Christophe Drouet
email: christophe.drouet@ensiacet.fr

Rjmat Carnot

Applied predictive thermodynamics (ThermAP). Part 2. Apatites containing Ni 2 +, Co 2 + , Mn 2 +, or Fe 2 + ions

Keywords:

Apatites are minerais encountered in many fields including geochemistry, nuclear and environmental sciences as well as medicine. This ubiquity is likely related to the diversity of ion substitutions that the apatite structure can accommoda te, maki ng of it an excellent "ion reservoi r" ei ther in natural settings or for the intentional production of doped systems with tailored properties. Despite this widespread interest for apatite compounds, however, only few studies are dedicated to study their thermodynamic properties. Yet, their knowledge becomes necessary for assessing stabi lity domains and understanding evolutionary trends in solution or upon heating, for example. Recently, the experi mental thermodynam ics of 33 phosphate apatite compounds (deriving from the composition M 1 o(POù;X 2 ) have been reviewed and their comparison allowed the development of the additive predictive mode! "ThennAP' (Applied Predictive Thermodynamic:s) capable of adequately predicting properties such as standard enthalpies (liHr), Gibbs free energies of formation (iiGr), or entropies (S 0 ) at T 298 K. for any composition involving ions among M 2 • ca2•. Ba 2 •• sr2•. Mg2•, Cd 2 •• Pb 2 •• cu 2 •• zn 2 • and X OH , f , Cl or Br . Although experirnental data for apatites involving other divalent cations such as Ni 2 •, Co 2 •. Mn 2 • or Fe 2 • do not seem to be available, the exploration of apatites doped with these ions is appealing from a practical and fundarnental viewpoint. for example for understanding geochemical events, or when using apatite precipitation for the elimination of metal cations from industrial wastewaters, or else for conferring magnetic properties to apatite systems in medicine. Based on multiple physico chemical correlations, the present contribution extends the additive predictive mode! ThermAP to Ni , Co , Mn(II) and Fe(II) doped apatites. It provides for the fi rst time estimations of enthalpies, Gibbs free energies of formation and entropies, unveiling the general stability ranking Mn(II) apatite > Fe(II) apatite > Co apatite � Ni apatite. This additive approach also allows one to estimate these properties for any composition in view of enabling thermodynamic calculations for applicative or fundamental purposes.

Apatites are ubiquitous minerais found in many domains, as was commented previously [1 ). In particular, phosphate apatites deriving from the generic formula M t <JPO 4 )6'( 2 represent an important subclass [2 4). Whether for understanding natural phenomena such as biomineralization and geological phosphate rock formationfweathering, for explaining minerai evolutions upon immersion or heating, or else for optimizing industrial processes involving such compounds, the knowledge of their thermodynamic properties becomes a primordial prerequisite (besides adequate physico chemical characterization).

Unfortunately, this thermodynamic aspect is sometimes disregarded, which can be concluded from the rather limited literature existing on this domain. In a previous recent work [1 ), the available reported experimental data on phosphate apatite end members, ranging from the 1950s to nowadays has been reviewed and compiled in a single reference document [1 ), giving an overview of standard formation energetics (LiHl and iiG;J, entropy s 0 and molar heat capacity c;,m• In that review work, the compilation of data for 33 phosphate apatite end members involving the ions ea 2 •, Ba 2 •, Sr2 ♦ • Mg2•, Cd 2 •, Pb 2 •, Cu 2 •, Zn 2 •, PO! , HPO� , OH , F , Cl and/or Br , then made it possible to quantify the "contributive" rote of each of these ions in the energetics of apatitic compounds, thanks for a contributional additive model [1,[START_REF] Drouet | Applied predictive thermodynamics, Therm'AP free calculation program[END_REF] that can be referred to as ThennAP for Applied Predictive Thennodynamics.

In the present work, the idea was to expand this list of ions in phosphate apatite systems to the cations Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ . Indeed, these ions may be found in some contaminated water from mines or in industrial waste waters, thus leading to pollution issues [START_REF] Thiruvenkatachari | Permeable reactive barrier for groundwater rernediation[END_REF], and their precipitation in the form of a stable solid phase can then be seen as a clever removal process for depollution pur poses [7,8]. In this field, the setup of permeable reactive barrier systems for groundwater remediation has in particular been devel oped [START_REF] Thiruvenkatachari | Permeable reactive barrier for groundwater rernediation[END_REF], including using apatite as solid remediator [7,8]. Also, these ions exhibit specific electron features that may prove of interest in the design of magnetic apatite systems (e.g. for the mag netic guidance of doped nanoparticles towards a diseased tissue) [9] or for exploiting cathodoluminescence features. In the geo chemical field, manganese is also known to be found in igneous, metamorphic and sedimentary apatites [10,11]. The possible pre cipitation of Mn , Ni , Co , or Fe containing apatites or sorption of these ions on apatite were also reported (see for example refer ences [12 19]).

In the present contribution, I will first summarize the main aspects concerning the ThermAP additive thermodynamic model that was recently developed to evaluate the energetics of phos phate apatite compounds (the initial paper, reference [1], being considered as ''Part 1" of the series of papers dedicated to the ThermAP approach). Then, I will unveil and exploit the existence of multiple correlations between the nature of doping ions and the energetics of apatitic phases, in view of deriving the energetic contribution of Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ ions in an apatitic matrix (thereby extending the ThermAP model to these ions). I finally use these findings to draw, for the first time, estimates of the heat of formation, Gibbs free energy of formation and entropy of Ni , Co , Mn(II) and Fe(II) doped apatites; and some illustrative exam ples of the use of such data are given.

Materials and methods

Summary of the ThermAP additive model applied to phosphate apatites

The aim of this section is to recall the main aspects of the ''pre dictive thermodynamic model" ThermAP that was developed recently for evaluating and understanding the energetics of phosphate apatite compounds: the full concept has been explained in detail elsewhere [1,[START_REF] Drouet | Applied predictive thermodynamics, Therm'AP free calculation program[END_REF] and will thus only be overviewed here.

When (reliable) experimental data are available, these should of course be considered in priority for solving any thermodynamics related question. The problem arises when no experimental data is accessible for a given composition of inter est. In this case, the best approach then lies in estimating this lacking information. In the case of apatites, several ''predictive" methods have been tested [1], including the Volume Based Thermodynamics approach (VBT) developed in the literature for investigating other solid systems [20,21] as well as some additive contributional models. The VBT approach was found to be interesting but limited, in the case of apatites, by the fact that compositions with alkali earth cations had a clearly distinct behavior as compared to other divalent cations. In contrast, an additive model proved to be more ''universal" and well adapted to estimate (within less than 1% of relative error, and often less than 0.5%) the energetics of phosphate bearing apatite com pounds involving not only alkali earth cations (i.e. in the ''s block" of the periodic chart of elements) but also cations from the ''d block" (Cu, Zn, Cd) and in the ''p block" (Pb).

This ''predictive model" was elaborated from a thorough com parative analysis of experimental datasets accessible in the litera ture about 33 phosphate apatite end members involving the following ions: Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ , Zn 2+ , PO 3 4 , HPO 2 4 , OH , F , Cl and/or Br . It is based on a simple ''additive" concept, stating that the energetics (e.g. DH f or DG f ) of a complex oxide are directly related to the energetics of the constituting bin ary oxides/compounds (e.g. CaO, P 2 O 5 and CaF 2 in the case of fluo rapatite Ca 10 (PO 4 ) 6 F 2 ) provided that a corrective factor is applied on each contribution. This type of additive model was previously used for evaluating the thermodynamics of other types of solids such as zeolites, simple phosphates or silicates [22 27]. Generally speak ing, by considering a family or sub family of compounds, it becomes possible to constrain the corrective factors to be applied in this additive model. In practice, this refinement is undergone by finding the best fit between the values of DG f , DH f or S°calcu lated by the model and experimental ones, using a solver approach.

The concept of such corrective factors is by itself easy to catch: the energetic contribution of a given ion, e.g. Ca 2+ for example, is likely to be rather similar for solids presenting the same crystalline structure (in the present case, the apatite structure). Using such an [1]). The relative uncertainty associated to each contribution is of the order of 1%.

Contributing sub-units in the apatite formula (T = 298 K, 1 bar) additive model, the contributions of various ions were recently evaluated for phosphate apatite compounds [1,[START_REF] Drouet | Applied predictive thermodynamics, Therm'AP free calculation program[END_REF], with an appre ciable fit to experimental data (typically with mean relative errors less than 1% and often less than 0.5%, as mentioned above). In this case of phosphate apatites, table 1 summarizes these contributions in terms of enthalpy of formation, Gibbs free energy of formation and entropy, as reported previously [1]. Therefore, the ThermAP model (for which additive contributions are given in table 1) was established by determining the energetic contributions (in enthalpy or free energy or entropy) of ions once contained in a solid apatitic matrix. These contributions were eval uated by comparison to experimentally available data on several substituted apatites using a solver type of approach (multifactorial resolution of equations). In a first stage, the contribution of each ion was initially set to the enthalpy (or free energy or entropy respectively) of the related binary oxide. For a given apatite com position, this led to a rough estimate of the DH f (or DG f or S°) of the apatite compound. In a second step, these rough estimates were then corrected (by a multiplying corrective factor of the bin ary oxide energetics) so as to minimize the difference between estimates and experimental values. This minimization was made possible by solving a multi equation mathematical system (an equation providing from each apatite composition available exper imentally). The final set of energetic contributions found for each ion tested so far is reported in table 1. For example, as shown in table 1, the Gibbs contribution g i of one calcium ion in a phosphate apatite compound, denoted g Ca , is 740 kJ/mol. This value corre sponds, in this additive concept, to the contribution of one calcium ion to the Gibbs free energy of formation DG f of the considered apatite compound. When compared to the free energy of formation of the binary oxide involving calcium, namely CaO (DG f (CaO) = 603.1 ± 0.9 kJ/mol [28]), one can see that the correcting factor a corr,G,Ca = 740/603.1 = 1.22 (7) has to be applied to the free energy of CaO to find the Gibbs contribution of one calcium ion in apatite: g Ca = a corr,G,Ca * DG f (CaO). Then, taking for example the case of fluo rapatite, Ca 10 (PO 4 ) 6 F 2 , the application of the model leads to DG f estimated = 10 * g Ca + 6 * g PO4 + 2 * g F = 12835.9 kJ/mol, and in a similar way to DH f estimated = 13,544 kJ/mol and S°e stimated = 770.3 J mol 1 K 1 , which are in good agreement (within 1%) with the recommended experimental based values: DG f = 12,899 kJ/mol, DH f = 13,658 and S°= 778 J mol 1 K 1 (see experimental values compiled and discussed in reference [1]).

g i /(kJ Á mol 1 ) h i /(kJ Á mol 1 ) s i /(J Á mol 1 Á K 1 )
In a more formal way, this additive model can thus be described via a linear combination of the energetics of the constituting binary oxides/compounds. Thereafter, the energetics of an apatite com pound M 10 (PO 4 ) 6 X 2 (or a more complex solid solution. . .) is treated as the linear combination of the energetics of the constituting bin ary compounds MO, P 2 O 5 , and MX 2 . For DH f and DG f , in particular, this gives:

DH f ðcomplex oxideÞ X i v i ða corr;H i DH f ðbinary oxideÞÞ; DG f ðcomplex oxideÞ X i v i ða corr;G i DG f ðbinary oxideÞÞ;
where v i represents the stoichiometry coefficient for the ion i found in a given chemical formula, and where a corr,H and a corr,G are the respective corrective factors. This model has the advantage to allow estimating, at least in a first good approximation, thermodynamic data not available in the literature, especially for solid solutions of any complex composi tion as well as for non stoichiometric systems, as was also shown previously [1].

Methodology used in this work for extending the ThermAP model to nickel, manganese, iron and cobalt divalent ions

The above recapitulated model was initially established thanks to comparisons between calculated values and experimental data reported for apatite end members involving the cations Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ and Zn 2+ . The present contribution focuses on the possible extension of this model to other divalent cations of potential applicative interest, namely Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ . With this objective in mind, a methodology was thus followed, in two steps:

First, the energetic contributions obtained with the model for Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ and Zn 2+ cations were re examined for underlining some correlations between physico chemical characteristics of these ions and the overall ener getics (DG f , DH f ) of the related apatite end members. Second, after discussion the plausible reliability of extending such correlations to the cases of nickel, cobalt, manganese and iron, these correlations were exploited to evaluate the energetic contri butions of the latter in an apatite structure.

Results and discussion

3.1. Evaluation of the energetic contributions h i and g i of Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ in phosphate apatites, based on the ThermAP additive model

The energetic contributions of Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ and Zn 2+ cations in an apatitic setting have been listed, using the ThermAP approach, in along with the energetics of the corre sponding binary oxides (typically from thermodynamic databases such as reference [28]) and the ''corrective" factors as defined in the previous section. One observation is that the contributions h i and g i for the divalent cations M 2+ tested are more negative than the value characteristic of the corresponding binary oxide MO. This indicates that the ion in the apatitic array is more stable than in the corresponding binary oxides. ln the same kind of idea, a neg ative value of experimental enthalpies of formation calculated "from the oxides" (and not from the elements in their standard state) has been notioed for apatitic compounds, whether in the case of non stoichiometric nanocrystalline apatites [29] or for lanthanide doped apatites [30]. lt may be noted that the ions considered in this table are located in very distinct places in the periodic chart of elements, in the "s", "d" or "p" blacks. At this point, it was thus interes ting to look for possible correlations versus specific features of these very diverse ions. ln particular, the plots of Oro rr,H or Ororr.c versus the effective ionic radius f;on (taken from Shannon and Prewitt reference work [31 ], data for coordination number 6), can reasonably well be described by linear fits, with R 2 = 0.9688 and 0.9248 respectively, as shown in figure 1. Interestingly, the ions are well distributed along the line, independently of their position in the periodic clas sification of elements. The variation of these corrective factors ver sus other parameters was also tentatively followed. In particular, a respectable quadratic correlation (R 2 $ 0.7 to 0.8) was found ver sus the first ionization energy of the elements (data from reference [START_REF] Dean | Lange' s Handbook of Chernistry[END_REF]). This datum, referring to the easiness of removing one elec tron from a given atom, is indeed related to the relative attraction of the nuclei for peripheral electrons and can thus be considered for understanding energetic phenomena. Another correlation was sought versus the parameters D H O 2 M or D G O 2 M, proposed by Tardy and collaborators [24,[START_REF] Tardy | [END_REF] to describe differences between the energetics of binary oxides and of the corresponding aqueous ions. Again, a quadratic relationship (with R 2 > 0.78) could be identified.

The existence of such clear correlations between physico chemical parameters of the cation M 2+ in phosphate apatite and their energetics (enthalpy and Gibbs free energies of formation) strongly supports, again, the use of this additive model for thermodynamic explorations of apatites, especially when no experimental data are available.

On the basis of these correlations, and keeping in mind that they were obtained for very diverse divalent ions (from either the ''s", ''p" or ''d" block) of the periodic chart of elements, it then became appealing to extrapolate the energetic contributions of other divalent ions. As mentioned above, particular interest is found for 4 additional divalent cations of transition metals: Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ . With this objective, the a corr = f(r ion ) relation ships determined here both for enthalpies and Gibbs free ener gies were used to evaluate the corrective factors a corr,H and a corr,G corresponding to these 4 new ions. The obtained values are tabulated in table 3 along with the energetics of the corresponding binary oxides (data from reference [28]), as well as the resulting contributions h i and g i for these 4 ions. It may be noted that the electronic structure of Mn(II), Fe(II) and Co divalent ions allows either high spin (HS) or low spin (LS) configurations, and the ionic radius is known to depend on the configuration. However, the cal cium Ca(1) and Ca(2) sites in apatite are typically coordinated with OH -, O 2 , F , Cl . . . which are known to be low field ligands; in these conditions only the high spin configuration is suspected to occur in apatites [10]. Therefore, the HS ionic radii values were considered for the construction of table 3.

Again, these estimated contributions are more negative than the energetics of the binary oxides, indicating that for these ions as well, chemical environments in the apatitic structure are energetically favorable as compared to simple oxides. It is impor tant to remark also that the two other correlations pointed out on figure 1 continue to be valid for the Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ cations (indicated by cross marks), thus corroborating these values. Moreover, contributions of the same order of magnitude although not identical due to different chemical environments were also found for Fe 2+ , Ni 2+ and/or Co 2+ in additive models reported by Chermak et al. [START_REF] Chermak | [END_REF] and La Iglesia [START_REF]A La lglesia[END_REF], respectively for silicates and simple phosphates.

This estimation of energetic contributions h i and g i for Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ ions in phosphate apatites, using the ThermAP approach, now theoretically allows one to evaluate the standard enthalpy of formation DH f and Gibbs free energy of for mation DG f (from the elements taken in their standard state), at T = 298 K, for any phosphate apatite composition containing cations among Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ , Zn 2+ , Ni 2+ , Co 2+ , Mn 2+ and Fe 2+ and anions among OH , F , Cl and Br . Consequently, the knowledge of both DH f and DG f could then also lead to an estimation of entropy. However, as was mentioned previously [1], the propagation of errors (relative errors of the order of 1%) from such estimated DH f and DG f to entropy could lead to non negligible uncertainties. Therefore, it is generally a more recommended approach to estimate first entropies by a par allel method, and to combine it to DH f (or DG f ) predicted by the ThermAP approach (unless known experimentally) to finally reach DG f (or DH f respectively), thus leading to a more reliable and self consistent set of DH f /DG f /S°values (i.e. obeying the law DG f DH f T Ã DS f ). For information, it can be added here that comparisons of calcu lated data using this model to existing experimentally available data on various compounds incorporating phosphate and Ni 2+ /Co 2+ /Mn 2+ /Fe 2+ has been performed, for stressing the validity of this ThermAP approach even beyond purely apatitic compounds. Although differences are expected due to this ''non apatitic" nat ure, a very good accord is found (typically with relative errors 62%), substantiating the validity of the model (see Supplementary Information, figure S1). Considering what was pre viously observed systematically for other cations using the ThermAp approach, and the comparative data of figure S1 for non apatitic compounds (for which a greater error is expected), then the estimated relative errors on the corrective factors a corr and on g i and h i energetic contributions are expected to be lower than 2%.

The following section will first examine the question of the entropy, and will then combine entropy and enthalpy evaluations to derive Gibbs free energies, with the objective to ultimately pro pose self consistent thermodynamic estimates for Ni , Mn(II) , Fe(II) or Co apatites.

Entropy evaluations and derived considerations

Entropies or heat capacities have only rarely been determined for apatitic systems (see a compilation in reference [1]); therefore methods are needed to approximate them for unknown systems. The ThermAP approach proved to allow such evaluations, but to this day only for apatite compositions involving the ions Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ , Zn 2 , OH , F , Cl or Br (including solid solutions or non stoichiometric samples) for which some entropy values had been reported in the literature. But in the case of Ni 2+ , Co 2+ , Mn 2+ or Fe 2+ ions, no experimental values seem to be accessible. In this context, other estimation methodologies have thus to be utilized to evaluate entropies.

In reference [1], Drouet tested two additional approaches: (1) a Volume Based Thermodynamic (VBT) method, and (2) a method based on a comparison to the sum of entropies of the constituting binary compounds. VBT is an approach developed by Glasser and Jenkins (see for example [20,21]) for application to various types of systems. In a few words, this approach uses the fact that many physical characteristics of a solid, including its entropy, can be rather straight forwardly related to the volume ''V m " occupied by one unit formula. When applied to apatites with V m ranging between 500 and 700 Å 3 , a linear S°(apatite) = f(V m ) correlation was established [1], in first approximation. However, if considering for example the nickel, manganese (II), cobalt, or iron (II) hydrox yapatite end member compositions, i.e. M 10 (PO 4 ) 6 (OH) 2 with M = Ni, Mn(II), Co or Fe(II), the evaluated values of V m go down to 382 Å 3 and thus fall significantly away from the 500 to 700 Å 3 range, making the previous correlation probably inaccurate for these 4 additional ions.

The second approach tested and discussed previously by Drouet [1] ,.

-

-8729 1017

Co '' ,.

•7636 •8313 939

Ni 1•

,.

-7$53 -8273 791

C u 2• F• Z n 2+ F• C d 2• F• •9Jl8 •9122 .9777
.9753 1041 1061

,,,, . a• Fe l• .,. -7220 -7224 . 7939 -7915 816 841 correlation was pointed out for values of D in the range 750 to 1150 J mol 1 K 1 • For more information on this result, please refer to the previously published paper (reference [11). In the case of Ni , Co , Mn(II), or Fe(II) hydroxyapatite end members, the entropy sum I:s for binary compounds lies for example between 770 and 980 J mol 1 K 1 (as calculated from reference [28]), thus contained in the previously considered range of 750 to 1150 J mol 1 K 1 • Consequently, the previously determined Jinear like 5 ° � f(I: 5 ) relationship may probably also be applied to the cases of nickel, cobalt, manganese (Il) and iron (Il). This Jeads to entropy values 5 ° ( at T = 298 K) Jisted in the first column of table 4. Interestingly, rather similar values (within Jess than 3.5% of dif ference) are also found if estimating entropy from reactions of the type:

C u 2' C u 2 • Cl• Br z n i+ Z n 2+ c1• Br• C d 2' C d 2• c1• Br• s• Pb 2• Pb 2' Pb 2• Pb 2' Pb OH• F• Cl•
Ca1o(P04) G {OH) i + 10Ni0-+ Ni1o(P04) G {OH) i + 10Ca0
and considering a change of entropy, à5 °,e action, close to zero (as often approximated for reactions in the solid state for which the global variation of volume is small). This good accord thus supports the orders of magnitude found here for standard entropies 5 ° of the theoretical end member compositions Ni 10 (P0 4 ) 6 (0H) i , Cot<{P04)s(0H}i, Mn10(P04)s(0H)i and Fe10(P04)s(0H)i.

The enthalpy contributions hi estimated in the previous section (see table 3) for Ni 2 •, Co 2 •. Mn 2 • and Fe 2 • allow the evaluation of the enthalpy of formation !!.Hl of the above 4 end member compo sitions. Theo, by combining 5 ° (in the form of the entropy of forma tion from the elements Li5;) and !!.Hl, the related Gibbs free energies of formation LiGl can also be expressed thanks to the rela tionship t!.Gl !!.Hl T * t!.5;. Ail these values have been adde d in table 4; they represent the best estimates to date for Ni , Co , Mn(II) or Fe(II) hydroxyapatite end members (whether theoreti cal or accessible by experiment).

If the same approach is carried out for the fluor , chlor or brom apatites series (considered at Jeast theoretically), then it can similarly Jead to the estimation of their standard thermody namic properties. This has been done, and figure 2 gathers this information in the form of an easy to use periodic table of recom mended thermodynamic data (shade d areas on this figure correspond to previously published data [1] and were thus not reproduced here for the sake of clarity). In ref [1 ], it was shown that the ThermAP contributive approach was valid not only for OH but also F, Cl and Br, independently of the cations present. The validity of these estimations for F , Cl and Br bearing apatites can thus be reasonably expected.

The knowledge of 5 ° for these end member compositions (whether theoretical or actually obtainable by experiment) also makes it possible to determine the entropy contributions "st rela tive to the ions Ni 2 •, Co 2 •. Mn 2 • and Fe 2 • ions when incorporated in phosphate apatites. They can easily be obtained with the ThennAP approach by simply subtracting the contributions of the phosphate ion PO! and of the related monovalent anion X = OH , F , Cl or Br contained in the M 10 (PO 4 )6'( 2 formulas (for M = Ni, Co, Mn or Fe). The obt ained s i contributions for these 4 ions have thus been calculated and added in table 3. No experimental data on the entropy of Ni /Co /Mn /Fe apatites being available, entropies 5 ° can however be tentatively estimated for other phosphate bearing compounds using the Si contributions given in the present study in order to validate the order of magnitude of these data. This was performed with some metal phosphates (see 5upplementary Information, figure 52). Despite the non apatitic nature of these compounds, the entropy estimated here for such metal phosphates have the same order of magnitude (with less than 20% of relative error) as experimental results. In the case of apatites, let's recall that maximal errors of up to 15% had been observed in some systems (see details on reference [11). This arises from greater uncertainties on reported entropies as opposed to enthalpies or Gibbs free energies. However, since in the case of apatites the energetics are mostly driven by enthalpy rather than entropy (the enthalpy contribution is always clearly predominant for apatites as pointed out in reference [1 ]), then these compar isons allow considering the entropy contributions Si estimated in the present manuscript as reasonably well fi tted. 5uch s i values can prove helpful for evaluating with the ThennAP mode! the thermodynamic properties (.iH;, .iG; and 5 °) of apatite solid solutions or nonstoichiometric samples involving some amount of nickel, cobalt, manganese (Il) or iron (Il). Again, it is reminded that the recommended approach in this case (for limiting propagated errors) consists infirst evaluating the entropy 5 ° and the enthalpy .iH;, by adding respectively the Si and hi contributions of the involved ions, and then to calculate the free energies .iG; through .iG; .iH; T * .i5;. For information, it can be remarked that the comparison of the .iG; thus evaluated to val ues calculated directly fr om the gi contributions (listed in tables 2 and 3) leads to an appreciable accord, with a mean relative error of 0.8%.

Sorne examples of uses and implications

As a graphical illustration of the above fi ndings, figure 3 shows for example the evolution of .iH; and .iG; ( at T = 298 K. 1 bar) for solid solutions between hydroxyapatite Ca 10 (PO 4 ) 6 (OH h and respectively the end members Ni t <JPO 4 ) 6 (OH) i , Co 10 (PO 4 �(OH) i , Mn 10 (PO 4 ) 6 (OH) i and Fe 1 0(PO 4 ) 6 (OH) i . It should be reminded here that these graphs do not imply "automatically" that ail the mem bers of these solid solutions can indeed be obtained by experiment; indeed a maximal limit of calcium substitution by nickel, cobalt, manganese (Il) or iron (Il) in hydroxyapatite may exist due to crys tallographic or other constrains. However, several literature reports (see introduction section) have envisioned the possible incorporation of these ions, at least to some extent, in the apatitic structure. Any such composition accessible by the experiment will then have a given set of thermodynamic data which may be found in fi gure 3 ( or in equivalent graphs for fluor , chlor apatites ... ). In any case, the results obtained in the present contribution indicate that the fol lowing stability order can be evidenced (for a given amount of metal ion): Mn(II)-apatite > Fe(/1)-apatite > Co-apatite ;;;,: Ni-apatite with Cu apatite (studied previously [1 ], with .iG; � 6524 kj/mol) being even less stable than nickel apatite. This stability ranking is not without remembering the order of the Mn, Fe, Co, Ni, Cu elements in Mendeleev's classification, however this cannot be "generalized" as, for example, Zn apatite stability (with a .iG; � 8623 kJ/mol [1]) remains between those of Mn(II) and Fe(II) apatites.

It may be remarked that in ail additive models, the energetics of mixing is, by defi nition, neglected. Previous comments [1] however pointed out that the order of magnitude of enthalpies of mixing were extremely small (only a few tens of kj/mol) compared to .iH; or .iG;, and in any case lower than the propagated uncertainty linked to the .iH; and .iG; estimates.

Another illustration of this work findings can for instance be found in the exploration of the energetics of precipitation reactions involving metal ions. Their precipitation in the form of apatitic solids has for example been proposed for the treatment of indus trial waste waters; but what is the driving force of such reactions? Is this remediation route appropriate and energetically favorable? The knowledge of the energetics of apatites loaded with various metal ions now allows such questions to be answered. Let us con sider, for example, the case of Ni 2+ ions and their possible The calculation of the change in standard Gibbs free energy of reac tion, ôG 0 re.1ction, leads to values ranging from 686 to 541 kj/mol for x increasing from O to 10. This shows the strong spontaneous character of this reaction at T= 298 K, in standard conditions. Theo, in any given situation, it becomes possible to calculate the actual change ôG,e.1c ti on from the relationship: ôG 0 ,eacrion = ôG 0 ,e.1ction + RT Ln(K) where K is the equilibrium con stant of the reaction Supposing for example a simple scenario with the following ionic activities in the medium (Caf, -: 0) = 0.1 M, (Nit�)= 0.2 M and (HP� ( a q) l = 0.1 M, and a pH at 8, this would lead to values of ôG,e.1ction between 321 and 215 kj/mol, thus still energetically favorable. By calculating this amount for various val ues of ionic activities, the most appropriate conditions can then be searched. An example is given for illustrating this point on figure 4. ln this case, one can see that increasing the amount of cal àum and/or phosphate in solution, the ôG reactio n value gets more and more negative (favorable situation). If the same plot is done for other nickel activities, on the other hand, (Nit�) can be shown to leave almost unchanged ôGre.1c ti on which remains between 316 for 0.2 M of nickel and 301 for 0.001 M. Therefore , it is pos sible to use such calculations to identify the most crucial experi mental parameters to modulate, and thus the best preàpitation scenario to be used in a given industrial setting.

Concluding statement

The above example illustrated one possible use of the thermo dynamic data reported in the present contribution; however the evaluation of energetic contributions of nickel , cobalt , man ganese (Il) or iron (Il) doped apatites should now allow more numerous and precise thermodynamic calculations involving these phases, whether in natural settings (geochemistry, bone diagene sis ... ) or for applicative purposes. The general stability ranking Mn(II) apatite > Fe(II) apatite > Co apatite ;;,: Ni apatite was established.

From a more general viewpoint, this work stresses again the interest of additive contributional approaches for the evaluation of unknown thermodynamic properties. Data for Ni 2 •, co 2 •, Mn 2 • and Fe 2 • ions obtained in the present work will therefore help enriching the apatite database of the ThennAP calculation sheet, freely accessible online (www.christophedrouet.com/thermAP).

6 RR 1 FIGURE 1 .

 611 FIGURE 1. Identification of correlations between the corrective factors Ocorr)1 and a,..,. . .: (for divalent cations in phosphate apatites) and the corresponding Ionie radius ( r 100), the first ionization energy, à tt c2 M and à,:O 2 M. Dark circles (e) represent divalent ions data treated in reference (1) and crosses ( x) represent the additional ions Ni�. Co 2 *, Mn� and Fe 2 * treated in the present work.

• 3 •

 3 Energetic contributions evaluated in the present work for Mn 2 •. Fe 2 *, Ni2* and Co� ions in phosphate apatites. Values estimated in this work from the a .,,, ,• f( r 100 ) linear correlations, see text. Effective Ionie radii (coord. 6) from Shannon and Prewitt, reference (31); high-spin (HS) radii were considered when applicable.

TABLE 4 •

 4 Recommended thermodynamic data (at T•298 K. 1 bar) for Ni 1 o{P0 4 )E(OH)2, Co,o(PO.J,,(OH}i, Mn 1o{ PO.)6(OH) 2 and Fe 1 o{PO.J,,(OHh theoretical end-members. Estimated by applying the previously determined relationship S°(M 1o{P04 )6X 2) � 1.0428 • I:s (where I:s is the entropy sum for binary compounds 9 MO+ 3 P2O5 + 1 MX2). • Estimated in this work by extension of the ThermAP additive model c Estimated from application of âGr âHr -T • 115r.Periodic Table of Phosphate Apatite end-members* gl+ M g2 +Mn� OH' Mn> F" Mn:lo Cl' Mn 1' . , , Co, Ni aparires: wirhin 2 % for LIH / and LIG/ and within 20% fors•

  o• Ni 1• .,.

FIGURE 3 .

 3 FIGURE 3. Craphical evolution of âH 1 {a) and âC 1 {b), for Ca-Ni-hydroxyapatite, Ca-Co-hydroxyapatite, Ca-Mn(U)-hydroxyapatite and Ca-Fe{U)-hydroxyapatite solid solutions {whether theoretical or obtainable by experiment}

  Example: at pH = 8, results for x N; = 5. s

FIGURE 4 .

 4 FIGURE 4. EX.lmple of evolution of .ô.Greact1on (precipitation of nickel-<loped hydroxyapatite Cil(,o xJl'li.{P04)s(OH)2) at pH• 8, for (Ni��)• 02 M, and consider ing an 00average• nickel amount of x • 0.5, versus calcium and hydrogenphosphate activities in solution (at 298 K} preàpitation in the form of a Ni doped calcium hydroxyapatite, at pH 8. This scenario can be described by the expected reaction: (10 x)Caf, "':o + xNif, "' :o + 6HPO! (aqJ + 8OH (aqJ -> Ca ( 10 x)Ni,(P04MOH) i + 6H 2 O ( liq)•

TABLE 1

 1 Contributive energetic values in the ThermAP predictive model, for the estimation of DG

	f , DH f or S°of phosphate apatites (at T = 298 K, 1 bar) from their composition (reproduced

TABLE 2

 2 Energetic contributions of Ca 2+ , Ba 2+ , Sr 2+ , Mg 2+ , Cd 2+ , Pb 2+ , Cu 2+ and Zn 2+ divalent cations in a phosphate apatite matrix (ThermAP model, from Drouet, JCT 2015[1]), and comparison to the corresponding binary oxides.

	Ion	Enthalpy

a Contribution to the energetics of the considered apatite phase. b Effective ionic radii (coord. 6) from Shannon and Prewitt, reference

[31]

.

TABLE 3

 3

  led to the approximated relationship: S°(M 10 (PO 4 ) 6 X 2 ) ffi 1.0428 * P

S , where

P

S represents the sum of entropies of the cor responding binary compounds: 9 MO + 3 P 2 O 5 + 1 MX 2 . This

Br• * whether theoretical or experimentally accessible FIGURE 2. Periodic table summarizing recommended thermodynamic values (at T • 298 K. 1 bar} for theoretiGII M 10 {P0 4 )6'( 2 end-members with M • Ni, Co, Mn(U) or Fe(U) and X• OH, F, Cl or Br. The shaded zones have already been reported previously (see Table 5 in reference (1)) and therefore were not detailed here.

Acknowledgments 1 deeply thank Dr. Pierre Alphonse for precious technical help conoerning the setup of the computerized ThennAP free calculation sheet (www.christophedrouet.com/thermAP) .

Appendix A. Supplementary data

Supplementary data associated with this article can be found, in the online version, at http:// dx.doi.org/10.1016/j.jct.2015.06.016.


J Roger

M Avenel

L Lapuyade

Characterization of SiC ceramics with complex porosity by capillary infiltration: Part A -Filling by hexadecane at 20°C

Keywords: SiC, Pore characterization, Infiltration kinetic, Infiltration mechanism, Washburn equation

The characterization of the porosity of ceramic materials is of prime importance for numerous applications. The work is being presented in two parts: A and B. In this part (Part A), several porous SiC samples exhibiting different sizes and porosity spectra were characterized by capillary infiltration of hexadecane at room temperature. The model materials contain a bimodal pore population generated by cracks and porogens. The monitoring of the weight increase in time for each sample evidence the occurrence of two successive filling stages with specific kinetics, except for the reference samples which contain only one kind of pore. The first kinetic is found clearly larger than the second. A new analytical function was used to identify the transition times and the kinetics of each stage. The various values of the kinetics revealed the occurrence of different infiltration mechanisms. The same model correlates silicon infiltration observations as well (Part B).

Introduction

Liquids rising in porous media plays a determinant role in numerous fields such as geology, concrete, oil recovery, catalysis, textile, paper industry, food processing, biotechnology, and in the production of structural materials by melt-infiltration [START_REF] Tsunazawa | An experimental study on the rate and mechanism of capillary rise in sandstone[END_REF][START_REF] Hanzic | Capillary absorption in concrete and the Lucas-Washburn equation[END_REF][START_REF] Perwuelz | Experimental Study of Capillary Flow in Yarns[END_REF][START_REF] Gruener | Capillary rise dynamics of liquid hydrocarbons in mesoporous silica as explored by gravimetry, optical and neutron imaging: Nanorheologyand determination of pore size distributions from the shape of imbibition fronts[END_REF][START_REF] Mirzajanzadeh | Water rise in a cellulose foam: By capillary or diffusional flow?[END_REF][START_REF] Forny | Wetting, disintegration and dissolution of agglomerated water soluble powders[END_REF][START_REF] Berthier | Microfluidics for biotechnology[END_REF][START_REF] Vafai | Porous Media: Applications in Biological Systems and Biotechnology[END_REF][START_REF] Aghajanian | The fabrication of metal matrix composites by a pressureless infiltration technique[END_REF]. Structural composites (CMC and MMC) are advanced materials widely used for aerospace, automotive, defense, and general engineering applications [START_REF] Rosso | Ceramic and metal matrix composites: Routes and properties[END_REF].

Pressureless infiltration of liquid metals into packed ceramic particulates or fiber preforms is obviously an efficient process for the fabrication of this kind of composite [START_REF] Ibrahim | Particulate reinforced metal matrix composites -a review[END_REF][START_REF] Rajesh | Infiltration of liquid metals in porous compacts: Modeling of permeabilities during reactive melt infiltration[END_REF][START_REF] Sree Manu | Liquid Metal Infiltration Processing of Metallic Composites: A critical review[END_REF][START_REF] Fu | Kinetics of molten metal capillary flow in non-reactive and reactive systems[END_REF]. It offers the potential to be an alternative, simple and low-cost processing route for the manufacture of ceramic reinforced metal matrix composites without remaining porosities. Recently, special interest has been focused on SiC ceramicreinforced metal-matrix composites. It was demonstrated that high SiC volume fraction composites can be produced by infiltrating liquid metals into packed SiC powder. This way was applied successfully at relative low temperatures with aluminum and copper [START_REF] Sharifi | Investigation into the kinetic behavior of molten aluminum pressureless infiltration into SiC performs[END_REF][START_REF] Zhang | SiC/Cu composite prepared by spontaneous infiltration of copper alloy into porous SiC ceramic[END_REF]. MMC can be tailored to have superior properties such as enhanced high-temperature performance, high specific strength and stiffness, increased wear resistance, better thermal and mechanical fatigue, and creep resistance than those of unreinforced alloys. This way is also applied for the fabrication of SiC/Si composites by capillary infiltration of molten silicon into porous silicon carbide [START_REF] Zollfrank | Microstructure Evolution and Reaction Mechanism of Biomorphous Si/SiC Ceramics[END_REF][START_REF] Zhu | A highstrength SiCw/SiC-Si composite derived from pyrolyzed rice husks by liquid silicon infiltration[END_REF][START_REF] Roger | SiC/Si composites elaboration by capillary infiltration of molten silicon[END_REF]. Macroporous silicon carbide is widely used in various industrial applications including filtration for gas and water, absorption, catalyst supports, concentrated solar power, thermoelectric conversion, etc [START_REF] Eom | Processing and properties of macroporous silicon carbide ceramics: A review[END_REF][START_REF] Bru | Macroporous Silicon: Technology and Applications[END_REF]. In broad terms, the mechanism of liquid transport involves capillary rise, but the precise details are lacking. In particular, it remains a lack of understanding about the filling of bimodal or even multimodal porous networks containing for example cracks and/or bubbles [START_REF] Lewandowska | Infiltration in a double-porosity medium: Experiments and comparison with a theoretical model[END_REF][START_REF] Suo | Modelling Imbibition Processes in Heterogeneous Porous Media[END_REF]. The Lucas-Washburn equation is a fundamental expression used to describe liquids rise in capillaries but also surprisingly in porous media [START_REF] Washburn | The dynamics of capillary flow[END_REF][START_REF] Abraham | Characterization of porous media by the kinetics of liquid penetration: vertical capillaries model[END_REF][START_REF] Trong | Characterization of porous materials by capillary rise method[END_REF]. In this case, this equation is based on the pore radius, the liquid viscosity, the surface tension, the contact angle, and the time. In a previous paper we confirm that the Washburn equation can be successfully applied to hexadecane filling of SiC compacts with narrow range porosities [START_REF] Marchais | Capillary infiltration of hexadecane in packed SiC powder and in SiC/SiC preforms: Pore description and calculation of molten Si infiltration[END_REF]. But these compacts are not always representative of real materials that commonly contain more complex porous networks with cracks and cavities. The aim of the present study is to investigate the infiltration kinetics of hexadecane at room temperature into SiC porous compacts with different pore distributions by the insitu measurement of the weight gain versus time. Then, the experimental measurements can be mathematically analyzed to determine the effect of larger porosities on the capillary infiltration mechanism and kinetic.

Background on the capillary infiltration into porous SiC

Several methods can be applied to measure porosity: buoyancy, gas expansion (He porosimetry), gas adsorption (BET), mercury intrusion porosimetry and saturation or imbibition [START_REF] Anovitz | Characterization and analysis of porosity and pore structures[END_REF]. All these methods are complementary as they reveal different characteristics of the porosities. Particularly, the imbibition by capillary rise is a useful and simple tool to characterize porous materials as it measures the progression of the liquid in the pore network [START_REF] Trong | Characterization of porous materials by capillary rise method[END_REF][START_REF] Karoglou | Capillary rise kinetics of some building materials[END_REF]. The Lucas-Washburn equation describes the front height evolution in porous materials during capillary infiltration [START_REF] Washburn | The dynamics of capillary flow[END_REF][START_REF] Abraham | Characterization of porous media by the kinetics of liquid penetration: vertical capillaries model[END_REF][START_REF] Trong | Characterization of porous materials by capillary rise method[END_REF] where 𝑡 is the time, 𝑟 is the radius of the medium pore, ℎ is the infiltration height, 𝜌 is the density, 𝜎 is the surface tension, and 𝜃 is the contact angle, which is assumed to be 0 due to the very low surface tension of hexadecane and 𝜇 its dynamic viscosity. This equation can be simplified in the case of small radii capillaries with creeping flow. Then,

𝑑 𝑑𝑡 . (ℎ. 𝑑ℎ 𝑑𝑡
) is close to 0. Equation 1 can be modified as follows [START_REF] Kumar | Capillary infiltration studies of liquids into 3D-stitched C-C preforms Part A: Internal pore characterization by solvent infiltration, mercury porosimetry, and permeability studies[END_REF]: where 𝑑 𝑒𝑓𝑓 is the effective diameter of the pores equal to 2. 𝑟. From the linear fit of the curves

h 2 = f(t),
it is possible to determine the effective diameter 𝑑 𝑒𝑓𝑓 of the pores. Maximum capillary height ℎ ∞ , called Jurin height, which may be attained by a liquid by capillary rise can be determined from Equation 5 [START_REF] Mirzajanzadeh | Water rise in a cellulose foam: By capillary or diffusional flow?[END_REF]:

ℎ ∞ = 4. 𝜇. 𝑐𝑜𝑠𝜃 𝜌. 𝜎. 𝑑 𝑒𝑓𝑓 Equation 5
The mass increase during the infiltration of the liquid can also be used to calculate the value of the effective diameter 𝑑 𝑒𝑓𝑓 by modifying Equation 4:

𝑚 2 𝑡 = 𝜎. 𝑐𝑜𝑠𝜃 4. 𝜇 . 𝜌 2 . 𝜀 2 . 𝑆 2 . 𝑑 𝑒𝑓𝑓 Equation 6
where 𝑚 is the infiltrated mass of the liquid, 𝜌 is the density of the liquid, 𝜀 is the porosity rate of the sample, and 𝑆 is the infiltration section of the sample. According to Equation 6, the weight increase varies linearly with time as the effective diameter 𝑑 𝑒𝑓𝑓 is constant for a narrow pore distribution as found in our previous work [START_REF] Marchais | Capillary infiltration of hexadecane in packed SiC powder and in SiC/SiC preforms: Pore description and calculation of molten Si infiltration[END_REF]. The Darcy'law describes the volumetric flow rate 𝑄 of a porous media of length 𝐿: 𝑄 = 𝑝 𝜇 .

𝑆. 𝑃

𝐿

Equation 7

where 𝑝 is the permeability, 𝑃 is the capillary pressure. If the pores are all of the same diameter 𝑑 𝑒𝑓𝑓 , the flow rate can be expressed with Equation 8 [START_REF] Benavente | Predicting water permeability in sedimentary rocks from capillary imbibition and pore structure[END_REF]:

𝑝 = 𝜀. 𝑑 𝑒𝑓𝑓 2 32
Equation 8

Then Equation 6 can be modified to:

𝑝 = 𝑚 4 𝑡 2
𝜇 2 2. 𝜎 2 . (𝑐𝑜𝑠𝜃) 2 . 𝜌 4 . 𝜀 3 . 𝑆 4

Equation 9

Few recent studies were focused on the effect of bimodal pore distribution and cracks on the capillary rise and the saturation time [START_REF] Li | Permeability of the porous Al2O3 ceramic with bimodal pore size distribution[END_REF][START_REF] Jiang | Pore network extraction for fractured porous media[END_REF][START_REF] Suo | Modelling imbibition processes in heterogeneous porous media[END_REF]. To our knowledge, no analysis method was developed to examine and understand the capillary infiltration of liquid in porous media containing a multimodal pore size distribution with or without cracks. Such a method could be useful to the production of porous or composite materials of SiC [START_REF] Esposito | Microstructure and properties of porous b-SiC templated from soft woods[END_REF].

Experimental

Starting materials

Commercial green β-SiC powder (99.8% purity, Alfa Aesar GmbH & Co KG, Germany) was used as raw material. This powder exhibits a bimodal distribution with a population centered at 0.2µm and a larger population centered between 1 and 2µm as shown in Figure 1-a,b). As reported elsewhere, this SiC powder was found to contain mainly cubic 3C-SiC polytype (β-SiC) [START_REF] Marchais | Capillary infiltration of hexadecane in packed SiC powder and in SiC/SiC preforms: Pore description and calculation of molten Si infiltration[END_REF]. The porous SiC compacts were obtained from mixtures composed of 3g of freshly prepared or 1-year-aged β-SiC powder and 1cm 3 of commercial ethanol as lubricant which were uniaxially pressed into rectangular specimens of about 50×5×5mm 3 at 100MPa using stainless steel die. The use of aged powder is known to create cracks because of the agglomeration of the grains; the height of these samples was limited to 36mm because of their mechanical weakness. The compacts from fresh, R#, and aged powders, C#, were labeled R1, R2, R3, C1, C2, C3, and C4. Porous SiC ceramics with bimodal and multimodal pore size distributions were also prepared from the mixture of fresh SiC powder with porogens: i) PMMA micro balls (Goodfellow) with a mean size of 131±46 µm (Figure 1-c)); or ii) with pyrolyzed LLDPE fibers (CANOE Technology Center, France) with a mean diameter of 32±3 µm and a length of 1.0±0.1 mm (Figure 1-d)). These pore agents and the SiC powder were slowly mixed by a planet ball miller at a low speed of 300rpm (Emax, Retsch, France) with zirconia balls as milling medium and ethanol as dispersion medium. The slurry was dried to remove the ethanol. The contents of PMMA balls (B#) or LLDPE fibers (F#) additives were 2.5%vol., 5.0%vol., 7.5%vol., and 10%vol., and the corresponding samples were labeled as B2.5, B5.0, B7.5, B10, F5.0, F7.5 and F10, the figures indicating the percentage of additive. Compacts were obtained from these powders following the method described herein. Then, the compacts were fired at 600°C for 2h in air atmosphere to decompose the PMMA balls and LLDPE fibers. All green β-SiC compacts were presintered at 1400°C for 3h under high vacuum (2×10 -4 mbar) to obtain strengthened and manipulable porous SiC ceramics. Industrial porous SiC compacts (Galtenco Solutions, France) were considered in this study, they are labeled as G1, G2, G3 and G4. Commercial foams of SiC (SiC4 grade, SICAT, France) were also used in this study, they are labeled as F1, F2, F3 and F4. These foams contain less than 500ppm of total impurities. Industrial compacts and commercial foams were cut into rectangular specimens of about 50×5×5mm 3 and 45×8×5mm 3 , respectively. All compacts were settled by polishing with a precision of at least 100 µm to homogenize their section. This stage is of importance as the capillary infiltration kinetic of a liquid depends on the cross section. Laboratory grade hexadecane provided by Sigma Aldrich Co. (USA) was used for the capillary infiltration experiments at 20°C. In Table 1 we give the contact angle on SiC, the dynamic viscosity, the density and the surface tension of hexadecane [START_REF] Onken | The Dortmund Data Bank: a computerized system for retrieval, correlation, and prediction of thermodynamic properties of mixtures[END_REF].

Characterization of the materials

Pore-size distributions of each type of sample were obtained from mercury intrusion porosimetry (Autopore IV, Micromeritics Instrument Corp., USA) for a sample volume between 0.6 and 1.2 cm 3 .

The porosity ε of the samples was calculated by applying the formula 𝜀 = 1 -𝑑 𝑎 𝑑 𝑡 with the apparent density 𝑑 𝑎 obtained from mercury intrusion porosimetry and the true density 𝑑 𝑡 measured from helium pycnometry (AccuPyc 1340, Micromeritics Instrument Corp., USA). The pore size distributions of the SiC consolidated compacts are shown in Figure 2. The pore size distribution of the porous SiC compacts presents a single peak with a narrow width, signifying a uniform pore size distribution centered on 0.22µm. The dimensions and the porosity percentage of all compacts are given in Table 2. The capillary infiltration experiments were performed on a commercial modified tensiometer (TensioCad W TE-5200 MO, CAD Instrument, France) operating at 20°C in a temperature-regulated room. Figure 3-a) shows the experimental setup of the infiltration experiments and the picture of a compact. The SiC compacts were hanged along their longest dimension to a modified analytical balance, which was connected to a computer. A crystallizer containing the solvent (hexadecane) was automatically raised until the contact with the solvent. The mass variations were recorded every second. Losses due to solvent evaporation were estimated to be less than 1% of the imbibed solvent and were hence considered as negligible.

Indeed, the saturatation vapor pressure of hexadecane at 20°C is very low since equal to 4.8810 -2 Pa [START_REF] Forziati | Physical Properties of n-Hexadecane, n-Decylcyclopentane, n-Decylcyclohexane, 1-Hexadecene and n-Decylbenzene[END_REF], so its evaporation can be neglected. Hexadecane was selected because of its low volatility. The time dependence of hexadecane infiltration in the SiC compacts was recorded by monitoring the weight gain increase. A typical curve of a compact is plotted in Figure 3-b). One can note that for the shortest durations, hexadecane rises very quickly by capillary action. Then, the infiltration kinetic is lower for longer durations and the weight of the meniscus can be deduced after the breakage of the contact with the liquid. The computed tomography X-ray scans of the compacts were performed using a 140kV Xray source with an intensity of 60μA (GE v|tome|x s research edition, General Electrics, USA). A 500ms exposure time and a spatial resolution of 8μm per voxel unit: this is the highest achievable resolution as dictated by the dimension of the specimens and the performance of the CT machine. These measurements give views of the size and distribution of the largest porosities (cracks, cavities), as shown in Figures 4 and5. The reconstructed images were obtained with Image Analysis System (ImageJ) [START_REF] Schneider | NIH Image to ImageJ: 25 years of image analysis[END_REF]. It can be seen in these figures that, compared to the reference compact, more or less large cracks and/or holes are present. It must be noted that the large pores do not connect directly to each other.

Results and discussion

Infiltration of hexadecane

It is assumed that the infiltration of the liquid through the porous SiC compacts is isothermal. No interfacial chemical reaction occurs in the process and the SiC compacts are considered as rigid bodies.

The curves plotted in Figures 6-a) show the time dependence of the hexadecane filling percentage of compacts representative of the six types of SiC compacts studied. Theses curves indicate that after 14400s, i.e. 4h, the complete filling is only obtained for Reference SiC (1), called R1, which is free of macro porosities. The squared weight gains corrected from the meniscus weight and normalized to the infiltration section and to the pore fraction, i.e. ( 𝑚-𝑚 𝑚𝑒𝑛𝑖𝑠𝑐𝑢𝑠 𝜀.𝑆 ) 2 = 𝑓(𝑡), are plotted in Figure 6-b). It can be seen that all the different types of samples exhibit similar behaviors with a first linear variation follows by a second linear variation but with a noticeable lower slope. For R1, the second slope is found null. It must be noted that the first slopes of the compacts made from aged SiC powder or with balls or with fibers are higher than the slope of the reference SiC compact. At the opposite, the first linear slopes of the foam and SiC from Galtenco are lower. The second slopes of all samples seem also practically linear, this second slope for SiC Galtenco appears later and is not visible in Figure 6-b). The transition between the two filling rates occurs at various durations without apparent relationship with the nature of the porosity and their quantity. On a general manner, after 4h, the compacts with cracks and voids are far from the complete filling. The first linear variations are in good agreement with the Washburn equation, but this law cannot describe the second linear variations. These behaviors are confirmed on the graphs in Figure 7 which show the ( 𝑚-𝑚 𝑚𝑒𝑛𝑖𝑠𝑐𝑢𝑠 𝜀.𝑆 ) 2 = 𝑓(𝑡) curves for all the samples of this study compared to R1 compact. Three reference SiC compacts are considered with the higher (R2) and lower (R3) slopes obtained with these compacts, R1 has an intermediate slope (Figure 7-a)). These three compacts have similar behaviors with a relative rapid saturation, the main difference lies in the values of their slopes which are more or less elevated. It must be noted that all these measurements are totally reproducible.

To our opinion, this variation is induced by the microstructure of the area in contact with the liquid. The shapes of the curves of the four compacts made from aged SiC powder, called C1, C2, C3 and C4, are relatively similar with two successive linear variations (Figure 7-b)). The main differences lie: i) in the values of their slopes, which are close (C1, C2) or higher (C3, C4) than the references, and ii) in the transition time between the two-filling kinetics. The curves of the four compacts obtained from PMMA balls addition are shown in Figure 7-c). They exhibit also two successive rates of filling and their first respective slopes are also close (B10) or higher (B2.5, B5.0, B7.5) than the reference compacts, and there is no evidence of any influence of the added balls content. Similar observations can be deduced from Figure 7-d) showing the curves of the compacts elaborated with fibers. The slopes are clearly larger than the references and no effect of the fraction of fibers added is evidenced. The curve shapes of the SiC foams (F1 to F4) are close during the two successive filling kinetics (Figure 7-e)). The first slope is within the range of the reference compacts but the second one is very low. The squared normalized weight gain curves of the SiC from Galtenco exhibit also two filling kinetics, the first being clearly slower than the reference compacts (Figure 7-f)). Two kinds of behaviors are identified for these compacts with a faster initial kinetic and a sooner transition time leading to a more important filling rate for G1 and G2 samples than for the G3 and G4 samples. On a general manner, these experimental results show variations of the filling kinetics induced by different filling mechanisms. The reference compacts (R1 to R3) are the only ones that strictly verify the Washburn equation as they exhibit a unique filling slope until saturation. The filling mechanisms of all the others compacts involve apparently two consecutive stages. The first kinetics are relatively rapid and close, higher or lower than the filling kinetic of the references. The second kinetics are blatantly sluggish with kinetics around ten times lower.

All these different observations confirm the determinant role plays by cracks and large porosities on the filling kinetic by capillary infiltration. But such complex behaviors cannot be described by the Washburn equation. That is why, an analytical analysis is proposed below.

Analytical analysis of hexadecane infiltration

As mentioned above, many interesting studies were realized about the capillary infiltration of liquid in capillaries and some works deal with the infiltration in porous materials. This last case is more complex as the porosities may be irregular and some large defects can easily occur which may have a significant effect on the kinetic of saturation. Thus, it was previously identified in porous materials that a difference of time occurs between the filling of quickly penetrated capillaries and the complete filling of slow penetrating capillaries which implies the presence of at least two characteristic pore sizes [START_REF] Trong | Characterization of porous materials by capillary rise method[END_REF]. Kumar et al. also consider the effect of two types of pore with different sizes on the infiltration of 3D-stitched C-C preforms [START_REF] Kumar | Capillary infiltration studies of liquids into 3D-stitched C-C preforms Part A: Internal pore characterization by solvent infiltration, mercury porosimetry, and permeability studies[END_REF]. In those two cases, a modification of the Washburn equation is made to explain the measurements with two pore sizes based on one infiltration kinetic. These considerations are not directly applicable to the results reported herein as two filling kinetic are clearly identified. This can be justified by the elevated difference in size between the small pores and the large defects in our model materials.

Others classical analytical theories can be considered depending on the limiting phenomenon.

Particularly, the asymptotic solution describing the saturation of porous materials can be expressed by a relaxation law which consists in an exponential variation of the front height in time [START_REF] Zhmud | Dynamics of capillary rise[END_REF]. Another interesting analytic solution to describe the saturation effect implies the Lambert equation [START_REF] Yang | Capillary Microflows in Porous Metals, 21 st Australasian Fluid Mechanics Conference[END_REF]. In our opinion, no one of these theories can accurately model the two-staged weight increases reported in the present work. Consequently, a new model is proposed here. As two successive linear variations of the squared weight gain with time are observed, it can be considered that two kinetics, called 𝑘 1 and 𝑘 2 , occur successively in time. For a given solvent i, the kinetic 𝑘 1 𝑖 relates to the fast filling of the smaller pores, and the kinetic 𝑘 2 𝑖 describes the slow filling of the larger and isolated pores. In a first step, the squared weight gain, ∆𝑚 2 , increases linearly in time according to the 𝑘 1 𝑖 kinetic (g 2 .mm -4 .s -1 ). Then, in a second step, it is supposed to increase slowly according to the 𝑘 2 𝑖 kinetic (g 2 .mm -4 .s -1 ). The overall squared-weight gain ∆𝑚 2 is normalized to the squared porous section (𝑆. 𝜀) and its evolution is calculated with Equation 6: Reference SiC (1). This confirms that the filling of the reference occurs only with one mechanism, contrary to all the others samples. The experimental and fitted curves giving the mass increase versus time for these six compacts are shown in Figures b,c,d,e,f). These figures show also the fitted variations of the mixing functions 𝑓(𝑡) and 1 -𝑓(𝑡) for the six samples. As expected from the elevated correlation percentages, all the fitted weight gain curves are very close to the experimental measurements. The fit of the reference compact (R1) describes pretty well the saturation of this sample which occurs after a duration of about 14100 seconds, accordingly with visual observations. This value is very close to the crossing time of the 𝑓(𝑡) and 1 -𝑓(𝑡) functions at 13500 seconds, i.e. 𝑓(𝑡) = 1 -𝑓(𝑡) = 0.5. The S-like shape of the 𝑓(𝑡) is only due to the sudden variation of the experimental slope which cannot easily be described mathematically, even with the function considered here. The fits could be improved by adding a weighting scheme to increase the weight of the points close to the transition.

∆𝑚 2 𝑆 2 . 𝜀 2 =
However, it seems that the fitting of the curves provides not only significant values of 𝑘 1 𝐻𝑒𝑥 and 𝑘 2 𝐻𝑒𝑥 kinetics, but also the duration for saturation of the smallest porosities. Indeed, for all the others samples, it is found that the crossing times of the 𝑓(𝑡) and 1 -𝑓(𝑡) functions are also in good agreement with the experimental transition times, as shown in Figures c,d,e,f). These values are not directly comparable as they depend on the height of each sample. The 𝑘 1 𝐻𝑒𝑥 and 𝑘 2 𝐻𝑒𝑥 kinetics are compared in Figure 9. It can be seen that the 𝑘 1 𝐻𝑒𝑥 values are relatively different. On a general manner, the 𝑘 1 𝐻𝑒𝑥 values of C1 and B10 compacts are close and slightly higher than for the references, as already seen in Figure 7. Much higher values are obtained with the compacts prepared by fibers addition, the calculated kinetic is almost twice as high as the kinetic of the reference R1. The 𝑘 1 𝐻𝑒𝑥 kinetics of the Galtenco samples are almost half of the reference R1, and the kinetics of the foams are close to the kinetic of the reference R1. The 𝑘 2 𝐻𝑒𝑥 kinetics are clearly lower than the 𝑘 1 𝐻𝑒𝑥 kinetics and are less than one-tenth of the corresponding 𝑘 1 𝐻𝑒𝑥 value. This second kinetic is clearly correlated to the presence of larger porosities (cracks, holes). A scrupulous examination of the fittings confirms their relevance. Indeed, the calculated values of the weight gains at the transition durations with the 𝑘 1 𝐻𝑒𝑥 kinetics differ from the experimental data by a maximum value of 1.3% (R1: 1.3%; C1: 1.0%; B10: 1.2%; F7.5: 1.0%; G1: 1.0% and F1: 1.0%). The deviations of the weight gains at the last measured points are equal to or less than 1.3% (R1: 1.3%; C1: 1.0%; B10: 0.2%; F7.5: 0.4%; G1: 0.9% and F1: 0.2%). It can be then considered that the fitted values as significant and can be relevant for a better understanding of the filling mechanisms. It is then clear that the largest porosities can have a positive and a negative effect on the filling kinetic and, consequently, on the filling time. In accordance with the porosimetry measurements (Figure 2-a)), the presence of cracks or holes in the SiC compact slightly increases the mean size of the smaller pores, so that the filling kinetics are increased as experimentally measured. At the opposite, a lower mean size of the smaller pores induces a lower kinetic of filling as measured for Galtenco compacts and SiC foams (Figure 2-b)). The very low values of the 𝑘 2 𝐻𝑒𝑥 kinetics are dramatically detrimental to the completeness of the filling. There are clearly correlated to the filling of the largest pores. But it depends also on the interconnectivity of these large pores as evidenced by the SiC+balls compacts. Indeed, the holes generated by the balls are large and poorly interconnected, consequently, the direct and relative rapid filling from the liquid source is almost impossible. The Galtenco compacts have a particular behavior induced by their heterogeneity: a part of the compacts contains no large pores and another part contains the large pores. This can probably justify the particular and unique shape of the 𝑓(𝑡) curve in Figure 8-e) with a very slow decrease to 0. It occurs probably a rapid filling of the small pores in the area containing the large pores and, in parallel, a slower progression in the area without these large pores. A slow decrease of the 𝑓(𝑡) curve below 0.5 could reveal the highly heterogeneous character of a material evidenced by the rounded shape of the squared weight gain in time. The low 𝑘 1 𝐻𝑒𝑥 kinetics obtained with the Galtenco compacts are correlated to the low size of their smallest pores, which confirms that the smallest pores are filled first. From the fitting equations, one can deduce also that: i) the time needed for the saturation of each samples, as shown in Figure 10-a), and ii) the percentage of filling associated to each kinetic, as shown in Figure 10-b). The comparison of these figures reveals that the duration for the complete filling increase in the presence of large porosities, as the reference samples need the shortest durations for full filling. The presence of large and practically isolated porosities induced a large increase of the saturation time as highlighted by the compacts prepared with balls. It should be noted that the isolated large pores cannot be completely filled because of the presence of air trapped in these pores. The estimation of the filling volumes during the two kinetic stages can be compared with mercury porosimetry analyses to determine the limiting size of pores filled during the first stage assuming that the smaller porosities are filled first. In this way, the transition occurs for pores close to 1µm for C1 and B10 compacts. The upper limit size of filled pores during the first stage is found close to 10µm for F7.5 and 3µm for F1 samples. A limit value of 0.2µm is identified for G1 sample. The findings reported herein could be highly useful for a better understanding of the liquid filling of porous materials with bimodal or multimodal pores distribution. From these results, a partial permeability value can be obtained for each stage. Indeed, for each kinetic a value of the filled porosity can be estimated and these values can be used in Equation 9to calculate the two partial permeabilities 𝑝 1 𝐻𝑒𝑥 and 𝑝 2 𝐻𝑒𝑥 . These values are reported in Table 3. The 𝑝 1 𝐻𝑒𝑥 values of R1, C1 and B10 samples are close, which confirms the filling of the smaller pores. While the values of F7.5 and F1 are respectively 4 and 16 times greater than the reference R1. This difference can be correlated to the complex pore network of F7.5 and F1, which promotes the flow of the liquid. The 𝑝 2 𝐻𝑒𝑥 values are clearly lower than the corresponding 𝑝 1 𝐻𝑒𝑥 values which demonstrates that the filling of the larger pores is slower. It can be considered that the ratios 𝑝 1 𝐻𝑒𝑥 /𝑝 2 𝐻𝑒𝑥 are all the more elevated than the larger pores are large. Further works are expected to improve the analysis of the capillary rise measurements and to establish precisely the filling mechanisms and their limitations. Others solvents can be used in this aim on the condition that their saturation vapor pressure remains sufficiently low to limit the losses by evaporation. A solution to work around this limitation consists in operating in a closed and temperate setup under saturatation pressure of the solvent. It was also recently demonstrated that the penetration velocity is sensitive to the geometry of the porous media, and can be tuned by varying the evaporation rate for a given geometry [START_REF] Liu | Tuning capillary penetration in porous media: Combining geometrical and evaporation effects[END_REF]. The inhibition of liquid penetration by evaporation can be offset by varying the geometry of the porous layer. Consequently, the application of the capillary infiltration to the characterization of a porous object with complex shape and relatively large surface is possible on the condition of operating with a solvent having a suitable saturation vapor pressure to limit the losses. This parameter must be considered carefully before any test under penalty of erroneous measurements.

Conclusion

The infiltration of several different porous SiC materials by hexadecane have been investigated. These

Figure 1 :Figure 3 :Figure 4 :Figure 5 :Figure 6 :Figure 7 :Figure 8 :

 1345678 Figure 1: Starting materials: a) grain size distribution of the SiC powder, b) secondary electrons image (magnification x30000) of the SiC powder; c) optical microscopy (magnification x5) of the PMMA balls and d) optical microscopy (magnification x5) of the LLDPE fibers.

Figure 9 :

 9 Figure 9: Graphic comparison of the 𝑘 1 𝐻𝑒𝑥 and 𝑘 2 𝐻𝑒𝑥 kinetics for the six main samples.

  (Equation1):

	𝑑 𝑑𝑡	(𝜋. 𝑟 2 ℎ. 𝜌.	𝑑ℎ 𝑑𝑡	) = (	2. 𝜎. 𝑐𝑜𝑠𝜃 𝑟	-𝜌. 𝑔. ℎ -	8. 𝜇. ℎ 𝑟 2 .	𝑑ℎ 𝑑𝑡	) . 𝜋. 𝑟 2
									Equation 1

  The expression of 𝑓(𝑡 𝑖 ) function is given in Equation7which varies between 0 and 1. It describes the more or less rapid transition between the two filling processes:

	𝑓(𝑡 𝑖 ) =	1 1 + (𝛼 + 𝛽. 𝑡 𝑖 + 𝛾. 𝑡 𝑖	3 ). 𝑒𝑥𝑝 ( 𝑡 𝑖 𝛿	)
				Equation 7
	where 𝛼, 𝛽, 𝛾 and 𝛿 are adjustment coefficients. The validity of Equation 6 was tested by fitting the
	experimental curves with this equation by adjusting the values of 𝑘 1 𝑖 , 𝑘 2 𝑖 , 𝛼, 𝛽, 𝛾 and 𝛿 with the help of
	a numerical solver. The obtained values of 𝑘 1 𝐻𝑒𝑥 and 𝑘 1 𝐻𝑒𝑥 kinetics with hexadecane (𝐻𝑒𝑥) are given is
	Table 2 for the six main compacts (R1, C1, B10, F7.5, G1 and F1). The correlation percentages are quite
	acceptable since they are greater than 99.9%. All the samples have a non-null 𝑘 2 𝐻𝑒𝑥 kinetic except
	𝑗
	𝑖=1 𝑡 0 =0
				Equation 6

∑ 𝑘 1 𝑖 . (𝑡 𝑖 -𝑡 𝑖-1 ). 𝑓(𝑡 𝑖 ) + 𝑘 2 𝑖 . (𝑡 𝑖 -𝑡 𝑖-1 ). (1 -𝑓(𝑡 𝑖 ))

Table 1 :

 1 Physical properties of hexadecane at 20°C[START_REF] Onken | The Dortmund Data Bank: a computerized system for retrieval, correlation, and prediction of thermodynamic properties of mixtures[END_REF]

	Physical properties	Hexadecane (20°C)
	Contact angle on β-SiC (°)	0
	Dynamic viscosity μ (g.cm -1 .s -1 )	3.04×10 -2
	Density ρ (g.cm -3 )	0.77
	Surface tension σ (g.s -2 )	27.42

Table 2 :

 2 Porosity, dimensions and filling kinetics of the six main samples

Table 3 :

 3 Experimental partial permeabilities 𝒑 𝟏 𝑯𝒆𝒙 and 𝒑 𝟐 𝑯𝒆𝒙 , and their ratio for the six main samples

	Samples	Reference	Dimensions (mm)	Porosity ε (%) (±1)	𝒌 𝟏 𝑯𝒆𝒙 (g 2 .mm -4 .s -1 ) 𝒌 𝟐 𝑯𝒆𝒙	Correlation percentages (%)
	Reference SiC (1)	R1	5065	47	8.77410 -8	-	99.973
	Cracked SiC (1)	C1	361110	53	1.00910 -7	4.38810 -9	99.968
	SiC + 10%vol. balls	B10	4775	59	1.03710 -7	2.07010 -9	99.941
	SiC + 7.5%vol fibers	F7.5	4655	57	1.82510 -7	1.11010 -8	99.914
	SiC Galtenco (1)	G1	5055	45	5.83110 -8	3.61410 -9	99.974
	SiC foam (1)	F1	4585	73	9.70010 -8	8.34610 -9	99.986
	Samples		Reference	𝒑 𝟏 𝑯𝒆𝒙	𝒑 𝟐 𝑯𝒆𝒙	𝒑 𝟏 𝑯𝒆𝒙 /𝒑 𝟐 𝑯𝒆𝒙
					(m 2 )
	Reference SiC (1)	R1	6.33110 -15	-	-
	Cracked SiC (1)	C1	6.86010 -15	4.87810 -18	1406
	SiC + 10%vol. balls	B10	7.78910 -15	5.68710 -19	13696
	SiC + 7.5%vol fibers	F7.5	2.82010 -14	1.85810 -19	1518
	SiC Galtenco (1)	G1	1.87510 -15	3.07910 -18	609
	SiC foam (1)	F1	1.03410 -14	1.23910 -17	835

Acknowledgement

The authors wish to thank Olivier Caty and Muriel Alrivie from the Laboratory of ThermoStructural Composites (LCTS -UMR 5801 -UB-CNRS-CEA-SAFRAN) for their kind assistance. The Galtenco company is thanked for the supply of custom samples. A part of this work has been performed within the frame of the training of Kylia Marcus-Cesaire from the University of Bordeaux, France.