S Féréol

Caluch A Fodil

Pelle G &

Daniel Isabey
email: daniel.isabey@inserm.fr

Force distribution on multiple bonds controls the kinetics of adhesion in stretched cells

Keywords: parallel bonds, zipper bonds, lifetime, dissociation rate, bond strength. 2

We show herein how mechanical forces at macro or micro scales may affect the biological response at the nanoscale. The reason resides in the intimate link between chemistry and mechanics at the molecular level. These interactions occur under dynamic conditions such as the shear stress induced by flowing blood or the intracellular tension. Thus, resisting removal by mechanical forces, e.g., shear stresses, is a general property of cells provided by cellular adhesion. Using classical models issued from theoretical physics, we review the force regulation phenomena of the single bond. However, to understand the force regulation of cellular adhesion sites, we need to consider the collective behavior of receptor-ligand bonds.

We discuss the applicability of single bond theories to describe collective bond behavior.

Depending on bond configuration, e.g., presently "parallel" and "zipper", the number of bonds and dissociation forces variably affect the kinetics of multiple bonds. We reveal a marked efficiency of the collective organization to stabilize multiple bonds by sharply increasing bond lifetime compared to single bond. These theoretical predictions are then compared to experimental results of the literature concerning the kinetic parameters of bonds measured by atomic force microscopy and by shear flow. These comparisons reveal that the force-control of bonds strongly depends on whether the force distribution on multiple bonds is homogeneous, e.g., in AFM experiments, or heterogeneous, e.g., in shear flow experiments. This reinforces the need of calculating the stress/strain fields exerted on living tissues or cells at various scales and certainly down to the molecular scale.

Introduction

It is widely accepted that mechanical forces control normal biological processes and their pathological alterations through cell-environment interactions [START_REF] Discher | Tissue cells feel and respond to the stiffness of their substrate[END_REF].

These interactions occur under dynamic conditions, e.g., the flowing blood; the internal tension, through adhesion sites which are protein assemblies through clusters which link the cytoskeleton to the extracellular environment [START_REF] Geiger | Assembly and mechanosensory function of focal contacts[END_REF]. In particular, receptor-ligand bond anchoring cells in the presence of flowing liquids are stretched by tensile mechanical forces that balance the drag force on the bound cell. Therefore, stabilization under force and resisting removal to shear stresses -and to intra/extracellular stretching forces -are fundamental properties of cellular adhesion. Moreover, weakening or strengthening adhesion by mechanical force remains an opened question if one wants to understand flow-controlled key cellular interactions such as: macrophage-epithelium, leukocyte-endothelium or bacteria-host cell [START_REF] Thomas | Catch bonds in adhesion[END_REF][START_REF] Thomas | Biophysics of catch bonds[END_REF].

A fundamental aspect to consider is the stochastic nature of molecular adhesion and the subsequent bond detachment under force which are extensively described in the literature [START_REF] Bongrand | Reproducible dissociation of cellular aggregates with a wide range of calibrated shear forces: application to cytolytic lymphocyte target cell conjugates[END_REF][START_REF] Evans | Energy landscapes of biomolecular adhesion and receptor anchoring at interfaces explored with dynamic force spectroscopy[END_REF]. A less considered aspect deals with the complexity of adhesion structures since cellular adhesion involves not only one but several bonds, which stabilize and resist mechanical forces through different strategies of association and reinforcement [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF][START_REF] Evans | Probing molecular attachments to cell surface receptors: image of stochastic bonding and rupture processes[END_REF][START_REF] Leckband | Long-range attraction and molecular rearrangements in receptor-ligand interactions[END_REF].

These strategies are most likely influenced by force transiting through intra/extracellular environment and a number of biochemical and biomechanical processes implicating signaling pathways [START_REF] Vogel | Cell fate regulation by coupling mechanical cycles to biochemical signaling pathways[END_REF]. From initial adhesion to focal complex and beyond the mature focal adhesion, the adhesion site configuration changes and the maximal force supported by adhesion sites varies by six orders of magnitude (pN-µN) [START_REF] Bruinsma | Non-specific binding by macrophages: evaluation of the influence of medium-range electrostatic repulsion and short-range hydrophobic interaction[END_REF]. Although these adhesion sites are made of similar weak bond units -i.e., a non covalent molecular link given to support forces up to the pN -their collective organization allows them supporting a huge range of stress [START_REF] Geiger | Assembly and mechanosensory function of focal contacts[END_REF]. The structural changes associated to stabilization and strengthening include: (i) increase in cell-substrate contact area (spreading) (Capo et al., 1981), (ii) recruitment of receptors to anchoring sites including their lateral association (clustering) [START_REF] Andre | Use of conjugates made between a cytolytic T cell clone and target cells to study the redistribution of membrane molecules in cell contact areas[END_REF][START_REF] Bell | Cell adhesion. Competition between nonspecific repulsion and specific bonding[END_REF][START_REF] Cluzel | The mechanisms and dynamics of (alpha)v(beta)3 integrin clustering in living cells[END_REF], (iii) interaction with cytoskeleton elements that lead to enhanced force distribution among bound receptors via local cytoskeleton stiffening (focal adhesion assembly) [START_REF] Pasternak | Lymphocyte mechanical response triggered by cross-linking surface receptors[END_REF]. These mechanisms of adhesion strengthening are supported by numerous observations from various cellular systems [START_REF] Balaban | Force and focal adhesion assembly: a close relationship studied using elastic micropatterned substrates[END_REF][START_REF] Choquet | Extracellular matrix rigidity causes strengthening of integrin-cytoskeleton linkages[END_REF][START_REF] Galbraith | The relationship between force and focal complex development[END_REF][START_REF] Giannone | Talin1 is critical for force-dependent reinforcement of initial integrin-cytoskeleton bonds but not tyrosine kinase activation[END_REF][START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF]. Trying to understand these mechanisms of adhesion stabilization or reinforcement, and how they are regulated by force, remains challenging questions which require new theoretical knowledge and pertinent analysis of data.

Recognizing the need to expose actions at the submicroscopic level, many research groups have employed ultrasensitive force techniques to probe extremely small regions of adhesive contact between surfaces functionalized by biological molecules [START_REF] Evans | Probing molecular attachments to cell surface receptors: image of stochastic bonding and rupture processes[END_REF][START_REF] Leckband | Long-range attraction and molecular rearrangements in receptor-ligand interactions[END_REF]. Most direct measurements of single bond strength have been performed with three types of ultrasensitive probes: the Atomic Force Microscope (AFM) where force is sensed by deflection of a thin silicon nitride cantilever [START_REF] Binnig | Atomic force microscope[END_REF], the Biomembrane Force Probe (BFP) where force is sensed by a glass microsphere glued to the pole of a micropipette-pressurized membrane capsule [START_REF] Evans | Sensitive force technique to probe molecular adhesion and structural linkages at biological interfaces[END_REF][START_REF] Simson | Micropipet-based pico force transducer: in depth analysis and experimental verification[END_REF]; and the laser optical tweezers where force is sensed by displacement of a microsphere trapped in a narrowly focused beam of laser light [START_REF] Ashkin | Forces of a single-beam gradient laser trap on a dielectric sphere in the ray optics regime[END_REF][START_REF] Ashkin | Force generation of organelle transport measured in vivo by an infrared laser trap[END_REF]. Common limitations associated with these techniques deal with thermal fluctuations, hydrodynamic interactions, and noteworthy the difficulty of controlling the number and the configuration of bonds implicated.

Since most of adhesive interactions are initiated and then sustained under flow conditions (in priority blood but no only), cell resistance to removing by fluid shear stress using a "laminar flow chamber" has become a relevant method to measure the kinetics of receptor-ligand adhesion [START_REF] Bongrand | Reproducible dissociation of cellular aggregates with a wide range of calibrated shear forces: application to cytolytic lymphocyte target cell conjugates[END_REF]Pierres et al., 1996a;Pierres et al., 1996b). Basically, particles or circulating cells covered by receptors are flowing near a surface coated by a ligand and arrested by the formation of limited number of adhesion bonds, e.g., one per cell in many experiments, as far as hydrodynamic forces do not exceed a few pNs [START_REF] Kaplanski | Granulocyte-endothelium initial adhesion. Analysis of transient binding events mediated by E-selectin in a laminar shear flow[END_REF]Pierres et al., 1996a;[START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF]. The relationship between the number of arrested particles and the duration of their arrest allows (i) determining detachment curves whose initial slope provides the dissociation rate [START_REF] Chang | The forward rate of binding of surface-tethered reactants: effect of relative motion between two surfaces[END_REF][START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF][START_REF] Tissot | Motion of cells sedimenting on a solid surface in a laminar shear flow[END_REF] and moreover the (ii) study the mechanisms of flow-enhanced cell adhesion [START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF].

The present paper concerns some under considered aspects of the regulation by mechanical force of collective bond organization implicated in the structuration of different adhesion sites. We purposely apply theories of bond dissociation -issued from stochastic approaches -to different types of bond configurations, e.g., parallel or zipper bonds, on which dissociation forces are exerted. The experimental data purposely used are issued from two different methods: the atomic force microscopy (AFM) and the laminar (viscous) shear flow. It appears that force-control of multiple bonds strongly depends on whether the force distribution on multiple bonds is homogeneous, such as during certain AFM experiments, or heterogeneous, such as during shear flow experiments. Thus, to understand the cell adhesion response to stretch, the knowledge of force distribution is needed down to the smallest scale possible. This reinforces the need of calculating the stress/strain fields exerted on living tissues or cells at various scales and certainly down to the molecular scale.

Methods

Kinetic theory of bond dissociation

Single bond

Weak noncovalent interactions between large molecules mediate many of life's functions in cells [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF]. The strength of interaction is the level of force to disrupt a bond on a particular time scale, knowing that these molecular interactions have limited lifetimes and thus fail under any level of force if pulled on for modest period of time [START_REF] Evans | Energy landscapes of biomolecular adhesion and receptor anchoring at interfaces explored with dynamic force spectroscopy[END_REF]. For instance, an isolated bond has no strength on times scales longer than its natural lifetime for spontaneous dissociation (entropy-driven): 𝑡 𝑜𝑓𝑓 0 . On the other hand, when a bond is pulled apart faster than needed for diffusive relaxation, it is classically thought to resist to detachment by a force given by the maximum gradient, 𝑓 𝑏𝑜𝑛𝑑 ~ 𝜕𝐸(𝑥) 𝜕𝑥 𝑚𝑎𝑥 , where 𝐸(𝑥) is the interaction potential or energy contour defined along the distance 𝑥 of separation. Thus, from the slow limit set by spontaneous transition, 𝑡 𝑜𝑓𝑓 0 , to the ultrafast limit for diffusive relaxation, 𝑡 𝐷 , ( 10  -10  s, i.e., the Brownian time), strength of non covalent bonds is governed by thermally activated kinetics under external force and depends on how force is applied over time. To describe bond kinetics, the simplest model of chemical binding is a single-sharp barrier model which takes into account the major impact of force on thermally activated kinetics that is an exponential amplification of the forward rate for dissociation characterized by a small force scale:

𝑘 𝑜𝑓𝑓 (𝑓)~ 𝑒 -𝐸 𝑏 (𝑓) 𝑘 𝐵 𝑇
(1)

𝐸 𝑏 𝑓 = 𝐸 𝑏 0 -𝑓 × 𝑥 𝛽 (2)
In this model, the energy barrier of a single bond is characterized in its simplest form by its height at zero force, 𝐸 𝑏 0 , and the distance in nanometer, or a fraction of it, at which this energy barrier is located, 𝑥 𝛽 . The impact of a force 𝑓 depends on how energy landscape is deformed by the force (given by Eq. 2). This model (Eq. 1 and 2) predicts that the rate of failure increases exponentially once force reaches a small scale 𝑘 𝐵 𝑇 𝑥 𝛽 , well below the adiabatic limit -𝐸 𝑏 (𝑓) 𝑘 𝐵 𝑇 . Thermal activation sets the scale for force:

𝑓 𝛽 = 𝑘 𝐵 𝑇 𝑥 𝛽 (3)
𝑓 𝛽 provides a reference force scale that remains rather small since the thermal energy per molecule: 𝑘 𝐵 𝑇pNnm, (with 𝑘 𝐵 = 1.381 × 10 -23 𝐽/𝐾, Boltzman constant and T absolute room temperature in °K), and 𝑥 𝛽 0.1-1 nm, (e.g., 𝑓 𝛽 ≈ 13 𝑝𝑁  for integrin-VCAM linkage [START_REF] Evans | Using force to probe single-molecule receptor-cytoskeletal anchoring beneath the surface of a living cell[END_REF]; 𝑓 𝛽 ≈ 3 -9 𝑝𝑁 for adhesion between platelets and subendothelial matrix [START_REF] Arya | Dynamic force spectroscopy of glycoprotein Ib-IX and von Willebrand factor[END_REF]). Using this scale, the rate of dissociation increases exponentially with force as proposed in the pioneering work by Bell [START_REF] Bell | Models for the specific adhesion of cells to cells[END_REF]:

𝑘 𝑜𝑓𝑓 𝑓 = 𝑘 𝑜𝑓𝑓 0 × 𝑒 𝑓 𝑓 𝛽 (4)
Eq. 4 provides the dissociation rate of a single bond on which a force 𝑓 is exerted. It describes a stochastic phenomena analogue to the disintegration of radioactive substance.

The dissociation rate 𝑘 𝑜𝑓𝑓 𝑓 is linked to the probability 𝑃 𝑡 that the bond survives by a simple differential equation which provides the expression of 𝑃 𝑡 :

𝑑𝑃 (𝑡) 𝑃(𝑡) = -𝑘 𝑜𝑓𝑓 𝑓(𝑡) × 𝑑𝑡 or equivalently 𝑃 𝑡 = 𝑒𝑥𝑝 -𝑘 𝑜𝑓𝑓 𝑓 𝑡 ′ × 𝑑𝑡 ′ 𝑡 0
(5)

The prefactor of the dissociation rate given by the Bell model (Eq. 4) corresponds to the zero-force dissociation rate which classically depends on initial barrier height 𝐸 𝑏 0 . (i) The response to mechanical force of an assembly of bonds obeys to the principle issued from the single bond theory which states that the response to force depends on the deformation of energy landscape which thereby resembles a mountainous terrain made of multiple energy barriers [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF][START_REF] Tsukasaki | Role of multiple bonds between the single cell adhesion molecules, nectin and cadherin, revealed by high sensitive force measurements[END_REF], (ii) Associated bonds behave independently (non cooperative bonds). This supposes sufficient linkage for assembly but insufficient linkage to force the bonds to behave identically to reach simultaneous dissociation (cooperative bonds). Indeed, preserving the stochastic nature of a given bond assembly is probably the most pertinent assumption for the weak non covalent interactions that characterize interfacial biological bonds, (iii) different bond configurations (e.g., parallel or zipper) result in different force regulations (Fig 1). In parallel bonds, the load is redistributed among the surviving bonds [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF].

𝑘 𝑜𝑓𝑓 0 = 1 𝑡 𝐷 𝑒 -𝐸 𝑏 0 𝑘 𝐵 𝑇 (
In "zipper" bond, only the leading bond experiences the entire loading force [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF].

Whatever the predominant bond structure, multiple bonds can break at random times from first to last. Thus the global dissociation rate is a function of the instantaneous force and the distribution of rupture times can merely be described by a first order Markov process with time-dependent rate constants [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF][START_REF] Williams | Analytical descriptions of dynamic force spectroscopy: behaviour of multiple connections[END_REF]. Hence, far from thermal equilibrium, it is still possible to generalize the single-dynamic barrier model (Eq. 1-4) to more complex barriers. This model predicts that the transition rate for escape past a cascade of N sharp barriers is governed by the sum of times 𝑘 𝑜𝑓𝑓 𝑛→𝑛-1 (𝑓) -1 needed to transit individual barriers, starting from the bound state [START_REF] Evans | Energy landscapes of biomolecular adhesion and receptor anchoring at interfaces explored with dynamic force spectroscopy[END_REF][START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF]:

𝐾 𝑜𝑓𝑓 𝑁→0 𝑓 = 1 𝑘 𝑜𝑓𝑓 𝑛 →𝑛 -1 𝑓 -1 𝑛 =𝑁 𝑛 =1 (7)
According to this approach, each individual bond acts as the barrier of a more complex bond made of N identical bonds. 𝑘 𝑜𝑓𝑓 𝑛→𝑛-1 (𝑓) represents the transition rate through each barrier of a multivalent bond. It depends on bond configuration as shown below.

Equal sharing of the force between parallel bonds has the effect of multiplying the force scale for these bonds [START_REF] Williams | Analytical descriptions of dynamic force spectroscopy: behaviour of multiple connections[END_REF]. Moreover, at each failure event, there are less and less possibilities of rebinding and therefore, the time scale for event is decreased and the dissociation rate is increased [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF]:

𝑘 𝑜𝑓𝑓 𝑛→𝑛-1 = 𝑛 𝑘 𝑜𝑓𝑓 0 × exp 𝑓 𝑛𝑓 𝛽 (8)
For zipper bonds, the transition rate is constant at given force:

𝑘 𝑜𝑓𝑓 𝑛→𝑛-1 = 𝑘 𝑜𝑓𝑓 0 × exp 𝑓 𝑓 𝛽 (9)
Combining Eq. 6 and 7 leads to following expressions of the overall dissociation rates across N bonds in parallel:

𝐾 𝑜𝑓𝑓 𝑁→0 = 𝑘 𝑜𝑓𝑓 0 1 𝑛 × exp - 𝑓 𝑛 𝑓 𝛽 𝑁 𝑛=1 -1 (10)
And a pure analytical expression for zipper bonds:

𝐾 𝑜𝑓𝑓 𝑁→0 = 𝑘 𝑜𝑓𝑓 0 × 1 𝑁 × exp 𝑓 𝑓 𝛽 (11)
The series given by Eq. 10 has no analytical expression. However, we presently propose an approximation of this series:

𝐾 𝑜𝑓𝑓 𝑁→0 ≈ 𝑘 𝑜𝑓𝑓 0 × 2 × 𝑓 𝑁𝑓 𝛽 0.5 × 𝑒𝑥𝑝 𝑓 𝑁𝑓 𝛽 (12)
Eq. 11 and 12 show that the global dissociation rate (noted with capital letters: 𝐾 𝑜𝑓𝑓 𝑓 ) for multiple bonds working uncooperatively can be written in the form of a generalized Bell's model:

𝐾 𝑜𝑓𝑓 𝑓 = 𝐾 𝑜𝑓𝑓 0 (𝑓) × exp 𝑓 𝐹 𝛽 (13)
The global dissociation rate, 𝐾 𝑜𝑓𝑓 𝑓 , its prefactor, 𝐾 𝑜𝑓𝑓 0 (𝑓), which may also depend on 𝑓, and the global force scale, 𝐹 𝛽 , show the kinetics of failure of multiple bonds with a necessarily composite energy barrier. Eq. 13 can be linked to an expression of probability 𝑃 𝑡 of bond survival having the form of Eq. 5. The new prefactor 𝐾 𝑜𝑓𝑓 0 (𝑓)

characterizes how the width of the overall chemical energy barrier is modified by force [START_REF] Evans | Dynamic strength of molecular adhesion bonds[END_REF]. The exponential term is a damping factor which takes into account In many experimental cases, the force experienced by an attachment is not constant but increases with time. Bonds under slow loading rates have long lifetimes but only withstand small forces, whereas bonds under fast loading have shorter lifetimes and withstand larger forces [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF][START_REF] Evans | Forces and bond dynamics in cell adhesion[END_REF][START_REF] Thomas | Biophysics of catch bonds[END_REF].

Hence, a classical way to explore the chemical energy landscape traversed in force-driven dissociation is to use adjustable force-ramp representing many loading rates from very slow to extremely fast. An alternative method has been recently proposed by Evans and Kinoshita [START_REF] Evans | Using force to probe single-molecule receptor-cytoskeletal anchoring beneath the surface of a living cell[END_REF] in which a force-clamp method is applied to a single bond. The force-clamp method can be seen as the easiest method to establish the lifetime profile as a constant dissociation force is applied. Depending on the contours of energy function with the distance of separation (one or more chemical energy barriers), some bonds weaken as force increases (as predicted by Bell model, Eq. 4) while some other bonds strengthen as force increases [START_REF] Sokurenko | Catch-bond mechanism of force-enhanced adhesion: counterintuitive, elusive, but[END_REF][START_REF] Thomas | Catch bonds in adhesion[END_REF][START_REF] Thomas | Biophysics of catch bonds[END_REF].

Results and Discussion

Theoretical aspects

Fig. 2 summarizes the theoretical predictions for lifetimes of multiple uncooperative bonds dissociating randomly (Eqs. 10-12). The lifetime of multiple bond (n=N) normalized by the bond lifetime of a single and isolated bond (n=1) is plotted versus the number N of identical bonds constituting the parallel and "zipper" configurations. The main message is that although the lifetime of the collective bonds still decreases with increasing force (as shown by Eq. 11 and 12), bonds organized collectively have a much higher capacity to stabilize the adhesive structure compared to isolated bonds. This stabilization effect increases with increasing the bond number for both parallel and zipper configurations. Note that by contrast, a series of uncorrelated bonds would dissociate faster than single bonds [START_REF] Evans | Energy landscapes of biomolecular adhesion and receptor anchoring at interfaces explored with dynamic force spectroscopy[END_REF]. In parallel and zipper bonds, this effect is particularly marked at the beginning of bond association: n= 2, 3…. Thus, the theory predicts that the early phase of lateral bond association is particularly efficient for stabilizing adhesion.

If both zipper and parallel configurations stabilize adhesion, the parallel configuration, where the dissociation force is equally divided between parallel bonds, reveals a marked effect of the force in the reduction of the global dissociation rate. This is basically due to the associated reduction in the exponential amplification classically described by Eqs 1 and 4 for single bonds, and generalized to collective bonds by Eqs. 11-13.

However, at first glance, the force-dependent reduction in the dissociation rate of bonds characterized by a collective organization may be seen as a counter intuitive result if we forget to consider that both single and collective bonds still dissociate more under a force increase (Eqs. 4,11 and 12). Indeed, the force-induced exponential decrease in bond lifetime is cut down by the collective bond organization, most remarkably in the parallel case. We underline here that the effect of increasing bond stabilization caused by a homogeneous force distribution has been under-considered before. This is nevertheless a fundamental mechanism for the control of cellular adhesion by force, including the force-dependence of dissociation or stabilization rates as well as the force-dependence of weakening or strengthening processes. In other words, lateral bond association (or clustering) is capable to promote an efficient reduction in the force-dependent increase in the dissociation rate of adhesion sites in a context where both collective and individual bonds all tend to fail at higher force.

By contrast, the zipper configuration provides a relevant example of the effect of a heterogeneous distribution of force on the global failure rate. For N uncorrelated bonds in the zipper configuration, the global failure rate is divided by N, meaning that the stabilization mechanism directly depends on the bond number but, contrary to the parallel configuration, not on the force applied. It is equivalent to say that a heterogeneous force distribution limits the beneficial effect -in terms of stabilization -of dividing the force through several bonds.

Another important message brought by the theoretical results plotted in Fig. 2 is that the analytical equation given by Eq. 10 provides a satisfactorily curve fitting while we noticed that previously proposed analytical forms such as in [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF]:

𝐾 𝑜𝑓𝑓 𝑓  𝑘 𝑜𝑓𝑓 0 × 𝑓 𝑁𝑓 𝛽 × exp 𝑓 𝑁𝑓 𝛽
, fail to provide a satisfactorily curve fitting. Note however that the form presently proposed (in Eq. 12) to fit the multiple bond lifetime (Eq. 10) remains consistent with the general form previously proposed to describe dissociation rates of bonds with necessarily complex energy barriers: [START_REF] Evans | Dynamic strength of molecular adhesion bonds[END_REF]. Present results (Fig 2 ) bring a new confirmation of the initial intuitions of Evans [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF][START_REF] Evans | Dynamic strength of molecular adhesion bonds[END_REF]

𝐾 𝑜𝑓𝑓 𝑓  𝑘 𝑜𝑓𝑓 0 × 𝑔(𝑓) × exp 𝑓 𝑁𝑓 𝛽 with 𝑘 𝑜𝑓𝑓 0 × 𝑔 𝑓 = 𝐾 𝑜𝑓𝑓 0 (Evans

Comparison with experimental results

 Parallel bonds (bonds equally stretched by AFM tip)

We attempted to analyze the feature of the parallel configuration described by the multiple bond theory above. The most straightforward example is provided by the reference data obtained by the group of A. Noy [START_REF] Friddle | Counting and breaking individual biological bonds: Force spectroscopy of tethered ligand-receptor pairs[END_REF][START_REF] Sulchek | Dynamic force spectroscopy of parallel individual Mucin1antibody bonds[END_REF] who used dynamic force spectroscopy via AFM to characterize binding of several individual singlechain variable fragments (scFv) antibody to the Mucin1 (MUC1) peptide (MUC1-Ab). This interaction is the main targeting mechanism for a family of experimental radioimmunotherapics for cancer treatment [START_REF] Albrecht | Production of soluble ScFvs with C-terminal-free thiol for site-specific conjugation or stable dimeric ScFvs on demand[END_REF]. Indeed, clinical researchers have used multivalency of protein adhesion as an affinity-enhancing approach in a variety of immunotherapies and imaging techniques to target specific tissues [START_REF] Souriau | Recombinant antibodies for cancer diagnosis and therapy[END_REF]. Noteworthy, authors [START_REF] Friddle | Counting and breaking individual biological bonds: Force spectroscopy of tethered ligand-receptor pairs[END_REF][START_REF] Sulchek | Dynamic force spectroscopy of parallel individual Mucin1antibody bonds[END_REF] found a direct experimental verification of the uncorrelated multiple bond rupture model because they were able to impose a uniform force distribution through a few parallel bonds and compared with the single bond. A remarkable confirmation of the theoretical predictions presented above was brought by these authors because the dissociation rate found for two (respectively three) bonds in parallel drop precipitously compared to the single (respectively two) bond. The present theoretical approach revealed that in addition to the bond number effect, the effect of the force magnitude above the reference force scale, was also a critical parameter. Indeed, in these AFM experiments, force increased sharply in parallel to the bond number (e.g., 𝑓 𝑓 𝛽 ≈ 10 for N=1, 𝑓 𝑓 𝛽 ≈ 20 for N=2, 𝑓 𝑓 𝛽 ≈ 30 for N=3). Theoretical results shown in Fig. 2 reveal that as force increases, the bond number effect on dissociation rate is dramatically increased in uncorrelated parallel bonds with homogeneous force distribution. Present data demonstrate that the main beneficial effect of multivalent interactions compared to single bond is the reduction on the kinetic off-rate and the corresponding increase in lifetime to the extent that the force distribution remains homogeneous. Such a parallel configuration would certainly be desired to improve the binding efficiency of receptor-targeted immunotherapic treatment [START_REF] Souriau | Recombinant antibodies for cancer diagnosis and therapy[END_REF].

 Zipper bonds (Bonds unequally stretched by shear flow)

The viscous shear flow allows studying the motion of receptor-coated cells or particles along a flat ligand-coated surface. At low hydrodynamic force (e.g.,  a few pN) and for low ligand concentrations (e.g., from 9 to 170 ligands/µm 2 in [START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF]), detectable cell/particle arrests reveal single molecular bond formation [START_REF] Kaplanski | Granulocyte-endothelium initial adhesion. Analysis of transient binding events mediated by E-selectin in a laminar shear flow[END_REF]Pierres et al., 1996a;[START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF]. Neglecting rebinding if receptor density is low, the dissociation rate could be deduced from the slope of the natural log of number of cell-surface interactions versus the duration of each event [START_REF] Kaplanski | Granulocyte-endothelium initial adhesion. Analysis of transient binding events mediated by E-selectin in a laminar shear flow[END_REF]Pierres et al., 1996b;[START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF][START_REF] Tissot | Motion of cells sedimenting on a solid surface in a laminar shear flow[END_REF]. In [START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF], the measured unstressed dissociation rates 𝑘 𝑜𝑓𝑓 0 of single bond remained in a close range of values: 2.4-3.8 s -1 , in spite of the wide variety of ligand tested (i.e., P-selectin, E-selectin, L-selectin and PNAd). This is consistent with data obtained by other authors on the same cells [START_REF] Alon | The kinetics and shear threshold of transient and rolling interactions of L-selectin with its ligand on leukocytes[END_REF][START_REF] Alon | Lifetime of the P-selectin-carbohydrate bond and its response to tensile force in hydrodynamic flow[END_REF][START_REF] Kaplanski | Granulocyte-endothelium initial adhesion. Analysis of transient binding events mediated by E-selectin in a laminar shear flow[END_REF][START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF]. The shear stressdependence on dissociation rates, 𝑘 𝑜𝑓𝑓 (𝑓), was affected by the change in receptor-ligand affinity to take values in a higher range: 5-36 s -1 as expected from Bell model (Eq. 4) [START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF]. Beyond, this type of experiment is important as it allows to study the flowenhanced cell adhesion mechanisms [START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF]. Concerning the limits of the single bond assumption in sheared cells, it has been suspected that, as ligand density or shear stress increases, an increased number of multiple bond events is involved in cell adhesion. However, their effect was not thought to deeply affect the estimated values of dissociation rate 𝑘 𝑜𝑓𝑓 (𝑓) [START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF].

In single bond models of cell arrest under shear, an oblique direction of the bond force reacting to the shear force is assumed [START_REF] Chang | The forward rate of binding of surface-tethered reactants: effect of relative motion between two surfaces[END_REF][START_REF] Smith | A direct comparison of selectin-mediated transient, adhesive events using high temporal resolution[END_REF][START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF]. Hydrodynamic theory predicts that once a bond is created with the substrate, it tends to re-equilibrate the translational and rotational velocity components at the cell wall, thus cancelling out the slip velocity [START_REF] Chang | The forward rate of binding of surface-tethered reactants: effect of relative motion between two surfaces[END_REF][START_REF] Hammer | Simulation of cell rolling and adhesion on surfaces in shear flow: general results and analysis of selectin-mediated neutrophil adhesion[END_REF]. This suggests that the direction of bond traction during shear flow is necessarily oblique, i.e., resulting from tangential and normal tensile (or compressive) force components at an adhesive interface [START_REF] Gallant | Cell adhesion strengthening: contributions of adhesive area, integrin binding, and focal adhesion assembly[END_REF]. Based on previous analysis [START_REF] Evans | Detailed mechanics of membrane-membrane adhesion and separation. II. Discrete kinetically trapped molecular cross-bridges[END_REF][START_REF] Ward | A theoretical analysis for the effect of focal contact formation on cell-substrate attachment strength[END_REF], detachment of a sheared cell is expected to occur via peeling of the leading edge of the cell. Accordingly, bond loading is highly non uniform along the contact area since bond forces are maximal at the periphery and decay rapidly toward the center of the cell [START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF]. Thus, in shear flow experiments where independent bonds start to respond collectively, i.e., at increasing shear [START_REF] Zhu | Mechanisms for flow-enhanced cell adhesion[END_REF] or at increasing ligand/receptor concentrations [START_REF] Gallant | Cell adhesion strengthening: contributions of adhesive area, integrin binding, and focal adhesion assembly[END_REF], the zipper configuration appears to be a satisfactorily model because force distribution on bonds is highly non uniform. It means that the huge effects on lifetimes of multiple parallel bonds (Fig 2) such as those reported from AFM experiments [START_REF] Sulchek | Strength of multiple parallel biological bonds[END_REF] could not be expected in shear flow experiments. Accordingly, studying cell adhesion strengthening in sheared fibroblasts, [START_REF] Gallant | Cell adhesion strengthening: contributions of adhesive area, integrin binding, and focal adhesion assembly[END_REF] were surprised to observe that parallel association of integrins by focal adhesion (respectively clustering) accounts for only 30% (respectively 60%) of the strengthening response to shear. We presently postulate that this result is linked to the highly unequal distribution of the resultant force during shear. These mechanical factors would clearly favor the zipper configuration which exhibits the lowest sensitivity to the bond number and no force-dependence for the prefactor of the dissociation rate. We conclude that the homogeneous or heterogeneous character of the force distribution, which is indeed associated to the resulting direction of the force applied on interfacial bonds, is a determinant factor of the kinetics of adhesion. relative decay in dissociation rate, (i.e., increase in lifetime), of collective bonds compared to single bond is enhanced by increasing stretching forces.

In the absence of precise knowledge about the force distribution in living tissues down to the nanoscale, it is easier to consider that the heterogeneous force distribution better fits the standard biological conditions. However, this is not necessarily true. We can say that the parallel organization of bonds constitutes an impressive means to cut down the straightforward process of force-induced reduction in lifetime. This remarkable property of parallel bonds -and to a lesser extent of zipper bonds -has not been so clearly pointed out before. Experiments by AFM constitute the most ideal conditions in which multiple parallel bonds can be stretched uniformly and indeed these AFM results seem consistent with the theory [START_REF] Sulchek | Dynamic force spectroscopy of parallel individual Mucin1antibody bonds[END_REF]. As soon as force distribution on multiple bonds is heterogeneous, the stabilization process will still be effective but much less marked than in the case of uniform distribution. This is because the totality of the force is exerted on the leading bond(s), i.e., the proper of the zipper configuration, and thereby the rupture force of the leading bond, is rapidly attained. We point out here that the force distribution applied on adhesion sites of sheared cells is likely highly heterogeneous. This situation most likely corresponds to many biological situations including the frequent case of cells sheared by the flowing blood and also the pathological cases in general. Taking again the example of the deep lung, the highly elevated stress sites found in diseased emphysematous lungs suggests that heterogeneous stress distribution could weaken adhesion sites in the deep lung and thereby contribute to the disease [START_REF] Gefen | Analysis of stress distribution in the alveolar septa of normal and simulated emphysematic lungs[END_REF]. Finally, another consequence of the present approach is to show that, as soon as bonds are multiple, with the same bond number, but stretched by different force/stress fields, one can expect very different data in terms of bond kinetics. Indeed, specific bond kinetics experiments performed with tightly controlled homogeneous and heterogeneous force distributions have to be designed to further test the presently developed concept.

Figure Legend

Figure 1: The two main configurations considered in the present study of the kinetics of multiple bonds are the parallel configuration (upper schema) and the zipper configuration (lower schema). In the parallel configuration, the dissociation force, f, is uniformly distributed amongst identical bonds. After breaking, force is redistributed amongst the remaining bonds. In the zipper configuration, the leading bond experiences the entire force until breaking while it is passed on to the next. In the two cases, it is assumed that bonds break at random times from first to last, (i.e., described by a Makovian sequence) because the lack of mechanical linkage between bonds ensures a sufficient independence between bonds. (presently in the range: 2-30). Lifetime being by definition the inverse of the dissociation rate, this ratio provides equivalently the evolution of the kinetic rate constant, (or dissociation rate), of the single bond comparatively to the collective bond on which the same force is exerted. In parallel bonds, the force amplifies -by orders of magnitude -the N-dependent increase in relative lifetime. In zipper bonds, the ratio of multiple to single bond lifetimes also increases, but remains proportional to N. The zipper-bond relationship can be superimposed to the

  Note that due to the exponential dependence shown by Eq. 6, bond lifetime is considerably higher than 𝑡 𝐷 . For instance, an integrin-VCAM bond with an energy barrier of around 20 k B T[START_REF] Evans | Using force to probe single-molecule receptor-cytoskeletal anchoring beneath the surface of a living cell[END_REF] has a bond lifetime 𝑡 𝑜𝑓𝑓 0  𝑠 . Bond lifetimes would reach a few hours for energy barriers of 30 k B T. Due to the enormous gap in time scale between ultrafast Brownian diffusion and characteristic times of laboratory experiments (10 -4 s -min), applied forces are essentially stationary on the Brownian scale. Multiple bonds Three assumptions are made to describe the dynamics of the failure of multiple uncorrelated bonds (Fig 1):

  the higher likelihood of reaching the top of the new energy barrier lowered by force in proportion to the global thermal force scale, 𝐹 𝛽 .For parallel bonds: and 15 which relate collective bond properties to properties of single bond show explicitly how the collective bond organization modifies the kinetic adhesion parameters in the case of parallel and zipper configurations respectively.

  about the relevance of the Bell type model (Eq. 13) to describe dissociation rates of multiple parallel and zipper bonds which are individually in agreement with the Bell type model (Eq 11). Advantageously, these analytical expressions allow to delineate the specific effects of bond number and dissociation force on the kinetic adhesion parameters in various bond configurations (Fig 2).

Figure 2 :

 2 Figure 2: The kinetics of uncorrelated multiple independent bonds is compared to the

Acknowledgement

For this work, we acknowledge receipt of a grant from Agence Nationale de la Recherche (ANR-09-PIRI-002-03).

Significance of present results

Theoretical approach and experimental results consistently show that the distribution of force amongst bonds organized collectively and working stochastically fundamentally affects the kinetics of adhesion. Depending on the size of adhesion site, the scale at which these mechanisms operate cover a wide range of scales: from the nanometer scale for the weak initial adhesion sites constituted by only a few molecules (e.g., integrins and talins) to ten micrometers for the mature focal adhesion which contains more than fifty adhesion proteins [START_REF] Zamir | Components of cell-matrix adhesions[END_REF]. The list of adhesion molecules continue to grow and key mechanisms that regulate signal transduction pathways and modulate adhesion assembly still need to be identified [START_REF] Lo | Focal adhesions: what's new inside[END_REF]. Future challenges include the determination of how these molecules work together as a physico-chemical complex and how their coupling contributes to an adapted biological response. This reinforces the need of calculating stress/strain fields exerted on living tissues at various scales and certainly up to the smallest molecular scales, e.g., in the lung [START_REF] Suki | Extracellular matrix mechanics in lung parenchymal diseases[END_REF]. In this perspective, computational fluid/solid mechanics appears for the future as an indispensible tool to get mechanical parameters at various scales, in various living structures and functional conditions, e.g., in the lung [START_REF] Elad | Air-conditioning in the human nasal cavity[END_REF]. More particularly, we find that, if the force exerted on a given adhesion site is distributed homogeneously amongst the bonds, one may expect a force-dependent stabilization process which dramatically depends on bond number. Noteworthy, the stochastic theory of parallel bonds stretched uniformly predicts a drastic reduction in the dissociation rate compared to single bonds. Unexpectedly, the parallel-bond relationship obtained for 𝑓 𝑓 𝛽 = 2. Diamond symbols represent the relationships predicted by the stochastic theory of uncorrelated multiple parallel bonds (see Eq. 10 and Williams [START_REF] Williams | Analytical descriptions of dynamic force spectroscopy: behaviour of multiple connections[END_REF]). Empty squares represent the relationships predicted by the theory of uncorrelated multiple zipper bonds (see Eq. 11 and Williams [START_REF] Williams | Analytical descriptions of dynamic force spectroscopy: behaviour of multiple connections[END_REF] and Evans [START_REF] Evans | Probing the relation between force--lifetime--and chemistry in single molecular bonds[END_REF]). The continuous curves correspond to the proposed adjustment of the stochastic theory using the explicit form given in Eq. 12 to describe the kinetics of uncorrelated parallel bonds.


L Romagnosi

N Gascoin

E El-Tabach

C Strozzi

Transpiration cooling of structures by means of hydrocarbon fuel is a way to sustain the large heat load encountered in combustion chambers. An European collaboration has been set up to investigate the fluid flow through porous media by means of numerical simulation in order to help analysing the experimental results obtained under pyrolysis regime (1200 K, 6 MPa). A parametric study is conducted to understand the effect of numerical modelling of velocity in the porous media, multi-component diffusion, the deployment of in-situ adaptive tabulation for chemistry and the use of perfect or real fluid properties. Due to the complexity of multi-physics phenomena, this work intends to ensure the reliability of the calculations before comparing them with experiments. Two semi-detailed kinetic schemes for dodecane pyrolysis (up to 15 species and 20 reactions) are tested.

Nomenclature

I. Introduction

Porous media are encountered in a wide variety of applications (geology, aerospace, biomass, separation process). In the framework of aerospace structures cooling, ceramic matrix composite materials are used to sustain the large heat load due to combustion [START_REF] Bouchez | PTAH-SOCAR Fuel-Cooled Composite Materials Structure: 2009 Status, 16th AIAA/DLR/DGLR International Space Planes and Hypersonic Systems and Technologies Conference[END_REF][START_REF] Nguyen | Innovative porous SiC-based materials: From nanoscopic understandings to tunable carriers serving catalytic needs[END_REF][START_REF] Arnold | Film Cooling of Accelerated Flow in a Subscale Combustion Chamber[END_REF][START_REF] Kanda | Effect of regenerative cooling on rocket engine specific impulse[END_REF]. Active transpiration cooling enables the internal cooling of the structures by convection between the filtrating fluid and the solid [START_REF] Langener | Experimental Investigation of Transpiration Cooling with Subsonic and Supersonic Flows at Moderate Temperature Levels[END_REF][START_REF] Langener | Experimental Investigations on transpiration Cooling for Scramjet Applications Using Different Coolants[END_REF][START_REF] Kuhn | Application of Transpiration Cooling for Hot Structures, RESPACE -Key Tech. for Reusable Space Sys[END_REF]. In case of using hydrocarbon fuel as a coolant, these heat and mass transfers imply also chemical reactions due to the fuel pyrolysis [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF]. Fluid such as methane are generally considered for rocket engines application [START_REF] Di Matteo | Modelling of Stratification in Cooling Channels and its Implementation in a Transient System Analysis Tool[END_REF]. Multi-species flow is obtained and multi-phase flow can even be generated due to the formation of light species and heavy ones (up to solid coke particles). The presence of solid carbon deposit inside the porosity of the solid modifies its permeability [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF]. Similarly ablation involves chemical reactions with surface effect, radiation and fluid flow in porosity caused by the formation of permeable char [START_REF] Vignoles | Modelling of carbon-carbon composite ablation in rocket nozzles[END_REF][START_REF] Lachaud | Multiscale Approach to Ablation Modeling of Phenolic Impregnated Carbon Ablators[END_REF]. The coal pyrolysis is driven by the same phenomena [START_REF] Patisson | Coal Pyrolysis In A Rotary Kiln: Part I. Model Of The Pyrolysis Of A Single Grain[END_REF].

The understanding of the multiphysics' coupling is necessary in order to predict the efficiency of the cooling technique. The micro-scale of the pores in permeable materials makes this work challenging. The characteristic timescales and length scales of the involved phenomena differ strongly. Either experimental or numerical studies are possible to investigate separately the fluid flow, the heat transfer and the chemistry aspects. Nevertheless experimentally investigating the coupling of these three aspects is quite complex. The formation of pyrolysis products inside the porous medium is both dependent on temporal and spatial dwelling. The residence time is different for each species (due to multi-component diffusion) [START_REF] Taylor | Multicomponent Mass Transfer[END_REF] while the thermal and hydraulic gradients depend on microscopic considerations such as fibers' layer, heterogeneity of the material, tortuosity, permeability [START_REF] Zhou | Analytic characterization of the permeability of dual-scale fibrous porous media[END_REF]. The coke formation impacts the permeability, and hence also the fluid flow and the associated heat transfer [START_REF] Zhang | Numerical Simulations of Fluid Flow and Convection Heat Transfer Through Fluid/Porous Layers[END_REF]. Moreover, catalycity might be another phenomenon to be considered depending on the solid material composition [START_REF] Gascoin | Thermal and hydraulic effects of coke deposit in hydrocarbon pyrolysis process[END_REF]. Consequently, to help analysing the experimental results a numerical approach is proposed in this work obtained previously with a specific bench of heavy hydrocarbons pyrolysis within porous media [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF].

A brief overview of existing studies from open literature is proposed to clearly state what is currently achieved in the framework of porous reacting flow simulation. Numerical modelling of pyrolysis in permeable structure is found to be taken into account by either a one-step Arrhenius law or with semi-detailed kinetic [START_REF] Moghtaderi | The state-of-the-art in pyrolysis modelling of lignocellulosic solid fuels[END_REF]. The work of Ratte et al. is interesting despite it is not specifically dedicated to porous materials because they investigate numerically the chemical composition profile in a wood particle [START_REF] Ratte | Mathematical modelling of slow pyrolysis of a particle of treated wood waste[END_REF]. Their chemistry description is limited (about 10 species) but it simulates the spatial profiles of coupled thermal and chemical parameters. Another modelling of reacting system with explicit simulation of the porous zone is proposed by Margolis and Telengator [START_REF] Margolis | Influence of Sublimation and Pyrolysis on Quasi-Steady Deflagrations in Confined Porous Energetic Materials[END_REF]. Despite they do not consider the chemical reaction pathway, their description enables investigating closely the effect of the pressure in case of exothermic degradation within the energetic materials. The work of Ravi et al. [START_REF] Ravi | Development of a semi-empirical model for pyrolysis of an annular sawdust bed[END_REF] is similar in terms of modelling despite the applications differ. A single chemical reaction (data derived from ThermoGravimetry measures) is considered and the void fraction produced during the degradation is modelled. This enables considering the thermal and fluid coupling. The thermal properties of the solid vary with the modification of the porosity, which is related to the variation of the chemical composition (char and void content). Through these three studies [START_REF] Ratte | Mathematical modelling of slow pyrolysis of a particle of treated wood waste[END_REF][START_REF] Margolis | Influence of Sublimation and Pyrolysis on Quasi-Steady Deflagrations in Confined Porous Energetic Materials[END_REF][START_REF] Ravi | Development of a semi-empirical model for pyrolysis of an annular sawdust bed[END_REF], it is seen that the modelling either focuses on the porous flow or on the chemistry but rarely on both key points, or only with simple kinetic schemes.

Concerning the porous zone, its modelling is still not obvious (internal 3-D heterogeneous geometry, radiation, physical properties, fluid-solid interactions). Jaojaruek and Kumar [START_REF] Jaojaruek | Numerical simulation of the pyrolysis zone in a downdraft gasification process[END_REF] consider the porous media with constant average properties in equilibrium with the fluid phase. The material is assumed to be a grey body. A permanent regime is considered and a single Arrhenius law is used for solving the chemistry. Their code allows the modelling of gasification process involving the conversion of hydrocarbon fuels. Similarly for wood thermal degradation, Sand et al. [START_REF] Sand | Numerical prediction of the transport and pyrolysis in the interior and surrounding of dry and wet wood log[END_REF] use Navier-Stokes equations and energy balance to compute the porous flow. The permeability varies with the solid conversion rate and a mixed series-and-parallel law takes the different constituents into account. A twosteps mechanism with the formation of six solid, liquid and gas species is used. The radiation is determined globally for a grey body. The emissivity varies with the production of char. Some experiments are performed to validate this code. The thermal and catalytic cracking of hydrocarbon feedstocks also conduct to study the formation of light hydrocarbons in porous catalytic materials [START_REF] Van Mao | The thermocatalytic cracking process for the production of light olefins and transportation fuels from gas oils[END_REF]. The numerical work of Mesalhy et al. [START_REF] Mesalhy | Numerical study for enhancing the thermal conductivity of phase change material (PCM) storage using high thermal conductivity porous matrix[END_REF] explicitly considers the Darcy's and Forchheimer's terms for the fluid flow permeation. The thermal and porosity properties of the solid are observed with respect to their effects on the change of phase of the inner fluid. Non-equilibrium conditions between solid and liquid are used, which is not so common in open literature. In these works [START_REF] Jaojaruek | Numerical simulation of the pyrolysis zone in a downdraft gasification process[END_REF][START_REF] Sand | Numerical prediction of the transport and pyrolysis in the interior and surrounding of dry and wet wood log[END_REF][START_REF] Van Mao | The thermocatalytic cracking process for the production of light olefins and transportation fuels from gas oils[END_REF][START_REF] Mesalhy | Numerical study for enhancing the thermal conductivity of phase change material (PCM) storage using high thermal conductivity porous matrix[END_REF], neither turbulence model nor highly detailed mechanism is used.

Numerous papers can be found on the theory of fluid filtration within porous media [START_REF] Zhang | A fractal model for gas permeation through porous membranes[END_REF]. The different analytical laws describing the fluid flow are now well-known [START_REF] Shahnazari | Permeability Prediction Of Porous Media With Variable Porosity By Investigation Of Stokes Flowover Multi-Particles[END_REF]. Analytical equations for transport equations of species in porous media are also found in studies related to permeation through membranes for separation process [START_REF] Dingemans | Mass transfer characteristics for VOC permeation through flat sheet porous and composite membranes: The impact of the different membrane layers on the overall membrane resistance[END_REF]. These equations are simpler than those which are necessary to describe multi-species reacting flow. Some analytical work can be mentioned considering the non-equilibrium thermal flow within porous media [START_REF] Mahmoudi | Analytical investigation of heat transfer enhancement in a channel partially filled with a porous material under local thermal non-equilibrium condition[END_REF][START_REF] Lee | Momentum boundary layer and its influence on the convective heat transfer in porous media[END_REF][START_REF] Shokouhmand | Simulation of laminar flow and convective heat transfer in conduits filled with porous media using Lattice Boltzmann Method[END_REF].

In addition to pyrolysis processes, combustion is also observed in porous media with an increased efficiency as the solid heat conditions enhance the heat transfer from the burned gases to the fresh gases. In case of solid fuel, its porosity plays a role on the combustion efficiency [START_REF] Madrigal | Effective Heat of Combustion for Flaming Combustion of Mediterranean Forest Fuels[END_REF]. The combustion studies within energetic porous materials are of interest for solid propulsion and explosives applications [START_REF] Asay | The Role Of Gas Permeation In Convective Burning[END_REF]. Reactive flow permeation can even be found under detonation regime, which shows the wide diversity of studies related to porous media. Al-Hamamre et al. [START_REF] Al-Hamamre | Combustion Of Low Calorific Gases From Landfills And Waste Pyrolysis Using Porous Medium Burner Technology[END_REF] present an experimental work with some numerical results to investigate the combustion of gases (CH 4 and CO 2 mainly) resulting from waste recycling. The detailed chemistry is split between heat and mass transfers which are solved separately. The 3-D code solving the fluid flow in the porous media considers homogeneous properties of the porous solid with the inner fluid.

The modelling of porous flow also refers to the solid itself. A large part of studies related to carbon based porous materials is related to the production of carbon films for their adsorbent and catalyst characteristics. Experiments on their pyrolysis enable studying the associated changes of porosity [START_REF] Ania | Pyrolysis of activated carbons exhausted with organic compounds[END_REF][START_REF] Krzesinska | The effect of pyrolysis temperature on the physical properties of monolithic carbons derived from solid iron bamboo[END_REF]. The determination of the pore size is conducted either by visualisations or by permeation method [START_REF] Kim | Evaluation of permeable pore sizes of macroporous materials using a modified gas permeation method[END_REF]. The experiments are generally conducted under steady-state conditions but some transient works also exist [START_REF] Capek | Permeation of gases in industrial porous catalysts[END_REF]. The oxidative pyrolysis is also investigated to determine the formation of pores (void fractions and pore diameter distribution) in the solid phase due to the heating and chemical reactions (function of the conversion rate) [START_REF] Shchipko | Catalytic pyrolysis of Kansk-Achinsk lignite for production of porous carbon materials[END_REF][START_REF] Ruiz | Study of porous development in pyrolysis chars obtained from a low-volatile coal[END_REF]. The pores' size and distribution in case of coal pyrolysis at 1023 K has been experimentally investigated by Yardim et al. in presence of water [START_REF] Yardim | Formation of porous structure of semicokes from pyrolysis of Turkish coals in different atmospheres[END_REF]. The repartition of solid, liquid and gas products has been observed. Takeichi et al. [START_REF] Takeichi | Preparation of porous carbon films by the pyrolysis of poly(urethane-imide) films and their pore characteristics[END_REF] propose SEM images of the porous carbon materials in case of pyrolysis of raw organic material up to 1173 K. The pore size distribution is driven by the initial composition of the material. The synthesis of carbon hollow spheres and tubes is also obtained by hydrocarbon raw material pyrolysis [START_REF] De | Hollow porous carbon microspheres obtained by the pyrolysis of TiO2/poly(furfuryl alcohol) composite precursors[END_REF][START_REF] Quan | Synthesis of carbon nanotubes and porous carbons from printed circuit board, waste pyrolysis oil[END_REF]. All of these works are related to the degradation of the porous material itself while no such a solid pyrolysis is expected in the present work. Nevertheless, these references are interesting because they show the link between porosity and chemistry in case of degradation. It has been shown in previous work that the permeability should also change in case of hydrocarbon pyrolysis within porous material because of coke formation [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF]. It is thus clear that the coupling between porous material and porosity changes are mainly studied by experiments when coupled to chemistry. Very few works are dedicated to the simulation of hydrocarbon pyrolysis without direct porosity modification (i.e. in case of solid degradation). In case of fuel degradation, the formation of coke modifies the permeability but this is quite indirect. It depends on the fluid flow and on the size of the particles.

It can be concluded from this literature survey that the numerical simulation involving chemistry, porous material and heat and mass transfers are quite rare. Only few consider explicitly the radiation through global law and none of them investigate the multi-species flow. The physical properties are generally computed as average values for the porous zone between the solid and fluid data. The microscopic description of the through flow is rarely encountered in case of reactive flow.

The present study aims at proposing a numerical simulation of an experimental bench made for fuel pyrolysis [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF]. To ensure the reliability of this simulation, numerous verifications and validations are performed since the modelling of reactive fluid flow in porous media is still scarce in the open literature. Using semi-detailed kinetic mechanism is required to propose a fine description of the chemical process. The maximum size of the kinetic scheme compatible with Computation Fluid Dynamics (CFD) code is investigated. This simulation is then used for comparison to improve the analysis and the understanding of the experimental results.

II. Description of the test cell geometry and of the numerical code

A steady-state and axisymmetric simulation is proposed using the commercial CFD code, Ansys/Fluent, to support the experiments performed previously [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF][START_REF] Travajno | Calcul multiphysique dans une paroi poreuse[END_REF]. The experimental apparatus to be simulated is referred in the following as the "test cell". Large details on the numerical methods are given in the present section since this is as much important as the results themselves to ensure their consistency.

A. Design of the test cell and boundary conditions

The test cell [START_REF] Gascoin | High temperature and pressure reactive flows through porous media[END_REF], made of stainless steel, is composed of a pipe of 5 mm diameter where the fluid is injected, of a central body (36 mm long and 16 mm of inner diameter) and of an exit pipe, which is identical to the inlet one (Figure 1). The porous sample is located slightly offset from the centreline of the cell. It has a radius of 15 mm and a thickness of 3 mm. For the present work, a non-reacting and non-catalytic medium made of stainless steel is selected (the so-called SS3 sample with porosity of 30.4 %, grain diameter of 14.1 m and pore diameter of 4.1 m) [START_REF] Gascoin | High temperature and pressure reactive flows through porous media[END_REF]. The grid consists in 14,242 structured mesh cells with typical characteristic size of 0.2 mm, distributed uniformly along the test cell. The reference position along the longitudinal axis is taken at the inlet of the porous zone and the one along the radial axis is on the symmetry line (Figure 1). The x-axis is oriented from the left to the right and the y-axis from the bottom to the top. Some elements of the real test cell (like the sensors, the drilling and the threading) have not been modelled to simplify the numerical case. These elements are supposed to have a negligible impact on the turbulence development and thus on the heat transfers.

The boundary conditions are the one of the experimental cases to be simulated. The walls temperature is constant and it is the one experimentally measured for stabilized conditions. The inlet fluid presents the same temperature because of a pre-heating zone in the experiments [START_REF] Fau | Fuel Pyrolysis through Porous Media: Coke Formation and Coupled effect on Permeability[END_REF]. The downstream pressure is adjusted to get the experimentally measured one at the inlet of the system. The mass flow rate is imposed at the inlet. Two fluids (nitrogen and dodecane) and several operating conditions have been tested. Dodecane has been preferred for comparisons purpose but methane should be considered in latter steps. The tests cases with nitrogen were formerly compared to the experimental work (Figure 2) to define the best computing methods for the inert gas flow (work of Fedioun et al. since 2008 [START_REF] Blondel | [END_REF][START_REF] Orear | [END_REF][START_REF] Aungier | A Fast, Accurate Real Gas Equation of State for Fluid Dynamic Analysis Applications[END_REF]). The results corresponding to the hydrocarbon only are presented in this paper to specifically focus on reactive flow, which is the aim of the present work.

B. Numerical modelling and assumptions

The porous zone is treated as a fluid zone where the pressure drop is simulated by adding a sink term S i (Eq. 1) in the momentum equation (Eq. 2). S i is defined as the sum of viscous loss term (the Darcy's term, the first one on the right hand side of Eq. 1) and of the inertial loss term (the second one on the right-hand side of Eq. 1). The permeability and the inertial resistance values of the SS3 porous sample are taken directly from the experimental results (α = 2.061•10 -13 m 2 and C 2 = 4,924,000 m -1 ) [START_REF] Gascoin | High temperature and pressure reactive flows through porous media[END_REF]. The modified governing equations for a single phase flow are given by Eq.2, Eq. 3 and Eq. 4.

i i i v v C v S    2 1 2    (1)       i f f f S g P v v t v                        (2)     0       v t f f    (3)               h f i i i eff f f s s f f S v J h T k P E v E E t                          1 (4)
In this paper, the fuel physical properties are determined as a function of temperature during the simulation. Concerning the computations of real fluid properties, the modified Redlich-Kwong equation developed by Aungier was chosen [START_REF] Aungier | A Fast, Accurate Real Gas Equation of State for Fluid Dynamic Analysis Applications[END_REF]. It is specifically developed for solving supercritical flows. It requires knowing the species critical properties such as the critical temperature, pressure and volume and the acentric factor. The mixture critical properties are evaluated by Aungier-Redlich-Kwong mixture law. This cubic law takes into account the reduction in pressure due to the intermolecular attractive forces and the lower available volume due to the repulsive forces considering also that each particle occupies its own volume.

In the porous zone, the physical properties are obtained combining properly the fluid characteristics with the solid ones. For example, by default, eff k is the thermal conductivity of the binary system fluid plus solid evaluated with a parallel law (Eq. 5). The thermal conductivity of the SS3 material is calculated as a function of the temperature with an averaged method between the Parallel and the Maxwell Upper Limit methods (k SS3 = -1.1873•10 -5 •T 2 + 2.6261•10 -2 •T + 6.7647) [START_REF] Travajno | Calcul multiphysique dans une paroi poreuse[END_REF]. The radiative contributions are not considered due to the low thermal gradients which are typically found in pyrolysis processes (in opposite to combustion cases). In addition, using radiation is difficult because of a numerical code limitation which models the porous zone as a fluid zone. Some computations were conducted by considering the radiation but no visible discrepancies have been observed. For this reason, no radiation is explicitly considered in the present work.

  s f eff k k k      1 (5)
Several tests have been performed to choose the best computing options. A pressure based solver is employed and the gradients are solved with a Least Squares Cell based method. A 2 nd order scheme is used for the momentum and the energy equations. When ambient temperature tests are simulated with dodecane, incompressible fluid hypothesis is preferred. When reactive dodecane flow is considered (roughly for temperature higher than 800 K), the species transport equations are added to the Navier-Stokes equations system (Eq. 2 to Eq. 4) and the case is slightly simplified switching from second to first order scheme. In this case, the transport mass balance is computed for all the species except for the major one to ensure a mass fraction sum equal to unity. The low mass flow rate ensures that the turbulence effects are negligible in the entire domain and low total pressure drop due to the dissipative effect at the wall. For this reason, the turbulence model is turned off and the wall friction is set to zero in this case.

All of the numerical test cases to be considered in this work correspond to laminar flows. Anyway, it is possible to include the small turbulence effect within the cell in order to obtain a more realistic solution and to consider the flow through the porous media as laminar. That means that the turbulence at the inlet of the SS3 sample will be transported downstream and its effect on the mixing and on the momentum will be ignored. The assumption, in considering the flow as laminar in the porous medium, has been verified as follows. Analytically, by reversing the Forchheimer's equation and by imposing the pressure drop obtained numerically, it is possible to calculate an equivalent velocity. For a nitrogen case with the highest mass flow rate value which has been simulated (pressure drop of 2.5 MPa and mass flow rate of 3.98 g.s -1 ), the fluid speed is about 0.63 m.s -1 . This mass flow rate is the 6/14 maximum one in the test cases which have been performed. It is much lower for the dodecane tests (lower than 0.1 g.s -1 ). For this equivalent velocity of 0.63 m.s -1 , the Reynolds number based on the pore diameter is Re dp = 3.63 and the one based on the grain diameter is Re dg = 12.49. The first value suggests that the flow inside the porous medium is laminar while the second one is slightly over the limit given by Sobieski and Trykosko for such a flow [START_REF] Sobieski | Darcy and Forchheimer laws in the experimental and simulation studies of flow through porous media[END_REF]. Indeed, they found that the ranges of validity for the Darcy's law is for a Reynolds number from 10 -5 to 10, the Forchheimer's range is from 10 to 80 and the Brinkman's equation (Eq. 6) is valid over 120. As a consequence, the flow is laminar in all the cases presented in this work. This is also probably valid for the above nitrogen case (Re dg = 12.49) because the limit at Re = 10 is not so accurate. The definition of the Reynolds number (either a function of pore or grain diameters) is even questionable.

v v K v K x P F D    2 2           (6)
In case of reactive flow, a semi-detailed kinetic mechanism can be used. The one proposed by You et al. [START_REF] You | Detailed and simplified kinetic models of n-dodecane oxidation: The role of fuel cracking in aliphatic hydrocarbon combustion[END_REF] is considered (56 species and 289 reactions). It has been manually reduced by the present authors to get different size mechanisms. First, a very simple scheme has been tested (Eq. 7 and Eq. 8). It involves four species (dodecane, ethylene, propyl and butyl) and two reactions. Other mechanisms have also been generated and tested. Among them, one with 15 species and 20 reactions will be presented in the results (section III).

C 12 H 26 => 3 C 2 H 4 +2 C 3 H 7 (7)
C 12 H 26 => 2 C 2 H 4 +2 C 4 H 9 (8)
The In-Situ Adaptive Tabulation (ISAT) method has been used for some reactive cases. It consists in tabulating the chemistry after having construct data tables during the simulation by using the detailed kinetic mechanism. The fundamental principle of the algorithm is to avoid the direct integration over time of the chemical source terms (one for each species) in updating the new species mass fractions. Indeed, it is possible to cover all the realizable initial states and time steps, saving these in tables and solving the chemistry by "look-up tables". Creating the tables a priori including also none encountered states is a waste of time and for this reason the tables must be built during the simulation (in-situ), based on realistic states. More details can be found in [START_REF] Pope | Computationally efficient implementation of combustion chemistry using in situ adaptive tabulation[END_REF] on the ISAT theory. When not mentioned in the paper, the tolerance for the ISAT convergence is 10 -6 to compute the residuals. More generally, all the residuals are below 10 -5 in order to ensure the right convergence of mass, momentum and energy. The authors noted during their work that the pressure based solver shows a faster convergence than the density based one, when ISAT algorithm is turned on.

The multi-species diffusion has been considered in the numerical simulation. The effects of the enthalpy transport due to the species diffusion are taken into account adding the following term in the energy equation:

    N j j j J h where j mix j j Y D J     , 
is the diffusion flux of the j-species and mix j D , is the diluted diffusion (i.e.

the diffusion coefficient of j-species in the mixture). Considering the Soret effect, the diffusion flux vector J j must be rewritten as (Fick's law):

T T D Y D J j T j mix j j       , ,



. It allows considering the mass diffusion due to the thermal gradient, which implies that heavy molecules diffuse less rapidly than the light ones. When the full multi-component diffusion is turned on, the Maxwell-Stefan equation is used to compute the diffusion flux vector obtaining a generalized Fick's law where mix j D , is replaced by the mutual diffusion of each species in other species j i D , [START_REF] Taylor | Multicomponent Mass Transfer[END_REF].

This model requires computing N 2 co-factor determinants of size    

1 1    N N
and one determinant of size

N N  [52],
where N is the number of species.

III. Parametric study on the modelling options of the numerical code

The aim of this work is to investigate the effects of several parameters on the numerical results to determine the best way to simulate the complex multiphysics compressible unsteady reactive process. Two velocity models are considered as well as the species diffusion, the In-Situ Adaptive Tabulation (ISAT) parameters, the kinetic scheme and the perfect or real fluid behaviours.

A. Fluid velocity modelling in the porous zone

1. Inert fuel at the ambient temperature For a mass flow rate of 65 mg.s -1 and an inlet pressure of 3.5 MPa, the cold kerosene surrogate enters the cell at ambient temperature with a velocity of 4.4 mm.s -1 and a density of 746 kg.m -3 . Crossing the porous zone, it undergoes a pressure drop of about 7600 Pa (Figure 3a). The high viscosity of the fluid (µ = 0.0014 Pa•s at the inlet) bounds the recirculation zone close to the upstream step and it nearly forms a uniform flow at the porous medium inlet. Observing the contour of axial velocity in Figure 3b, no velocity variation within the porous medium appears. It is due to the superficial velocity model used in the numerical code. The superficial velocity is defined as the ratio between the volume flow rate [m 3 /s] and the geometrical cross section [m 2 ]. This value is constant from the inlet to the outlet. The pressure drop and its impact on the density are not considered. The porosity is not used in the transport equations, which means that the cells of the porous medium are considered to be entirely open. As a consequence, this model cannot predict the increase in velocity due to reduced cross section. For this reason, this model is judged to be not adapted despite it is interesting for computation time and stability reasons. The velocity gradients inside the porous medium should not be neglected because they become important when chemical reactions occur. Thus, the physical velocity model is tested. The effective cross section is about 65 mm 2 when considering the porosity value (instead of 201 mm 2 which is the real cross section exposed to the flow). Thus, the axial velocity increases by 233 % along the axis of symmetry from 0.37 mm.s -1 to 1.24 mm.s -1 (Figure 4a, for non-reactive isothermal dodecane at 300 K). This sudden jump of velocity at the porous medium inlet is due to the mass conservation. The reduction in cross section at the interface fluid-porous zone (upstream side) imposes that the product ρU x increases for the conservation of mass flow (A eff < A exp ). Since the density is set constant, the flow speed increases and it remains constant along each streamline. At the outlet, the opposite of the inlet appears (sudden decrease of the fluid speed). By calculating varying density, the effect of the pressure, which decreases through the permeable sample, could have been seen. This would correspond to a velocity increase. Considering the velocity change in the porous zone modifies the Reynolds number (Re dp = 2.78.10 -3 for the physical model) and moreover the residence time (t r_phy = 2.42 s instead of t r_sup = 8.27 s). This is of major importance because the residence time directly impacts the heat transfers between the porous solid and the fluid. As a consequence, it plays a role on the pyrolysis process. For this reason, the way the fluid filtrates through the porous zone needed to be investigated.

Reactive fuel at temperature over the pyrolysis one

In case of hot flow (inlet fluid: 3.5 MPa and 1200 K; walls: 1230 K; 2 reactions scheme), the velocity increases through the porous zone due to the heating process (decrease of the density) (Figure 4b). This trend is less visible for the superficial model than for the physical one. The rise is reinforced for the physical model because of the pressure drop effect on the density (due to pyrolysis). In addition, the fluid flow model in the porous zone has a major impact on the chemistry. To investigate this point, the species mass fractions have been averaged for each x-position over the entire y-position and the difference of values between both cases is observed (Figure 5a). At the test cell outlet, the superficial model presents a higher pyrolysis rate than the physical one because the remaining mass fraction of dodecane is lower (Figure 5b). This is attributed to a longer residence time of the fluid in the porous medium with the superficial model as presented in previous section. Nevertheless, the differences remain quite small.

B. Effect of the chemical reactions on the heat and mass transfers

To identify the effect of the chemistry (fuel cracking), a similar computation to the one above at high temperature (1200 K) is achieved by turning the chemical reactions off in the porous zone. The thermal effect is investigated (without chemistry). Along the porous medium in the x-direction, the axial velocity clearly increases when the chemical reactions are solved (Figure 6a). This is due to the decrease of the molecular weight during the pyrolysis process. The fluid constant (r=R/M) increases resulting into an expansion of the flow. Because of the boundary conditions, the outlet pressure in both reacting and none reacting cases are similar but the inlet one slightly differs (Figure 6b). It is higher for the reactive case than for the other one, which is in agreement with the higher fluid velocity. At the porous zone inlet, a temperature rise is observed (Figure 6c). This is due to the conductive heat transfer between the porous medium and the fluid. The high thermal conductivity of the solid enhances the conductive heat transfer and the temperature increases of about 50 K in less than 0.1 mm for the case without reactions. Then, it continues to increase along the SS3 material (22 K more over the 3 mm). Around x = 3 mm, close to the outlet of the porous zone, the temperature decreases sharply because a part of this accumulated heat is transferred to the downstream system. By comparing both test cases, the heat transfers are clearly impacted by the chemistry. The temperature gradient between points 0 mm and 3 mm is divided by a factor 2 when considering the chemistry (34 K for the reactive case and 72 K for the inert one). Due to the endothermic effect of the pyrolysis, the inlet temperature rise is the lowest for the reactive case, only 20 K instead of 50 K (Figure 6c). The maximum fluid temperature reached in the reactive case is 1002 K (starting from 968 K at the inlet) instead of 1078 K for the inert case (starting from 1006 K at the inlet). The combined effect on the density decrease of the thermal heating of the fluid and of the pressure drop is different between both cases (Figure 6d). For the reactive case, the chemistry tends to smooth the sharp gradients observed at the porous zone inlet because of the difference in terms of cross section. For the inert case, the density is almost constant around 30 kg.m -3 ; which shows that the slight pressure drop (around 7000 Pa) has a minor impact on the density. A similar conclusion for the thermal effect on the density can be proposed. As a consequence, the density changes in the porous zone are mainly due to the chemistry since the reacting case presents much density variations. This imposes to investigate the way the multi-species flow is treated. The computation method of the density will be investigated in section III.E.

C. Tabulation of the chemistry during kinetic computations

The stiff problems are characterized by different chemical timescales. The ISAT algorithm offers a faster convergence but the accuracy of the results strictly depends on the configuration parameters (table size, number of trees and tolerance). The first two parameters mainly affect the time required by the solver to build the tables. Therefore a sensitivity analysis has been performed on the ISAT error tolerance. As an example, the results of two calculations with 10 -4 and 10 -5 tolerances have been compared and the resulting differences in terms of species mass fraction are presented in Figure 7a,b. Although all residuals are converged, the reduction of the ISAT tolerance from 10 -4 to 10 -5 causes some changes in the species mass fraction to the third decimal. This may be considered as negligible but when using more complex kinetic schemes, this difference could be of importance. Indeed, some important radicals are formed with much lower amount (generally less than 10 -3 ). The higher discrepancies are found for the dodecane in the porous zone (up to 0.009 in terms of mass fraction) (Figure 7a). The differences become less significant when reducing further the ISAT value to 10 -6 (Figure 7c,d). When decreasing the tolerance by one order of magnitude (10 -5 to 10 -6 ), the difference is divided by a factor 4. As a consequence, this process should be repeated until the species variation is lower than a fixed limit. However, a lower tolerance requires a higher computation cost and a compromise between speed and accuracy should be found. For the simplest kinetic scheme (4 species and 2 reactions), the ISAT tolerance of 10 -5 is judged to be acceptable while a tolerance of 10 -6 is preferred for the more detailed schemes, where some fluctuations of species could appear for coarse tolerances. The ISAT tolerance also impacts the heat and mass transfers. As an example, a discrepancy of 10 % is found between two cases on the fluid temperature in the recirculation zone, before the porous zone inlet. Indeed, a temperature around 940 K is calculated for a tolerance of 10 -4 and it is about 1040 K for a tolerance of 10 -6 .

D. Effect of the multi-component diffusion

To evaluate the contribution of multi-component diffusion in the coupling between chemistry and mass transfer, a more complex chemical kinetic mechanism has been used (15 species and 20 reactions). It includes lighter species such as radicals H and CH 3 . This mechanism is a simplified version of the mechanism developed by You et al. [START_REF] You | Detailed and simplified kinetic models of n-dodecane oxidation: The role of fuel cracking in aliphatic hydrocarbon combustion[END_REF]. A maximum around hundred reactions was found by the authors to run the numerical simulation within a reasonable time (less than one month for steady-state calculation) with a desktop computer (8 cores). Using a 20 reactions scheme enables faster computations. Two computations have been done with and without diffusion. The effect of the diffusion is observed on dodecane and on hydrogen (Figure 8). For a number of y-positions, the x-axis evolutions of the mass fractions are slightly shifted downstream when considering the diffusion. As a consequence, the pyrolysis process is slightly delayed due to the faster diffusion of light species such as the radicals. The explanation is the following. Upstream to the porous medium, the thermal gradient along the x-coordinate is negative. This means that when radical species are produced their axial diffusion velocity increases due to the Soret effect. Consequently, the reaction zone slightly shifts downstream despite the effect of the diffusion is clearly less than the other heat and mass transfer effects. This result is still interesting because in other reactive process like combustion system, the temperature gradients between hot products and cold reactants are so high that the hot species upstream diffusion ( 0

   x T
) enhances the reactions. It is now demonstrated to be not the case in the present pyrolysis application.

E. Effect of the fluid flow properties

Previous results were obtained with the restrictive assumption of ideal behaviour of kerosene under high pressure and high temperature conditions. These operating conditions are not within the range of validity of the ideal gas law. The kerosene enters the cell under supercritical state and the real fluid model must be used (the so-called real gas law). The computation method has been detailed in section II.B. Two test cases with the same inlet fluid conditions (890 K and 3.5 MPa) are conducted with ideal and real gas laws. The effects of the real gas model are found on the heat and mass transfers (Figure 9). A lower velocity is encountered at the cell inlet for the real gas model (32 mm.s -1 vs. 38 mm.s -1 for the ideal gas law) (Figure 9a). This is due to the boundary conditions because the velocity depends on the density for a constant mass flow rate. The ideal gas approach underestimates the density: ≈ 80 kg.m -3 instead of ≈ 94 kg.m -3 for the real fluid law (Figure 9b). The two models mostly match for 0 mm < x < 15 mm and then they differ slightly again downstream. This is clear on the temperature profile (Figure 9c). The small difference found on the outlet velocity between the two laws, about 1.5 %, may explain the thermal discrepancies (about 2% on the temperature field). This is related to the lower fluid velocity (in the case of real-gas model), which enhances the residence time and decreases the heat transfers from the walls to the fluid. The effect on the pyrolysis rate is quite negligible (Figure 9d). No significant impact of the properties computation method is found on the pressure distribution.

A useful parameter to evaluate where the real fluid law is more appropriate is the compressibility factor Z. It is defined as the ratio P/RT. Its value is less than unity (≈ 0.86) at the inlet of the porous medium and it increases along the cell (Figure 10). It is close to unity downstream to the porous material which explains why the density does not differ so much between both calculations. The numerical simulation is demonstrated to present important advantages for the reactive permeation studies despite experiments are mandatory to ensure the validity of the computations, where all the microscopic phenomena (heterogeneous reactions, pores size and distribution) cannot be accurately simulated but only fairly modelled.

IV. Conclusion

The numerical investigations of hydrocarbon fuel pyrolysis in porous media are generally found in extremely limited number in the open literature. Most of the related works fall in the ablation topics and wood pyrolysis fields. Experimentally, fuel decomposition is slightly more investigated despite the existing studies remain rare. To enhance the comprehension of such coupled heat and mass transfers with chemical reactions, numerical simulation is of interest to grant access to the residence time distribution, to the flow heterogeneities, to the not uniform permeation profiles and to the thermal and chemical gradients. The present work, conducted within a European collaboration, enables to settle a correct representation of an experimental reactive permeation bench. Numerous numerical points (radiation, multi-component diffusion, fluid properties, fluid velocity in the porous zone and computation way of the chemistry) have been investigated to ensure the validity of the results. Thanks to this simulation development, reactive fluid will be handled with comparison to experiments in order to provide a deeper analysis of these experimental results. A comparison to the one dimensional analytical Forchheimer's equation will enable observing the effect of two dimensional heterogeneous flow profile on the inert fluid permeation. Methane should be later investigated for rocket engines application.

2 C

 2 medium directly exposed to the flow (m²) eff A = effective cross section of the porous material related to the porosity (m²) = inertial resistance factor (m -1 ) general quantity  upstream and downstream to the porous medium ()

Figure 1 .Figure 2 .

 12 Figure 1. Geometry of the numerical test cell

Figure 3 .Figure 4 .

 34 Figure 3. Streamlines, absolute pressure (a) and axial velocity (b) for dodecane at T = 300 K and P = 3.5 MPa (ρ constant).

Figure 5 .

 5 Figure 5. Spatial difference of averaged species mass fraction, along y-axis for each x-position, between the physical model and the superficial one (a) and dodecane mass fraction along the x-axis for the two velocity formulations (b).

Figure 6 .

 6 Figure 6. Velocity (a), pressure (b), temperature (c) and density (d) distribution along the axis of symmetry of the porous medium with (green line) and without reactions (red line) in the porous zone

Figure 7 .

 7 Figure 7. Differences of C 12 H 26 (a and c) and C 2 H 4 (b and d) mass fractions computed between ISAT tolerances of 10 -5 and 10 -4 (a and b) and of 10 -6 and 10 -5 (c and d).

Figure 8 .

 8 Figure 8. Comparison of dodecane (a) and of monatomic hydrogen (b) mass fraction distribution with (solid line) and without (dash-dot line) multi-component diffusion.

Figure 9 .

 9 Figure 9. Effect of the equation of state on the velocity (a), density (b), temperature (c) and dodecane mass fraction (d).

Figure 10 .

 10 Figure 10. Compressibility factor distribution.

Acknowledgments

This work was supported by the ESA-ESTEC, Contract no.: 3-12861/09/NL/PA. The authors would like to sincerely thank F. Travajno and B. Le Naour from MBDA-France, G. Fau and P. Gillard from PRISME, D. Blanc and D. Courilleau from IUT Bourges, D. Davidenko from CNRS-ICARE Orléans, A. Ingenito from the University of Rome La Sapienza (Italy) and all the students involved in the projects head by I. Fedioun for their technical and scientific support involving the computations and the experiments.

Pyrolysis in Porous Media:

Numerical model and parametric study.