Celine Barreteau

Jean-Claude Crivello

Jean-Marc Joubert

Eric Alleno

Optimization of criteria for an efficient screening of new thermoelectric compounds: the TiNiSi structure-type as a case-study

Keywords: Thermoelectricity, DFT calculations, Semiconductors, Intermetallic

High-throughput calculations can be applied to a large number of compounds, in order to discover new useful materials. In the present work, ternary intermetallic compounds are investigated, to find new potentially interesting materials for thermoelectric applications. The screening of stable non-metallic compounds required for such applications is performed by calculating their electronic structure, using DFT methods. In the first section, the study of the density of states at the Fermi level, of pure elements, binary and ternary compounds, leads to empirically chose the selection criterion to distinguish metals from non-metals. In the second section, the TiNiSi structure-type is used as a case-study application, through the investigation of 570 possible compositions. The screening leads to the selection of 12 possible semiconductors. The Seebeck coefficient and the lattice thermal conductivity of the selected compounds are calculated in order to identify the most promising ones. Among them, TiNiSi, TaNiP or HfCoP are shown to be worth a detailed experimental investigation.

INTRODUCTION

Thermoelectricity may be used as an interesting technology able to reduce the energy consumption through the recycling of waste heat. The conversion efficiency of thermoelectric devices is related to the dimensionless figure of merit ZT of their constituting materials, which is defined as

𝑍𝑇 = 𝑆 2 𝜎𝑇 𝜅 ⁄ (1)
where: S represents the Seebeck coefficient; σ the electrical conductivity; κ the thermal conductivity and T the absolute temperature. Even though efficient thermoelectric devices exist, they are limited to niche markets because of their quite high cost compared to their conversion yield. One way to improve this technology is to find new materials based on more abundant, less toxic and less expensive elements. These limitations led us to define a set of elements, we could investigate to search for new interesting compounds [START_REF] Barreteau | Looking for new thermoelectric materials among TMX intermetallics using high-throughput calculations[END_REF] . For this set of elements, hundreds of compounds have already been reported in crystallographic databases, and even if some of them are known to be good thermoelectric compounds, there is still a very large number of possibly unknown compounds. In order to screen more efficiently the possible combinations, we use high-throughput calculations. We will present an efficient method to screen ternary intermetallic compounds in order to find new promising thermoelectric compounds. In this article, we will focus on the equimolar compounds and use the orthorhombic TiNiSi structure type as a case-study. We choose the TiNiSi structure-type because even though it is a large family with more than 70 compounds already reported in the crystallographic databases [START_REF] Villars | Pearson's Crystal Data: Crystal structure database for inorganic compounds[END_REF][START_REF] Belsky | New developments in the Inorganic Crystal Structure Database (ICSD): accessibility in support materials research and design[END_REF] , it has been less subject of investigations compared with other equimolar structure-types. In contrast, numerous high-throughput calculations and screenings have already successfully shown very promising results for the equimolar half-Heusler family [START_REF] Feng | Electronic structure as a guide in screening for potential thermoelectrics: demonstration for half-Heusler compounds[END_REF][START_REF] Guo | Electrical properties dominated promising half-Heusler thermoelectrics through high-throughput materials computations[END_REF][START_REF] Legrain | Materials screening for the discovery of new half-Heuslers: Machine learning versus ab initio methods[END_REF][START_REF] Ma | Computational investigation of half-Heusler compounds for spintronics applications[END_REF][START_REF] Zeeshan | Ab initio design of new cobalt-based half-Heusler materials for thermoelectric applications[END_REF] .

Before estimating the thermoelectric properties of the compounds pertaining to the TiNiSi family, the main goal of the present work is to define an efficient method to discriminate among the investigated compounds those, which might be both stable and semiconducting. In our previous investigation [START_REF] Barreteau | Looking for new thermoelectric materials among TMX intermetallics using high-throughput calculations[END_REF] , we compared the enthalpy of formation of the same combination of elements in four structure types to determine the most stable one. In the present work, for the compounds stable in the TiNiSi structure-type, investigation of the electronic band structure will be carried out, to exclude the metallic compounds and keep only the possible semiconductors and semi-metals. A definition and the application of our exclusion criterion will be presented in this article. For the most remarkable compounds, the Seebeck coefficient and the thermal conductivity will be estimated using calculations based on the BoltzTrap package [START_REF] Madsen | A code for calculating band-structure dependent quantities[END_REF][START_REF] Madsen | BoltzTrap2, a program for interpolating ban structures and calculating semi-classical transport coefficients[END_REF] and the Slack equation [START_REF] Slack | Nonmetallic crystals with high thermal conductivity[END_REF][START_REF] Slack | The thermal conductivity of nonmetallic crystals[END_REF] respectively. All the results will be compared with the literature and the selected compounds will be highlighted.

COMPUTATIONAL DETAILS AND METHODOLOGY

In our previous work 1 , a combinatorial approach has been implemented to screen all possible T-M-X combinations with the 1:1:1 stoichiometry generated within a restrained set of chemical elements. In this set of elements, T is a transition metal from the Ti, V, and Cr columns or Sr, Ba, Y and La. M is an element from the first line of transition metals and X a sp element (Al, P, Si, Sn, and Sb). All T-M-X combinations have been calculated in four different structure-types:

the orthorhombic TiNiSi (Pnma), the cubic MgAgAs (F-43m, half-Heusler), the two hexagonal ZrNiAl (P-62m) and BeZrSi (P63/mmc). The 570 possible configurations have been calculated in these 4 structure-types by systematically assigning an unique element from our T-M-X nomenclature to each crystallographic site, discarding "inverse" configurations or site mixing, yielding 2280 different compounds.

The calculations of the stability and ground-state properties are based on the density functional theory (DFT). This technique allows us to obtain the enthalpy of formation at 0 K for each compound (by total energy difference with the elemental reference state) as well as its electronic structure. The calculations were conducted using the projector augmented wave (PAW) method, implemented in the Vienna Ab initio Simulation Package (VASP) [START_REF] Kresse | Ab-initio molecular dynamics for liquid metals[END_REF][START_REF] Kresse | From ultrasoft pseudopotentials to the projector augmented-wave method[END_REF][START_REF] Kresse | Efficient iterative schemes for an ab-initio total-energy calculations using plane-wave basis set[END_REF][START_REF] Hafner | The Vienne Ab-Initio Simulation Package[END_REF] ). These parameters ensured robust results in the total energy within convergence tolerance for the calculations set within 10 -6 eV of difference. We performed volume and ionic relaxation steps for each structure-type, and we considered spin polarization for all the configurations. Blöchl correction was applied in a final step calculation [START_REF] Blöchl | Improved tetrahedron method for Brillouin-zone integrations[END_REF] . Phonon dispersion bands were obtained by computing the atomic forces for different, finite atomic displacements [START_REF] Parlinski | First-principles determination of the soft mode in cubic ZrO2[END_REF] and a subsequent calculation and integration over the corresponding phonon frequencies, within the harmonic approximation. i.e. without any volume dependence, using the Phonopy code [START_REF] Togo | First-principles phonon calculations in materials science[END_REF][START_REF] Togo | First-principles calculations of the ferroelastic transition between rutile-type and CaCl2-type SiO2 at high pressures[END_REF] . The phonon calculations were carried out within a 1 × 2 × 1 supercell (24 atoms) with 12 atomic displacements for the TiNiSi structure-type.

The Seebeck coefficient was calculated from the Boltzmann transport theory [START_REF] Ashcroft | Solid State Physics[END_REF] as implemented in the BoltzTrap simulation package (v2) [START_REF] Madsen | A code for calculating band-structure dependent quantities[END_REF][START_REF] Madsen | BoltzTrap2, a program for interpolating ban structures and calculating semi-classical transport coefficients[END_REF] . The BoltzTrap code uses the electronic structure calculated with a dense k-mesh (12 x 20 x 10) from VASP as an input data and assumes a constant relaxation time for the electrons.

The first step of our selection concerns the stability of the calculated compounds. A positive enthalpy of formation leads to the exclusion of the combination, as it is obvious that the ternary compound is unstable toward the decomposition into its pure elements. In the second step, the most negative enthalpy of formation of the four-competing structure-types at the same composition, allows the selection of the most probable stable structure. Comparison with the literature and the crystallographic database has confirmed the robustness of our model 1 .

Definition of the non-metallic criterion

In addition to the stability criterion, it is important to find a simple and robust way to exclude metallic compounds for thermoelectric application. For this purpose, we decided to focus our screening method on the density of states (DOS) at the Fermi level (EF). By definition, when the DOS is zero at the EF, a semiconducting or insulating ground-state is expected. However, DFT-PBE is known [START_REF] Tran | Importance of the kinetic energy density for band gap calculations in solids with density functional theory[END_REF][START_REF] Yakovkin | The problem of the band gap in LDA calculations[END_REF] to not only systematically underestimate the gap magnitude but also to provide a non-zero DOS for known semiconductors. To compensate for this spurious effect, a non-zero discrimination criterion to select semiconductors is thus required to assess an unknown ground state (not too strict so as to miss interesting compounds and not too high for being efficiently selective). To estimate this value, preliminary calculations have been performed on the chemical elements in their stable structure (table S1) at 0 K and on known semiconducting binary compounds. Most of the elements of the periodic table have been investigated, except the 4f elements and the noble gas, whereas the investigated binary compounds have been chosen among known semiconductors or semi-metals.

The calculated DOS are presented in Fig. 1. As expected, all the metallic elements (s, d blocks)

have a non-zero DOS at EF, consistent with their metallic ground-state. In group III, (B, Al, Ga, In), the last three elements present a non-zero DOS at EF in agreement with their metallic ground-state whereas B has a zero DOS consistent with its semiconducting ground-state [START_REF] Greiner | Electrical Resistivity of Boron[END_REF] .

When looking at the first elements from group IV (C, Si, Ge), it is obvious that all the DOS at EF exhibit a zero value. The last two elements of this group, Sn and Pb, have a non-zero DOS at EF, 0.05 state eV -1 by atom and 0.5 state eV -1 by atom, respectively. This very low DOS at EF for Sn agrees with previous results [START_REF] Yu | Ab initio calculation of the properties and pressure induced transition of Sn[END_REF] which have shown that "gray tin", the α form of Sn (A4), which is the stable allotrope at low temperature [START_REF] Busch | Seminducting properties of gray tin[END_REF] , is semiconducting. Elemental Pb displays a DOS at EF similar to the metals of the block d. In group V (P, As, Sb, Bi), P has a zero DOS at EF whereas As, Sb, and Bi have a DOS smaller than 0.1 state eV -1 by atom. This confirms the semi-metallic ground-state of these last three elements, while elemental P is a semiconductor as reported in reference [START_REF] Goodman | Valence-band structures of phosphorus allotropes[END_REF] . In group VI with Se and Te, the DOS at EF is zero as these elements are reported as semiconductors [START_REF] Reitz | Electronic band strucuture of Selenium and Tellurium[END_REF][START_REF] Gaspar | Electronic structure of semi-conducting selenium and tellurium[END_REF] . In addition, all the binary compounds presented in Figure 1 exhibit a zero DOS or a very low DOS at EF, smaller than 0.005 state eV - reported semiconductors and semi-metals are in yellow and metals are in blue.

Reported results in Figure 1 allow to show an interval between the DOS at EF for known semiconductor or semi-metallic elements such as Sn or As and these of Ga and Zn, two metallic elements. Thus, with a selection limit set at 0.1 state eV -1 by atom for the DOS at EF, semiconductors, as well as elements with a semi-metallic ground state such as As or Sb, are kept whereas metals are excluded. Only one exception to this criterion is found: with this limit, Be, which is known as a metal is also selected. However, Be is known to exhibit a very low DOS at EF, despite its metallic behaviour [START_REF] Loucks | Band structure and Fermi surface of Beryllium[END_REF][START_REF] Mitchell | Electrical resistivity of beryllium[END_REF][START_REF] Chou | Ab initio study of structural and electronic properties of beryllium[END_REF][START_REF] Papaconstantopoulo | Handbook of the band structure of elemental solids[END_REF] . Overall, this limit ensures a judicious selection with a very small number of ambiguous results.

These calculations on elemental and known binary semiconductors allowed to establish an empirical discrimination criterion. However, before using it to screen unknown compounds, an additional test has been made on reported compounds from the Heusler family (full and half-Heusler). This family has been chosen because it has already been the object of several studies 4-5, 35 , Heusler compounds are very similar chemically to our screened compounds and metallic, semi-metallic and semiconductors compounds are reported in this family. The results are presented in the figure 2 with the DOS at EF for the considered ternary intermetallic compounds. As presented in Figure 2, our discrimination criterion is efficient to sort metals for semiconductors. Most of the reported semiconductors in the Heusler family exhibit a zero DOS at EF. However, some reported semiconductors or semi-metals do not and would be excluded without our criterion, even though, according to the literature, they are promising thermoelectric compounds. This is the case for ZrNiSn, TaCoSn or Fe2VAl [START_REF] Li | n-type TaCoSn-based half-Heusler as promising thermoelectric materials[END_REF][START_REF] Maier | Order-disorder transitions in the Fe2VAl Hesuler alloy[END_REF][START_REF] Kuentzler | Gap at the Fermi level and magnetism in RMSn ternary compounds (R= Ti, Zr, Hf and M= Fe, Co, Ni)[END_REF][START_REF] Hu | Interpreting the combustion process for high-performance ZrNiSn thermoelectric materials[END_REF] . Indeed, ZrNiSn is a very good example of both the necessity and the relevance of our discrimination criterion. As presented Twelve selected compounds in the TiNiSi structure type Among the 570 unique configurations calculated in the 4 structures type, 317 compounds are more stable in the TiNiSi structure-type than in any other structure type. But only 12 compounds are found to be non-metallic, according to the electronic criterion set at 0.1 state eV -1 by atom for the DOS at EF. They are reported in Table 1. All of these compounds are mechanically stable (no imaginary frequency is observed in their phonon bands, figure S2).

Half-metals [START_REF] De Groot | New class of materials: halfmetallic ferromagnets[END_REF] , which can be interesting for spintronic applications have also been searched. For these compounds, the majority-spin electrons are metallic whereas the minority-spin electrons are semiconducting. The criterion was applied to each spin directions. Eleven compounds complied with this criterion but only four are stable in the TiNiSi structure-type: CrNiSi, SrTiSn, BaVSb and BaCrSn (the DOS of these compounds are presented in figure S3 and phonon bands in figure S4). Nonetheless, these compounds are not interesting for thermoelectricity since they display a metallic total DOS and will not be investigated in detail in the following. We will first review the existence and structure of each selected compound. Then, for each candidate, a more specific study of the electronic structure will be presented. Finally, to assess their thermoelectric properties, the calculated values of the Seebeck coefficient and the thermal conductivity, will be examined.

Seven compounds among the 12 have not been reported in crystal structure databases. However, a careful study of phase diagrams may provide additional information. For the Ti-Cu-Sb, Zr-Zn-Si, Hf-Cu-Sb, and Ta-Cu-Si systems, no equimolar compound has been reported in the ternary phase diagram [START_REF] Dolotko | Interaction of Components in the Zr-Zn-Si System at 600 °C[END_REF][START_REF] Kobluyk | Interaction between the components in the Ti-Cu-Sb system at 770 K[END_REF][START_REF] Melnychenko-Koblyuk | Interaction between components in Hf-Cu-Sb ternary system at 770K[END_REF][START_REF] Savitskii | Copper effect on structure and superconductive properties of transition metal silicides [V(Nb, Ta, Cr, Mo, W)-Si-Cu alloy systems[END_REF][START_REF] Reid | Thermodynamics of (Cr, Mo, Nb, Ta, V or W)-Si-Cu ternary systems[END_REF] . It is therefore unlikely that an equimolar compound is stable in the ternary orthorhombic structure, but it cannot be formally discarded. For the 3 other unreported compounds, HfMnSb, MoZnSi, and HfZnSi, no phase diagram is available for the corresponding systems, it is therefore difficult to conclude on their existence.

Among the 5 reported compounds, TiMnSb has been reported in the literature to crystallize in the BeZrSi structure type [START_REF] Noda | High-pressure synthesis and magnetic properties of MMnSb (M= Ti, V, Cr) with Ni2In type structures[END_REF] in contrast with our calculations of a stable TiNiSi prototype.

However, as discussed in our previous work [START_REF] Barreteau | Looking for new thermoelectric materials among TMX intermetallics using high-throughput calculations[END_REF] , this discrepancy can be explained by the experimental conditions for the synthesis. Indeed, Noda et al. [START_REF] Noda | High-pressure synthesis and magnetic properties of MMnSb (M= Ti, V, Cr) with Ni2In type structures[END_REF] have obtained TiMnSb under extreme conditions (high pressure and high temperature) possibly giving rise to a nonequilibrium phase.

Finally, four compounds are reported in the crystallographic databases in the orthorhombic TiNiSi structure. In the cases of TiCoP and TaNiP, the literature only deals with their crystal structure and there is no investigation of their electronic properties [START_REF] Rundqvist | The Crystal Structure of ZrFeP and Related Compounds[END_REF][START_REF] Lomnytska | Component interaction in the system titanium-cobalt-phosphorus[END_REF] . More information is provided for the other two compounds: TiNiSi and HfCoP. Recently, Huang et al. [START_REF] Huang | Preparation and thermoelectric properties of pseudogap intermetallic (Ti1-xVx)NiSi solid solutions[END_REF] investigated the thermoelectric properties of TiNiSi-based solid solutions. Theoretical results, based on DFT calculations with PBE-GGA, confirmed the presence of a pseudo-gap in the electronic band structure of TiNiSi. Measurements of the electric conductivity versus the temperature confirmed that TiNiSi is a semi-metal while thermoelectric properties were improved by the substitution of Ti by V [START_REF] Huang | Preparation and thermoelectric properties of pseudogap intermetallic (Ti1-xVx)NiSi solid solutions[END_REF] . These results confirm a previous report by Kong et al. [START_REF] Kong | Structural and elastic properties of typical equiatomic ternary intermetallic compound TiNiSi under pressure[END_REF] who have also highlighted the semi-metallic ground state of TiNiSi, thanks to their calculations. HfCoP, based on an extended Hückel calculation by Kleinke et al. [START_REF] Kleinke | Synthesis Crystal Structure and Properties of HfM'P (M'= Fe, Co, Ni) in Comparison to ZrNiP[END_REF] , was considered as a metal. However, neither the electronic density of states nor the electrical resistivity was reported in this reference. The ground state of HfCoP is hence currently unknown.

Seebeck coefficient and lattice thermal conductivity at finite temperature

The Seebeck coefficient was calculated at finite temperatures within the semi-classical Boltzmann transport theory with the BoltzTrap package, assuming a constant relaxation time.

The possibility of calculating the electronic thermal conductivity and electrical conductivity offered by the package was not implemented because these quantities depend on the electronic relaxation time, which then plays the role of an adjustable parameter. This is not the case for the Seebeck coefficient, which is independent on the relaxation time when it is assumed to be constant. Only, TiNiSi and TiCuSb exhibit a different behaviour with Seebeck coefficient at 300 K higher than at 800 K. Even though this trend is maintained for the negative Seebeck coefficient, with the absolute value at 800 K higher than at 300 K, four compounds exhibit an opposite behaviour, TiMnSb, TiCuSb, MoZnSi and HfMnSb.

Several compounds like ZrZnSi or HfZnSi do not show values of the Seebeck coefficient larger than ± 50 μV K -1 and therefore do not seem to be very promising in contrast with TiMnSb or MoZnSi for example. First extremum of the calculated Seebeck coefficient close to the Fermi level are listed in Table 2. In this table, values at 700 K are reported, to compare the compounds in the same way as it will be done for the lattice thermal conductivity (values at 300 K are reported in table S2). We choose to focus on the maximal values of the Seebeck coefficient when μ is near the Fermi level for an easy comparison. As mentioned before, ZrNiSi, HfZnSi, HfCuSb, and TaCuSi present quite low Seebeck coefficient near the Fermi energy, regardless their conduction type. Considering both positive and negative values, TiMnSb and HfMnSb can be considered as the best compounds when other compounds could be promising at most in one of the two conduction types. For example, TiNiSi has a negative Seebeck coefficient twice as large as its positive counterpart whereas it is the opposite for MoZnSi. These results are comparable at other temperatures where similar trend can be noticed (table S2). No clear tendencies can be noticed in the evolution of the Seebeck coefficient regarding the nature of each T, M or X element in these compounds.

Lattice thermal conductivity, which is essential for the determination of ZT, can be estimated by the Slack's equation [START_REF] Slack | Nonmetallic crystals with high thermal conductivity[END_REF][START_REF] Slack | The thermal conductivity of nonmetallic crystals[END_REF] by assuming that heat is conducted only by acoustic phonons [START_REF] Hou | Machine-Learning-assisted development and theoritical consideration for the Al2Fe3Si3 thermoelectric material[END_REF][START_REF] Shiota | A new semiconductor Al2Fe3Si3 with complex crystal structure[END_REF] , through the equation 2:

𝐾 𝐿 = 3.1 × 10 -6 𝑀 𝜃 𝐷 3 𝛿 𝛾 2 𝑛 2 3𝑇
(2

)
where θD is the Debye temperature in K, δ 3 is the volume per atom in Å [START_REF] Belsky | New developments in the Inorganic Crystal Structure Database (ICSD): accessibility in support materials research and design[END_REF] , n is the number of atoms in the primitive unit cell, Mav is the average atomic mass in the atomic unit and γ is the average Grüneisen parameter. The Debye temperature and the Grüneisen parameter are evaluated from the transverse and longitudinal sound velocity, which are obtained from the calculated elastic modulus. The calculated elastic properties (table S3) have been obtained from the phonon calculations for each compound (see figure S2). Slack's model neglects the transport of heat by optical phonons and this could lead to underestimated values of lattice thermal conductivity. However, it also assumes that the acoustic phonons are scattered only by acoustic phonons, neglecting all the other scattering mechanisms (isotopes, impurities, vacancies, charge carriers, dislocations, grain boundaries, crystal boundaries, etc.) This last approximation has the strongest impact on Slack's values [START_REF] Slack | The thermal conductivity of nonmetallic crystals[END_REF] , which are systematically overestimated (up to an order of magnitude) as soon as the relative atomic mass difference of the constituting elements is larger than roughly ~ 100% 12 : this condition is indeed fulfilled for the 12 screened compounds.

Despite these drawbacks, Slack's model provides reliable trends at temperatures larger than the Debye temperature, when phonon -phonon scattering dominates over the other scattering processes. To fulfil this condition for all the compounds, the values presented in Table 2 were calculated at 700 K.

The range of lattice thermal conductivity values is quite broad with estimated values between 3 and 15 W m -1 K -1 which cannot be compared to experimental values due to lack of data. These values are quite large compared to the best thermoelectric materials such as Bi2Te3 or Sb2Te3 which exhibit low thermal lattice thermal conductivity (~ 2 W m -1 K -1 ) [START_REF] Elsheikh | A review on thermoelectric renewable energy: Principle parameters that affect their performance[END_REF] . They are, nonetheless, still promising because they were calculated for a perfect crystal and thus as above discussed, overestimated. Moreover, it is known that defects or impurities can strongly decrease the lattice thermal conductivity. For instance, the half-Heusler HfNiSn displays a thermal conductivity of 10 W m -1 K -1 at 300 K [START_REF] Uher | Transport properties of pure and doped MNiSn (M=Zr, Hf)[END_REF] while the solid solution Hf0.75Zr0.25NiSn0.99Sb0.02displays a thermal conductivity of 5 W m -1 K -1 at 300 K and ZT = 0.9 at 800 K [START_REF] Yan | Enhanced Thermoelectric Figure of Merit of p-Type Half-Heuslers[END_REF] .

Keeping the X-element constant and considering the T element, as expected, the lattice thermal conductivity is, in general, lower for "heavier" T in the same columns of the periodic table as for TCuSb or TCoP. In the case of TMnSb, the trend is reversed as κL(TiMnSb) is lower than κL(HfMnSb). Keeping the T element constant and considering the X element, we noticed that the antimonides present the lowest thermal conductivity compared to the silicides or the phosphides, in accordance with the larger atomic mass of antimony.

DISCUSSION

Keeping in mind that the values presented here are from calculations based on DFT-PBE and constant relaxation time approximation implemented in a perfect crystal, our study offers the possibility to identify the most promising compounds among the 12 selected TiNiSi-type compounds from our screening of 570 possible compositions. For example according to Table 2, even though, lattice thermal conductivity is quite low for the two TCuSb compounds, their Seebeck coefficients are also systematically low. ZrNiSi and TaCuSi exhibit also low Seebeck coefficient, which are for these compounds accompanied by high lattice thermal conductivity values. Thus, if one adds to these low values the fact that, according to the reported ternary phase diagram, these four compounds are not stable, it is probably not worth trying their synthesis.

On the other hand, HfMnSb and MoZnSi combine large Seebeck coefficients and a small lattice thermal conductivity that could lead to interesting properties. Nonetheless, their existence is not established.

Regarding the four compounds which have already been reported in the literature, a ranking is also possible. The highest value of the Seebeck coefficient is reached in TiNiSi, which has already been reported as interesting for its thermoelectric properties in the solid solution Ti1-xVxNiSi [START_REF] Huang | Preparation and thermoelectric properties of pseudogap intermetallic (Ti1-xVx)NiSi solid solutions[END_REF] . The presently calculated values of the Seebeck coefficient seem to confirm the calculated values calculated by Huang et al. [START_REF] Huang | Preparation and thermoelectric properties of pseudogap intermetallic (Ti1-xVx)NiSi solid solutions[END_REF] , with S ≈ ± 65 µV. K -1 at 300 K for both type of charge carrier whereas our calculation gives S = 64 µV. K -1 and S = -69 µV. K -1 at 300 K.

However, the measured values (-16 µV. K -1 at 300K) in Ref. [START_REF] Huang | Preparation and thermoelectric properties of pseudogap intermetallic (Ti1-xVx)NiSi solid solutions[END_REF] differ from the calculated one, but this difference can be explained by the charge carrier concentration of the pristine TiNiSi, which is around 2 × 10 21 cm -3 . This relatively large value indicates the presence of intrinsic defects in the sample. Finally, HfCoP seems to be more interesting than TiCoP as its Seebeck coefficient is larger and associated with a lower lattice thermal conductivity. TaNiP exhibits the lowest lattice thermal conductivity of the four reported compounds and negative Seebeck coefficient of the same order of magnitude as in HfCoP.

CONCLUSION

After determining the stability of thousands of TMX compositions in 4 structure types, a more focused investigation has been undertaken on the compounds which crystallize in the TiNiSi structure type. In a second stage, the electronic band structure has been examined and an effective criterion of the DOS at the Fermi level has been established by comparing calculated DOS of known semiconducting elements and compounds. Then, this empirical criterion of 0.1 state eV -1 by atom has been applied to the TMX compounds, stable in the orthorhombic structure. It has led to the selection of 12 non-metallic compounds among the 317 possible compounds while it was commonly admitted that all the TiNiSi are metallic [START_REF] Landrum | The TiNiSi Family of Compounds: Structure and Bonding[END_REF] . Additional phonon and BoltzTrap calculations allowed to rank the 12 selected compounds to identify those which are the most interesting to be synthesized from the point of view of their anticipated thermoelectric properties. Although MoZnSi could be an interesting compound as it presents the highest Seebeck coefficient associated with a small lattice thermal conductivity, there is no information about its possible stability. So, it seems more relevant to focus on HfCoP and TaNiP as starting points. Both compounds are already reported in the literature as crystallizing in the TiNiSi structure type [START_REF] Lomnytska | Component interaction in the system titanium-cobalt-phosphorus[END_REF][START_REF] Kleinke | Synthesis Crystal Structure and Properties of HfM'P (M'= Fe, Co, Ni) in Comparison to ZrNiP[END_REF] and present the best compromise between high Seebeck coefficient and a low lattice thermal conductivity.

Figure S3

: Electronic density of states for the four (half-metals) TMX compounds calculated by DFT using GGA-PBE approximation. The Fermi level is chosen as the origin of energies.

Both spin directions are presented. Table S2: Values of the first extremum of the calculated Seebeck coefficient close to the Fermi level at 300 K for the TiNiSi compounds.

Table S3: The bulk modulus (B), shear modulus (G), Young's modulus (E) and Poisson ratio (ν) of the screened TMX compounds.

TABLE OF CONTENTS

Figure 1 .

 1 Figure 1. Density of states (DOS) at EF for a selection of pure elements and binary compounds;

Figure 2 :

 2 Figure 2 : Density of states (DOS) at EF for a selection of half-Heusler and Heusler compounds; reported semiconductors and semi-metals are in yellow and metals are in blue.

Figure 2 ,

 2 Figure2, the DOS(EF) is equal to 0.07 states/eV by atom, which is not strictly zero, however

Figure 3 .

 3 Figure 3. Calculated Seebeck coefficient as a function of the chemical potential for each compound. Each curve corresponds to a temperature varying by steps of 100K in the range 300 K (yellow) to 800 K (black). The Seebeck coefficient has been plotted as a function of the chemical potential (μ) for each

Figure S4 :

 S4 Figure S4: Phonon curves dispersion of the four (half-metals) TMX compounds calculated

Figure 1 .

 1 Figure 1. Density of states (DOS) at EF for a selection of pure elements and binary compounds;reported semiconductors and semi-metals are in yellow and metals are in blue.

Figure 2 :

 2 Figure 2 : Density of states (DOS) at EF for a selection of half-Heusler and Heusler compounds; reported semiconductors and semi-metals are in yellow and metals are in blue.

Figure 3 .

 3 Figure 3. Calculated Seebeck coefficient as a function of the chemical potential for each compound. Each curve corresponds to a temperature varying by steps of 100K in the range 300 K (yellow) to 800 K (black).

  Figure 1

Table 1

 1 Non-metallic compounds in the TiNiSi structure type.

	TMX	Crystal	Enthalpy of	DOS at the Fermi
	Compounds	databases	formation	level (state. eV -1 by
			(kJ mol -1 )	atom)
	TiCuSb	Unreported	-20.90	0.092

Table 2 .

 2 Values of the first extremum of the calculated Seebeck coefficient close to the Fermi level and Slack's lattice thermal conductivity at 700 K for the TiNiSi compounds.

	TMX compounds	Smax (μV K -1 ) > 0	Smax (μV K -1 ) < 0	Slack's lattice Thermal
		(700K)	(700K)	conductivity at 700K (W m -1
				K -1 )
	TiMnSb	85	-78	5.4
	TiCoP	52	-55	15.7
	TiNiSi	54	-87	13.7
	TiCuSb	52	-28	3.9
	ZrZnSi	30	-36	9.9
	HfMnSb	70	-98	6.1
	HfCoP	71	-72	11.6
	HfCuSb	44	-47	3.7
	HfZnSi	33	-41	7.4
	TaNiP	41	-76	7.6
	TaCuSi	58	-54	9.7
	MoZnSi	95	-48	6.5

Table 1

 1 Non-metallic compounds in the TiNiSi structure type.

Table 2 .

 2 Values of the first extremum of the calculated Seebeck coefficient close to the Fermi level and Slack's lattice thermal conductivity at 700 K for the TiNiSi compounds.

Corresponding author: barreteau@icmpe.cnrs.fr

by atom, which is consistent with their semiconducting ground states.

ACKNOWLEDGEMENTS

DFT calculations were performed using HPC resources from GENCI-CINES (Grant 2019-A0060906175) and financial support from the Agence Nationale de la Recherche is acknowledged (ANR-18-CE05-0010-01).

SUPPORTING INFORMATION AVAILABLE

Table S1: Elemental, binaries and ternaries compounds information's


J Gavinet

J-M Ourcival

J Gauzere

L García De Jalón

J-M Limousin

Jean-Marc Limousin
email: jean-marc.limousin@cefe.cnrs.fr

Drought mitigation by thinning: benefits from the stem to the stand along 15 years of experimental rainfall exclusion in a holm oak coppice

Keywords: Quercus ilex, drought resistance, tree mortality, stool density, climate change, adaptive management, stem growth, resprouting, ecosystem manipulation

    Drought mitigation by thinning: Benefits from the stem to the stand along 15 years of experimental rainfall exclusion in a holm oak coppice

Introduction

Forest ecosystems are currently facing fast changes, such as unprecedented climatic conditions and landuse changes, that already affect tree growth, mortality and reproduction [START_REF] Allen | A global overview of drought and heat-induced tree mortality reveals emerging climate change risks for forests[END_REF][START_REF] Doblas-Miranda | A review of the combination among global change factors in forests, shrublands and pastures of the Mediterranean Region: Beyond drought effects[END_REF][START_REF] Linares | Interacting effects of changes in climate and forest cover on mortality and growth of the southernmost European fir forests[END_REF][START_REF] Mcdowell | Darcy's law predicts widespread forest mortality under climate warming[END_REF]. Such changes are expected to further accelerate in the future and strongly challenge current forest management practices [START_REF] Keenan | Climate change impacts and adaptation in forest management: a review[END_REF][START_REF] Millar | Climate Change and Forests of the Future: Managing in the Face of Uncertainty[END_REF][START_REF] Stephens | Operational approaches to managing forests of the future in Mediterranean regions within a context of changing climates[END_REF]. Decreasing stand density by thinning is traditionally used to improve residual tree growth and dimensions, but it has gained a renewed interest as an adaptive management practice aiming at improving forest resistance or resilience to drought [START_REF] Bottero | Density-dependent vulnerability of forest ecosystems to drought[END_REF][START_REF] Navarro-Cerrillo | Is thinning an alternative when trees could die in response to drought? The case of planted Pinus nigra and P. Sylvestris stands in southern Spain[END_REF][START_REF] Vilà-Cabrera | Forest management for adaptation to climate change in the Mediterranean basin: A synthesis of evidence[END_REF]. By reducing stand leaf area, thinning is expected to reduce rainfall interception and stand transpiration, thereby improving water availability for the remaining trees [START_REF] Bréda | Effects of thinning on soil and tree water relations, transpiration and growth in an oak forest (Quercus petraea (Matt.) Liebl.)[END_REF][START_REF] Del Campo | Hydrology-oriented (adaptive) silviculture in a semiarid pine plantation: How much can be modified the water cycle through forest management?[END_REF][START_REF] Giuggiola | Improvement of water and light availability after thinning at a xeric site: which matters more? A dual isotope approach[END_REF]. Reduced competition for water after thinning has been shown to increase mean tree growth and to reduce growth vulnerability to extreme drought events [START_REF] Aldea | Thinning enhances the species-specific radial increment response to drought in Mediterranean pine-oak stands[END_REF][START_REF] Bottero | Density-dependent vulnerability of forest ecosystems to drought[END_REF], although with a high variability depending on species, climate and thinning intensity (Sohn et al., 2016b). How thinning regimes, i.e. their frequency, intensity and the characteristics of the trees harvested affect the forest resistance to drought is still poorly known. Thinning effects on tree growth and mortality strongly depend on individual characteristics such as tree size and local competitive environment [START_REF] Bose | Tree-level growth and survival following commercial thinning of four major softwood species in North America[END_REF], which also affect the drought responses [START_REF] Bennett | Larger trees suffer most during drought in forests worldwide[END_REF][START_REF] Colangelo | Size Matters a Lot: Drought-Affected Italian Oaks Are Smaller and Show Lower Growth Prior to Tree Death[END_REF][START_REF] Trouvé | Growth partitioning in forest stands is affected by stand density and summer drought in sessile oak and Douglas-fir[END_REF]. Furthermore, positive thinning effects may vanish after a few years [START_REF] Bréda | Effects of thinning on soil and tree water relations, transpiration and growth in an oak forest (Quercus petraea (Matt.) Liebl.)[END_REF][START_REF] Cotillas | Growth response of mixed mediterranean oak coppices to rainfall reduction: Could selective thinning have any influence on it? For[END_REF][START_REF] Del Campo | The impact of adaptive forest management on water fluxes and growth dynamics in a water-limited lowbiomass oak coppice[END_REF][START_REF] Sánchez-Humanes | Increased drought reduces acorn production in Quercus ilex coppices: thinning mitigates this effect but only in the short term[END_REF], and even reverse in the long-term [START_REF] D'amato | Effects of thinning on drought vulnerability and climate response in north temperate forest ecosystems[END_REF] because stand leaf area can quickly recover to its pre-thinning values either by an increase in tree leaf area [START_REF] Bréda | Effects of thinning on soil and tree water relations, transpiration and growth in an oak forest (Quercus petraea (Matt.) Liebl.)[END_REF][START_REF] Giuggiola | Reduction of stand density increases drought resistance in xeric Scots pine forests[END_REF] or by the development of understory vegetation following increased light availability [START_REF] Gebhardt | The more, the better? Water relations of Norway spruce stands after progressive thinning[END_REF][START_REF] Simonin | The influence of thinning on components of stand water balance in a ponderosa pine forest stand during and after extreme drought[END_REF]. Thinning is also accompanied by an increase in light, temperature and evaporative demand in the canopy and understory [START_REF] Gavinet | Hardwood seedling establishment below Aleppo pine depends on thinning intensity in two Mediterranean sites[END_REF][START_REF] Lechuga | Managing drought-sensitive forests under global change. Low competition enhances long-term growth and water uptake in Abies pinsapo[END_REF], possibly resulting in increased water losses by evapotranspiration that can reduce or offset thinning benefits for water availability [START_REF] Del Campo | The impact of adaptive forest management on water fluxes and growth dynamics in a water-limited lowbiomass oak coppice[END_REF][START_REF] Gebhardt | The more, the better? Water relations of Norway spruce stands after progressive thinning[END_REF]. Moreover, when thinning reduces population sizes too strongly, genetic drift may lead to the loss of rare alleles and restrict local adaptation processes [START_REF] Lefèvre | Considering evolutionary processes in adaptive forestry[END_REF]. Thinning is thus generally considered as a short-term adaptation option [START_REF] Vilà-Cabrera | Forest management for adaptation to climate change in the Mediterranean basin: A synthesis of evidence[END_REF], and thinning regimes are to be determined depending on species, stand type and site conditions [START_REF] Ameztegui | Managing stand density to enhance the adaptability of Scots pine stands to climate change: A modelling approach[END_REF]Skov et al. 2009;Sohn et al. 2016b).

Whether thinning is a suitable strategy to improve tree vitality and resistance to stress in resprouting stands is not straightforward because an important resprouting dynamic can take place after thinning and quickly offset the benefits of competition reduction [START_REF] Cotillas | Growth response of mixed mediterranean oak coppices to rainfall reduction: Could selective thinning have any influence on it? For[END_REF][START_REF] Ducrey | Recherches et expérimentations sur la conduite sylvicole des peuplements de chêne vert[END_REF][START_REF] Ducrey | Effect of cleaning and thinning on height growth and girth increment in holm oak coppices (Quercus ilex L.)[END_REF]. The resprouting ability allows trees to persist under disturbances or drought [START_REF] Bond | Ecology of sprouting in woody plants: the persistence niche[END_REF][START_REF] Clarke | Resprouting as a key functional trait in woody plantschallenges to developing new organizing principles[END_REF][START_REF] Zeppel | Drought and resprouting plants[END_REF] and is widespread in broadleaved species such as Mediterranean oaks. These oak forests represent the main late-successional stands around the Mediterranean basin (Quézel & Médail 2003), and because of their strong resprouting ability, they have been generally managed as coppices for the provision of fuel wood, charcoal or tannins. With rural abandonment and the decreasing use of firewood and charcoal as energy sources, oak coppices have been progressively abandoned since the middle of the XXth century and form now ageing stands with frequent signs of dieback, lack of seed regeneration and low levels of productivity [START_REF] Camarero | Prior height, growth, and wood anatomy differently predispose to drought-induced dieback in two Mediterranean oak speciesk[END_REF][START_REF] Gentilesca | Drought-induced oak decline in the western Mediterranean region: an overview on current evidences, mechanisms and management options to improve forest resilience[END_REF]. Besides, in coppice stands, the long history of coppicing has led to complex forest structures comprising both unique genets and multi-stemmed clumps [START_REF] Valbuena-Carabaña | Coppice forests and genetic diversity: A case study in Quercus pyrenaica Willd. from Central Spain[END_REF], where individual stem growth and vigour are also influenced by stool characteristics [START_REF] Salomón | Clonal structure influences stem growth in Quercus pyrenaica Willd. coppices: Bigger is less vigorous[END_REF]. Whether management by thinning could improve productivity and drought resistance of aged oak coppices in the long-term remains to be determined, especially in a context of climate change pushing drought stress beyond its historical range. Thinning from below at a moderate intensity, i.e. by removing the smallest suppressed stems and by keeping at least one stem on multi-stemmed stools, may limit resprouting and maintain stand-level productivity (Cabon et al., 2018b;[START_REF] Ducrey | Recherches et expérimentations sur la conduite sylvicole des peuplements de chêne vert[END_REF] as well as a large population size of unique genets. This type of thinning from below can additionally decrease fire risk by limiting the vertical continuity of the vegetation cover and limit the impact on understory temperature and evaporative demand. Thinning from below has potentially two types of effects: the selected trees are the biggest, which are often more vigorous and more efficient in resource acquisition and transformation [START_REF] Binkley | Thinking about efficiency of resource use in forests[END_REF], while the removal of small trees can release the competition for belowground resources (Cabon et al., 2018b;[START_REF] Giuggiola | Competition for water in a xeric forest ecosystem -Effects of understory removal on soil micro-climate, growth and physiology of dominant Scots pine trees[END_REF].

However, the relative importance of size-selection and competition release, the influence of stool-level variables, and the impact of thinning on genetic diversity on the long-term remain to be described in a context of increasing drought.

In this study, we investigate the effects of thinning on holm oak growth, mortality, resprouting and stool diversity under ambient and increased drought severity, by using a 15-year experiment in a mature coppice combining rainfall exclusion (-27% of incoming precipitations) and thinning from below (-30% basal area). A previous study after 5 years of treatment showed that thinning improved tree growth and survival but also increased tree resprouting, suggesting that thinning effect could be quickly offset by resprouting [START_REF] Rodríguez-Calcerrada | Is selective thinning an adequate practice for adapting Quercus ilex coppices to climate change?[END_REF]. Ten years later, we aimed to i) analyse the long-term evolution of thinning effects and resprouting dynamics, ii) understand the individual characteristics driving tree demographic rates (growth, survival, resprouting) and their response to drought and iii) explore the stand-level evolution of basal area, stem density and stool densitywhich relates to the number of unique genotypesunder the different treatments.

Material & method

Experimental site

The experiment was conducted in southern France (35 km northwest of Montpellier), on a flat area in the Puéchabon State Forest (43°44'29'' N;3°35'46'' E,270 m a.s.l.). This forest has been coppiced for centuries with clear cuts approximately every 25 years, until the last clear cut in 1942. The evergreen Quercus ilex L. forms a dense canopy with a height of c. 5.5 m, a mean basal area of 30 m 2 .ha -1 and a density of c. 6000 stems. ha -1 , representing c. 4000 genetically different individuals (stools).ha -1 . The evergreen species Buxus sempervirens, Phyllirea latifolia, Pistacia terebinthus and Juniperus oxycedrus, compose a sparse understory layer with c. 20% cover and represent less than 3% of stand basal area. The bedrock is a hard Jurassic limestone and the soil is extremely rocky with c. 75% of stones and rocks in the top 0-50 cm and 90% below. The stone-free fraction of the soil within the 0-50 cm layer is a homogeneous silty clay loam (38.8% clay, 35.2% silt and 26% sand). The area has a Mediterranean-type climate with a mean annual temperature of 13.2°C (on-site meteorological station, 1984-2017), the coldest month being January (5.5°C) and the hottest month July (22.9°C). The mean annual precipitation is 910 mm with a range of 550 -1549 mm (1984-2017). Rainfall mainly occurs during autumn and winter with about 80% between September and April.

Experimental design: rainfall exclusion and thinning experiments

In March 2003, a factorial combination of throughfall exclusion and thinning treatments was set up on three 20 × 20 m replicated blocks located 200m apart one from the other. Each block comprises four treatments applied on a 10 × 10 m plot: control, throughfall exclusion (further "dry"), thinned, and thinned with throughfall exclusion (further "thinned dry"). For the throughfall exclusion treatment, half of the block was equipped with parallel 14m long and 0.19m wide PVC gutters hung below the canopy with a slope, between 1m and 2m in height, so as to cover 33% of the ground area under the tree canopy.

Taking into account rainfall interception and stemflow, the net input of precipitation was reduced by 27% compared with the control plots [START_REF] Limousin | Modelling rainfall interception in a mediterranean Quercus ilex ecosystem: Lesson from a throughfall exclusion experiment[END_REF]. On the other half of the blocks, identical gutters were installed upside down to homogenize albedo and understorey micro-climate without reducing precipitation inputs. Thinning from below was applied once, in winter 2003, on half of the plots by removing 30% of the initial plot basal area (27% of Q. ilex basal area, Table 1), thereby reducing the stem density by 49 % and stool density by 36% (Table 1). Dead, diseased and suppressed stems were felled and the understory cleared. One to four stems were cut on multi-stemmed stumps, in order to leave the larger stems. Thinning resulted in an increase of the mean Q. ilex stem diameter by 25% (Table 1). Measures of pre-dawn water potential on a subsample of trees during the first 5 years of the experiment confirmed that rainfall exclusion increased tree water stress while thinning decreased it [START_REF] Rodríguez-Calcerrada | Is selective thinning an adequate practice for adapting Quercus ilex coppices to climate change?[END_REF].

Demographic data

In each plot, stems with a diameter at breast height (DBH) larger than 2 cm were individually tagged, mapped and their DBH was measured before and after treatment application in March 2003.

Neighboring trees located within a 3m band around the blocks were also measured to assess the competitive environment of the trees inside the plots. Annual inventories of stem DBH were then conducted every winter from 2004 to 2019 with a diameter tape. Annual stem basal area increments were calculated for each year n from 2003 to 2018 as BAIn = 𝜋 4 (𝐷𝐵𝐻 𝑛+1 2 -𝐷𝐵𝐻 𝑛 2 ). Trees were recorded as dead when only brown leaves remained on the tree without crown regreening during the following years. Resprouts were inventoried twice, in 2008 and 2018. For each stump, the number of resprouts, the diameter of all resprouts and the height of the dominant resprout were measured. We considered the resprouts taller than 1.30m as established resprouts. None of them reached the threshold of 2cm DBH to be considered as recruits.

Stool identification and clonal structure

The stool clonal structure of the coppice was visually assessed during the first inventory by grouping as clones the stems that were visibly connected to the same stump and root system. This assessment was facilitated by the shallow and rocky soil of the site where the stumps are generally clearly visible above the ground. This visual census of genetically different individuals was then confirmed by genetic analyses performed in 2018 on a subsample of 166 individuals selected in the denser unthinned plots.

Leaves were collected on one stem of every identified stump and frozen on dry ice. The DNA was extracted using DNeasy 96 Plant kit (QIAGEN GmbH, Hilden, Germany) and genotyping was performed using 70 validated SNPs markers [START_REF] Bonal | 70 highly validated SNP markers for Quercus ilex[END_REF] developed from ddRAD-Seq data by the Genome Transcriptome Facility of Bordeaux (PGTB) following the methods described in [START_REF] García | A comparative analysis between SNPs and SSRs to investigate genetic variation in a juniper species (Juniperus phoenicea ssp. turbinata)[END_REF]. To be considered clones, individuals had to be compared at more than 50 loci with exact matching of their SNP markers. The results confirmed the relative reliability of the visual assessment as only 25 individuals among the 166 sampled were unidentified clones (15% of undetected clones). The mean distance between unidentified clones was 1.05 m and never exceeded 3 m (see a map of one block in Fig. S1). These results show that stools extent was spatially restricted and that it could reliably be described from visual inventories. We therefore decided to retain the visual stool identification as a proxy for the number of genetically different individuals within the coppice.

Competition indexes

We calculated the stem-level competition index as the sum of the neighbors' basal area in a 3m radius around each stem (BAn). BAn is a distance-independent, size-symmetric competition index. Size symmetric competition occurs when each tree exerts a competitive effect proportional to its size. This is typically assumed to be the case when belowground resources such as water are the main limitation.

We also computed other widely used competition indexes: two size and distance-dependent competition indexes, the Hegyi competition index [START_REF] Hegyi | A simulation model for managing jack-pine stands[END_REF] and the Weiner competition index [START_REF] Weiner | Neighbourhood Interference Amongst Pinus Rigida Individuals[END_REF] and the size asymmetric index Basal Area of Larger trees within a 3m radius (BAL). The models including these indexes were however always less performant than with BAn (lower or similar R 2 , higher AIC, with particularly no effect of BAL, data not shown). We thus retained BAn to characterize the local competition intensity as this index is simpler and not related to tree size, allowing us to separate mechanisms of competition and size-selection. We tested the effect of competition absolute values (BAn) and competition release, expressed as the percentage of BAn removed by thinning or neighbor's mortality.

Data analysis

Annual stem basal area increment (BAI) were averaged over the 15 years of the experiment for each tree. This mean growth variable was first analysed as a function of thinning, rainfall exclusion and their interaction using a linear model. Then, we constructed a stem growth model as a function of individual characteristics using stem size in interaction with indexes of competition intensity (BAn), competition release (percentage difference of BAn before and after thinning, see above) and rainfall exclusion as explanatory variables. Thinning was not included in these latter models because it was redundant with competition indices and strongly inflated models VIF [START_REF] Zuur | A protocol for data exploration to avoid common statistical problems[END_REF]. We selected the best model based on the Akaike criterion (AIC). When AIC difference between two models was lower than 2, we considered the models equivalent and selected the most parsimonious one. The variance explained by each of the predictors was assessed using the R package relaimpo [START_REF] Grömping | Relative Importance for Linear Regression in R: The Package relaimpo[END_REF]. Finally, the temporal evolution of treatment effects on tree growth was analysed by separating the dataset in three 5-years periods and analysing stem growth as a function of treatments, period and their interactions using a linear mixed model with stem identity as a random factor in order to account for repeated measures. In all models, stem growth was log-transformed to satisfy conditions of normality and homoscedasticity.

Holm oak stem and stool mortality was analysed similarly as a function of i) treatments and ii) individual characteristics. Mortality probability was modelled using a logistic regression with a logit link. Rainfall exclusion effect on the diameter of dead trees in unthinned plots was tested with a t-test.

The year 2006 had the driest spring season over the course of the experiment, with only 52mm of precipitations between April and June (average 2003 -2018 = 211mm). Spring water limitation has been shown to be the most important environmental factor for stem growth and stand productivity in our site [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF][START_REF] Lempereur | Growth duration is a better predictor of stem increment than carbon supply in a Mediterranean oak forest: implications for assessing forest productivity under climate change[END_REF]. We thus used data from the year 2006 to quantify the impact of an extreme drought event on individual stem growth. We calculated the drought resistance for each stem as the ratio of stem growth in 2006 to mean growth in the three preceding years [START_REF] Lloret | Components of tree resilience: effects of successive low-growth episodes in old ponderosa pine forests[END_REF]. Drought resistance was analyzed as described for stem BAI and mortality as a function of rainfall exclusion and thinning treatments, then as a function of individual characteristics. Drought resistance was log-transformed to obtain normal residuals and 13 outliers were discarded, corresponding to stems with either a negative BAI in 2006 (typically due to a measurement error) or a mean BAI during the preceding years almost null, 8 of which died some years later.

Resprouts number, maximum height and total basal area per stool were analysed as a function of rainfall exclusion, thinning treatments, inventory date (5 and 15 years after thinning) and their interaction, using stool identity as a random factor to account for repeated measures. A negative binomial distribution was used for the number of resprouts to account for overdispersion (function glmer.nb of the R package lme4), and linear mixed models were used for log-transformed resprout height and basal area. When interactions between factors were significant, we performed Tukey post-hoc tests. Then, we modelled resprout basal area (the most integrative index of resprout dynamics) as a function of stool size and the percentage of stool basal area that had been cut by thinning (stool competition release).

To examine treatment effects at the stand-level, we calculated for each year and plot the relative stand basal area, stand basal area increment, stem density and stool density as the ratio of these variables to their initial values (pre-thinning). We used a linear mixed model to test the evolution of stand characteristics as a function of years since treatment application, thinning, rainfall exclusion and their interactions with plot as a random factor to account for repeated measures.

Results

Stem growth and resistance to drought

Mean stem growth over the 15-year period was strongly affected by thinning (P < 0.001, Figure 1A, explained variance of 12.2%) but not by rainfall exclusion (P = 0.77). Stem growth was positively related to stem DBH in interaction with stem-level competition release, so that biggest stems responded more to competition release (Table 2; Figure 1B). Stem DBH explained 53% of stem growth variance, competition release an additional 7% and the interaction between these factors 0.5%. Absolute competition (BAn) had an additive negative effect and explained 1.9% of stem growth variance. (see Table 2).

Thinning effect on stem growth decreased over time (Thinning x Period interaction, P < 0.001; Figure 2) but remained significant (P < 0.001, post-hoc Tukey tests) during all the periods. Rainfall exclusion had no effect on stem growth in any of the periods (Rainfall exclusion main effect P = 0.43, interaction with Period P = 0.32). The extreme spring drought in 2006 induced a strong reduction of stem growth to about 50% of the average of the preceding years in the unthinned plots (Figure 3A). Thinning strongly increased stem resistance to drought (P < 0.001; Figure 3) as tree growth was maintained in thinned stands. In contrast, stem growth resistance to drought was not affected by rainfall exclusion (P = 0.77) nor by the interaction between rainfall exclusion and thinning (P = 0.25). Thinning explained 8.1% of the variance in stem growth resistance to drought. Stem resistance to drought was positively related to stem DBH (explained variance = 2.5%) and competition release (explained variance = 4.5%; Figure 3B; Table 2). 2).

Stem and stool mortality

The average stem mortality rate was 2.2%.year -1 in control stands and was not significantly increased by rainfall exclusion (P = 0.14). The mortality rate was, however, strongly decreased to only 0.2%.year - 1 in thinned stands (P < 0.001, explained variance of 18%): only 3 and 4 stem deaths were recorded in thinned and thinned dry stands, respectively. Similarly, at the stool-level, the mortality rate was 2.1%.year -1 in the control stands and increased to 2.5%.year -1 in the dry treatment, although this difference was not significant (P = 0.22). Thinning strongly reduced stool mortality (P < 0.001) to 0.03%.year -1 , regardless of the rainfall exclusion treatment (Figure 4B). Different letters indicate differences between treatments (P< 0.05).

Stem mortality probability was influenced by stem DBH, competition release and rainfall exclusion (Table 2, Figure 5). Stem mortality risk was higher for smaller stems: all trees under 3 cm of DBH died during the 15 years of the experiment while the mortality probability was almost null for DBH higher than 10 cm. Stem DBH explained 40% of the variance of the mortality risk, rainfall exclusion 4% and competition release 4%. Taking into account stem diameter, the mortality risk was higher in the dry treatment (P=0.003, Table 2) and decreased with competition release. The mean diameter of dead stems was higher in the rainfall exclusion treatment in unthinned plots (t = -3.16, P = 0.002, Figure 5A). The stool mortality risk was also higher for smaller stools, and rainfall exclusion shifted the mortality risk toward bigger stools in the unthinned treatment (Rainfall exclusion x Stool BA, P = 0.01, Supp. Figure S2). (variance explained by stem DBH 40%, Rainfall exclusion 4% and Competition release 4%).

Resprouting dynamics

None of the resprouts characteristics was influenced by rainfall exclusion (Table S1). Resprouts number was higher in thinned plots at the two inventory dates (Table S1, Figure 6A). Resprouts maximal height was higher in thinned plots after 5 years but did not differ among treatments after 15 years (Table S1, Fig 6B). The established resprouts (height > 1.30m) after 15 years represented 22% of the resprouts in thinned plots and 13% in unthinned plots. Resprouts basal area was higher in thinned plots in both inventory dates, but it strongly decreased between the two inventory dates in thinned stands resulting in a lower difference between treatments after 15 years (Table S1, Figure 6C). Resprouts basal area increased with stool size (P < 0.001) and with the proportion of stool basal area cut (P < 0.001), with a positive interaction between these two factors (P < 0.001) so that stool size influenced less the resprouting when it had not been thinned (data not shown).

Evolution of the stand basal area and density

Stand basal area increment was 0.16 ± 0.05 m 2 •ha -1 •year -1 in the control plots and decreased to 0.06 ± 0.07 m 2 •ha -1 •year -1 in the dry plots, although the difference was not significant (P = 0.3). Stand basal area increment was more than doubled in thinned stands, up to 0.41 ± 0.06 m 2 •ha -1 •year -1 (difference from the control : P<0.001), with no significant decrease under drier conditions (0.39 ± 0.06 m 2 •ha - 1 •year -1 , Table 3). Since the start of the experiment, a net annual loss of basal area (when mortality exceeds the growth of the remaining stems) occurred in 3 years in the control plots and in 7 years in the dry plots but never occurred in thinned plots.

The relative basal area (as the percentage of plot initial basal area) was strongly influenced by the interaction between time, rainfall exclusion and thinning (P = 0.001, Table 3). In the unthinned treatments, rainfall exclusion led to a slower progression of stand basal area (0.3 vs 0.8%.year -1 ; Figure 7A). After the initial reduction of basal area due to thinning, the increase in basal area was much faster in thinned stand, reaching 1.4%.year -1 with no difference between the thinned and thinned dry plots. The differences between thinned and unthinned plots thus progressively decreased over time: after 15 years thinned stands had recovered about 95% of their initial basal area.

Stem and stool relative density were both influenced by the 3-way interaction between time, thinning and rainfall exclusion (Table 3). Stem and stool densities decreased faster in the dry (slopes of -1.8 and -2.1 %.year -1 , respectively) than in the control treatment (slopes of -1.3 and -1.5%.year -1 ). In the thinned treatments, after an initial reduction of stem and stool densities by about 50% and 36% due to thinning, they remained stable in time (slope not different from 0) in both thinned and thinned dry plots. At the end of the experiment, stem density was still higher in the unthinned plots (Figure 7B), but from 2010 (after 7 years of treatment) stool density was not significantly different between thinned and unthinned plots (Figure 7C).

Discussion

Thinning and rainfall exclusion effects from the stem to the stand

Thinning strongly stimulated stem growth and increased the resistance to an extreme drought event, thereby confirming the interest of this management practice to alleviate drought effects on stem growth (Cabon et al., 2018b;[START_REF] Navarro-Cerrillo | Is thinning an alternative when trees could die in response to drought? The case of planted Pinus nigra and P. Sylvestris stands in southern Spain[END_REF]Sohn et al., 2016b). Partial rainfall exclusion, in contrast, had no effect on individual stem growth, probably because growth occurs mainly during spring rainy periods when the difference in tree water stress between rainfall exclusion treatments is low [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF][START_REF] Rodríguez-Calcerrada | Is selective thinning an adequate practice for adapting Quercus ilex coppices to climate change?[END_REF], and also because of a change in tree allocation in favor of wood compared to leaves and fruits [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF]. The rainfall exclusion effect appeared, however, as a shift of the sizemortality risk relationship toward the bigger individuals, which led to more frequent losses of basal area and consequently a slower basal area increment at the stand-level. In thinned stands, mortality was almost completely cancelled which, combined with the faster growth, resulted in a stand basal area increment more than doubled compared to the control, similarly to what was observed by Ogaya et al. ( 2019) in another Q. ilex stand. In our experiment, thinning proved to have larger effects on stem growth and mortality than rainfall exclusion, in line with recent studies at larger scales showing that competition is more important than climate in determining forest demography [START_REF] Fernández-De-Uña | Stand Competition Determines How Different Tree Species Will Cope with a Warming Climate[END_REF][START_REF] Ruiz-Benito | Patterns and Drivers of Tree Mortality in Iberian Forests: Climatic Effects Are Modified by Competition[END_REF][START_REF] Zhang | Half-century evidence from western Canada shows forest dynamics are primarily driven by competition followed by climate[END_REF].

The loss of stools from mortality appeared to happen faster than the loss of stems in the unthinned plots, and this tendency was reinforced by the rainfall exclusion. The consequence of this phenomenon is a progressive loss of unique genotypes within the coppice stand and an increasing level of clonality (ratio of stems overs stools), especially under more stressful dry conditions. This result suggests that multistemmed stools survive better to the increased competition for water under rainfall exclusion (Fig. S2), either because they grow in more favourable micro-habitats or because these individuals are better genetically adapted to drought. Thinning, on the other hand, reduced the initial stool density in the plots by 36%, but our experiment revealed that thinning maintained population sizes and did not induce further loss of unique genotypes in the long-term. Interestingly, all treatments converged toward a similar stool density by the end of our 15-year study. Drastic reductions in the number of unique genotypes (by thinning and/or mortality) can lead to the loss of rare alleles [START_REF] Danusevicius | The effect of tending and commercial thinning on the genetic diversity of Scots pine stands[END_REF].

However, previous studies have shown that moderate thinning intensities have rather limited consequences on genetic diversity in the long-term [START_REF] Danusevicius | The effect of tending and commercial thinning on the genetic diversity of Scots pine stands[END_REF][START_REF] Koskela | Translating conservation genetics into management: Pan-European minimum requirements for dynamic conservation units of forest tree genetic diversity[END_REF][START_REF] Valbuena-Carabaña | Coppice forests and genetic diversity: A case study in Quercus pyrenaica Willd. from Central Spain[END_REF]. On the contrary, thinning has been shown to increase acorn production in our experimental site [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF], so it may result in positive effects on sexual regeneration and seedling diversity (Lefevre et al., 2014).

Stem performances and thinning effects are linked to stem size and competition release

Stem size had a strong positive effect on stem growth, probability of survival and, to a lesser extent, growth resistance to drought. The strong decrease in mortality risk with tree size has been often evidenced in Mediterranean forests [START_REF] Colangelo | Size Matters a Lot: Drought-Affected Italian Oaks Are Smaller and Show Lower Growth Prior to Tree Death[END_REF][START_REF] Hülsmann | How to predict tree death from inventory datalessons from a systematic assessment of European tree mortality models[END_REF][START_REF] Ruiz-Benito | Patterns and Drivers of Tree Mortality in Iberian Forests: Climatic Effects Are Modified by Competition[END_REF] but contrasts with results from more humid forests [START_REF] Bennett | Larger trees suffer most during drought in forests worldwide[END_REF]. This positive effect of stem size on all stem performances may be linked to a better access of big trees to deep water [START_REF] Kerhoulas | Tree size, stand density, and the source of water used across seasons by ponderosa pine in northern Arizona[END_REF] and suggests that this holm oak stand, even 70 years after the last clear-cut, is still in a selfthinning dynamic where suppressed stems are more likely to die while dominant trees are still vigorous.

By selecting the biggest, most vigorous stems, thinning from below directly improved the mean stem growth, survival and resistance to drought. Thinning thus anticipated the self-thinning: we can estimate that about 60% of the stems and all the stools that have been cut would have died during the following 15 years (see Table 1; Figure 6b; Appendix. S2). This size-selection effect accounted for about half of the thinning effect on stem growth. Conversely, the fact that rainfall exclusion limited the positive influence of stem or stool size on the mortality risk could indicate a deeper depletion of water resources in the dry treatment.

Thinning effects on stem growth were stronger for bigger stems, which were more able to respond to competition release. A similar disproportionate effect of thinning on the bigger stems has previously been found for Q. ilex (Cabon et al., 2018b;[START_REF] Mayor | Growth response of holm oak (Quercus ilex L) to commercial thinning in the Montseny mountains (NE Spain)[END_REF], and for Q. faginea [START_REF] Cañellas | Growth response to thinning in Quercus pyrenaica Willd. coppice stands in Spanish central mountain[END_REF], while on boreal species [START_REF] Bose | Tree-level growth and survival following commercial thinning of four major softwood species in North America[END_REF] found that tree size had no influence on the thinning response and [START_REF] Jones | Canopy tree growth responses following selection harvest in seven species varying in shade tolerance[END_REF] that small trees responded more. In our experiment, thinning was made from below by removing only the smallest stems so the increase in light availability was very limited for the bigger, dominant trees. Competition release, that explained roughly half of the thinning effects on stem growth and resistance to drought, is thus probably related to a release of competition for water.

Competition for water is indeed more symmetric than competition for light [START_REF] Schwinning | Mechanisms determining the degree of size asymmetry in competition among plants[END_REF], so that the suppression of small trees or even of understory shrubs can increase water availability for big trees [START_REF] Giuggiola | Competition for water in a xeric forest ecosystem -Effects of understory removal on soil micro-climate, growth and physiology of dominant Scots pine trees[END_REF]. The stronger effect of competition release on bigger trees suggests that they had a greater capacity to expand their root system to take up this additional water. [START_REF] Kerhoulas | Tree size, stand density, and the source of water used across seasons by ponderosa pine in northern Arizona[END_REF] observed that bigger trees, but also trees growing in thinned stands, rely on deeper water sources than small trees or trees in denser stands, suggesting that thinning indeed foster root expansion.

In addition, the higher relative importance of competition release compared to stem size for resistance to an extreme drought event also points to the improvement of water availability as the main mechanism for thinning effects during an extreme drought. This improvement of water availability may arise both from a deeper rooting and from a decrease in stand-level transpiration [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF], delaying the depletion of water resources in thinned stands (Cabon et al., 2018b).

Remarkably, absolute values of local competition index had less influence on stem growth than competition release. Given that the stand history is homogeneous, the low explanatory power of neighbor basal area may reflect the fine-scale spatial heterogeneity of the soil rock fraction and water holding capacity. Indeed, in this very rocky soil, a low value of neighbor basal area probably reflects bad local conditions because the vegetation is expected to be in equilibrium with the local soil carrying capacity (Cabon et al., 2018a;[START_REF] Eagleson | Ecological optimality in water-limited natural soil-vegetation systems: 1. Theory and hypothesis[END_REF]. Absolute values of neighbor basal area after thinning are thus a poor predictor of competition intensity, as these values already differed before thinning due to soil heterogeneity. Interestingly, this may suggest that the same percentage of basal area removal could be applied successfully across variable local soil conditions, although this remain to be assessed in a more variable context.

Temporal evolution of treatment effects

Thinning effect on tree growth decreased with time as expected, but it remained positive 15 years after thinning even though the stand almost recovered its initial basal area. This long lasting thinning effect contrasts with studies in temperate forests where thinning effect can vanish in a few years (e.g. [START_REF] Bréda | Effects of thinning on soil and tree water relations, transpiration and growth in an oak forest (Quercus petraea (Matt.) Liebl.)[END_REF]. Slow growing species such as Q. ilex take longer to recover to their pre-thinning basal and leaf area. In line with this idea, [START_REF] Bose | Tree-level growth and survival following commercial thinning of four major softwood species in North America[END_REF] also observed that thinning effects lasted longer for shade-tolerant species. In Mediterranean forests, thinning effects were observed to be maintained to a similar level after 8 years in a Q. faginea stand [START_REF] Cañellas | Growth response to thinning in Quercus pyrenaica Willd. coppice stands in Spanish central mountain[END_REF] and to remain significant after 12 years in a Q. ilex stand [START_REF] Mayor | Growth response of holm oak (Quercus ilex L) to commercial thinning in the Montseny mountains (NE Spain)[END_REF] and after 20 years in Pinus halepensis stands [START_REF] Manrique-Alba | Long-term thinning effects on tree growth, drought response and water use efficiency at two Aleppo pine plantations in Spain[END_REF]. In contrast, in a mixed stand including Q. ilex and Q. faginea, [START_REF] Cotillas | Growth response of mixed mediterranean oak coppices to rainfall reduction: Could selective thinning have any influence on it? For[END_REF] found that growth stimulation after thinning vanished after only 2 years and attributed this to a fast and intense resprouting in their low-density stand. Here, resprouting was rather intense in the first years after thinning (Figure 5; [START_REF] Rodríguez-Calcerrada | Is selective thinning an adequate practice for adapting Quercus ilex coppices to climate change?[END_REF], particularly for big stools that had been intensively cut as resprouting is known to be stimulated by previous stool size and thinning intensity (e.g. [START_REF] Adamec | Development of Sessile Oak and European Hornbeam Sprouts after Thinning[END_REF]. We then observed a decrease over time in resprouts number and summed basal area in the thinned plots, meaning that resprout mortality was not compensated by resprout growth probably because the thinning from below did not increase sufficiently light availability in the understory. This contrasts with results by [START_REF] Retana | Regeneration by sprouting of holm-oak (Quercus ilex) stands exploited by selection thinning[END_REF] who observed that the decreasing number of resprouts was compensated by their increasing size, so that the total biomass of resprouts remained constant over time. Besides, no vegetative recruits were produced in the 15-year period following thinning in our study. This confirms that our low-intensity thinning from below was effective at controlling the resprouting dynamic in thinned stands, which probably participates in explaining the long-lasting effects of thinning in our site.

As a consequence of this absence of recruitment and of the self-thinning mortality in unthinned stands, stool densities are now similar in thinned and unthinned stands, and stand basal area is expected to become similar in the coming years. In the future, both stands will thus converge to similar structures, but generated either by management or by natural mortality. If drought sensitivity is partly genetically determined, then natural mortality may select for the best adapted trees. On the contrary, thinning may relax this selection and lead to maladaptive responses on the long-term, questioning the relevance of such management for forest adaptation to climate change (Lefevre et al. 2014). Once competition release has been offset by the growth of neighbors, frequent thinning may be needed to maintain positive effects on drought responses (Sohn et al., 2016a). Here, the positive thinning effect on growth resistance to drought has been highlighted only 3 years after thinning application and remains to be confirmed for further extreme drought events. In our experiment, thinning removed the small trees that proved to be more sensitive to drought in terms of growth and mortality, and has probably fostered root expansion, which may also induce long-term positive effects on drought resistance. Whether this will be enough to sustain an improved resistance to drought on the long-term, or whether natural mortality does a better job in selecting for drought resistance, remain to be assessed in the future decades.

This ageing oak coppice was still in a self-thinning phase 70 years after the last clear-cut, with a progressively decreasing number of individuals and a low stem and stand-level growth. Long-term increased drought by experimental rainfall exclusion led to the death of bigger trees than in the control treatment, which slowed the stand basal area increment. Thinning from below can cancel the drought effects by stimulating stem growth and avoiding mortality, even under chronic or extreme drought. The basal area increment at the stand level was twice faster in the thinned stand, while the loss of unique genotypes was not larger in the long-term than the one induced by natural mortality. Size-selection of the most vigorous stems and competition release for belowground resources both participate to explain thinning effects. Positive thinning effects diminished over time but remain significant after 15 years. By applying a moderate thinning intensity (-30%), keeping the dominant tree layer and keeping at least one stem on multi-stemmed stool, resprouting dynamic was rather limited, which may have participated to the long-lasting thinning effects and is important for fire prevention. A test of this thinning regime in other stands with differing climatic and edaphic conditions is required to confirm its benefits under a wide range of situations. Thinning has also been shown to improve stand-level above-ground productivity and fruit production while reducing water consumption [START_REF] Gavinet | Rainfall exclusion and thinning can alter the relationships between forest functioning and drought[END_REF], thus confirming that it can be an efficient way to improve stand and tree vitality and to adapt the Mediterranean coppices to climate change. Thinning oak coppices will also result in forests with bigger trees, which may be used for other products than traditional fuel wood.

Figure 1 :

 1 Figure 1 : Quercus ilex stem mean annual basal area increment (BAI) as a function of : A) treatments (means ± standard error of 99 to 173 stems per treatment, different letters indicate differences between treatments) and B) model projections and experimental data of stem BAI as a function of stem DBH and competition release (CR), confidence intervals of the models are represented in grey. CR=0 represent trees unthinned (represented by grey points), CR=0.25 is the median (white points) and CR=0.4 the upper quartile (crosses) of CR values in the thinned treatments. Point symbols are determined by the closest value of CR. Competition release is calculated as the change of neighbors' basal area in a 3m radius compared to its initial value before treatment application in 2003. Full model R 2 is 62%, stem DBH explains 53% of the variance and CR 7%

Figure 2 :

 2 Figure 2: Temporal evolution of thinning effects on Quercus ilex stem growth (mean ± standard error). Thinning effect decrease with time but is always significant (P < 0.001). Different letters indicate significant differences between treatments and periods (Tukey post-hoc tests).

Figure 3 :

 3 Figure 3 : Quercus ilex stem growth resistance to the extreme 2006 drought as a function of A) treatments (means ± standard error of 85 to 123 trees per treatment) and B) model projections of drought resistance as a function of stem DBH and competition release (CR) with confidence intervals in grey and ticks representing data distribution. Full model R 2 is 7%, stem DBH explains 2.5% of the variance and CR 4.5% (see Table2).

Figure 4 :

 4 Figure 4 : Quercus ilex stem (A) and stool (B) annual mortality rates as a function of treatment.

Figure 5 :

 5 Figure 5 : Quercus ilex stem mortality probability as a function of A) stem DBH and rainfall exclusion treatment -the inset shows dead stems DBH as a function of rainfall exclusion treatment in unthinned stands, asterisks denote the significance of a t-test (**, P < 0.01) -and B) competition release, for a median DBH value (7.3 cm). Grey areas represent the confidence intervals of the models. Ticks represent dead (top) and alive (bottom) trees according to their DBH and rainfall exclusion treatment (black = Control, grey = Dry). Full model R 2 = 48%

Figure 6 :

 6 Figure 6: Quercus ilex resprouts number per stool (A), height of the dominant resprout (B) and total basal area of the resprouts (C) as a function of thinning treatment (black: unthinned, grey: thinned) and period. Data are means ± standard error of 246 to 284 stools per thinning treatment (pooled over the rainfall exclusion treatment which had no effect). Different letters indicate differences between treatments and inventory date (P < 0.005, Tukey post-hoc tests)

Figure 7 :

 7 Figure 7: Evolution of Quercus ilex relative stand basal area (A), stem density (B) and stool density (C) as a percentage of pre-thinning values in the different treatments. Data are means ± standard error of 3 plots per treatment.

Figure S1 :Figure S2 :

 S1S2 Figure S1 : Map of one experimental block showing the spatial aggregation patterns of stems and stools. Each stem is a point (or a triangle for secondary species) and the first number represents the stool number followed by the stem number within the stool (#stool_#stem). Stems cut by thinning and dead stems in 2018 are also indicated

Table 1 :

 1 Stem, stool and stand characteristics before and after thinning. Means and range (within brackets) of the different metrics are shown. BA = Basal Area. 495 Neighbors basal area is calculated within a 3m-radius circle around each stem.

	496
	497
	498
		DBH	Neighbors BA	Stool BA	Plot BA	Stem density	Stool density
		(cm)	(cm 2 /m 2 )	(cm 2 )	(m 2 /ha)	(nb/ha)	(nb/ha)
	Control	6.6 [1.3 -15.3]	25.9 [6.1 -53.6]	69.9 [1.6 -531.4] 24.2 [23.4 -25.1] 5767 [4000 -6900] 3470 [1700 -4500]
	Dry	7.2 [1.7 -15.6]	26.8 [9.6 -57.0]	76.2 [4.1 -701.7]	27.6 [22.4 -31.0] 6133 [5100 -7000] 3630 [3300 -4200]
	Thinned	7.2 [1.0 -15.7]	26.3 [4.5 -46.5]	70.6 [0.7 -706.2] 27.9 [19.6 -32.6] 5867 [4000 -8000] 4030 [2600 -5900]
	Thinned Dry 6.6 [1.2 -13.5]	26.3 [6.2 -41.3]	70.5 [1.1 -453.7]	29.7 [22.5 -36.4] 6933 [6000 -8100] 4230 [3500 -5000]
	After thinning
	Thinned	9.0 [4.0 -15.7]	20.1 [3 -34.6]	86.8 [12.9 -486.0] 20.8 [15.9 -23.5] 3133 [2400-3900] 2470 [1700 -3200]
	Thinned Dry 8.3 [2.2 -13.5]	19.0 [4.3 -31.7]	77.4 [7.1 -393.3]	21.9 [16.9 -25.1] 3433 [3200 -3600] 2830 [2600 -3100]
	After 15 years
	Control	8.5 [3.4 -16.4]	32.0 [6.3 -64.5]	112 [11.3 -622]	27 [26.2 -28.2] 4167 [3400 -4900] 2400 [1500 -2900]
	Dry	9.2 [3.7 -17.1]	30.7 [9.1 -51.0]	128 [13.9 -802]	28.9 [23.6 -31.7] 4133 [3000 -4900] 2230 [1900 -2600]
	Thinned	10.5 [4.8 -18.9] 28.4 [3.7 -47.2]	117 [18.1 -654]	26.6 [21.1 -29.7] 3033 [2300 -3600] 2267 [1600 -2800]
	Thinned Dry 9.7 [3.6 -15.2]	25.9 [7.3 -45.1]	102 [24.4 -541]	27.9 [22.7 -32.2] 3300 [3100 -3400] 2733 [2500 -3000]
	499

Table 2 :

 2 Parameters of the best individual model for stem growth, growth resistance to drought and mortality risk. See model predictions on Figure 1B, 3B and 5.

	Stem growth -log(BAI+1), cm 2 .year -1
		Estimate	Standard error	P -value
	Intercept	-0.21	0.05	<0.001
	DBH	0.10	0.005	<0.001
	Competition Release (CR)	-0.38	0.25	0.14
	DBH:CR	0.08	0.03	0.01
	Basal Area of the neighbors (BAn)	-0.004	0.002	0.02
	r 2	0.62
	Stem resistance to drought -log(BAI2006/mean(BAI2003,BAI2004,BAI2005) +1)
	Intercept	0.09	0.08	0.2
	DBH	0.02	0.01	0.02
	Competition Release (CR)	0.64	0.15	<0.001
	r 2	0.07
	Mortality probability -log odds (mortality)
	Intercept	4.9	0.6	<0.001
	DBH	-1.04	0.11	<0.001
	Competition Release (CR)	-4.87	1.55	0.1
	Rainfall exclusion	0.8	0.3	0.003
	Pseudo-r 2	0.48

Table 3 :

 3 Results of the linear mixed model analysis for stand-level evolution of basal area increment, relative basal area, relative stem density, and relative stool density (as a fraction of initial situation at the start of the experiment in 2003) according to time and treatments. Significant P-values are highlighted in bold.

			Basal Area Increment Relative Basal Area Relative Stem Density	Relative Stool Density
		DF	F	P -value	F	P -value	F	P -value	F	P -value
	Rainfall Exclusion (E) 1	1.3	0.28	2.43	0.16	1.5	0.25	0.03	0.86
	Thinning (Th)	1	31.1	<0.001	58.8	<0.001	72.8	<0.001	6.3	0.03
	Time (Ti)	15	0.9	0.35	1415.3	<0.001	528.9	<0.001	426.3	<0.001
	E × Th	1	0.6	0.44	0.01	0.94	0.01	0.92	0.6	0.47
	E × Ti	15	0.1	0.76	31.1	<0.001	11.49	<0.001	11.3	<0.001
	Th × Ti	15	0.6	0.43	263.1	<0.001	440.7	<0.001	370.5	<0.001
	E × Th × Ti	15	0.1	0.92	12.125	0.001	11.14	0.001	10.7	0.001

Table S1 :

 S1 Resprouts characteristics according to treatments and inventory date : ANOVA Table.Significant P-values are highlighted in bold.

		D.F.	Number of resprouts per	Height of the dominant	Basal Area of resprouts
			stool		resprout	per stool
			LR  2	P	LR  2	P	LR  2	LR  2
	Rainfall exclusion (E)	1	1.8	0.17	0.1	0.7	1.4	0.23
	Thinning (T)	1	50.1	<0.001	23.1	<0.001	84.4	<0.001
	Inventory date (D)	1	1.4	0.24	19.2	<0.001	1.5	0.22
	E × T	1	1.2	0.28	0.7	0.38	0.4	0.51
	E × D	1	0.8	0.37	0.5	0.48	0.6	0.44
	T × D	1	0.9	0.33	9.4	0.002	4	0.04
	E × T × D	1	0.01	0.95	2.3	0.12	0.03	0.86

Acknowledgements

The Puéchabon experimental site belongs to the OSU OREME (UMS 3282) and the SOERE F-ORE-T, supported annually by Allenvi and the French national research infrastructure ANAEE-France (ANR-11-INBS-0001). This research was supported by the "Fond Stratégique Forêt Bois" (FSFB) funded by the French ministry of Agriculture and Food through the INNOV'ILEX project (MAA 2017). The authors are indebted to Serge Rambal and Richard Joffre who designed and initiated the long-term rainfall exclusion experiment in Puéchabon, and to Alain Rocheteau, David Degueldre and Karim Piquemal whose efforts greatly contributed in maintaining it throughout the years. SNP genotyping was performed at the Genome Transcriptome Facility of Bordeaux (Grants from Investissements d'Avenir, EquipEx Xyloforest ANR-10-EQPX-16-01).

Supplementary Material