Noureddine Lebaz

Arnaud Cockx

Mathieu Spérandio

Jérôme Morchain
email: jerome.morchain@insa-toulouse.fr

POPULATION BALANCE APPROACH FOR THE MODELING OF ENZYMATIC HYDROLYSIS OF CELLULOSE

In this numerical work, a population balance based model is proposed in order to describe the cellulose particles size evolution during the enzymatic hydrolysis. Two kinds of actions are considered: endoglucanase activity that cleaves randomly β-1,4-glycosidic linkages of cellulose and exoglucanase activity which reduces the particles size with chain-end-cleaving producing cellobiose (a dimer of two glucoses linked by a β-1,4-glycosidic bond). A discretization method with a fixed pivot technique is used for the endoglucanase action and a moving pivot technique for exoglucanase attack. The numerical resolution is then validated by analytical solutions available in literature. Afterwards, the combination of the two actions is investigated for different enzyme ratios in order to reproduce the endo-exo synergism numerically. Since the biodegradation of cellulose releases D-glucose as a final product due to β-glucosidase which hydrolyzes cellobiose into two molecules of glucose, numerical kinetic model predicting the fractional conversion of cellulose is derived from the population balance developed model. The enzymes activity is strongly affected by the accumulation of the end-products (cellobiose and glucose) during the hydrolysis, the inhibition effect is thereby incorporated in the model. The numerical model prediction is compared to experimental data in the case of combined activity and shows a promising approach for the modeling of cellulose-cellulase systems.

Introduction

The use of lignocellulosic materials, known as the largest renewable carbohydrate source on earth, [START_REF] Jørgensen | Biofuels Bioprod[END_REF] for biofuels and chemicals production has been widely studied since the global supply of fossil fuels is limited. [2][3][4] The focus has been on the bioconversion of biomass by means of hydrolytic enzymes capable of degrading these insoluble materials to soluble and fermentable sugars. However, both the enzyme systems and the understanding of cell wall structure are still subject to improvement. [START_REF] Kristensen | Enzymatic Hydrolysis of Lignocellulose[END_REF] Indeed, the hydrolysis reaction, considered as one of the main limiting steps of the bioconversion process, [START_REF] Chauve | [END_REF] has a multidimensional mechanism. It depends on the chemical and structural substrate features (chemical composition especially lignin content, accessibility, reactive area …) [7,8] as well as on the biocatalyst (nature and enzymatic cocktail composition) [9,10] and on the operating conditions (temperature, pH and mixing). [11,12] In addition to that, one has to consider complex substrate-biocatalyst interactions which remain to be elucidated. [13] This complexity makes the kinetic modeling of such systems delicate and the experimental determination of the role of each parameter on the whole reaction unobvious.

Nevertheless, several kinetic models based on different assumptions have been developed for the enzymatic hydrolysis of pure cellulose and lignocellulosic materials, reviewed by Bansal et al. [14] Michaelis-Menten based models, developed initially for homogenous reaction (such as the action of -glucosidases on soluble cellobiose) have been used by many authors to describe the whole hydrolysis process. [15,16] Since hydrolysis is a heterogeneous reaction, the overall rate depends on the amount of adsorbed cellulases. Adsorption models based generally on Langmuir adsorption isotherm are, thereby, used for fitting experimental data. [17,18] In these models it is assumed that the limiting step is related to the enzyme adsorption. Irreversible adsorption on lignin has been reported for lignocellulosic materials [19,20] as well as nonreactive adsorption on crystalline cellulose. [21,22] It is important to observe that Michaelis-Menten and adsorption models share the same mathematical formalism. The parameters have different meaning but they are only apparent and global constants. For example, the fact that the specific area changes during the reaction because of the substrate liquefaction cannot be taken into account by these models. [23] More comprehensive models have been proposed integrating the rate limiting factors such as the synergistic interaction of cellulases, the activity inhibition by products accumulation and the enzyme accessibility to the hydrolysable fractions. [18,24] But they are based on the same approach as those described before. Movagarnejad et al. combined Langmuir isotherm concept and shrinking particle theory [25] assuming that the particles are cylindrical and preserve their shape during the shrinkage process with an invariant diameter to height ratio.

The weakness of such models lies in their inability to integrate heterogeneous particles size distribution with different shapes.

In the light of these considerations, population balance modeling can be an alternative for developing mechanistic models for systems as complex as lignocellulosic materials hydrolysis. In fact, population balance equations, used for diverse applications, draw strength from their ability to describe dynamic systems evolution taking into consideration particular level processes (breakage, aggregation ...). Except Hosseini and Shah, [26,27] no literature could be found using this approach for enzymatic lignocellulosic materials hydrolysis. In their discrete models, Hosseini et al. assimilated the cellulose hydrolysis to a linear chain depolymerisation. But for cellulose hydrolysis by endoglucanase, the discretization scans all possible polymer lengths leading to an onerous resolution since the number of equations is proportional to the maximum initial chain length.

In this paper, we develop models for cellulose particles size evolution during their hydrolysis by endoglucanase and exoglucanase enzymes. For the resolution of the population balance equations, we use the discretization method developed by Kumar and Ramkrishna [29][30][31] with a fixed pivot technique for the endoglucanase activity and a moving pivot technique for the exoglucanase attack. Once these two models are validated, a combined model for both activities is proposed and solved using the moving pivot technique. This model predicts the particles size distribution during hydrolysis as well the conversion yield of the substrate. The originality of the work proposed here is to incorporate enzyme activity mechanisms into the robust framework of population balance modeling. Finally, the extension of the developed model in order to take into account the end-products inhibition effect is explored.

Model framework and numerical method

Population balance general equations describe particulate systems behavior when submitted to a dynamic evolution of their physical or chemical properties. In our case, we consider pure cellulose particles assimilated to polymer chains with a given initial size distribution. Enzymatic hydrolysis of this substrate leads to a progressive decrease in the chains length due to breakage and chain-end-cleaving processes. As a first step, each enzymatic activity is treated separately and the associated numerical method is presented. Then, the full model combining the two activities is described.

For a continuous and homogeneous system, the evolution of the particles size distribution undergoing binary break-up and constant growth is described by the kinetic equation: [30]  

' ) ' ( ) ' ( ) ' , ( ) ( ) ( ) ( ) ( ) ( dx x n x x x x x n x x n x G t x n x              (1)
The first term on the right hand side corresponds to the net accumulation related to the growth of particles which is a kinetic phenomenon. The second term represents the loss (death) of particles with a size x due to their break-up (with a frequency 𝛤(x)). The last term describes the production of new particles with a size x (birth) due to the breakage of larger particles.

The birth and death terms are discrete events. The last integral term represents the main issue in solving the equation (1). Except in some particular cases for which analytical solutions exist, population balance equations resolution is carried out by numerical techniques. Here, we chose the discretization method (method of classes) developed by Kumar and Ramkrishna. [28]

Discretization method

The discretization method consists of dividing the entire property range (size, volume, specific area …) into a finite number of intervals, noted M. All particles for which the property lies in the interval [i ; i+di] are attributed to the i th interval represented by its center xi called pivot. Thus, the initial problem is transformed into a set of M equations, regarding the number of particles, Ni(t), in the i th interval. The integral term in each equation is approximated numerically using the same discretization.

Modeling the endoglucanase activity

The hydrolysis due to endoglucanases solely is assimilated to a pure breakage process (G=0 in equation ( 1)). Endoglucanases are known for internally cleaving cellulose chains (figure 1).

Their random actions on β-1,4-glycosidic bonds lead to binary breakages without any position preference but as admitted by Hosseini and Shah, [26] the breakup frequency is a function of chain length since longer chains are more susceptible to be attacked. Thus, as a first approach, we consider the breakage frequency 𝛤 proportional to cellulose chains length.

As the polymers breakage is a discrete process, it has a probabilistic foundation. Considering a cellulose chain with an initial length x', we look for the probability (x,x') to form a new chain with a length x after the breakup of the chain x'. Since the breakage is binary, there are only two possible locations along the initial chain length x' on which the enzymes can act in order to form the chain x (figure 1). Thus, the probability expression is 2/x'. [31] During the population particles evolution, new particles are formed by breakup. The size of the newly formed particle may not correspond to an existing pivot. So this particle has to be affected to the neighboring pivots. Considering a new particle with a size  in the range [xi, xi+1], one has to calculate the fractions a(, xi) and b(, xi+1) that will be affected to the pivots xi and xi+1. Following the approach proposed by Kumar and Ramkrishna, we chose to preserve the particle number and mass, noted f1 and f2, which lead to equations ( 2) and (3). [28] )

( ) ( ) , ( ) ( ) , ( 1 1 1 1 1    f x f x b x f x a i i i i     (2) ) ( ) ( ) , ( ) ( ) , ( 2 1 2 1 2    f x f x b x f x a i i i i     (3)
The breakage process is defined as a continuous event over a given cellulosic chain length. It could be thought that this is in contradiction with the fact that breakup actually occur at the glycosidic linkages. In fact, the new formed chains appearing in the system are distributed on two adjacent pivots, as described previously. This solves the apparent contradiction since the pivots correspond to polymer lengths which are physically possible. given nature (binary in this scheme) and breakage frequency . The arrows show the two possible locations to have a particle of length x after the breakage of x'.

After several mathematical transformations, Kumar and Ramkrishna establish the discrete equation for the pure breakup process as the linear sum of the birth and death particles terms given below: [28] 

      M i k k k ik i i i t N n t N dt t dN ) ( ) ( ) ( (4)
The left hand side term of the equation ( 4) corresponds to the variation of the total number of the particles in the i th interval (1 ≤ i ≤ M) per unit time. This variation is the consequence of two phenomena appearing on the right hand side of the equation. The first term refer to the particles in the i th interval disappearing at the frequency i because of hydrolysis. It corresponds to the number of the particles of length xi undergoing breakage process per unit time, i.e. the second term on the right hand side of equation (1). The second phenomenon is related to the particles of length xk ≥ xi which may give, after breakage, at least one particle in the i th interval. The number of the new particles appearing in the i th interval per unit time is given by the second term of the equation ( 4) corresponding to the discretization of the integral term of the equation (1). For further details regarding the calculation of nik, please refer to the original work cited. It is worth noting that these nik values can be calculated a priori, i.e. before starting the resolution of the set of equations ( 4). Indeed, the nik values depend only on the properties f1 and f2 and the value of the pivots.

The numerical resolution of the discretized population balance equation was performed using the fixed pivot technique and was validated by analytical solutions available in literature. Ziff and McGrady proposed analytical solutions in the cases where the breakage frequency is (𝛤=1) or (𝛤=x α ). [31] Once the model is validated, other breakage frequency mathematical expressions can be easily implemented.

Modeling the exoglucanase activity

The second major hydrolytic enzyme is exoglucanase, called also cellobiohydrolase. It acts on the ends of cellulose chains and releases cellobiose (figure 2). This action causes a progressive decrease of the cellulose chains length. Unlike endoglucanases requiring probabilities for predicting daughter particles length, exoglucanase activity leads invariably to a chain with a degree of polymerization decreased by two after each attack. An equal reactivity, irrespective of the chain length, is assumed. We propose to model this hydrolytic activity as a pure growth [30] with a negative rate defined as the rate of decrease of the chains length (G(x)=dx/dt). The hydrolysis by exoglucanases is thus described as a continuous dissolution of the cellulose polymers. The equation to solve is obtained by setting =0 in equation (1). In that case, the moving pivot technique [29] is more appropriate for resolving the discretized population balance equation. The gradual decrease of the cellulose chains length can be related to the pivot moving rate, thus the technique reproduces indirectly the enzymatic hydrolysis rate. Generally, the moving pivot rate is either constant (G(x)= σ0) or depending on the particles size (G(x)=σ0x). In the two cases, analytical solutions exist [30] and are used in this work to validate the numerical results.

Compared with similar growth systems available in literature, [33][34] it should be noted that in our case, in view of negative moving rate (rate of the chains dissolution by the exoglucanase activity), we must avoid switching to negative pivots which is nonsense. Note that the first pivot corresponds to a chain length of 2 monomers of glucose (cellobiose). Considering this, we set the first pivot equal to the smallest chain length in the system (cellobiose). This pivot is not moving. The other pivots move to the left at a constant rate, this means that all the cellulose chains (except cellobiose) are attacked at a same frequency. The first pivot being fixed, knowing the initial chain lengths and the moving rate, the time necessary for the second pivot to reach the first one is calculated as:

) ( ) ( 1 2 x G x t x    (5)
The PBE is integrated over this time interval. At the end of the time interval, the second pivot is eliminated and its final number density or total particles number is assigned to the first pivot. The procedure is repeated until the final simulation time is reached.

Since the pivots level is a chain length, their moving to smaller chain lengths induces weight loss corresponding to the produced cellobiose which can be easily evaluated. The total mass of the polymers is conserved through the process. The first moment of the distribution multiplied by the molar mass of the monomer (glucose in the present case) corresponds to the total mass of all chains. Given the initial mass, mT, is known and that the PBM leads to the chain length distribution at any time, the difference between the initial mass and the mass of all chains longer than 2 corresponds to the mass, mc, that was transformed into cellobiose.

Then it comes that

T x Glu T T c m dx t x xN M m m t m     2 ) , ( ) ( N (6)
Where the integral term corresponds to the total mass of the chains larger than cellobiose.

Being capable of assessing the amount of the cellobiose produced during the process by means of the equation ( 6), one can deduce at any given time the fraction of converted substrate.

Combined model and synergistic interaction of cellulases

The discretized population balance equation ( 1) for a combined endoglucanase-exoglucanase activity is resolved with the moving pivot technique. This technique is chosen because it is judicious to treat the pure growth (exoglucanase activity) part of the equation (1). However since the pivot are now moving, the nik values involved in the breakage part of the equation are no longer constant and their values have to be updated throughout the solution procedure as proposed by Kumar and Ramkrishna. [30] The resolution of the discretized population balance equation in the different considered cases was performed with Matlab® software.

Results and discussion

Endoglucanase activity

As discussed previously, the discretized population balance equation in the case of the endoglucanase activity is resolved with a fixed pivot technique. The initial distribution used in order to validate the numerical approach is monodisperse (all the chains have the same initial length). The breakage frequency depends on the particles size (𝛤=x) and the analytical solution has been discussed by Ziff and McGrady. [31] The Figure 3a represents a comparison between analytical and numerical number distibutions during the hydrolysis reaction by endoglucanase. The Figure 3b represents the effectiveness of the breakage in terms of new particles formed in the system during the reaction. The numerical results are very close to the analytical solutions in the two cases, the error is insignificant at this level. The numerical approach is, thereby, validated.

Exoglucanase activity

The initial pivots grid is defined according to an arithmetic sequence where the first pivot corresponding to the elementary particles produced by the exoglucanase chain-end-cleaving is fixed, when the others move at a constant rate. The cellobiose as the smallest and indestructible particle in such system is chosen for the fixed pivot.

As expected in the case of pure growth, the Figure 4a shows a constant translation without any deformation of the initial distribution. The number of chains represented by the pivots weight remains unchanged while all the cellobiose produced during the reaction accumulates on the first pivot as we can see in the Figure 4b. It can be noticed that the cellobiose production rate is constant until the left side of the distribution reaches the first pivot. From that point onward, the number of chains decreases meaning that the rate of cellobiose production becomes limited by the amount of available hydrolysable chains. At the end of the process, it can be seen that the mass conservation is satisfied since a 100% conversion is reached. This mass conservation after total hydrolysis validates the accuracy of the numerical procedure. As we can see in the Figure 5a, the combined endo/exoglucanase activity results in a progressive deformation of the initial number distribution. In fact, the particles size is reduced as a consequence of the exoglucanase action while new particles are produced by the chains cleavage leading to an increase in the number of small particles. Here again the mass balance is preserved along the conversion as illustrated in Figure 5b. In addition to that, if we assume that all the produced cellobiose is converted into glucose (disregarding the inhibition effect and the concentration of the -glucosidases), the Figure 5b is also representative of the cellulosic fraction converted into glucose. Knowing the initial amount of substrate, the produced glucose can be assessed at any given time. Generally, the elementary mechanisms of both cellulolytic activities can be recognized. In order to test the robustness of the model, we have compared the results obtained with various values for  and G which is analogous to testing different enzymes ratios (see the Table 1). In the Figure 6a, we report the results obtained when only the endoglucanase activity increases (which can be experimentally correlated with the enzyme concentration). In the Figure 6b, the endoglucanase activity remains constant while that of exoglucanase changes. The increase on the activity of each enzyme leads to an enhancement of the cellulose hydrolysis rate, the maximum conversion is reached more rapidly. These results are consistent and show that the model offers an assessment of the contribution of each enzyme activity.

Furthermore, the results obtained by the combined model have been compared to those given by the two models developed initialy for the two cellulolytic activities as shown by the Figure 7. The comparison between the sum of conversions obtained by the endoglucanase and the exoglucanase activities separately and that given by the combined model respecting the same rates for the two activities shows a difference due to the endo-exo synergism. These results confirm that the global model can predict and assess the synergy effect numerically according to the initialy imposed enzymes activity. In order to test the applicability of the developed model to fit experimental data in term of fractional substrate conversion, we refer to the experimental work of Woodward et al. as an example. [34] In their experiments, they have studied the enzymatic hydrolysis of microcrystalline cellulose (Avicel PH-105) using a mixture of endoglucanase/exoglucanase from Trichoderma reesei and -glucosidase from Aspergillus niger. The initial average degree of polymerisation of Avicel PH-105 is 230. [35] We assume a normal initial size distribution of the polymers with a peak at 230.

Model parameters (breakage frequency and moving rate) have been estimated assuming a maximal conversion of 40% of initial substrate mass. As shown in the Figure 8, the model prediction is very close to the experimental results. It must be highlighted that, due to the recalcitrant nature of microcrystalline cellulose (high crystallinity index) and other parameters affecting the reaction rate, incomplete conversion is observed.

Figure 8: Substrate conversion of Avicel PH-105: (-) model prediction, (*) experimental results [35] (=0.5, G=-48).

Discussion

The developed models for each of the two main cellulolytic activities show predictable results, successfully validated against analytical solutions. This section mainly focused on the results obtained when the two enzymatic activities are considered. Regarding the substrate conversion, typical enzymatic hydrolysis rate behavior characterized by a rapid initial rate and a progressive slowdown is reproduced by the model. In our case, the decrease of the overall hydrolysis rate is due to the substrate limitation (the number of hydrolysable chains finally goes to zero), and reaction stops when the converted fraction reaches 100%.

As far as the number of chains is concerned, the total particles number increases because of the appearance of smaller chains in the system due to the endoglucanase activity that cleaves randomly the cellulose chains. The exoglucanase activity leads to a progressive dissolution of the chains. The associated particles mass loss produces cellobiose as a final product. The model gives access to the changes in the chains length distribution. The evolution of this distribution is markedly different depending on the magnitude of each enzymatic activity considered (Table 1). The increase in chain number is related to endoglucanases whereas the production of cellobiose is mostly due to the action of exoglucanases that leave the number of chains unchanged.

It is remarkable that the substrate conversion curves (figures 6a and 6b) are very similar to that observed in most experiments. In that sense the numerical results could be fitted with a kinetic model (Michaelis-Menten model) to identify the apparent constants. However, we can claim that the apparent values of the parameters Vmax and Km would be dependent on the actual values of  and G. It is significant the Michaelis-Menten like shape of the conversion curve does not come from the combination of multiple enzymatic reaction models: and G are true constants.

The numerical parameters  and G are considered as constants in this first step of the model development. It is proposed here to relate them to the enzymatic activities. On one hand, the break-up frequency (Γ) expressed as the number of chains breakages by time unit is easily turned into an enzymatic activity as defined in literature by a simple units conversion (breakages number to molar amount of reduced sugars).

N   U (6)
Where N is the Avogadro number and U the activity in mole of product (reduced sugars in our case) formed by unit time.

On the other hand, the amount of produced cellobiose can be deduced from the constant rate (G) used for modeling the exoglucanase activity, as follows.

During the hydrolysis by exoglucanases, the number of the chains is constant. Assuming that the chains are polymers of a known monomer (ex. glucose) with an initial total mass mT and that the degree of polymerization (DP) of the substrate is known, the initial number of chains is given by:

  DP M m N Glu T  0 (7)
If one assumes that the amount of enzyme is sufficient to attack all chains, the amount of cellobiose formed by time unit due to the action of exoglucanases is:

  Cello c M N G dt dm N 0 2  (8)
Recall that G=dx/dt, that is to say the rate of chain length reduction expressed in monomers.

The coefficient 2 indicates that two monomers are needed to form one cellobiose molecule.

Then G/2 times N(0) represents the number of molecules of cellobiose formed by unit time.

  N 0 2 1 N G dt dm M U c Cello   (9)
Knowing the cellobiose releasing rate, the exoglucanase activity (U) can be easily deduced.

When the three major hydrolytic enzyme activities (endoglucanase, exoglucanase and glucosidase) are combined, it is demonstrated (see Figure 7) that the hydrolysis rate exceeds the sum of the rates of the single enzymes showing a synergistic interaction of the cellulases system. [36] Different synergisms depending on the nature of the cellulase system have been reported, [9,38] the main one is between endoglucanase and exoglucanase called endo-exo synergism. [10,39] The endoglucanase internal chain-cleaving makes available additional chain ends for the action of exoglucanases, that's why, when the reaction is not limited by the concentration of exoglucanases, the hydrolysis rate increases.

One of the critically important barriers affecting the enzymatic hydrolysis rate is the endproducts inhibition. The accumulation of cellobiose and glucose is partially responsible for the deceleration of the cellulases activity even at low degree of conversion of the substrate. This means that the effectiveness of the hydrolytic enzymes falls as the end-products concentration increases in the reactor. Incorporating this phenomenon in the population balance model is quite straightforward. The model parameters characterizing the two main activities (breakage frequency and growth rate) have to be related instantaneously to the endproducts concentration. We give in the equations (11) and (12) the kinetic rate expressions used in the case of non-competitive inhibition. The degree of inhibition in the system is controlled numerically by the constants K and KG. The specific constants K and KG can be experimentally determined. Here, they are set to arbitrary values.  and G are no longer constants, they are inversely proportional to the end-products concentration.

      K C K C S s 0 ) ( (10) G S G s K C K G C G   0 ) ( (11)
In the Figure 9, we compare the degree of conversion of the same initial chains population with and without taking into account the inhibition effect. These numerical results are obtained by the population balance model. the availability of the substrate leading finally to a total conversion, the introduction of the inhibition effect reflects more the reality, the substrate conversion rate falls drastically thus, the substrate is partially converted.

The results obtained in this work in terms of the substrate particles size distribution evolution during enzymatic hydrolysis are qualitatively similar to those obtained experimentally by Mooney et al. [39] who measured the fiber length distribution of Douglas-fir kraft pulp during enzymatic hydrolysis (figure 10). A direct comparison between our model and these data was not attempted since the data correspond to fibers. Indeed, in the latter case, one has to think about how to express the effectiveness of each enzyme activity knowing that during the reaction, the two soluble/insoluble substrate systems coexist. Nevertheless, the similarity is encouraging and this tends to demonstrate that the numerical approach proposed here can be generalized to insoluble particles. from hydrolysis [40] .

If one assumes that the -glucosidase activity is such that the conversion of cellobiose into glucose is instantaneous, then the glucose concentration should evolve as the cellobiose concentration. But if one considers also the inhibitory effect on -glucosidase, the whole degradation kinetic will become dependent on the three enzymatic activities. These inhibitory effects have been identified [40] and can be further incorporated in the model through the introduction of a kinetic rate expression for the -glucosidase.

Experimental data, [34] as those of major works, [42,43] show partial substrate conversion because of its decreasing reactivity due to several parameters: changes in cristallinity index, [43] accessible surface area [8] and chemical composition. [START_REF] Zhao | Biofuels Bioprod[END_REF] Incorporating these factors and their interactive effects in the model requires the development of the appropriate parameters expressions (this is done in the case of end-products inhibition effect) or the addition of other dimensions to the existing model.

Conclusion

A population balance based kinetic model has been developed for the enzymatic hydrolysis of cellulosic materials. The pure cellulose has been assimilated to polymer chains with given initial length distribution leading, in the case of enzymatic hydrolysis, to a depolymerisation process. Firstly, the two main activities have been modeled separately in order to validate the numerical approach. A discretization method has been used for the resolution with a fixed pivot technique for the endoglucanase activity and a moving pivot technique for the exoglucanase action. Afterwards, a combined model for endoglucanase-exoglucanase activity has been proposed. Under some assumption, the kinetics of final product release are also readily accessible through this model. This model has been tested for different enzyme ratios and shows coherent results. In addition, the comparison between the conversion obtained by the global model and those given by the two initial models separately in the same conditions shows a difference attributed to the endo-exo synergism. Since the enzymatic activity is inhibited by the products, this critically important phenomenon is incorporated in the model.

Finally, the prediction of the population balance model compares favorably with the experimental results available in the literature.

In this study, the formulation of a one-dimensional model is based only on the hydrolytic enzyme elementary mechanisms known in the literature disregarding other factors affecting cellulose hydrolysis rate. The principal aim was the validation of the population balance approach using simple mathematical breakage frequency and growth rate expressions.

According to the results, this approach presents a real potential for a better description of the enzymatic hydrolysis of such complex systems.

Figure 1 :

 1 Figure 1 : Schematic representation of the endoglucanases activity. Random scission with a given nature (binary in this scheme) and breakage frequency . The arrows show the two possible locations to have a particle of length x after the breakage of x'.

Figure 2 :

 2 Figure 2 : Schematic representation of the exoglucanases activity. Each attack produces one cellobiose and reduces the length of the initial chain by two monomers.

Figure 3a :

 3a Figure 3a: Number distribution evolution: comparison between analytical (-) and numerical results (*). The initial distribution is monodisperse and undertakes pure breakup with a breakage frequency =x.

Figure 3b :

 3b Figure 3b: Total particle number evolution under the endoglucanase activity, (-) analytical solution (*) model simulation.

Figure 4a :

 4a Figure 4a: Number distribution evolution under exoglucanase activity (G=-5)

Figure 5a :

 5a Figure 5a: Number distribution evolution under combined endo/exo activities (=0.25, G=-10).

Figure 5b :

 5b Figure 5b: Simulated substrate conversion during the hydrolysis considering endo and exo activities. The mass balance is conserved.

Figure 6a :

 6a Figure 6a: Substrate conversion evolution under different endoglucanase rates Γ.

Figure 6b :

 6b Figure 6b: Substrate conversion evolution under different exoglucanase rates G.

Figure 7 :

 7 Figure 7: Substrate conversion during hydrolysis ( … ) endo+exo separately, (-) endo and exo in cocktail (=0.5, G=-10). The arrow quantifies the gain due to the synergistic effect of the two cellulases.

Figure 9 :

 9 Figure 9: The inhibition effect. The cellulases activities slow down with conversion leading to a partial substrate conversion. Two cases with different inhibition constants are shown (=0.25, G(0)=-10). The implementation of the new parameters expressions in the population balance model gives interesting results. Unlike the previous cases where the enzymes activities are only limited by

Figure 10 :

 10 Figure 10: Changes in population fiber-length distribution of Douglas fir kraft pulp resulting

mT

  Number density of chains with a length x at time t G(x)Growth rate for chains of length x [monomer.t-1 ](x)Breakage frequency for a chain of length x [t -1 ] k Defined as (xk) (x, x')dx' Number of daughter chains formed in length range x to x+dx due to the breakage of chains of length x'  Property of the population a(x, xi) Number of chains of length xi assigned when a chain of length x, such that xi ≤ x ≤ xi+1 , is formed b(x, xi+1) Number of chains of length xi+1 assigned when a chain of length x, such that xi ≤ x ≤ xi+1 , is formed f1, f2 Properties (number and mass) which have to be preserved N(t) Total number of chains at time t Ni(t) Total number of chains in the i th interval at time t M Total number of classes (chains length ranges) used to represent the total chains population xi Representative length for the i th length range (pivot) nik Total number of chains assigned to representative length xi when a chain of length xk breaks  Characteristic time defined as the necessary time for the second pivot to reach Total chains mass at time t=0

Table 1 :

 1 Reaction time prediction (hours) for a conversion of 99% of the substrate under endo/exoglucanase activities.

	G	Γ	0,01	0,5	1	2
	0,5		168	30	16	11
	5		18,4	12	8,4	5,6
	10		9,1	7,2	5,8	4,6
	15		6	5,2	4,5	3,5


F Olivier

N Rochet

S Delpeux-Ouldriane

J Chancolon

V Sarou-Kanian

F Fayon

S Bonnamy

Strontium incorporation into biomimetic carbonated calcium-deficient hydroxyapatite coated carbon cloth: biocompatibility with human primary osteoblasts

Keywords: Activated carbon fiber cloth, Calcium-deficient hydroxyapatite, Strontium substitution, Sono-electrodeposition, Human primary osteoblast, Cell proliferation

It has already been shown that sono-electrodeposition can be used to coat activated carbon fiber cloth (ACC) with calcium phosphates (CaP) and we recently demonstrated that cathodic polarization at -1 V / Hg/Hg2SO4 was the best parameter to obtain a carbonated calcium deficient hydroxyapatite (CDA) coating with optimal uniformity and homogeneity. In the present study, we investigated whether this technique was suitable to dope this carbonated CDA coating by partial substitution with another bivalent cation such as strontium. We show here that a strontium-substituted carbonated CDA coating can be produced and quantitatively controlled up to at least 10 at.%. In this range we demonstrate that the presence of strontium does not modify either the textural or the structural properties of the carbonated CDA. Owing to the well-known effect of both carbonated CDA and strontium in bone formation, the biocompatibility of ACC coated or not with carbonated CDA or with strontium substituted carbonated CDA was tested using primary human osteoblasts. Our data revealed a positive and dose-dependent effect of strontium addition on osteoblast activity and proliferation. In conclusion, we show here that electrodeposition at -1 V is a suitable and easy process to incorporate cations of biological interest into CaP coating.

Introduction

The structure of the bone mineral phase is similar to that of nanocrystalline calcium phosphate (CaP) apatites [START_REF] Gómez-Morales | Progress on the preparation of nanocrystalline apatites and surface characterization: Overview of fundamental and applied aspects[END_REF], partially enriched with a large variety of substitution ions such as Na, K, Mg, Sr, Cl, F, Zn. These trace elements, and specifically the strontium (Sr), play a role in bone mineral phase remodeling [START_REF] Quarles | Cation Sensing Receptors in Bone: A Novel Paradigm for Regulating Bone Remodeling?[END_REF]. Therefore, it has been proposed to dope synthetic apatitic bone substitutes used in orthopedic surgery with small amounts of such ions.

Moreover, it has been shown that incorporation of Sr 2+ , Mg 2+ , F -, Cl -or CO3 2-in the crystal structure of apatites improves their biological behavior [START_REF] Laurencin | Magnesium incorporation into hydroxyapatite[END_REF][START_REF] Lijuan | Effect of different synthesis conditions on the microstructure, crystallinity and solubility of Mg-substituted hydroxyapatite nanopowder[END_REF][START_REF] Cox | Preparation and characterisation of nanophase Sr, Mg, and Zn substituted hydroxyapatite by aqueous precipitation[END_REF][START_REF] Abert | Wet chemical synthesis of strontium-substituted hydroxyapatite and its influence on the mechanical and biological properties[END_REF]. For examples, these ionic incorporation have a potential role in cell adhesion, biomineralization, bone metabolism and resorption processes [START_REF] Šupová | Substituted hydroxyapatites for biomedical applications: A review[END_REF].

Interestingly, calcium ion (Ca), with an ionic radius of 1.00 Å, can be easily substituted in the apatitic structure by strontium, which has an ionic radius of 1.18 Å [START_REF] Frasnelli | Synthesis and characterization of strontium-substituted hydroxyapatite nanoparticles for bone regeneration[END_REF][START_REF] Shannon | Revised Effective Ionic Radii and Systematic Studies of Interatomie Distances in Halides and Chaleogenides[END_REF][START_REF] Li | Chemical composition, crystal size and lattice structural changes after incorporation of strontium into biomimetic apatite[END_REF][START_REF] Drouet | Surface enrichment of biomimetic apatites with biologically-active ions Mg2+ and Sr2+: A preamble to the activation of bone repair materials[END_REF]. Strontium appears as an ion of choice to impart new properties to CaP since many studies have shown that this cation promotes bone tissue growth and decreases bone resorption without altering neither the mineralization process nor its chemical composition [START_REF] Boanini | Ionic substitutions in calcium phosphates synthesized at low temperature[END_REF][START_REF] Leilei | Double-layer TC4/Sr substituted hydroxyapatite bioactive coating for carbon/carbon composites[END_REF][START_REF] Querido | The effects of strontium on bone mineral: A review on current knowledge and microanalytical approaches[END_REF]. The degree of substitution depends mainly on the initial concentration of strontium (precursor solution) from 0 to 100 % [START_REF] Boanini | Ionic substitutions in calcium phosphates synthesized at low temperature[END_REF]. Moreover, strontium incorporation into bioceramics between 3 to 7 at.% stimulates preosteoblast activity and differentiation, inhibits osteoclast activity and differentiation and is used in osteoporosis treatment [START_REF] Capuccini | Interaction of Srdoped hydroxyapatite nanocrystals with osteoclast and osteoblast-like cells[END_REF][START_REF] Capuccini | Strontium-substituted hydroxyapatite coatings synthesized by pulsed-laser deposition: In vitro osteoblast and osteoclast response[END_REF][START_REF] Panzavolta | Strontium-Substituted Hydroxyapatite-Gelatin Biomimetic Scaffolds Modulate Bone Cell Response[END_REF]. However, it has being shown that strontium substitution has a dose-dependent impact on bone formation [START_REF] Boivin | Strontium distribution and interactions with bone mineral in monkey iliac bone after strontium salt (S 12911) administration[END_REF][START_REF] Verberckmoes | Dose-dependent effects of strontium on osteoblast function and mineralization[END_REF]. A recent study has showed that 10 at.% of strontium corresponds to the critical substitution rate into the biomimetic apatite [START_REF] Bussola Tovani | Formation of stable strontium-rich amorphous calcium phosphate: Possible effects on bone mineral[END_REF]. These authors showed that using higher amount leads to the formation of a stable strontium-apatite phase that may induce high-dose related pathologies.

Furthermore, it has been shown that incorporation of strontium has an influence on the mechanical properties of apatite. Increases in compressive strength and hardness were measured for a strontium substitution of less than 10 at.% [START_REF] Landi | Sr-substituted hydroxyapatites for osteoporotic bone replacement[END_REF][START_REF] Shen | Synthesis of strontium substituted hydroxyapatite whiskers used as bioactive and mechanical reinforcement material[END_REF].

As for the synthesis of doped CaP, several methods can be used such as sol-gel [START_REF] Renaudin | Effect of strontium substitution on the composition and microstructure of sol-gel derived calcium phosphates[END_REF][START_REF] Qaisar | Sol-gel synthesis and TEM-EDX characterisation of hydroxyapatite nanoscale powders modified by Mg, Sr or Ti[END_REF], precipitation [START_REF] Cox | Preparation and characterisation of nanophase Sr, Mg, and Zn substituted hydroxyapatite by aqueous precipitation[END_REF][START_REF] O'donnell | Structural analysis of a series of strontium-substituted apatites[END_REF], hydrothermal [START_REF] Lin | Strontium substituted hydroxyapatite porous microspheres: Surfactant-free hydrothermal synthesis, enhanced biological response and sustained drug release[END_REF], micro-arc [START_REF] Chung | Systematic strontium substitution in hydroxyapatite coatings on titanium via micro-arc treatment and their osteoblast/osteoclast responses[END_REF], pulsed-laser [START_REF] Capuccini | Strontium-substituted hydroxyapatite coatings synthesized by pulsed-laser deposition: In vitro osteoblast and osteoclast response[END_REF][START_REF] Pereiro | Pulsed laser deposition of strontium-substituted hydroxyapatite coatings[END_REF] and pulsed electro-deposition [START_REF] Huang | Corrosion behavior and biocompatibility of strontium and fluorine co-doped electrodeposited hydroxyapatite coatings[END_REF][START_REF] Drevet | Pulsed electrodeposition for the synthesis of strontiumsubstituted calcium phosphate coatings with improved dissolution properties[END_REF]. Several of these processes have thus been proposed for prosthesis coating with a doped CaP layer. The major use of CaP materials in the orthopedic field consists in prosthesis coating to improve their integration in the host body [START_REF] Røkkum | Total hip replacement with an entirely hydroxyapatite-coated prosthesis: 5 years' follow-up of 94 consecutive hips[END_REF][START_REF] Rogers | The ABG hydroxyapatite-coated hip prosthesis: One hundred consecutive operations with average 6-year follow-up1[END_REF]. Among these techniques, electrodeposition has been used to perform CaP coatings or doped CaP coatings on conducting metallic alloys [START_REF] Drevet | Pulsed electrodeposition for the synthesis of strontiumsubstituted calcium phosphate coatings with improved dissolution properties[END_REF]. Electrodeposition is a soft chemical process, performed at room temperature, in a water-based electrolyte. It consists in a controlled co-precipitation of calcium ions and/or other divalent cations and phosphate ions on the surface of a conducting electrode. In most cases, electrodeposition is used to perform a CaP coating on rigid metallic substrates such as stainless steel [START_REF] Gopi | A comparative study on the direct and pulsed current electrodeposition of hydroxyapatite coatings on surgical grade stainless steel[END_REF], titanium [START_REF] Dumelie | In vitro precipitation of electrodeposited calcium-deficient hydroxyapatite coatings on Ti6Al4V substrate[END_REF][START_REF] Eliaz | Electrochemical processes of nucleation and growth of hydroxyapatite on titanium supported by real-time electrochemical atomic force microscopy[END_REF], magnesium [START_REF] Walter | Influence of the cathodic activity of magnesium alloys on the electrochemical deposition of calcium phosphate[END_REF], and cobalt [START_REF] Yan | Effect of metallic nanoparticles on the biotribocorrosion behaviour of Metal-on-Metal hip prostheses, 17th[END_REF] alloys or flexible substrates such as carbon mesh [START_REF] Sandeman | The in vitro corneal biocompatibility of hydroxyapatite coated carbon mesh[END_REF].

Carbon mesh, carbon fibers and activated carbon fiber cloth (ACC) are biocompatible and bioinert materials [START_REF] Adams | Carbon fiber-reinforced carbon as a potential implant material[END_REF][START_REF] Minns | The repair of osteochondral defects in osteoarthritic rabbit knees by the use of carbon fibre[END_REF][START_REF] Minns | Mechanical and histological response of carbon fibre pads implanted in the rabbit patella[END_REF][START_REF] Mäkisalo | Tissue reactions around two alloplastic ligament substitute materials: experimental study on rats with carbon fibres and polypropylene[END_REF][START_REF] Carranza-Bencano | Carbon "ber implants in osteochondral defects of the rabbit patella[END_REF] that are often proposed not as a support but as a material designed to improve the mechanical properties of apatites for bone regeneration [START_REF] Slosarczyk | Hot pressed hydroxyapatite-carbon fibre composites[END_REF][START_REF] Wu | Biomimetic synthesis and characterization of carbon nanofiber/hydroxyapatite composite scaffolds[END_REF][START_REF] Wang | Design and fabrication of carbon fibers with needle-like nano-HA coating to reinforce granular nano-HA composites[END_REF][START_REF] Zhao | Controllable preparation of SiC coating protecting carbon fiber from oxidation damage during sintering process and SiC coated carbon fiber reinforced hydroxyapatite composites[END_REF][START_REF] Wang | Controllable preparation of a nano-hydroxyapatite coating on carbon fibers by electrochemical deposition and chemical treatment[END_REF][START_REF] Hong-Mei | Calcium phosphate sonoelectrodeposition on carbon fabrics and its effect on osteoblast cell viability in vitro[END_REF].

In a previous work, the sono-electrodeposition process was optimized for the coating of activated carbon cloth with calcium-deficient hydroxyapatite [START_REF] Olivier | Influence of electrochemical parameters on the characteristics of sono-electrodeposited calcium phosphate-coated carbon fiber cloth[END_REF]. We showed that using this process the thickness, morphology, microtexture, chemical composition and structure of the CaP coatings could be tailored. After optimization of the experimental conditions, carbonated calcium-deficient hydroxyapatite (CDA), which displays a homogeneous plate-like microtexture, was deposited on ACC. The aim of the present study was, firstly, to explore whether combining ACC and Sr-doped CaP coatings using sono-electrodeposition yielded a hybrid material. To achieve this, we investigated the effect of strontium addition (< 10 at.%) on the chemical, morphological and structural characteristics of the hydroxyapatite phases synthetized. Then, we investigated the biocompatibility and cell activity of different CaP coated materials (with and without strontium incorporation), in vitro, on human primary osteoblasts.

Materials and methods

Activated carbon fiber cloth

The substrate consisted of an activated carbon fiber cloth (ACC) referred to as FM50K. The FM50K prepared from a viscose precursor was supplied by the Zorflex® Company. Prior to its use, the ACC substrate was washed with boiling distilled water using a Soxhlet® extractor to remove any traces of contaminants coming from the activation and carbonization steps leading to its manufacturing. After 12 hours of washing, it was dried at 70 °C under vacuum.

Synthesis of the xSr-CaP coatings on activated carbon fiber cloth

The xSr-CaP coating (x being the strontium ratio in the electrolyte) on ACC was obtained by sono-electrodeposition under cathodic polarization. Different Sr/(Sr+Ca) ratios in the electrolyte were considered: 0, 5 and 10 at.% (samples respectively named 0Sr-CaP, 5Sr-CaP and 10Sr-CaP). An electrochemical cell was homemade to use ACC (FM50K) as the working electrode, a platinum basket surrounding the ACC was used as the counter electrode and Hg/Hg2SO4 (Radiometer analytical, France) was used as the reference electrode. The cathodic polarization of the carbon electrode was performed using a potentiostat/galvanostat (Biologic Science Instruments VMP-2, France)

with a constant potential of -1 V/Hg/Hg2SO4, which is the theoretical barrier potential of water reduction at initial pH = 4.8. Sonication (35 kHz, 30 %) was applied during polarization (ultrasonic bath, Transsonic Ti-H-10, Fisher Scientific France). Using this procedure, three materials were synthetized (materials 2 to 4), which are listed in Table 1. analyzer (Thermo Scientific France).

The porosity was characterized by N2 and CO2 adsorption respectively at 77 K and 273 K, using an Autosorb-1 (Quantachrome). The total pore volume was estimated from the amount adsorbed at P/P0 = 0.95. The micropore (pore size < 2nm) and mesopore (2 nm < pore size < 50 nm) volumes were determined by applying the density functional theory (DFT) method to the N2 isotherms assuming a slit pore shape. The Dubinin-Radushkevich (DR) theory was applied to the CO2 adsorption isotherm at 0 °C to estimate the volume of ultramicropores (pore size < 0.7 nm). The pKA distribution of functional groups grafted on the ACC surface and the pHPZC (pH value corresponding to an ACC net charge of zero) were measured by potentiometric titration with NaOH (0.1 mol/L), using a very low incremental volume (0.01 ml) in a wide range of pH = 3-11. The pore size distribution was measured by mercury intrusion porosimetry using a Poresizer 9320 porosimeter (Micromeritics). The mesopore and macropore sizes ranged from 4 nm to 350 µm.

In vitro biological procedures

The human primary osteoblasts were purchased from PromoCell (Heidelberg, Germany).

They were expanded on plastic culture flasks in complete medium constituted of in α-MEM medium (LONZA, Emerainville, France) supplemented with 10 % fetal calf serum (EUROBIO, Courtaboeuf, France), 100 U/ml penicillin, 100 U/ml streptomycin (Gibco Thermo Fisher Scientific, France), 2 % ultraglutamine (LONZA) and incubated at 37 °C under 5 % CO2. Four materials were tested, namely FM50K, FM50k+0Sr-CaP, FM50K+5Sr-CaP and FM50K+10Sr-CaP. Material patches of 8 mm diameter were pre-incubated in 48-well culture plates in completed medium for 24 h. Human primary osteoblasts were then seeded at a density of 7000 cells/well. After 2 and 4 days, cell viability and proliferation were analyzed.

Each material was tested in triplicate.

Results

Activated carbon fiber cloth characterization

The ACC material is made of twisted carbon fiber threads and the threads are knitted to form a flexible material (Fig. 1 (a)). The carbon fibers are 8-12 µm in diameter and have a lobed shape, as presented in the SEM images (Fig. 1 (b) and (c)). A heterogeneous sizing (coated during the industrial process) on the carbon surface is observed. The EDS analysis of the sizing shows the presence of aluminum and zinc on the carbon fiber surface (Fig. 1 (d)).

FM50K exhibits a narrow pore size distribution centered around 1 nm, with mainly ultramicropores of less than 0.7 nm and micropores of 0.7 to 2 nm. A negligible amount of mesopores (2 to 50 nm) is observed, representing less than 20 % of the total pore volume (Fig. 1 (e)). The mercury intrusion porosimetry curve reveals a distribution of macropores between 0.8 µm and 450 µm. The maximum of pore distribution at 95 µm corresponds to the average size of spaces between carbon fiber twists (Fig. 1 (f)). A pore size distribution between 0.8 and 60 µm is characteristic of different spaces between fibers of the same twist.

Data show that FM50K is a material with a multi-scale porosity.

The chemical and textural characteristics of FM50K are reported in Table 2. The surface chemistry characterization shows that the ACC substrate is slightly hydrophilic with 80.6 wt% C and a rather high oxygen content (4.8 wt%). Its pHPZC (pH value corresponding to an ACC net charge of zero) is almost neutral (7.8) and its pKa distribution of functional groups presents mainly acidic functional groups (2.93). The composition of the coatings was investigated by both SEM and TEM-EDS analyses. For all samples, the measured (Ca+Sr)/P molar ratio is 1.4, independent of the Sr content and significantly lower than that of the electrolyte solution (1.67). For the 5Sr-CaP coating, SEM and TEM-EDS analyses indicate a strontium concentration of 4.1 ± 0.6 and 3.0 ± 0.5 at.%, respectively. For the 10Sr-CaP coating, SEM and TEM-EDS measurements give a Sr content of 6.8 ± 0.8 and 6.6 ± 0.8 at.%, respectively. These values are in good agreement with those obtained by MP-AES (2.9 ± 0.5 and 5.5 ± 0.4 for the 5Sr-CaP and 10Sr-CaP coating respectively). It can be concluded that Sr concentrations of 5 and 10 at.% in the electrolyte solution lead to Sr contents of about 3 and 6 at.% in the xSr-CaP coating, respectively (Table 3). Additional P, Ca and Sr EDS-SEM elemental maps (Fig. 2

(e)-(f))

reveal that the strontium is uniformly distributed in the obtained coatings. These results show that electrolyte strontium concentration until 10 at.% leads to the formation of xSr-CaP coating with a plate-like morphology constituted by lamellar nanocrystals and a constant Ca+Sr/P ratio of 1.4. These features are very similar to that reported for Sr-doped stoichiometric HA (Ca/P = 1.67) obtained from precipitation in aqueous medium for which the Sr incorporation was found smaller than that expected from initial concentrations [START_REF] Bigi | Strontium-substituted hydroxyapatite nanocrystals[END_REF]. In order to better understand the effect of strontium incorporation on the chemical composition of the hydroxyapatite phase, FTIR analysis of non-doped and strontium-doped CaP coatings (samples respectively named 0Sr-CaP, 5Sr-CaP and 10Sr-CaP) was performed (Fig. 3.) The spectra of the three coatings obtained at constant potential are almost similar and are characteristic of AB-type calcium carbonated deficient hydroxyapatite (CDA), as described in previous papers [START_REF] Olivier | Influence of electrochemical parameters on the characteristics of sono-electrodeposited calcium phosphate-coated carbon fiber cloth[END_REF]. However, the presence of strontium in the coatings leads to a slight shift in the 5Sr-CaP and 10Sr-CaP spectra. In the three spectra, four main absorption domains are identified and are detailed for the 0Sr-CaP coating (Fig. 3 (b) and (c)). The first domain, between 450-800 cm -1 , consists of a double adsorption band located at 566 and 603 cm -1 with a shoulder at 580 cm -1 ; it is attributed to the ʋ4 bending mode of apatitic PO4 3-. The second shoulder observed at 633 cm -1 is assigned to apatitic OH - [START_REF] Rey | Physico-chemical properties of nanocrystalline apatites: Implications for biominerals and biomaterials[END_REF][START_REF] Vandecandelaere | Biomimetic apatite-based biomaterials: on the critical impact of synthesis and post-synthesis parameters[END_REF]. Other bands located at 535, 545 and 617 cm -1 are attributed to non-apatitic HPO4 2-, apatitic PO4 3-and non-apatitic PO4 3-[53]. The second domain, between 800 and 1300 cm -1 , displays the most intense and broadest band centered at 1036 cm -1 with the shoulder at 1094 cm -1 . It is attributed to the ʋ3 bending mode of apatitic PO4 3-. The weak band located at 962 cm -1 is assigned to the ʋ1 mode of PO4 3-[54]. The very weak absorption peak at 873 cm -1 can be ascribed to the presence of HPO4 2-and/or CO3 2-ions [START_REF] Rey | Bioactive Ceramics: Physical Chemistry[END_REF]. The third domain, between 1300 and 1580 cm -1 , shows a weak adsorption band at 1385 cm -1 that reveals the presence of a small amount of NO3 -ions in the CaP phase coming from the calcium nitrate precursor.

The absorption bands between 1400 and 1580 cm -1 indicate the existence of structural CO3 2- and these CO3 2-bands do not change with the evolution of the strontium rate [START_REF] Rey | Physico-chemical properties of nanocrystalline apatites: Implications for biominerals and biomaterials[END_REF][START_REF] Rey | The carbonate environment in bone mineral: A resolution-enhanced fourier transform infrared spectroscopy study[END_REF]. The last domain displays broad absorption bands between 1600 and 1700 cm -1 and between 2800 and 3600 cm -1 . They are assigned to structural water and adsorbed water bands, respectively [START_REF] Weber | Enzyme-functionalized biomimetic apatites: concept and perspectives in view of innovative medical approaches[END_REF]. and 2 ~ 50° (using K Cu = 0.154056 nm), corresponding respectively to (002), ( 10) and

(004) reflections of a carbon with a turbostratic structure [START_REF] Oberlin | Landmarks for graphitization[END_REF]. The other peaks observed on this diffractogram are characteristic of the sizing and correspond to the zinc aluminate (JCPDS n o : 39-473) present in the sizing composition. It is used as a surface treatment to enhance the fiber processability by preventing carbon fibers from breaking.

For the undoped CaP coating, the XRD diffractogram shows that the main diffraction peaks are located at 2 = 26.08°, 32.00° and 33.15° corresponding respectively to (002), ( 211) and

(300) reflections of the hexagonal hydroxyapatite (HA) structure (JCPDS no: 74-0566).

Indexation of all diffraction peaks confirm the occurrence of only calcium-deficient hydroxyapatite (CDA) as unique crystalline phase. For the two strontium-doped coatings (5Sr-CaP and 10Sr-CaP), the XRD diffractograms are similar to that of the 0Sr-CaP coating, excepted that a shift of the main reflections towards lower 2θ° values is observed when increasing the amount of strontium incorporated in the coating (Table 4). The lattice parameters (a and c) of the undoped and Sr-doped CDA deposited phases were determined using the Debye refinement of the XRD diffractograms (after removing the ACC substrate contribution). The observed increase of these parameters, a from 9.366 Å to 9.379 Å and c from 6.820 Å to 6.834 Å, is consistent with a calcium/strontium substitution in the apatite structure (the ionic radii of Ca and Sr being 1.00 and 1.18 Å, respectively). It could be noted that the relative increase of these parameters is slightly smaller than that reported for stoichiometric HA with the same Sr-doping obtained by precipitation [START_REF] Frasnelli | Synthesis and characterization of strontium-substituted hydroxyapatite nanoparticles for bone regeneration[END_REF] or high-temperature solid-state synthesis [START_REF] Bigi | Structural Refinements of Strontium Substituted Hydroxylapatites[END_REF]. This may suggest that the amount of Sr incorporated in the coherent domains of the HA nanocrystals is weaker than that expected from chemical analyses. Further structural investigation of the xSr-CaP coatings was performed using solid-state 31 P magic angle spinning (MAS) NMR spectroscopy. As shown in Fig. 5, the quantitative 31 P MAS spectra of the three samples display an intense resonance located at 2.7-2.9 ppm typical of a carbonated nanocrystalline CDA [START_REF] Panda | SEM and solid state NMR investigations of carbonate-containing hydroxyapatite nano-particles synthesized by hydroxide-gel technique[END_REF][START_REF] Aue | Solid-state phosphorus-31 nuclear magnetic resonance studies of synthetic solid phases of calcium phosphate: potential models of bone mineral[END_REF][START_REF] Roussière | Hybrid materials applied to biotechnologies: coating of calcium phosphates for the design of implants active against bone resorption disorders[END_REF][START_REF] Mellier | Design and properties of novel gallium-doped injectable apatitic cements[END_REF] in agreement with FTIR and XRD results. 31 P-1 H cross-polarization (CP) MAS experiments performed at different contact times reveal that the observed lineshape contains two overlapping peaks, a main one (at 2.8 ppm for 0Sr-CaP coating) and a broader one of weaker intensity (at about 2.0 ppm for the 0Sr-CaP sample), as previously observed for nanocrystalline CDA and HA [START_REF] Roussière | Hybrid materials applied to biotechnologies: coating of calcium phosphates for the design of implants active against bone resorption disorders[END_REF][START_REF] Mellier | Design and properties of novel gallium-doped injectable apatitic cements[END_REF][START_REF] Jäger | A solid-state NMR investigation of the structure of nanocrystalline hydroxyapatite[END_REF]. The 31 P isotropic chemical shifts and linewidths of the two peaks, assigned to the ordered core of the nanoparticle and to a disordered non-apatitic surface layer respectively [START_REF] Mellier | Design and properties of novel gallium-doped injectable apatitic cements[END_REF][START_REF] Jäger | A solid-state NMR investigation of the structure of nanocrystalline hydroxyapatite[END_REF][START_REF] Wang | Water-mediated structuring of bone apatite[END_REF][START_REF] Von Euw | Amorphous surface layer versus transient amorphous precursor phase in bone -A case study investigated by solid-state NMR spectroscopy[END_REF], were determined from simultaneous fits of CP-MAS spectra at different contact times and the obtained parameters were used to reconstruct the quantitative MAS spectra (see Fig. 5 (a), (b) and (c) where the thin green line is associated to apatitic PO4 3-groups of the particle core, while the broad purple line corresponds to PO4 3-and HPO4 2-surface groups). When increasing the Sr-content, slight shifts of the narrow (from 2.8 to 2.9 ppm) and of the broad (from 2.0 to 2.2 ppm) peaks are observed, in agreement with previous studies [START_REF] Bigi | Strontium-substituted hydroxyapatite nanocrystals[END_REF]. A significant increase of the relative intensity of the broad peak is also observed for the Srdoped coatings. This could suggest a preferential incorporation of Sr atoms in the disordered surface layer of the CDA nanoparticles and explain the moderate increase of cell parameters observed upon Sr-doping.

(a) (b) (c)
This assignment is supported by 2D 31 P-1 H heteronuclear correlation (HetCor) MAS NMR spectra, which allow deciphering the short-range proximities between P and H atoms. As shown in Fig. 6, the 31 P-1 H 2D HetCor MAS spectra of the three coatings (cross-polarization contact time of 1 ms) exhibit an intense correlation peak between the main 31 P resonance at 2.8-2.9 ppm and a narrow 1 H resonance at 0 ppm characteristic of the OH -hydroxyl ions of the apatitic structure. In addition, they show a broader correlation peak between the wider 31 P contribution at a slightly smaller chemical shift (2.0-2.2 ppm) and broad 1 H signals between 5 and 15 ppm corresponding to water molecules and HPO4 hydrogenophosphate ions. These features undoubtedly reveal that the undoped and Sr-doped CDA nanoparticles exhibit an ordered apatitic core and a disordered surface layer involving water molecules and nonapatitic PO4 3-and HPO4 2-groups [START_REF] Mellier | Design and properties of novel gallium-doped injectable apatitic cements[END_REF][START_REF] Jäger | A solid-state NMR investigation of the structure of nanocrystalline hydroxyapatite[END_REF][START_REF] Wang | Water-mediated structuring of bone apatite[END_REF][START_REF] Von Euw | Amorphous surface layer versus transient amorphous precursor phase in bone -A case study investigated by solid-state NMR spectroscopy[END_REF], thereby showing similarities with bone mineral nanoparticles [START_REF] Wang | Water-mediated structuring of bone apatite[END_REF][START_REF] Von Euw | Amorphous surface layer versus transient amorphous precursor phase in bone -A case study investigated by solid-state NMR spectroscopy[END_REF].

Biological characterization

Fluorescence analysis showed that human primary osteoblasts exhibited a fusiform shape along the ACC fibers (Fig. 7 (a) to (d)). After 2 and 4 days of culture of living cells (green fluorescence), very few dead cells (red fluorescence) were observed at the surface of the four materials. Cell enumeration (Fig. 7e) revealed that the cell density after 2 days was similar and statistically higher on FM50K (46 ± 10 cell/mm 2 ), FM 50K+5Sr-CaP (35 ± 10 cell/mm 2 ) and FM 50K+10Sr-CaP (35 ± 10 cell/mm 2 ) compared to FM 50K+0Sr-CaP (19 ± 3 cell/mm 2 ).

After 4 days, cell proliferation was observed on all materials, except on FM50K (from 46 ± 10 to 32 ± 11 cell/mm², i.e. a decrease of 30%) and was the highest on strontium-doped materials, especially FM50K+10Sr-CaP.

Discussion

In the present study, FM50K activated carbon fiber cloth was chosen as a working electrode for its specific characteristics and properties, such as conductivity, hydrophilicity due to the presence of acidic groups on the ACC surface, and multi-scale porosity. We hypothesized that microporosity and hydrophilic properties would allow high ACC impregnation in aqueous electrolyte and that electrical conductivity would be optimal for CaP coating by sono-electrodeposition. Our results established that FM50K-ACC was an appropriate substrate to obtain a uniform CaP coating with possible strontium substitution. The synthesis of strontiumsubstituted CaP coatings led to the expected final strontium range, i.e. between 3 and 7 at. % as based on previous studies. Interestingly, it was observed that the materials (xSr-CaP coated ACC substrate) were easy to handle allowing to perform all characterizations without damaging them. This also means that chemical bonds between xSr-CaP coatings and ACC fibers are tough enough.

It was relevant to analyze and then to compare the morphological, chemical and structural characteristics of the coatings, with and without strontium. Firstly, it is shown that strontium incorporation (< 10 at. %) does not modify the morphology and microtexture of xSr-CaP coatings. This leads to an uniform thickness (estimated around 1.5 µm) and a plate-like morphology. In TEM, whatever the materials, the coatings consist in lamellar nanocrystals which size is ranging from a few tens to a few hundreds of nanometers. Secondly, FTIR, XRD and solid-state NMR analyses reveal that the xSr-CaP coatings consist of carbonated CDA nanoparticles with an ordered apatitic core and a disordered hydrated surface layer.

XRD data highlighted variations of cell parameters of the hexagonal crystallographic structure: a parameter increases from 9.366 Å to 9.379 Å and c parameter increases from 6.820 Å to 6.834 Å. This is an evidence of a calcium/strontium substitution in the coherent domains of the HA nanocrystals. It is confirmed by the 31 P MAS spectra of the xSr-CaP coatings which suggest that the Sr atoms are incorporated into both the core and the disordered surface layer of the CDA nanoparticles.

This incorporation mechanism is expected to lead to different Sr release rates in solution or in vivo since the Sr cations located in the hydrated disordered layer are loosely bound and can be more easily exchanged in solution than those incorporated into the more ordered core of the nanoparticles [67]. Biological experiments revealed that ACC and CaP coated ACC are biocompatible. However, ACC substrate without any coating does not favor osteoblast proliferation, hence the interest of CaP coating. It is worth noting that the best proliferation rate is obtained with the highest Sr-substituted CaP coating. This result is consistent with previous studies showing the benefit of using Sr to increase bone regeneration [START_REF] Capuccini | Strontium-substituted hydroxyapatite coatings synthesized by pulsed-laser deposition: In vitro osteoblast and osteoclast response[END_REF].

Conclusion

The originality of the present study concerns the incorporation of a strontium cation into a biomimetic apatite coating deposited on conductive ACC material by sono-electrodeposition.

This process is an interesting and efficient technique to deposit a strontium-substituted carbonated CDA coating on a carbon fiber cloth. A single crystalline Sr-CaP phase is obtained, in which calcium and strontium cations occupy similar positions in the apatitic lattice. The substitution of calcium by strontium ions, having a different ionic radius and electronegativity, leads to slight structural modifications of the CDA phase but does not modify the nanocrystal organization or the microtexture of the CDA coatings. The coatings display a plate-like morphology consisting of lamellar nanocrystals of carbonated CDA. The structure of the nanocrystals comprises an apatitic ordered core and a non-apatitic disordered and hydrated surface layer, typical of a biomimetic apatite. After hybrid material synthesis and characterization, cell viability assays were performed on human primary osteoblast cells, cultured on the carbon fiber cloth substrate without and with xSr-CaP coatings. As these materials show a strong affinity with human primary osteoblasts, the incorporation of strontium in the CDA phase structure has a beneficial role on cell proliferation.

  The xSr-CaP coating morphology and thickness were investigated by scanning electron microscopy (SEM-Hitachi S4500) operating at 5 kV. To determine the microtexture of the xSr-CaP coatings, transmission electron microscopy (TEM) images were acquired with a TEM-PHILIPS CM20 operating at 200 kV. For sample preparation, 1 mg of coating was scraped from the carbon surface, dispersed in ethanol and sonicated for 30 sec; a droplet of the suspension was deposited on an amorphous carbon coated copper grid (Electron Microscopy Sciences). Ca/P and Ca+Sr/P atomic ratios were determined by Energy Dispersive X-ray Spectroscopy (EDS) analysis performed in SEM at 10 kV and in TEM at 200 kV. The atomic ratios were calculated from the intensities of the Ca peak (K Ca = 3.691 keV), P peak (K P = 2.015 keV) and Sr peak (K Sr = 14.164 keV) present in the SEM-EDS or TEM-EDS spectra.The strontium concentration present in xSr-CaP coatings was measured by microwave plasma atomic emission spectroscopy (MP-AES) using a 4200 MP-AES Agilent Technologies. For that, the material (ACC/xSr-CaP) was immersed in a 35 % nitric acid solution (Sigma Aldrich®) for 3 h. Supernatant containing the xSr-CaP was diluted to a concentration of 0.5 % nitric acid to prevent the formation of oxides. The solution obtained was sprayed into a nitrogen plasma which reached the temperature of 5000 °C. At this temperature, the atoms in the nebulized sample were dissociated, atomized, ionized and excited. The emitted photons were measured using optical detection at the characteristic wavelengths of the elements of interest (λCa = 393.366 nm and λSr = 407.771 nm) in order to quantify the elements of interest in the liquid solution by comparing the signal measured with the calibration curve.The chemical composition of the xSr-CaP coatings was analyzed by Fourier transform infrared spectroscopy (FTIR). Spectra were acquired in transmission mode using a Nicolet 6700 spectrometer on KBr pellets in the range of 4000-400 cm -1 (resolution = 4 cm -1 , 256 scans). Pellets were prepared by admixing 1 mg of ACC/xSr-CaP hybrid materials with 99 mg of KBr.The crystallographic structure of the xSr-CaP coatings was analyzed by X-Ray diffraction (XRD), using an INEL diffractometer (CSP120), working at 40 kV and 30 mA (Cu Kα = 1.5418978 Å). Data were collected in the range 2θ = 10°-80°, with a curve detector (angular instrumental resolution of 0.05°), in transmission mode, using a Si(111) monochromator. To determine the lattice parameters of xSr-CaP crystalline phases, refinements were performed by the Lebail method using the Jana2006 software. Solid-state nuclear magnetic resonance (NMR) experiments were performed on a Bruker Avance I spectrometer operating at a magnetic field of 7.0 T ( 1 H and 31 P Larmor frequencies of 300.13 MHz and 121.50 MHz, respectively), equipped with a 4 mm magic-angle spinning (MAS) NMR probehead. The 31 P 1D quantitative MAS, { 1 H}-31 P cross-polarization (CP) MAS and 2D { 1 H}-31 P heteronuclear correlation spectra were recorded at a sample spinning frequency of 10 kHz.1 H and 31 P NMR chemical shifts were referenced to TMS and H3PO4 in solution (85 %), respectively.

  Cell viability was assessed using the fluorescent assay LIVE/DEAD TM Viability/Cytotoxicity Kit (Life Technologies SAS, France). At each endpoint ACC patches were gently washed with phosphate-buffered saline (PBS) and incubated at room temperature with PBS containing 2 µmol/L calcein-AM and 4 µmol/L ethidium-homodimer-1. Calcein-AM is a nonfluorescent cellpermeant fluorescein derivative, which is converted by cellular esterase activity into cellimpermeant and highly fluorescent calcein. Calcein accumulates inside live cells having intact membranes, which results in a green fluorescent signal. Ethidium-homodimer-1 enters dead cells with damaged membranes and undergoes a 40-fold enhancement of fluorescence upon binding to their DNA, leading to a red fluorescent signal. After 20 min incubation, the ACC patches were put onto a glass slide, covered with a coverslip and viewed under a fluorescence microscope (Zeiss Observer Z1 with FITC filters and Texas Red) equipped with FITC and Texas Red filters. The number of live and dead cells was estimated using the ImageJ Software. 2.5. Statistical analysis Quantification of strontium rates in in vitro results are expressed as a mean ± standard deviation of independent experiments. Differences between in vitro results were statistically analyzed using Student test and differences were considered statistically significant at *p < 0.05, **p < 0.01 and ***p < 0.001.

Fig. 1 .

 1 Fig. 1. SEM images of pristine activated carbon fiber cloth (FM50K) at three magnifications: a) x100, b) x1000, c) x3000 and d) the EDS analysis of carbon fiber surface e) Pore size distribution of FM50K measured by N2 adsorption at 77 K f) Pore size distribution of FM50K measured by mercury intrusion porosimetry. (2-column fitting image)

3. 2 .

 2 Calcium phosphate coating characterization SEM micrographs show the occurrence of similar morphologies for xSr-CaP coatings deposited by sono-electrodeposition, independently of the strontium ratio present in the electrolyte (Fig.2 (a)-(d')). A uniform coating is observed for the three strontium ratios, with a plate-like morphology and a uniform thickness estimated around 1.5 µm. Different observations suggest that the incorporation of strontium does not influence the thickness, uniformity or morphology of the coating. TEM characterization did not reveal any difference in microtexture between coatings without (Fig.2 (b'')) and with strontium (Fig.2 (c'') and (d'')). The three coatings are composed of lamellar nanocrystals varying in size from a few tens to a few hundreds of nanometers. Due to the thinnest of the nanocrystallites, a folding of their edges was often observed (white arrows).

Fig. 2 .Table 3 .

 23 Fig. 2. SEM images of ACC with and without coatings: a and a') FM50K (ACC pristine), b and b') FM50K+0Sr-CaP, c and c') FM50K+5Sr-CaP and d and d') FM50K+10Sr-CaP. TEM images of corresponding coatings: b'') 0Sr-CaP, c'') 5Sr-CaP and d'') 10Sr-CaP; e), f) EDS-SEM elemental maps of P, Ca and Sr performed on 5Sr-CaP and 10Sr-CaP coatings. (2-column fitting image)

Fig. 3 .

 3 Fig. 3. (a) FTIR spectra of the xSr-CaP materials (bottom-up: FM50K+0Sr-CaP; FM50K+5Sr-CaP; FM50K+10Sr-CaP, respectively); (b) detail of ʋ4 mode of PO4 3-domain with peak deconvolutions and (c) detail of ʋ1 and ʋ3 modes of the PO4 3-domain. (2-column fitting image)

Fig. 4 .Table 4 .

 44 Fig. 4. XRD diffractograms of the pristine ACC and the ACC/xSr-CaP materials (bottom-up: FM50K; FM50K+0sr-CaP; FM50K+5sr-CaP; FM50K+10sr-CaP, respectively). (2-column fitting image)

Fig. 5 .

 5 Fig. 5. 31 P quantitative MAS NMR spectra of the ACC/xSr-CaP materials: a) 0Sr-CaP; b) 5Sr-CaP and c) 10Sr-CaP. Fits of experimental spectra are shown in red, the narrow peak assigned to the ordered core contribution and the broad one corresponding to the disordered surface layer are shown in green and purple. (1-column fitting image + color)

Fig. 6 .

 6 Fig. 6. Two-dimensional { H} 31 P HetCor MAS NMR spectra of the sono-electrodeposition coatings: a) 0Sr-CaP; b) 5Sr-CaP and c) 10Sr-CaP. (2-column fitting image)

Fig. 7 .

 7 Fig. 7. (a)-(d) Typical fluorescence micrographs obtained after 2 and 4 days of culture of human primary osteoblasts on the four materials, namely FM 50K (a); FM 50K+0Sr-CaP (b); FM 50K+5Sr-CaP (c) and FM 50K+10Sr-CaP (d). Cell staining: live=green and dead=red. Mean ± standard deviation, *p < 0.05, **p < 0.01 and ***p < 0.001. e) Quantification of cell number by mm² on the four materials after 2 and 4 days of culture using the ImageJ software. (2-column fitting image + color)

Table 1 .

 1 The investigated materials as a function of the electrolyte composition. (

2-column fitting image)

	[Ca] mmol/L	-	5	4.75	4.5
	[Sr] mmol/L	-	0	0.25	0.5
	[P] mmol/L	-	3	3	3
				FM50K	FM50K	FM50K
			FM50K	+	+	+
				0Sr-CaP	5Sr-CaP	10Sr-CaP
	2.3.	Characterization techniques of materials
	The characteristics of the ACC substrate and the xSr-CaP coatings were analyzed using
	several techniques.
	Elemental analysis (weight %) of C and O were determined using a Flash 2000 (CHON-S)
	Electrolyte
			Material 1	Material 2	Material 3	Material 4
	composition

Table 2 .

 2 Chemical and textural characteristics of FM50K. (

2-column fitting image) Elemental analysis (mass.%) pKa distribution (mmol/g) N2 and CO2 adsorption

						Acidic	Vmicro 1	Vmeso 1	Vultramicro 2
					7<pKa<1					Vtotal
	% C	% O	pHPZC	3<pKa<7		functional	N2, DFT	N2, DFT	CO2, DR
					1					cm 3 /g
						groups	cm 3 /g	cm 3 /g	cm 3 /g
	80.6	4.8	7.8	2.13	0.80	2.93	0.41	0.07	0.55	0.51

The microporous and mesoporous volumes were determined by applying the density functional theory method (DFT) to the N2 adsorption isotherms.

The ultramicroporous volume was determined by applying the Dubinin-Radushkevich theory (DR) to the CO2 adsorption isotherm.

Acknowledgements

This work was financially supported by the FP7-IRSES European program "ABREM" 2012-2015, the Region Centre-Val de Loire projects "MatBioReOs" 2014-2017 and "MatCCaPBio" 2018-20210.

The authors are grateful to K. Beck for mercury porosity, R. Nehmé for assistance in MP-AES experiments, A. Richard for the SEM images, F. Warmont for the TEM images and E.

Verron for calculating the hydroxyapatite lattice parameters.