J.-M Commandré

H Lahmidi

S Salvador

N Dupassieux

Pyrolysis of wood at high temperature: The influence of experimental parameters on gaseous products

Keywords:

The pyrolysis of wood was carried out in an Entrained Flow Reactor at high temperature (650 to 950 °C) and under rapid heating conditions (N 10 3 K s -1 ). The influence of the diameter and initial moisture of the particle, reactor temperature, residence time and the nature of the gaseous atmosphere on the composition of the gaseous products has been characterised. Particle size, between 80-125 and 160-200 μm, did not show any impact. Pyrolysis and tar cracking essentially happen in very short time period: less than 0.6 s; the products yields are only slightly modified after 0.6 s in the short residence times (several seconds) of our experiments. Higher temperatures improve hydrogen yield in the gaseous product while CO yield decreases. Under nitrogen atmosphere, after 2 s at 950 °C, 76% (daf) of the mass of wood is recovered as gases: CO, CO 2 , H 2 , CH 4 , C 2 H 2 , C 2 H 4 and H 2 O. Tests performed under steam partial pressure showed that hydrogen production is slightly enhanced.

Introduction

The necessity for alternative energy sources has caused renewed interest in the use of biomass, since this renewable fuel does not contribute to increased CO 2 concentration in the atmosphere [START_REF] Guerrero | Pyrolysis of eucalyptus at different heating rates: studies of char characterization and oxidative reactivity[END_REF][START_REF] Belosludov | Combinatorial computational chemistry approach as a promising method for design of Fischer-Tropsch catalysts based on Fe and Co[END_REF]. Because of tension on petroleum cost, second generation biofuels produced by Fischer-Tropsch synthesis after biomass gasification show an increasing interest. The physical and chemical mechanisms involved in biomass steam gasification are still poorly understood [START_REF] Dupont | Biomass pyrolysis experiments in an analytical entrained flow reactor between 1073 and 1273 K[END_REF], especially under high temperatures (1073 K b T b 1273 K) and high heating rate (N10 3 K s -1 ). Our work focuses on processes occurring at the particle scale during the first seconds of the transformation. This first stage-named "pyrolysis" in our study-gathers both biomass decomposition, often called "primary pyrolysis" or "devolatilisation", and the secondary cracking reactions of the vapours produced during primary pyrolysis.

In literature, investigations on pyrolysis characterization have been mainly carried on slow pyrolysis [START_REF] Vamvuka | Pyrolysis characteristics and kinetics of biomass residuals mixtures with lignite[END_REF]; high heating rate pyrolysis of biomass has not been studied as closely, particularly in the high temperature domain [START_REF] Bridgwater | Fast pyrolysis processes for biomass[END_REF][START_REF] Goyal | Bio-fuels from thermochemical conversion of renewable resources: a review[END_REF]. Published high temperature pyrolysis experiments on biomass involve reactors such as fluid bed reactors [START_REF] Gray | Pyrolysis of a wood-derived material. Effects of moisture and ash content[END_REF][START_REF] Rapagna | Steam-gasification of biomass in a fluidised-bed of olivine particles[END_REF], radiant heating techniques [START_REF] Di Blasi | Weight loss dynamics of wood chips under fast radiative heating[END_REF], drop tube furnaces and entrained flow reactors [START_REF] Van De Steene | Rice husk, straw and bark behaviour during pyrolysis, combustion, and gasification: fundamental study, 12[END_REF][START_REF] Shuangning | Devolatilisation characteristics of biomass at flash heating rate[END_REF][START_REF] Wagenaar | Co-combustion of miscanthus in a pulverised coal combustor: experiments in a droptube furnace[END_REF]. In the present paper, experiments are performed in a free fall reactor, also called entrained flow reactor (EFR). This reactor is an interesting analytical tool to study the global reaction of flash pyrolysis [START_REF] Dupont | Biomass pyrolysis experiments in an analytical entrained flow reactor between 1073 and 1273 K[END_REF][START_REF] Van De Steene | Rice husk, straw and bark behaviour during pyrolysis, combustion, and gasification: fundamental study, 12[END_REF][START_REF] Li | Fast pyrolysis of biomass in free-fall reactor for hydrogen-rich gas[END_REF][START_REF] Visentin | Experimental study of cellulose fast pyrolysis in a flow reactor[END_REF][START_REF] Cetin | Influence of pyrolysis conditions on the structure and gasification reactivity of biomass chars[END_REF][START_REF] Brink | A flow reactor technique for the study of wood pyrolysis. I. Experimental[END_REF][START_REF] Commandré | Formation des oxides d'azote lors de la combustion de cokes de pétrole dans des conditions de précalcinateur de cimenterie[END_REF][START_REF] Zanzi | Rapid high-temperature pyrolysis of biomass in a free-fall reactor[END_REF][START_REF] Zanzi | Rapid pyrolysis of agricultural residues at high temperature[END_REF], since very high heating rates of the particles can be obtained. Moreover, the residence time can be varied over a range of a few seconds, and the particles can be assumed to be isolated from each other, which makes interpretation easier [START_REF] Visentin | Experimental study of cellulose fast pyrolysis in a flow reactor[END_REF][START_REF] Commandré | Formation des oxides d'azote lors de la combustion de cokes de pétrole dans des conditions de précalcinateur de cimenterie[END_REF]. Up to now, few experiments have been performed under the operating conditions of interest, the most complete set of data being given by Zanzi [START_REF] Zanzi | Rapid high-temperature pyrolysis of biomass in a free-fall reactor[END_REF][START_REF] Zanzi | Rapid pyrolysis of agricultural residues at high temperature[END_REF]. Zanzi et al. have studied the effect of the process conditions, such as temperature and particle size, on the product of high heating rate pyrolysis of agricultural residues, including gas composition and char reactivity, in a free fall reactor. They showed that higher reaction temperature (800 to 1000 °C) led to lower yields of tars and higher yields of gaseous products. In addition, pyrolysis of smaller particles (0.5 to 1 mm) improved the char reactivity for combustion or gasification [START_REF] Di Blasi | Kinetics of primary products formation from wood pyrolysis[END_REF]. Results found by Guerrero et al. [START_REF] Guerrero | Pyrolysis of eucalyptus at different heating rates: studies of char characterization and oxidative reactivity[END_REF] showed that the reactivity to oxygen of the char produced by eucalyptus pyrolysis increases with the heating rate of biomass particles.

Several works [START_REF] Guerrero | Pyrolysis of eucalyptus at different heating rates: studies of char characterization and oxidative reactivity[END_REF][START_REF] Shuangning | Devolatilisation characteristics of biomass at flash heating rate[END_REF][START_REF] Li | Fast pyrolysis of biomass in free-fall reactor for hydrogen-rich gas[END_REF][START_REF] Dupont | Vapogazéification de la biomasse: contribution à[END_REF] have focused on gaseous products. They showed that gaseous amounts increase with temperature, residence time, and heating rate. Li et al. [START_REF] Li | Fast pyrolysis of biomass in free-fall reactor for hydrogen-rich gas[END_REF] showed that solid and liquid yields increase with particle diameters, whereas gaseous amounts decrease. The main gaseous products during flash pyrolysis of biomass are carbon monoxide (CO), hydrogen (H 2 ), carbon dioxide (CO 2 ), methane (CH 4 ) and other light hydrocarbons. Studies interested in these gaseous products [START_REF] Meesri | Lack of synergetic effects in the pyrolytic characteristics of woody biomass/coal blends under low and high heating rate regimes[END_REF] have shown that the ratio CO/CO 2 increases with temperature, but decreases when the heating rate increases. For a temperature range of 500-800 °C, Li et al. confirmed that H 2 and CO concentrations increase with temperature, and the H 2 /CO ratio increases with particle diameters. Moreover, Couhert et al. [START_REF] Couhert | Is it possible to predict gas yields of any biomass after rapid pyrolysis at high temperature from its composition in cellulose, hemicellulose and lignin?[END_REF] showed that gaseous products of high heating rate pyrolysis of wood particles are slightly influenced by wood nature, and are similar for soft and hard wood.

The parameters which can influence biomass gasification, and particularly pyrolysis as the first reaction step, are particle diameter, particle initial moisture, reactor temperature and residence time, and the nature of the gaseous atmosphere. In our study, flash pyrolysis of wood particles was characterised over a range for each parameter: The purpose of our study was to evaluate the influences of each one of these parameters on the nature of gaseous products and on their amounts during high temperature flash pyrolysis of wood. This reaction is the first step in the conversion of biomass to syngas (H 2 and CO) in a fluidised bed and/or entrained flow reactors under development at this time.

Experiments

The experiments were carried out with pine wood. Table 1 gives the results of the proximate and the ultimate analysis, as well as the natural moisture content of pine wood. Levels of ash and volatile matters are measured following the standard methods. Results are in agreement with those found by ECN in their Phyllis Database [23]: this wood contains 10.7% of natural moisture, low ash content and very low N and S content.

Experiments were carried out at laboratory scale in an Entrained Flow Reactor (EFR) under high temperatures (650 to 950 °C) and high heating rate (N10 3 K s -1 ). This EFR is schematically described in Fig. 1. The complete description of the apparatus and its thermal and hydrodynamic characterization were given in previous works from Van de Steene [START_REF] Van De Steene | Controlling powdered fuel combustion at low temperature in a new entrained flow reactor: laboratory experiments and modelling[END_REF] and Commandré [START_REF] Commandré | Formation des oxides d'azote lors de la combustion de cokes de pétrole dans des conditions de précalcinateur de cimenterie[END_REF]. A steam generator and a steam over-heater were added to conduct experiments under different steam partial pressure.

It is a 75 mm i.d. and 2 m long tube into which an isothermal reaction zone is fed with an electrically preheated laminar gas flow. Besides the nitrogen, a controlled steam flow could be added thanks to a water mass flow meter/controller and a steam generator that produced steam at 150 °C.

The powdered solid is injected through a hot-oil cooled feeding probe avoiding steam condensation, and dispersed over the cross-section of the reactor using a dispersion dome. It is pyrolysed for a controlled residence time before being sampled by a hot-oil cooled probe.

The solid injection consists of a low velocity conveyor belt that ensures an accurate solid mass flow rate. The belt is fed from a V-shaped rail along which a precisely weighted quantity of solid is regularly spread using a calibrated wedge. The mass flow of the biomass was 1 g/min, which ensures that particles may be assumed to be isolated from each other. In these conditions of gas and solid flow, the particles density is low in order to keep isothermal the flow even exo or endothermal reactions are involved. Temperature measurements with thermocouples on radial and axial position along the reactor have shown that the reactor is isothermal after 0.30 m away from the solid injection [START_REF] Commandré | Formation des oxides d'azote lors de la combustion de cokes de pétrole dans des conditions de précalcinateur de cimenterie[END_REF][START_REF] Dupont | Vapogazéification de la biomasse: contribution à[END_REF].

The gas and solid phases can be collected at different heights (z) using a hot-oil cooled sampling probe, allowing precise control of the residence time. Residence time of gas was calculated for the sampled section around the axis of the reactor by considering the average gas velocity in the section. Residence time for particles was calculated for each sample taking into account the slip velocity of the solid in the gas.

Downstream the sampling probe, a condenser was used to remove water from gases prior to analysis. The water remaining in gases after the condenser was quantified using a cooled-mirror hygrometer. The sum of the mass of condensed water and of the mass of water quantified by the hygrometer in gases represents the water contained in the sampled gases.

After separating gas from particles and water, the gases are forwarded to the analysers, via a heated line:

a Fourier Transform InfraRed (FTIR) analyser to quantify continuously CO, CO The experimental conditions for all experiments are summarized in Table 2. As illustrated in Fig. 2, the repeatability-checked for experiments at the two positions 50 cm and 90 cm and at 900 °C-was satisfactory (variation b6%).

Two size fractions (80-125 μm and 160-200 μm) were selected by pneumatic sieving to remove the fine particles under 80 μm, then by manual sieving to prepare the two size fractions.

In the following figures where gas species are plotted, the gas yields are expressed through the ratio of the mass of gas on the initial wood mass (wet, unless specified) injected.

Since the amount of ash is very low in wood, the normalised ashtracer method usually used to measure solid conversion could not be applied for this test.

Results and discussion

Observations of wood particles

A view of the initial pine wood particles under scanning electron microscope is presented in Fig. 3-a. One can observe different shapes of particles with various lengths; the majority is stick-shaped. Their dimensions are very different and the ratio between the length and the width of a particle varies from 2 to 10.

After pyrolysis at 900 °C for 2 s under N 2 , two categories of particles can be noticed (Fig. 3-b). The first ones have preserved their initial form while contracting. The second ones have undergone a drastic modification in morphology: they have become hollow spheres. The diameter of the final char particles approximately corresponds to the width of the initial wood particles. These phenomena may be an effect of temperature [START_REF] Bonelli | Effect of pyrolysis temperature on composition, surface properties and thermal degradation rates of Brazil nut shells[END_REF]. These observations let us think that a fusion phenomenon occurred and that the particles pass through a plastic phase with gas generation inside. Considering that biomass is a polymer, such behaviour of wood during pyrolysis seems logical and is in agreement with the "fusion-like phenomenon" of wood stated by Lédé et al. [START_REF] Lédé | Fusion-like behaviour of wood pyrolysis[END_REF] and the existence of an intermediate liquid compound (ILC).

Mass balance

From the results of experiments, it is possible to propose a mass balance in terms of yields of solid residue, gas and tars. Such results are particularly relevant for considerations of eventual steam reforming of the pyrolysis products to supply a syngas (CO and H 2 ). In some processes, the solid residue might also be separated from the gases and tars, to be burned separately and supply the energy required by the process. The experiment was carried out at 900 °C for around 2 s under nitrogen. For this specific experiment, a high gassampling flow rate was used to collect all the products (gas, tars and solid residues) in the reactor.

Char and tar mass quantification was carried out using the following procedure. During the experiment, the major part of the char residue was recovered in a settling box and in a filter at the outlet of the sampling probe and weighed once the experiment was over. Then, the sampling probe was cleaned using acetone. The liquid obtained contained tars, as well as chars that had been stuck on the probe by these tars. This liquid was filtered and solids were washed with acetone in order to separate the solid residue from the liquid phase; the mass of chars weighed was added to the previous one.

The remaining liquid was left under a hood for 48 h. The viscous residue after liquid evaporation was considered to be the mass of tars. The mass difference-after and before experiment-of the filter located downstream of the sampling probe and settling box was also quantified as tars.

The experiment was repeated three times; the results obtained are presented in Fig. 4, as a % of the mass of the initial biomass. The term "gas" represents the sum of the yields of the gases that were analysed, i.e. H 2 , CO, CO 2 , CH 4 , C 2 H 4 , C 2 H 2 and H 2 O. The quantification of gases appears to be very repeatable, and the mass released by particles is as high as 76%. The accuracy of char quantification is not as good: the char yields vary between 4.4% and 7.7%. This can be attributed to the difficulty in cleaning and collecting all particles stuck with tars in the sampling probe.

The quantification of tars is repeatable (1.4 to 1.5%). The measured mass of tars is in accordance with the value of 1% found by Zanzi [START_REF] Zanzi | Rapid high-temperature pyrolysis of biomass in a free-fall reactor[END_REF]. Nevertheless, doubts remains about the quantity of hydrocarbon species that are evaporated with acetone during the protocol.

Fig. 4 also shows that approximately 16%w of the wood initial mass injected in the reactor is non-quantified. It is difficult at this stage to attribute this mass to non-quantified gases, char or tar. In order to specify the nature of the missing 16%w, we carried out mass balances for the elements C, H and O. A balance for each element-detailed in Fig. 5-indicates that 10%w of carbon, 1%w of hydrogen and 6%w of

Influence of particle size

As far as particle diameter is concerned, experiments were carried out at 900 °C with two classes of particle size of pine: [80-125] and [160-200] μm. At this temperature, comparing the different gas species in the product gas at the reactor outlet, very similar results are obtained for both sizes (Fig. 6). Since the pyrolysis process is very fast, no conclusion can be drawn immediately in terms of the apparent reaction kinetics. One could expect the particles with different diameters to have undergone different heating rates. Nevertheless, as demonstrated previously [START_REF] Van De Steene | Controlling powdered fuel combustion at low temperature in a new entrained flow reactor: laboratory experiments and modelling[END_REF], the heating rate in such conditions is controlled essentially by the mixing of the flows of the cold transport gas and the hot atmosphere gas. It is likely that it was similar for both particle sizes, which could explain the similarity in the gases produced. In these conditions, this result shows that for both particle sizes there is no heat and mass transfer limitation for the solid-gas reaction, which is in accordance to [START_REF] Dupont | Biomass pyrolysis experiments in an analytical entrained flow reactor between 1073 and 1273 K[END_REF].

We should reiterate here that the observed gases are the result of primary pyrolysis and also of tar cracking. One might also imagine that the primary products of pyrolysis were different for the two particle sizes, but that tar cracking and gas phase reactions brought them to similar compositions.

Influence of residence time

Pine wood particles were pyrolysed at 850 °C in the EFR. Results of gas yields during pyrolysis are plotted in Fig. 7. Except for hydrogen, the yields of all species globally appear to be virtually unchanged after the first measurement point at 0.6 s. In the pyrolysis process, primary volatile matters are produced first. They consist mainly of some permanent gases and heavy hydrocarbons, which are then cracked to form additional permanent gases. Therefore, it is likely that both pyrolysis and tar cracking reactions only slowly evolve during the short residence time involved in our experiments.

This figure also shows that CO is the main gas released by wood pyrolysis. Hydrogen is important in the molar fraction, but it has the lowest yield in the mass fraction. The species C 2 H 4 shows increasing yields with time while C 2 H 2 and H 2 decrease. Based on these observations, the homogeneous reactions of hydrogenation

C 2 H 2 + H 2 ⇄ C 2 H 4
can be assumed to occur inside the reactor. However, the mass balance of this reaction can only explain 25% of the decrease observed on H 2 with residence time.

Influence of reactor temperature

As discussed in the Introduction, increasing the temperature and/ or the heating rates supports the gas production while minimising the production of solid residue. However, the nature of gases also varies with the reaction temperature. In our experiments, pine wood was pyrolysed under nitrogen between 650 and 950 °C and sampled at 80 cm from the injection point, i.e. when the pyrolysis had terminated.

Fig. 8 shows the gas yields versus temperature. Note that in this figure the yield of H 2 was multiplied by 10 and that of CO was divided by 10. Hydrogen production is clearly improved by an increase in temperature: its yield is zero at 650 °C, it appears at 700 °C and its mass percentage rises to over 1% at 950 °C. These results lead us to think that higher pyrolysis temperature would lead to an even higher percentage of hydrogen.

The increase in hydrogen yield is accompanied by an increase in C 2 H 2 from 0.3%w at 650 °C to more than 2%w at 950 °C. Contrary to C 2 H 2 , the mass percentage of C 2 H 4 tends to decrease at temperatures above 850 °C according to results found by Dupont et al. [START_REF] Klose | On the intrinsic reaction rate of biomass char gasification with carbon dioxide and steam[END_REF]. The reaction of C 2 H 4 consumption to produce C 2 H 2 and H 2 is thermodynamically favoured when temperature increases. This could explain why C 2 H 4 amounts decrease when temperature increases in favour of C 2 H 2 .

The concentrations of CO 2 increase slightly until 900 °C, and seem to stabilise above this temperature. The quantity of CH 4 increases from 2% at 650 °C to 6% at 950 °C. If CH 4 is regarded as a stable species that is difficult to convert to H 2 by further steam reforming, these results could indicate that pyrolysis at low temperature is recommended.

The CO concentration seems to increase with temperature from 22% at 650 °C to 47% at 800 °C, and is constant until 950 °C.

At this point, we compare the experimental results to the composition of the gas phase at thermodynamic equilibrium. More precisely, we compare the values of H 2 /CO obtained experimentally with pine wood to results found using the software FactSage. At each temperature, the C-H-O-N-S composition of the gas phase at the exit of the reactor was introduced into FactSage as initial conditions; the equilibrium composition was then computed by minimising the Gibbs energy. The presence of tars and of solid carbon residue was neglected in this approach. Results are plotted in Fig. 9. The computed H 2 /CO ratio is higher than that obtained experimentally, whatever the temperature. Temperature increase is accompanied by an increase in the experimental H 2 /CO ratio. Experimental gas phase composition is thus very far from the thermodynamic equilibrium at the lower temperature, but tends towards it at the higher temperature, despite not reaching equilibrium, even at 950 °C. At this temperature under 1000 °C, it is clear that reactions are under kinetic control.

Influence of particle initial moisture

As flash pyrolysis is a very rapid process, one can wonder if removing the moisture initially present in particles, a very endothermic reaction, would accelerate or slow down the process and/or impact on the amounts of the gas species formed. The H 2 O coming from wood natural moisture can also modify the apparent kinetics by reacting with biomass following the heterogeneous gasification reaction. Two tests were run, the first one with initial pine wood moisture of 10.7%, the second one with pine wood dried at 105 °C for 24 h (no moisture). These two woods (wet and dry) were pyrolysed at 850 °C for 0.9 to 1.5 s under nitrogen. The results obtained are presented in Fig. 10.

The results are very similar for all the gas species quantified. It can be concluded that there is no influence of natural moisture on the nature and the amounts of the gases finally resulting from pyrolysis.

Influence of the nature of the atmosphere

Until now, experiments have been carried out under N 2 . However, in a gasification process, steam-either resulting from wood drying or injected in the reactor-is present in large amounts. The presence of steam in the gas phase can affect the process through different mechanisms:

it can modify the heating rate of particles since the radiation is by far more efficient in the presence of steam than in nitrogen alone; -it can favour chemical reactions such as tar cracking or water-gasshift reaction.

The effect of steam was here characterised by performing two experiments in the EFR according to conditions described in Table 2. Pine wood was thus pyrolysed at 900 °C for 0.8 to 1.5 s in the reactor under atmosphere containing nitrogen with alternatively 0 or 10 vol.% of steam. One can note that the initial moisture of wood, after drying, led to the formation of steam in the reducing atmosphere. In the present work, it was equivalent to 0.86 vol.% steam only in the atmosphere.

In order to compare the effect of the nature of the atmosphere on pyrolysis, we have converted results to show the difference between gas yield under 10 vol.% of steam and under nitrogen. Values are positive when gas yields under 10 vol.% of steam are higher than gas yields under nitrogen, and vice versa. Results plotted in Fig. 11 show that the higher the initial steam concentration (10%), the higher the H 2 concentration and the lower the CO concentration. The occurrence of the water-gas-shift reaction might be an explanation for this result. From 0 to 10% of initial steam concentration, CO concentration decreases by 2.8 10 -3 mol/g wood while hydrogen increases by 1.2 10 -3 mol/g wood . However, CO 2 concentration increases by 4.9 10 -4 mol/g wood when it should have increased of circa 1.2 10 -3 mol/g wood (equal to hydrogen increase) if water-gasshift was the only reaction involved. It can be concluded that the water-gas-shift reaction alone cannot explain the changes in the gas phase composition when steam is added to the atmosphere.

No significant influence of the steam partial pressure in the atmosphere is observed on CH 4 , C 2 H 4 and C 2 H 2 concentrations.

As regards the possibility of an attack on the char residue by steam, this reaction is not likely to occur significantly at such a low temperature. It has indeed been shown [START_REF] Klose | On the intrinsic reaction rate of biomass char gasification with carbon dioxide and steam[END_REF][START_REF] Mermoud | Experimental and numerical study of steam gasification of a single charcoal particle[END_REF] that even at a temperature of 1000 °C, it takes several tenths of minutes to gasify char.

Conclusion

The experimental work carried out made it possible to characterise quantitatively the influence of the main process parameters of flash pyrolysis of thin wood particles seen here as the first step of the conversion of wood to syngas (H 2 and CO). Experiments in this field remain delicate and have to be repeated for a better accuracy. The quantification of gas species is highly repeatable, while that of char and tars needs to be improved.

At 850 °C, regardless of particle size between 80-125 and 160-200 μm, both the pyrolysis and tar cracking essentially happen in very short residence time period: less than 0.6 s; the products yields are only slightly modified after 0.6 s in the short residence times (several seconds) of our experiments, and each gas yields remain unchanged.

The stick-shaped wood particles sometimes keep their initial shape and shrink, and sometimes change to spherical hollow spheres with a diameter close to the initial diameter of the stick.

The amount of each gas species largely depends on the temperature level in the range 650 to 950 °C. Hydrogen concentration is largely improved by the temperature increase while CO decreases.

From the mass balance established after 2 s at 950 °C, 76% of the mass of wood is recovered as gases, 7% as char and 1.5% as tars. In our experimental conditions, the gas species taken alone (excluding char residue and tars) never reach the thermodynamic equilibrium, but tend towards it at the highest temperature.

The increase in the steam partial pressure favours the production of hydrogen, although yields remain low. Between 0 and 10% of steam, the H 2 mass fraction increases from 1.3 to 1.6%w and CO decreases from 50.0 down to 42.2%w. The presence of water induces a number of chemical reactions in the gas phase. The water-gas-shift reaction alone cannot explain the differences observed.

▪

  temperature: 650 to 950 °C; ▪ residence time: 0.5 to 2 s; ▪ size of wood particles: 80-125 μm and 160-200 μm; ▪ particle initial moisture: dried wood and 10.7% (natural moisture) ▪ nature of the atmosphere: N 2 and N 2 containing 10% H 2 O

Fig. 1 .

 1 Fig. 1. Schematic view of the Entrained Flow Reactor.

Fig. 2 .

 2 Fig. 2. Repeatability of pyrolysis experiments of wood at 900 °C at 50 (a) and 90 cm (b) from the injection point.

Fig. 3 .

 3 Fig. 3. Observation of initial and pyrolysed wood under electron microscope. (a) Initial wood; (b) pyrolysed wood under N 2 ; scale is 500 μm.

Fig. 5 .

 5 Fig. 5. Mass balance by element: C, H, O in % of the mass of initial biomass.

Fig. 6 .

 6 Fig.6. Gas species mass normalised to the initial mass of wood released after 2 s pyrolysis of two sizes of wood particles at 900 °C under N 2 .

Fig. 4 .

 4 Fig. 4. Mass balance (three experiments a, b, c) realised after wood pyrolysis at 900 °C under N 2 ; H 2 O is quantified with gas; NQ stands for Non Quantified mass.

Fig. 7 .

 7 Fig. 7. Time evolution of the gas composition (ppmv) during pyrolysis of wood at 850 °C under N 2 .

Fig. 8 .

 8 Fig. 8. Evolution of the gas species mass (normalised to the initial mass of biomass) versus temperature during pyrolysis of wood under N 2 and sampled at 80 cm downstream of the injection point (residence time ranged from 1.4 to 1.6 s).

Fig. 9 .

 9 Fig. 9. Comparison of H 2 /CO mole ratio measured in experiments of wood pyrolysis and computed by equilibrium calculations.

Fig. 10 .

 10 Fig. 10. Gas species mass (normalised to the initial mass of dry biomass) of dry (lines) and wet (dotted lines) wood pyrolysed at 850 °C under N 2 .

Fig. 11 .

 11 Fig.11. Difference of gas yield formed during wood pyrolysis at 900 °C between 10% and 0 % of steam concentration in nitrogen.

Table 1

 1 Proximate, ultimate, and moisture analysis of pine wood.

	Proximate analysis	wd%	Pine wood
	Ash		0.4
	Volatile matter		82.3
	Fixed carbon (by difference)		17.3
	Ultimate analysis	waf%
	C		45.46
	H		5.89
	O (by difference)		48.65
	N		b 0.1
	S		b 0.1
	Moisture	w%	10.7
	w%: mass percent.

wd%: mass percent dry. waf%: mass percent ash free.

Table 2

 2 The characteristics of experiments carried out under nitrogen and under steam.

	Influence of	Particle size μm	Flow of nitrogen l/min at STP	Flow of steam l/min at STP	Ratio Steam/N 2 vol.%	Temperature °C	Residence time s
	Residence time	160-200	18	0	0	900	0.5 to 1.5
	Temperature	160-200	23	0	0	650 to 950	1.4 to 1.6
	Particle size	80-125	18	0	0	900	1.9
		160-200	18	0	0	900	1.9
	Natural wood moisture	160-200	23	0	0	850	0.9; 1.3; 1.45
	Steam concentration	160-200	18	0	0	900	1.9
		160-200	14.2	3.8	10	900	1.9

Acknowledgement

We would like to thank the Institut Français du Pétrole, which has supported this research.


Frédéric Lebon
email: lebon@lma.cnrs-mrs.fr

Serge Dumont
email: serge.dumont@unimes.fr

Raffaella Rizzoni
email: raffaella.rizzoni@unife.it

Juan Carlos Lopez-Realpozo

Raul Guinovart-Diaz

Reinaldo Rodriguez-Ramos
email: reinaldo@matcom.uh.cu

Julian Bravo-Castillero

Frederico J Sabina

J C López-Realpozo

R Guinovart-Díaz

R Rodríguez-Ramos

Soft and hard anisotropic interface in composite materials

Keywords: Fibres, Interface/interphase, Parallelogram cell, Fibrous composites, imperfect contact, asymptotic homogenization method

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

In this work, micromechanical analysis method is applied to a periodic composite with unidirectional fibers and parallelogram cells. The analytical expressions of the homogenized elastic properties are calculated for two phase composite with hard and soft interfaces. The Asymptotic Homogenization Method (AHM), for two-phase fibrous periodic composites with imperfect adhesion and oblique cell is used for the calculation of the plane elastic effective coefficients. This contribution is an extension of previous works by the authors (Rodriguez-Ramos et al. 2011 [START_REF] Rodríguez-Ramos | Two analytical models for the study of periodic fibrous elastic composite with different unit cells[END_REF], Guinovart-Diaz et al. 2011 [START_REF] Guinovart-Díaz | Influence of parallelogram cells in the axial behaviour of fibrous composite[END_REF]), where only the perfect contact was considered for the antiplane problem. Besides, the present investigation is different of those of Lopez-Realpozo et al. 2011 [START_REF] López-Realpozo | Transport properties in fibrous elastic rhombic composite with imperfect contact condition[END_REF] and Rodriguez-Ramos et al. 2013 [START_REF] Rodríguez-Ramos | Different approaches for calculating the effective elastic properties in composite materials under imperfect contact adherence[END_REF] since the plane problem is solved for the calculation of the effective coefficients for composites with parallelogram cell. The novelty of the present work is that the imperfection of the interface in the composite with parallelogram cell is taken into account introducing two spring-type stiffnesses   n t K ,K for plane problems.

Using a classical approach [START_REF] Lebon | Modelling adhesion by asymptotic techniques[END_REF], the spring parameters can be identified from a three phase problem where the interphase coating the fiber is very thin. The paper is organized as follows. In the first part of the paper the derivation of the contact law mechanically equivalent to the interphase coating the fiber is reviewed on the basis of an energy method [START_REF] Rizzoni | Higher order model for soft and hard elastic interfaces[END_REF]. The method allows obtaining the spring-type interface law for a general anisotropic behavior of the interphase and the interface parameters   n t K ,K are explicitly given in terms of the elastic constants of the interphase.

The results of the micromechanical analysis presented in the second part of this paper are mainly focused on the impact of the arrangement of the fibers and the mechanic imperfection at the interface on the plane properties in the composites. Moreover, the theoretical approach is validated with some theoretical models.

Modeling of imperfect contact

assumed to be symmetric, with the minor and major symmetries, and positive definite.

The adhesive is assumed to be soft, i.e. ijkl ijkl C C    with ijkl C independent of  .

The adherents are subjected to a body force density The approach used in [START_REF] Rizzoni | Higher order model for soft and hard elastic interfaces[END_REF] to obtain the contact law is based on the fact that stable equilibrium configurations of the composite assemblage minimize the total energy:

        1 2 g 1 2 ijkl i,j k,l i i x ijkl i,j k,l i i x i i x S S ijkl i, j k,l x B 1 1 E u C u u f u dV C u u f u dV g u dA 2 2 1 C u u dV , 2                 
in the space of kinematically admissible displacements:

    3 u V u H S ; : u 0 on         ,
Where  

H S ;

  is the space of the vector-valued functions on the set S  which are continuous and differentiable as many times as necessary. Under suitable regularity assumptions, the existence of a unique minimizer u V   is ensured [START_REF] Ciarlet | Mathematical Elasticity. Volume I: Three-dimensional Elasticity[END_REF].

For the asymptotic analysis, it is convenient to introduce the following change of variables

    1 2 3 1 2 3 p : x , x , x z ,z ,z 
in the adhesive:

1 1 z x  , 2 2 z x  , 3 3 z x   , which gives 1 1 z x      , 2 2 z x      , 3 3 z x       . A change of variable     1 2 3 1 2 3 p : x , x , x z ,z ,z 
is also introduced in the adherents:

1 1 z x  , 2 2 z x  ,   3 3 z x 1 2 1    ,
where the plus (minus) sign applies whenever   1 2

x S x S     and one has

1 1 z x      , 2 2 z x      , 3 3 z x      .
After the change of variables p , the interphase occupies the domain

      3 1 2 3 1 2 B z ,z ,z : z ,z A, z 1 2      ,
and the adherents occupy the domains  

1,2 1,2 3 S S 1 2 1 i     . The sets       3 1,2 1 2 3 1 2 3 z ,z ,z : z ,z A,z 1 2       
are taken to denote the interfaces between B, and 1,2 S and Using the changes of variables given by p , p and denoting

1 u u p     and 1 û u p    
the displacement fields from the rescaled adhesive and adherents, respectively, the total energy takes the rescaled form:

        1 2 1 2 ijkl i i, , j j i k,,jl i i i, j z ijkl i i, , j j i k,,jl i i i, j z i i i z S S 1 33 3 3 ki k,3 i,3 ki i,3 k, ki i, k, 1 1 Ê u , u C u u f u dV C u u f u dV g u dA 2 2 1 ˆˆˆˆˆK u u 2K u u K u u 2                                      z B dV ,  
where the matrices jl K (with j,l 1, 2,3  ) are defined by the relation:

jl ki ijkl K : C  .
In view of the symmetry properties of the elasticity tensor ijkl

C it results that   T jl lj K K  , j,l 1, 2,3  .
Next, the existence of asymptotic expansions of the displacement fields with respect to the small parameter  is assumed:

  0 1 2 2 2 u u u u o       ,   0 1 2 2 2 ˆˆˆû u u u o       ,   0 1 2 2 2 u u u u o         .
Substituting these expansions into the rescaled energy, one obtains

      0 0 0 ˆÊ u ,u E u ,u o 1      , with         1 2 0 0 0 1 0 0 0 2 0 0 0 0 ijkl i i, , j j i k,,jl i i i, j z ijkl i i, , j j i k,,jl i i i, j z i i i z S S 1 33 0 0 ki k,3 i,3 z B 1 1 Ê u , u C u u f u dV C u u f u dV g u dA 2 2 1 ˆK u u dV . 2                The energy   0 0 0 Ê u , u is now minimized with respect to couples   0 0
û , u in the set of the rescaled admissible displacements:

        ^0 3 3 1,2 u 1,2 ˆV u, u H S ; B; : u 0 on , u u on .         
Using standard arguments, the following equilibrium equations are obtained:

  1,2 0 ijkl k,l i 1,2 , j C u f 0 in S   ,   1,2 0 ijkl k,l j i g C u n g on   ,   1,2 0 ijkl k,l j 1,2 C u n 0 on S     g u 1 , 2      ,   33 0 jk k,l ,3 K u 0 in B  ,   1,2 0 33 0 ijkl k,l ir r,3 1,2 Ĉ u K u on   ,
where n is taken to denote the outward normal. The first three equations are the equilibrium equations of the adherents at the order zero, with the suitable boundary conditions.

The remaining two equations imply the continuity of the traction vector, say

  0 1 , 2 0 3 3 i : C u i  
at the order zero across the rescaled interphase, i.e.

0 3 i 0   ,
where

    1 2
     indicates the jump in the quantity  across the rescaled interphase B . After integration of equation  

33 0 jk k,3 ,3 K u 0  with respect to 3 z , using the natural condition 1,2 0 33 0 i3kl k,l ir k,3 Ĉ u K u  on 1,2
 and the continuity of the displacement vector fields at the surfaces 1,2  , one obtains

0 3 3 0 3 i K u   ,
which represents the classical law for a soft interface. In [START_REF] Rizzoni | Higher order model for soft and hard elastic interfaces[END_REF] it is shown that this law can be rephrased in the composite's limit configuration obtained as 0   . In this case one obtains:

  0 3 3 0 3 i K u   ,
where    indicates the jump of  across the surface  defined on the adherents' limit configuration as 0   and 0 3 i  is the traction vector field defined on the same limit configuration.

Spring-like interface law for a plane anisotropic adhesive

The elastic coefficient of the adhesive's material enter the interface law through the matrix 33 K .

For the general case of anisotropic (triclinic) material, the matrix 33 K has components

C 0 C K 0 C 0 C 0 C            .
If the adhesive is monoclinic with symmetry plane orthogonal to 3 i , then 33 K is diagonal:

3 1313 33 monoclinic,i 2323 3333 C 0 0 K 0 C 0 0 0 C            .
The same holds if the adhesive is orthotropic and the coordinated planes are the symmetry planes of the material.

For an orthotropic adhesive, the non-vanishing coefficient in 3 33 monoclinic,i K can be expressed in terms of the engineering constants as follows [START_REF] Jones | Mechanics of Composite Materials[END_REF]:

1313 31 C G  , 2323 23 C G  , 12 21 3333 1 2 1 C E E      ,
where 31 G , 23 G are the shear moduli in the (3,1) and (2,3) planes, respectively

; i E , i 1,2,3  is the Young's modulus in the i-direction; ij  , i, j 1, 2,3  are the Poisson's ration,
and they satisfy the three reciprocal relations

ij ji i j E E    , i, j 1, 2,3  , i j  ,
1 2 3 1 2 E E E                .
If the material is transversely isotropic and 3 i is the direction of the axis of transverse symmetry, then 33 K is diagonal and 1313 2323 C C  . Finally, for an isotropic material, one has

33 isotropic 0 0 K 0 0 0 0 2                 , with ,   the Lamé constants. In terms of the engineering constants, E, ,G  one has G   ,      E 1 2 1 1 2           .

Contact law modeling a curved anisotropic thin adhesive

The spring-like interface laws obtained in the previous sections are local and thus they can be locally applied to curved thin adhesives. In particular, it is introduced a local reference system with axes parallel to the tangential t and the normal n direction of the adhesive midline let s denote the third direction, orthogonal to the adhesive's plane. Based on the previous analysis, one could propose a contact law of the form

      n n t t s s T u T K u T u                     
, with u the displacement, 0 3 T i   the traction vector and 33 K K  the spring stiffness, given as shown in the previous section. In particular, the case of K with constants coefficients would correspond to an adhesive material whose symmetry properties are independent of the directions, as represented in Figure 1.a.

However, one could also imagine an adhesive whose properties are "globally" constant, i.e.

they are fixed with respect to a global frame, a case schematized in Figure 1.b. In this second case, the matrix K entering the above interface law, which is written in the local frame, can be obtained by applying a transformation from the global reference frame properties "globally" constant.

General considerations for heterogeneous structures with imperfect contact

A two phase uniaxial reinforced material is considered here in which fibers and matrix have transversely isotropic elastic properties; the axis of transverse symmetry coincides with the fiber direction, which is taken as the 3 Ox axis. The fibers cross-section is circular.

Moreover, the fibers are periodically distributed without overlapping in directions parallel homogeneous medium (Fig. 2).

The overall properties of the above periodic medium are sought using the well-known AHM [START_REF] Bakhvalov | Homogenization Averaging Processes in Periodic Media[END_REF][START_REF] Pobedria | Mechanics of Composite Materials[END_REF], and the following considerations are assumed. Two variables are introduced, i.e. x and   y x . They are referred as the slow or macroscopic and fast or microscopic variables, respectively, where l L   is a small dimensionless parameter, L is a linear dimension of the body and l is the diameter of the unit cell.

Then, it follows that in terms of the fast variable y , the appropriate periodic unit cell S is taken as a regular parallelogram cell in the 1 2 y y  plane so that S S    , where the domain 2 S is occupied by the matrix and its complement 1 S (fiber) is considered by a circle of radius R and center at the origin O (Fig. 2). The common interface between the fiber and the matrix is denoted by  .

      ij ijkl k,l , C ,    x y y x y ,
where ij  , ij  are the stress and strain tensors respectively, related to the small parameter  . Assuming zero body forces the elastic equilibrium equation is

ij, j 0   in  ,
where the subscript comma denotes partial differentiation. The gradient equations, which are the strain-displacement equations

k l kl l k 1 u u 2 x x              ,
where

  i i u u , 
x y are the components of the mechanic displacement related to the small parameter  . Replacing and into a system of partial differential equations with rapidly oscillating coefficients can be obtained

      ijkl k,l , j C u , 0  y x y in  .
Equation represents a system of equations for finding i u  . For a complete solution, it is necessary to assign suitable boundary and interface conditions, for instance

0 i i u u  on u  ; 0 ij j i n S   on T  ,
where  imperfect interface. This may be due to the presence of an interphase but also due to interface bond deterioration caused for different reasons, such as fatigue damage or environmental and chemical effects. Within this approach, the composite is modeled as a two-phase material with imperfect interface conditions.

The vanishing value of n K  , t K  and s K  corresponds to pure debonding (normal perfect debonding), in-plane pure sliding, and out-of-plane pure sliding, respectively. The status of the mechanical bonding is completely determined by appropriate values of these constants.

For large enough values of the constants, the perfect bonding interface is achieved.

Using the vector notation and defining the mechanical displacement u , the traction vector T and the spring stiffness matrix K in the following manner,

n n n t t t s s s u T K 0 0 = u , = T , = 0 K 0 , u T 0 0 K                                  u T K
the mechanical imperfect condition ( [START_REF] Pobedria | Mechanics of Composite Materials[END_REF], [START_REF] Rodriguez-Ramos | Closed-form expressions for the effective coefficients of fibre-reinforced composite with transversely isotropic constituents-I. Elastic and square symmetry[END_REF]), in general, may be expressed as:

          1 1 2 0 , 1       T T T K u on  .
In the relation, the symbol  

    1 2 g g g  
indicates the jump in the quantity g at the common interface denoted by  between the fiber and the matrix. The superscripts in brackets   , 1,2    denote the matrix (1) and fiber (2), respectively. The magnitudes n u , t u , s u are the tangential and normal components of the mechanical displacement vector; n T , t T , s T are the tangential and normal components of the traction vector T  

I i j j T u n 
and n is the outward unit normal on  .

In order to study the imperfect contact conditions, the relations between the mechanic displacement and traction vectors are related to Cartesian mechanic displacement   i u and the traction   i T vectors respectively (Fig. 3) by the following expressions,

                                                                       
Fig. 3. The coordinate system in each point P   of the interface.

Two scale asymptotic homogenization method

By means of the asymptotic homogenization method, the original constitutive relations with rapidly oscillating material coefficients -are transformed into equivalent system

  * 0 ijkl l,kj C u 0  on  ,
with the boundary conditions

  0 0 i i u u  on u  ;   0 0 ij j i n S   on T  ,
with constant coefficients * C , which represent the elastic properties of an equivalent homogeneous medium. They are called as effective coefficients of the composite  .

In order to obtain average coefficients * C the periodic solutions of six pq L local problems on S in terms of the fast variable y , where p,q 1, 2,3  are found. Each local problem uncouples into sets of equations, i.e. plane-strain and antiplane-strain systems.

The pq L problem consists to find the displacements

    pq N  y in S  , 1, 2  
(double periodic functions with periods 1 w , 2 w ) as solution of the following system of partial differential equations,

  pq i , 0      in S  , (1)
where,

      pq i , i k pq k, C N          , (2)
the comma notation denotes a partial derivative relative to the y  component, i.e., , U U y

 
   the summation convention is also understood for Greek indices, which run from 1 to 2; no summation is carried out over upper case indices, whether Latin on Greek.

Thus, the eq. on  for the pq L problem can be expressed in the following indicial form,

    1 2 pq pq 0   T T , (3)
    1 pq n n pq n T 1 K N          ,     1 pq t t pq t T 1 K N          ,     1 pq s s pq s T 1 K N          , (4)
where pq n N , pq t N , pq s N and pq n T , pq t T , pq s T have the same meaning as , but adequate to pq L problems. To assure the only one solution of the pq L problems, the functions also satisfy that pq k N 0  . In the calculations the following relations are used:

  1 n 1212 n K C K R   ,   1 t 1212 t K C K R   and   1 s 1212 s K C K R  
where n K  , t K  and s K  are dimensionless imperfect parameters. The symmetry between the indices p and q shows right away that at most six problems needs to be considered. Once the local problems are solved, the homogenized moduli * ijpq C may be determined by using the following formulae

* ijpq ijpq ijkl pq k,l C C C N   .
The potential method of complex variables N  in S  of the problem defined by equations ( 1)-( 4) must be found among doubly periodic functions of half periods 1 w , 2 w (Fig. 2). Each local problem (1)-( 4) uncouples into sets of equations. An in-plane strain system for

  pq N   , 1, 2  
and an out-of-plane strain Laplace's equation for

  pq 3
N  has to be solved. Then the solution of the in-plane (outof-plane) strain problems involves the determination of the in-plane (out-of-plane) displacements, strains and stresses over each phase S  of the composite. Due to the non- vanishing components of the elastic tensors   ijpq C  , the only non-homogeneous problems, that have a non-zero solution, correspond to the four in-plane strain problems L  ,   1, 2,3   and 12 L , and the two out-of-plane strain ones 23 L and 13 L . Therefore, the solutions of both (in-plane and out-of-plane) local problems ( [START_REF] Rodriguez-Ramos | Closed-form expressions for the effective coefficients of fibre-reinforced composite with transversely isotropic constituents-I. Elastic and square symmetry[END_REF], [START_REF] Guinovart-Diaz | Closed-form expressions for the effective coefficients of fibre-reinforced composite with transversely isotropic constituents -I. Elastic and hexagonal symmetry[END_REF]) (six pq L problems need to be considered due to the symmetry between p and q ) lead to obtain the average coefficients of the composite given in Fig. 2. Imperfect contact for out-of-plane problem is well described in [START_REF] López-Realpozo | Transport properties in fibrous elastic rhombic composite with imperfect contact condition[END_REF]- [START_REF] Rodríguez-Ramos | Different approaches for calculating the effective elastic properties in composite materials under imperfect contact adherence[END_REF] where the imperfect parameter

  1 s 1212 s K C K R   is considered.
Therefore, from now on, only the solution of the plane local problems L  ,   1, 2,3   and 12 L are studied where the imperfect parameters Therefore, the methods of a complex variable z in terms of two harmonic functions

  1 n 1212 n K C K R   ,   1 t 1212 t K C K R   are involved..

Solution of plane problems L

  z   and   z  
, the Kolosov-Muskhelishvili complex potentials, are applicable. The potentials are related to the displacement and stress components by means of the formulae

                                    / 1212
2C N z iN z z z z z , 2 z z , 2i 2 z z z ,                                                    where     3 4       ,         1122 1111 1122 C . C C        A relevant representation of the complex potentials   z   and   z   of periods 1 w and 2 w is       p p k * * * * 1 0 p k pk 2 k p 1 k 1 p 1 k 1 p p k * * * * 1 0 p k pk 2 k p 1 k 1 p 1 k 1 z R z z z a a a , z c , R z R R z R z z (z) b b b , z d . R z R R                                                                                 
After some algebraic manipulations, a system of infinite equations is obtained for

                                    a a a a a a a a a a where   2 a R      a ,   12 12 a R m  a ,     2 * 1212 1 1212 C C   , 1 2 2 1 1 1 2 1 2 w w A w w        , 1 2 2 1 2 1 2 1 2 w w A w w w w      , , 1 2 2 1 3 1 2 1 2 w w A w w w w      , 1 2 2 1 4 1 2 1 2 w w A w w w w      ,   (1) 2 2 2 1p 3 1p 1 1p 1 1p Ĥ BA R A A R CR        ,   2 2 2 2p 4 1p 1p 1 1p 1 Ĥ A BR A R ACR       ,   2 0 kp 1p 1 1p k1 k kp ˆŴ AR C B Cr       ,       0 n p 2 0 kp 1p 1 1p k1 kp p k 2p k p 2 0 p K B ˆM AR C BG D A             ,             * * 2 3 n 1 3p 1p 1p 1 1p 0 p 2 2 K 2 p ˆĤ R A R P A                  ,   * 1 0 4p p 1p n 2 Ĥ 1 A K             ,   * 1 * 0 p n 2 p B 1 A K            , * (1)
* (2) * n 0 1 4 K ˆĈ B 2            ,         (2) * 4 0 n 0 2 ˆB 1 2B P 1 ˆ2 K                      ,     p k p kp k 2p p k p k k p 2 G p 2 k kR C T            ,   * kp i 2 p k i 2 i 1 r         ,     1 * n t k 0 k k K K 1 D D E     ,           1 * ( 2 ) * * n t n p 0 * (1) * p n p K K 1 K 2 p 2 D A K 1 2 p E                     ,         (2) * * n t n p 1 0 * p n p K K K 2 p 2 B B 1 K E                        ,     1 * n t k 0 k k K K 1 B C 1 E     ,         * * * p n t n t n t A K K p 2 p 1 K K p K K p 2                     ,     * * p t n n t B K K 1 K K p      ,         2 * * p n t n t n t C K K K K p K K p 2              ,   * p n t D K K p    , p p p p p A B E C D   ,         2 2 3 * 2 0 3 (1) Re A R 1 1 i I mA R 1 2                ,     2 3 * 2 0 3 (1) Re A R 1 1 ˆ1 Re A R ˆ1 2 F                   , (2)
* n 1 F 1 4 K      ,           1 2 2 1 1 2 2 2 2 1 1 1 1 2 C C C C           ,           1 2 1 2 2 1 1 1 1 2 2 2 2 2 C C C C              ,       2 2 3 2 2 1 1 2 C C       y       2 2 4 2 2 1 1 2 C C        .
We can write from the expression,

C C C C C C N N , C C 2iC C C C C N N 2i C                                2,1 1,2 * 33 33 1133 1,1 2233 2,2 N N , C C C N C N .                     * *
C C N C N , C C N C N , C C N C N , C C C N C N ,              
and the following effective coefficients can be calculated using the average operator, applying the Green theorem and considering

        2 1 1 , 1,2,3            ,     1111 1122 2k C C      and   1212 m C    we obtain                   * ( 1 ) * 2 1 2 1111 1111 2 11 1 1 1 n 1 k 1 m 2V k k k C C V k Re 1 2 Re 1 , m k K m k                            a                   * ( 1 ) * 2 1 2 2211 2211 2 11 1 1 1 n 1 k 1 m 2V k k k C C V k Re 1 2 Re 1 , m k K m k                            a               1 * * 2 1133 2 3311 1133 2 1133 1 1 n 1 1 2V k C k C C V C k Re 2 Re 1 m K m k                        ,     * ( 1 ) 1211 2 11 1 C V k Im 1        a ,                   * ( 1 ) 1122 1122 2 22 1 1 * 2 1 2 1 n 1 k 2 m C C V k Re 1 m k 2V k k k 2 Re 2 K m k                             a ,                   * ( 1 ) 2222 2222 2 22 1 1 * 2 1 2 1 n 1 k 2 m C C V k Re 1 m k 2V k k k 2 Re 2 K m k                             a ,               1 * * 2 1 1 3 3 2 3322 3322 2 2233 1 1 n 1 2 2V k C k C C V C k Re 2 Re 2 m K m k                        ,       * ( 1 ) 1222 2 22 1 C V k 1 Im    a ,                 * ( 1 ) 1133 1133 2 1133 33 1 1 2 2 1 1 3 3 1 * 1133 1 n 1 1133 k 3 C C V C Re 1 m 2V C k C 2 Re 3 , K m C                                a                 * ( 1 ) 2233 2233 2 2233 33 1 1 2 2 2 2 3 3 1 * 2233 1 n 1 2233 k 3 C C V C Re 1 m 2V C k C 2 Re 3 , K m C                                a             1 2 2 * * 2 1133 1133 3333 3333 2 1133 1 1 n 1 1133 3 2VC C C C V C Re Re 2 3 m K m C                     ,      * ( 1 ) 1233 2 1133 33 1 C V C 1 Im    a ,             (2) * 2 2 * ( 1 ) 2 2 1112 2 12 1 2 1 n1 k 4k m V C m V Re 1 Re m Km                    a ,             (2) * 2 2 * ( 1 ) 2 2 2212 2 12 1 2 1 n1 k 4k m V C m V Re 1 Re m Km                     a ,        (2) * * 2 2 2 2 3312 2 2 1 n1 4 m l V C m l V Re Re m Km               ,       1 * ( 1 ) 1212 1212 2 12 1 C C m V Im 1         a , where                 2 1 * * * 2 0 k1 k 1 1 k 1 1 2 * 2 1 2 * * * 2 0 k1 k 1 1 0 k 1 2 C 1 A R B 1 , ˆ2 1 2 C P 1 A R 1 2 ˆB B                                                                                          a a a a             2 1 * * 2 * 0 1 1 2 1 k 11 2 k k 1 2 2 * 2 1 * * 2 * 0 1 1 2 1 k 11 2 k k 1 2 C 1 A R B 1 2 1 2 C 1 A R B                                                                                a a a a
.

Numerical results

As a validation of the present model in Fig. 4 The material parameters used in the calculations are given in Table 1 2010 [START_REF] Otero | Semi-analytical method for computing effective properties in elastic composite under imperfect contact[END_REF].

In Tables 2345we consider the effective coefficients of the composite. Table 2 shows the behavior of the composite with perfect contact (

Conclusions

Two-scale asymptotic homogenization technique was used for effective coefficients calculation of an circular elastic fiber reinforced composite with imperfect bonding



  are assumed to be located far from the interphase and the fields of the external forces are sufficient regularity to ensure the existence of equilibrium configuration.

    is called the rescaled configuration of the composite body. Lastly, u  and   indicate the images of u

Fig. 1 .

 1 Fig. 1. (a) Adhesive properties independent of the directions n and t . (b) Adhesive

  complex numbers, which define the parallelogram periodic cell of the two phase composite. Therefore, the composite  consists of a parallelogram array of identical circular cylinders embedded in a

Fig. 2 .

 2 Fig. 2. Fibrous composite material and regular parallelogram cell with circular cell

  prescribed displacement, force on the boundary of the composite. An adhesive with monoclinic symmetry and constant elastic moduli is considered in the homogenization analysis via a layer of mechanical springs of zero thickness. The spring constants 33 K K  are the measures for the magnitude of the associated continuities on  , where n K  , t K  and s K  are the diagonal spring constant parameters, which have dimension of stress divided by length, as seen in previous sections. These constants are called the interface parameters and they are written with another notation for the monoclinic case in the diagonal matrix 33 K given in page. 7. It is seen that infinite values of the parameters imply vanishing of displacement jumps and therefore perfect interface conditions. At the other extremity, zero values of the parameters imply vanishing of interface tractions and therefore disbond. Any finite position values of the interface parameters define an

 and 12 L

 12 . Effective coefficientsThe aim of this work is to obtain the coefficients L for composites with parallelogram arrays.

Fig. 4 .

 4 Fig. 4. Comparison between the AHM-spring model and differential scheme reported in[START_REF] Guinovart-Díaz | Effective elastic properties of a periodic fiber reinforced composite with parallelogram-like arrangement of fibers and imperfect contact between matrix and fibers[END_REF].

Table 1 .

 1 . Material constants of the constituents used in the composite taken from Otero et al.

		1111 C	1122 C	1133 C	1212 C
	Matrix -Al	94.23	40.38	40.38	26.92
	Fibra -SiC	483.68 99.07	99.07	192.31

Table 2 .

 2  show the performance of the composite for different values of in-plane imperfect parameters for diverse volume fractions and different angle of the cell. Similar effect of the angle of the cell in the value of the effective coefficients is observed. Behavior of the composite with perfect contact (

	n K 10 	6	and	t K 10 	6	), different volume

Table 3 .

 3 Comportment of composite for different values of in-plane imperfect parameters for diverse volume fractions and different angle of the cell.

				n K 10 	6	,	t K 50 
	Vf		0.05	0.2		0.35	0.55	0.75
	90	0	1.054961 1.256744 1.533481 2.090600 3.088266
	* 1111 C
	0 75	1.054739 1.252523 1.518770 2.048378 2.985805
	90	0	1.042007 1.165717 1.283256 1.446467 1.837369
	* 1122 C
	0 75	1.042570 1.176434 1.320523 1.544647 1.951860
	90	0	1.022294 1.100148 1.200096 1.391695 1.747719
	* 1133 C
	0 75	1.022300 1.100262 1.200482 1.391651 1.731405

Table 4 .

 4 Comportment of composite for different values of in-plane imperfect parameters for diverse volume fractions and different angle of the cell.

					n K 10 	,	t K 50 
	Vf			0.05	0.2		0.35	0.55	0.75
	1111 C *	90	0	1.026128 1.118309 1.233205 1.418272 1.628926
		0 75	1.025962 1.115538 1.225215 1.404864 1.625590
	1122 C *	90	0	0.988111 0.926652 0.822606 0.628191 0.400022
		0 75	0.988509 0.933313 0.841944 0.661499 0.411230
	1133 C *	90	0	0.991006 0.962850 0.932785 0.889275 0.841009
		0 75	0.991004 0.962814 0.932658 0.888904 0.840381
	3333 * C	90	0	1.194071 1.776466 2.359162 3.136628 3.914840
		0 75	1.194071 1.776466 2.359161 3.136618 3.914793
		90	0	1.051650 1.215785 1.412449 1.803788 2.699235
	* 1212 C
		0 75	1.052204 1.223664 1.433339 1.833875 2.647238

Table 5 .

 5 Comportment of composite for different values of in-plane imperfect parameters for diverse volume fractions and different angle of the cell.

Acknowledgment

The authors gratefully acknowledge to the project supported by the Cooperation Scientifique Franco-Cubaine (FSP) program "Soft and Hard Interfaces in Composites Homogenized Analytically and Numerically" project N° 29935XH. Also the authors wish to thank the CNRS, Ecole Centrale de Marseille, Aix-Marseille University and to the project Composite Materials from University of Havana, Cuba. The funding of CONACYT project number 129658 is gratefully acknowledged. RR thanks the financial support from ENDIF through the "Fondo per necessità di base della ricerca-2014". This work was partially written while R. Rodríguez-Ramos was visiting the Computer Institute of Universidade Federal Fluminese (UFF) (Project Edital PI-UFF 2014) supported by the PROPPI and Department of International Relations -DRI -of this university.

.

between constituents and parallelogram periodic cell. Analytical expressions are obtained for all effective coefficients considering two different n K and t K (normal and tangential) values of imperfect mechanic parameters. The numerical computations for all effective coefficients are shown for different values of the aforementioned imperfect parameters and inclination angle of the cell. The influence in the effective coefficient behavior of the imperfect adhesion and the inclination angle of the cell is significant and this is reported. In order to generalize the present work to more general fibers shape, the study of fibers with elliptic shape could be the object of a future work.


Serge Dumont
email: serge.dumont@u-picardie.fr

Frédéric Lebon

Effective properties of linear random materials: application to Al/SiC and resin/glass composites

Keywords: Random materials, Wavelet-element method, Homogenization, Elasticity

Here we present a statistical study on the effective linear properties of random materials, i.e., microstructures which are random lattices described by a stochastic process. The local numerical procedure associated with the homogenization methods used here was based on a wavelet-element method. The resulting numerical results are compared with those obtained using classical theories. A new approach was developed in order to determine the effective properties in cases where the characteristics of the microstructure are not known.

Introduction

The aim of this study was to determine and compute the effective properties of linear random heterogeneous media. Although natural geomaterials are heterogeneous, they are generally regarded as homogeneous materials with effective properties which determine the overall behavior. Grounds, for example, are usually modeled as a body consisting of homogeneous materials with these effective properties. The behavior of these grounds cannot be predicted if the effective properties are not exactly determined. The need to incorporate more detailed information about the microstructure is now widely recognized. During the last few years, many studies have been published presenting measurements and describing empirical and theoretical techniques for determining these effective properties. The numerical cost of solving this problem directly for each material sample would obviously be prohibitive. We previously presented a method involving a numerical scheme combined with homogenization techniques which can be used to account for the details of microstructures [START_REF] Dumont | Representation of plane elastostatics operators in Daubechies wavelets[END_REF][START_REF] Dumont | Wavelet-Galerkin method for periodic heterogeneous media[END_REF][START_REF] Dumont | Wavelet-Galerkin method for plane elasticity[END_REF][START_REF] Dumont | A numerical tool for periodic heterogeneous media. Application to interface in Al/SiC composites[END_REF][START_REF] Lebon | Small parameter, homogenization and wavelets[END_REF]. These algorithms based on wavelet analysis make it possible to predict the effective properties from microstructural data and thus to optimally design linear composite materials from real images of the sample.

Unfortunately, in most cases, the microstructural details are not exactly known, and the effective properties are therefore assessed from statistical studies. In a well-known paper [START_REF] Torquato | Random heterogeneous media: microstructure and improved bounds on effective properties[END_REF], Torquato (see also [START_REF] Torquato | Random heterogeneous materials[END_REF]) presented results obtained on a case where some statistical informations, such as the correlation functions, was available. This improves the bounds of the effective properties of two-phase random heterogeneous media. Results of other kinds have been published in [START_REF] Kanit | Apparent and effective physical properties of heterogeneous materials: representativity of samples of two materials from food industry[END_REF][START_REF] Sab | On the homogenization and the simulation of random materials[END_REF]. Here it is assumed that only the volume fraction of each phase is known, and it is proposed to find the average effective properties of the composite. In addition, the microstructure is assumed to be given pixel by pixel, i.e., the characteristic shape of the microstructural components is taken to be a small square and its characteristic length to be the length of the pixel.

We will not discuss here the size of the representative volume, which is often a real problem with materials of this kind. We will assume the Hill's condition to be valid (see [START_REF] Hill | Elastic properties of reinforced solid: some theoretical principles[END_REF][START_REF] Wriggers | Aspects of the computational testing of the mechanical properties of microheterogeneous material samples[END_REF][START_REF] Wriggers | Introduction to computational micromechanics[END_REF]) and that exists a representative volume on which we will work. For further details on the problem of choosing a sample when there exists no representative volume, we refer Fig. 1 The homogenization process to [START_REF] Hazanov | Order relationships for boundary condition effect in hetegeneous bodies smaller than representative volume[END_REF][START_REF] Huet | Application of variational concepts to size effects in elastic heterogenous bodies[END_REF][START_REF] Huet | A method of substructuring large-scale computational micromechanical problems[END_REF][START_REF] Kanit | Determination of the size of the representative volume element for random composite: statistical and numerical approach[END_REF]. In addition, the method proposed in previous papers is improved here. In [START_REF] Dumont | Representation of plane elastostatics operators in Daubechies wavelets[END_REF] only Daubechies wavelets are used. In this paper B-spline wavelets are introduced. It will be shown below how these functions improve the numerical procedure.

In Sect. 2, the notations used are presented and the mechanical problem is defined. Section 3 is devoted to the numerical algorithm. In Sect. 4, a statistical study on AlSiC and resin-glass composites is presented. The results obtained are compared with the Hashin and Shtrikman bounds [START_REF] Hashin | A variational approach to the theory of the elastic behaviour of multiphase materials[END_REF] and with the self-consistent schema [START_REF] Hill | A self-consistent mechanics of composite materials[END_REF][START_REF] Peyroux | Caractéristiques thermoélastiques de matériaux composites à fibres courtes[END_REF]. At the end of the paper, the future perspectives of this approach are discussed.

The mechanical problem

Here it is proposed to study the behavior of a heterogeneous medium consisting of a two-phase isotropic elastic composite. A homogenization step is first performed, consisting in determining the effective properties of the equivalent homogeneous medium. Let us take a plane periodic composite Ω (macro-scale, variable x) and a rectangular periodicity cell Y (micro-scale, variable y). The bonds at all the interfaces of the composite are assumed to be perfect.

The elasticity coefficients C ε i jkl are assumed to be periodic, with period Y . We take the boundary ∂Ω of Ω to be split into two parts:

∂Ω = ∂ 1 Ω ∪ ∂ 2 Ω with ∂ 1 Ω ∩ ∂ 2 Ω = ∅.
A traction force equal to F is imposed on ∂ 1 Ω, and the displacement vanishes at ∂ 2 Ω (see Fig. 1).

In the elastostatic framework, the problem (P ε ) can be written

Problem(P ε ) Find u ε ∈ V such that a ε (u ε , v) = L(v), ∀ v ∈ V , where a ε (u, v) = Ω C ε i jkl e kl (u)e i j (v)dx, e kl (u) = 1 2 (u k,l + u l,k ), L(v) = Ω f vdx + ∂ 1 Ω
Fvdl,

and V = {v ∈ (H 1 (Ω)) 2 , v = 0 on ∂ 2 Ω}.
When ε → 0, it has been established (see for example [START_REF] Bensoussan | Asymptotic analysis for periodic structures, 1st edn[END_REF][START_REF] Sanchez-Palencia | Homogenization techniques for composite media[END_REF]) that the solution of the previous problem u ε tends to u in V , the space of the admissible global displacements, where u is the solution of problem P Problem(P)

Find u ∈ V such that a(u, v) = L(v), ∀ v ∈ V , where a(u, v) = Ω A i jkl e kl (u)e i j (v)dx and A i jkl = 1 meas(Y ) Y (C i jkl + C i j pq e pq (u kl ))dy.
To determine the displacements u kl , it is necessary to solve three local problems on Y (Problem P Y ) in the space H = (H 1 per (Y )) 2 , the set of admissible local displacements:

Problem(P Y ) E kl being given, find u kl ∈ H such that a Y (u kl , v) = l(v), ∀v ∈ H, ( 1
)
where

a Y (u, v) = Y C klst (y)e kl (u)e st (v)dy and l(v) = - Y C klst (y)E kl e st (v)dy.
Comment 1 Even if each component of the composite is isotropic so that the elasticity coefficients can be written using only the bulk and shear moduli K and µ, the homogeneous equivalent material is not necessarily isotropic and the elasticity tensor cannot be written using bulk and shear moduli.

Numerical procedure

In order to avoid having to remesh the cell Y after each draw, we decided to discretize problem P Y with a meshless method, the wavelet element method [START_REF] Daubechies | Orthonormal bases of compactly supported wavelets[END_REF][START_REF] Dumont | Representation of plane elastostatics operators in Daubechies wavelets[END_REF]. Moreover, as we will see below, the elements of the stiffness matrix can be computed with this method once and for all, without having to use quadrature formulae. Thanks to the localization of the wavelets in both the physical and spectral spaces, the stiffness matrix is not ill conditioned. This gives a good rate of convergence for solving the linear systems [START_REF] Dumont | Wavelet-Galerkin method for periodic heterogeneous media[END_REF]. This method can also be said to be a compromise between finite elements [START_REF] Garboczi | An algorithm for computing the effective linear elastic properties of heterogeneous materials; threedimensional results for composites with equal phase Poisson ratio[END_REF] and Fourier analysis [START_REF] Michel | A computational method for linear and nonlinear composites with arbitrary phase contrast[END_REF] methods.

Periodic conditions were applied here to the boundary of the cell Y , but boundary conditions of other kinds can be used with wavelet transforms, such as vanishing displacements [START_REF] Dumont | Meshfree methods and boundary conditions[END_REF].

Daubechies wavelets

In this section, we will briefly describe Daubechies wavelets [START_REF] Daubechies | Orthonormal bases of compactly supported wavelets[END_REF].

Let M be a positive integer. There exist a sequence of real numbers (h n ) n for n = 0, . . . , 2M -1, a function ϕ called the scaled function and a function ψ called the associated wavelet, such that

ϕ(x) = √ 2 2M-1 n=0 h n ϕ(2x -n), ∀x ∈ R, ψ(x) = √ 2 2M-1 n=0 g n ϕ(2x -n), ∀x ∈ R, (2)
where

g n = (-1) n h 2M-1-n , n = 0, . . . , 2M -1.
The functions ϕ and ψ have the following properties:

-ϕ and ψ are orthogonal, since

(ϕ, ψ) = R ϕ(x)ψ(x)dx = 0.
-The supports of ϕ and ψ lie in the closed interval [0, 2M -1].

Let

ϕ jk (x) = 2 j/2 ϕ(2 j x -k); V j = span{ϕ jk , k ∈ Z}; ψ jk (x) = 2 j/2 ψ(2 j x -k); W j = span{ψ jk , k ∈ Z}, (3)
the functions ϕ jk and ψ jk fulfill:

R ϕ jk (x)ϕ j (x)dx = δ k , ( 4
) R ψ jk (x)ψi (x)dx = δ k δ i j . ( 5
)
Relations (4) mean that (ϕ jk ) k is an orthonormal basis of V j for each j. Relations (5) mean that (ψ jk ) k is an orthonormal basis of W j for each j, and spaces W j and W i are orthonormal if j = i. We therefore have The set of spaces {V j , j ∈ Z} is called a multiresolution analysis of L 2 (R).

V j ⊂ V j+1 , V j+1 = V j W j , j∈Z = {0}, L 2 (R) = j∈Z V j .
Figure 2 gives plots of the scale function and the associated Daubechies wavelet of order M = 3.

We can then obtain a multiresolution of L 2 per (0, 1), the space of the periodic functions with a period equal to one, by periodizing the scaled function and the associated wavelet:

φ jk (x) = 2 j/2 n∈Z ϕ(2 j (x + n) -k), j ≥ 0, k = 0, . . . , 2 j -1; (6)
and likewise for the periodic associated wavelet.

We then perform a multiresolution of L 2 per (Ω), where Ω = (0, 1) × (0, 1) with a tensor product. The periodic two-dimensional scale function becomes [START_REF] Daubechies | Orthonormal bases of compactly supported wavelets[END_REF] and the three associated two-dimensional wavelets are:

Φ jk (x) = φ jk 1 (x 1 )φ jk 2 (x 2 ), ∀k = (k 1 , k 2 ) ∈ Z 2 , x ∈ Ω,
Ψ 1 jk (x) = φ jk 1 (x 1 )ψ jk 2 (x 2 ), ∀k = (k 1 , k 2 ) ∈ Z 2 , x ∈ Ω, Ψ 2 jk (x) = ψ jk 1 (x 1 )φ jk 2 (x 2 ), ∀k = (k 1 , k 2 ) ∈ Z 2 , x ∈ Ω, Ψ 3 jk (x) = ψ jk 1 (x 1 )ψ jk 2 (x 2 ), ∀k = (k 1 , k 2 ) ∈ Z 2 , x ∈ Ω. ( 8
)
If we take Λ j to denote the subset of Z 2

Λ j = k ∈ Z 2 , 0 ≤ k 1 ≤ 2 j -1, 0 ≤ k 2 ≤ 2 j -1 , ( 9
)
and {V j , j ≥ 0} the multiresolution of (L 2 per (Ω)) 2 , an element u in V j can be written

u(x) = u 1 (x) u 2 (x) with u d (x) = k∈Λ j u k d Φ jk (x) ∀x ∈ Ω.
Comment 2

-The Daubechies wavelets are orthogonal for the scalar product of (L 2 per (Ω)) 2 , and they are only implicitly given by the sequence {h n } n=0,...,2M-1 (see 2).

-Unfortunately, the length of the support is very large.

For example, if M = 3, the size of the support is in the [0, 5] range. These values can be compared with those obtained using P2 finite elements, where the support is in the [0, 2] interval. The number of non vanishing elements in the rigidity matrix will therefore be larger with the Daubechies wavelet method. -These functions are not very regular and they are not symmetric (see Fig. 2). The consequences of these properties will be seen in Sect. 4.

B-spline wavelets

In order to avoid the problems due to the orthogonality of wavelets, bi-orthogonal wavelets have been introduced [START_REF] Cohen | Biorthogonal bases of compactly supported wavelets[END_REF]. This procedure was based on a multiresolution analysis. It can be explained as follows:

-Define two sequences of embedded subspaces:

. . . ⊂ V -2 ⊂ V -1 ⊂ V 0 ⊂ V 1 ⊂ V 2 ⊂ . . . . . . ⊂ Ṽ-2 ⊂ Ṽ-1 ⊂ Ṽ0 ⊂ Ṽ1 ⊂ Ṽ2 ⊂ . . .
-Define W j ( W j resp.) the complementary space of V j into V j+1 ( Ṽj into Ṽj+1 resp.), such that

V j ⊥ W j , Ṽj ⊥ W j .
-We therefore have two scale functions (ϕ and φ) and two wavelets (ψ and ψ) such that

R ϕ(x) φ(x -k)dx = δ 0k , ∀k ∈ Z, R ψ(x) ψ(x -k)dx = δ 0k , ∀k ∈ Z, R ϕ(x) ψ(x -k)dx = 0, ∀k ∈ Z, R φ(x)ψ(x -k)dx = 0, ∀k ∈ Z, (10)
-A function f is decomposed into

s j 0 k = R f (x)ϕ j 0 k (x)dx, ∀k ∈ Z, d jk = R f (x) ψ jk (x)dx, ∀ j ≥ j 0 , k ∈ Z and f (x) = k∈Z s j 0 k φ j 0 k (x) + j≥ j 0 k∈Z d jk ψ jk (x).
-Therefore, four filters are necessary:

ϕ(x) = √ 2 N n=0 h n ϕ(2x -n); φ(x) = √ 2 Ñ n=0 hn φ(2x -n); ψ(x) = √ 2 N n=0 g n ϕ(2x -n); ψ(x) = √ 2 Ñ n=0 gn φ(2x -n). ( 11
)
Comment 3

-The scale function ϕ is often a B-Spline [START_REF] Cohen | Biorthogonal bases of compactly supported wavelets[END_REF][START_REF] Deslaurier | Une famille d'ondelettes bi-orthogonales sur l'intervalle obtenue par un schéma d'interpolation itérative[END_REF]. In this case, we define the B-spline scale function of order N as N convolutions of the characteristic function on [0, 1]. The B-spline of order 2 is the P1 finite element on R, and the B-spline of order 3 is plotted in Fig. 3. This function is differentiable, piecewise polynomial and of order two. The support of this function lies in the [0, 3] range, and this function is symmetric. -The choice of Ñ , the size of the filter h, is much more difficult. The larger Ñ is, the more regular φ and ψ will be. But the support of these functions will also increase.

In order to give all four functions ϕ, φ, ψ and ψ the

Operator discretization

In this section, we will focus on solving problem P Y . This problem is solved here using a wavelet-element method [36] with Daubechies wavelets [START_REF] Daubechies | Orthonormal bases of compactly supported wavelets[END_REF] or B-Spline wavelets [START_REF] Cohen | Biorthogonal bases of compactly supported wavelets[END_REF][START_REF] Deslaurier | Une famille d'ondelettes bi-orthogonales sur l'intervalle obtenue par un schéma d'interpolation itérative[END_REF]. This method is based on the standard form of the discretized elastostatics operator [START_REF] Beylkin | On the representation of operators in bases of compactly supported wavelets[END_REF]. Details of this method can be found in our previous papers [START_REF] Dumont | Wavelet-Galerkin method for periodic heterogeneous media[END_REF][START_REF] Dumont | A numerical tool for periodic heterogeneous media. Application to interface in Al/SiC composites[END_REF][START_REF] Dumont | Wavelets and PDE's: on two applications in solids mechanics[END_REF].

In what follows, without any loss of generality, we will consider a two dimensional problem using the plane deformations hypothesis. The stress and the deformation tensors will therefore be taken to be vectors of R 3 :

σ = ⎛ ⎝ σ 11 σ 22 σ 12 ⎞ ⎠
and so on for the deformation tensor ε. We will also identify the elasticity operator with a 3 × 3 matrix denoted C on the representative volume and A on the homogeneous equivalent material, such that

σ = Cε.
More specifically, the solution u of P Y is approximated by

u(y) = k∈Z 2 u jk Ψ jk (y) y ∈ Y, ( 12
)
where Ψ jk (y) = 2 j Ψ (2 j y -k) is a localized function around the point y = 2 -j k. Since the function Ψ is compactly supported, the sum can be written in a finite dimensional subspace of Z 2 , denoted Λ j . Using test functions of this kind, the projection of the linear elasticity tensor

a(u, v) = Y Ce(u) : e(v)dy is a matrix K defined by [K] = [K i ] i ∈Λ j ,
where K i is a 2 × 2 elementary matrix

K i = F 11
C(y) = n∈Z 2 C n Θ jn (y),
where Using the scale relation ( 2), these coefficients are computed once and for all without any approximation, after diagonalizing a small dimensional matrix [START_REF] Beylkin | On the representation of operators in bases of compactly supported wavelets[END_REF][START_REF] Dumont | Representation of plane elastostatics operators in Daubechies wavelets[END_REF][START_REF] Dumont | Wavelet-Galerkin method for periodic heterogeneous media[END_REF]. Due to the particular type of wavelet basis involved (see formula 5), the discrete stiffness matrix takes the form plotted in Fig. 4. When it is necessary to solve a linear system with a matrix of this kind, a multigrid method gives excellent results.

Θ jn (y) = 2 j θ(2 j y 1 -n 1 )θ (2 j y 2 -n 2 ).

Comment 4

In what follows, the image of the microstructure will be composed of square pixels (as in Fig. 1), and this explains the use of Haar wavelets for the discretization of the elasticity coefficients C αβ . One can also imagine that these coefficients could be given after performing digitization and filtering, using other types of wavelet based on a camera image. Note that the Jpeg2000 norm is based on the use of Daubechies wavelets [START_REF] Christopoulos | Jpeg2000: the new still picture compression standard[END_REF][START_REF] Taubman | JPEG2000 Image compression fundamentals. Standards and practices, series: the international series in engineering and computers science[END_REF].

In this case, the approach developed here can be extended only by replacing the Haar basis by the wavelet used in the filtering step.

Numerical results

In this part of the paper, numerical tests on Al/SiC and resin/ glass composites are presented. The effects of two parameters are studied: the volume fraction of each component (Aluminum and Silicium or Resin and Glass) and the distribution of each material in the microstructure. With each volume fraction, approximately 5000 draws are performed and results are studied statistically.

The elementary representative volume will be represented by an image with 32 × 32 pixels as proposed by Wriggers et al. [START_REF] Wriggers | Aspects of the computational testing of the mechanical properties of microheterogeneous material samples[END_REF]. The discretization of this domain will be carried out with 256 × 256 wavelet elements. The numerical results obtained in this study were not significantly different from those obtained using finer discretizations with 512 × 512 or 1024 × 1024 wavelet elements, for example.

The pixels corresponding to the second material are chosen randomly, using a uniform law.

The aim of this paragraph is to show that, even if after each draw s the homogeneous equivalent material with elasticity tensor equal to A s is not isotropic, the mean elasticity tensor on the draws à = A s s becomes isotropic and one can define "mean" bulk and shear moduli K and μ from Ã.

In [START_REF] Wriggers | Aspects of the computational testing of the mechanical properties of microheterogeneous material samples[END_REF], a different approach was used. For each sample s, the authors computed the effective bulk and shear moduli:

3K * s = trσ s /3 Ω trε s /3 Ω and 2µ * s = σ s Ω : σ s Ω ε s Ω : ε s Ω ,
where σ = σ -(trσ/3)I d and ε = ε -(trε/3)I d. They then studied the properties of the mean values of these quantities in the samples: K * = K * s s and μ * = µ * s s .

AlSiC composite

First we focus on AlSiC composites having the following properties: the Young's modulus is equal to E = 0.41 • 10 +6 MPa and the Poisson's ratio is equal to ν = 0.19 in the case of silicium, and the Young's modulus is equal to E = 0.72 • 10 +5 MPa and the Poisson's ratio is equal to ν = 0.32 in that of aluminum. We consider a two-dimensional material, and we assume the plane stress hypothesis to be valid.

The changes in each elasticity coefficient are presented in Fig. 5 in the case of about 5000 microstructures at a volume fraction of 70%. The fluctuations can be seen to be very large.

On the other hand, if we look at the histogram of each component of the elasticity tensor, we can see that the values show a Gaussian distribution (Fig. 6).

The mean values on the nth first draws :

A i j (n) = 1 n n s=1 A s i j
versus the number of draws are given in Fig. 7, for each coefficient i j of the elasticity tensor. Coefficients A i j (n) in the case of all i, j = 1, 2, 3 can be seen here to converge on a mean value when the number n of draws increases. If we study in detail the changes in the mean coefficients, we can see that the mean is obtained in quite a small number of draws (Fig. 8) and that the standard deviation is small in comparison with the coefficients. In the case presented in Fig. 8, i.e., at a volume fraction equal to 70%, the mean is equal to 12358 MPa and the standard deviation is equal to 768 MPa. The values are in the [12083,12682] range. In addition, after 2000 draws, the maximum standard deviation amounted to less than 0.01% of the coefficient. Therefore, 5000 draws suffice to obtain an accurate mean coefficient. The histogram in terms of polar coordinates (Fig. 9) shows the good distribution of the results in terms of the mean values. The values given in Fig. 7 show that the numerical limit is that of an isotropic material. In particular, we have A 13 = A 23 = 0 and A 11 = A 22 . It is therefore possible to compute the bulk and shear moduli of the limit material (Mean bulk and shear moduli).

In what follows, we take L-HS and U-HS to denote the Lower and Upper Hashin and Shtrikman bounds [START_REF] Hashin | A variational approach to the theory of the elastic behaviour of multiphase materials[END_REF]. We take WAV to denote the Mean Value obtained using the Wavelet-element method. Figure 10 shows three phases: in the first phase, at volume fractions lower than 0.25, the L-HS and WAV values are very similar. The second phase corresponds to volume fractions of 0.25 to 0.95. In this phase, the WAV value is mid-way between the two HS bounds. In the third phase, at volume fractions above 0.95, the WAV values are similar to the values of U-HS bounds. If we compare the results obtained with the self-consistent schema on the shear modulus and at volume fractions lower than 0.5, we can see from Fig. 11 that at the lower volume fraction values (lower than 0.1), the L-HS value, the results based on the self-consistant schema (denoted SCS) and the WAV value are fairly similar. At higher volume fraction values ranging from 0.1 to 0.25, the L-HS and WAV values are practically identical and the SCS value is lower than the other two. At volume fractions above 0.25 and below 0.4, the L-HS and SCS values are similar, and lower than the WAV value.

Resin-glass composite

In the resin-glass composites studied here, the Young's modulus is equal to E = 0.73 • 10 +5 MPa and the Poisson's ratio is equal to ν = 0.25 in the case of the glass, and the Young's modulus is equal to E = 0.30 • 10 +4 MPa and the Poisson's ratio is equal to ν = 0.40 in that of the resin. We assume the existence of a two-dimensional material, and the plane stress hypothesis is assumed to be valid. In this case, the difference between the properties of the two composites is greater than previously. The histograms plotted in polar coordinates in Figs. 11 and12 show the good distribution of the mean results obtained at various glass and resin volume fractions (20% in Fig. 12 and 80% in Fig. 13).

Figure 14 shows that the behavior of this composite was similar to that observed in the case of the AlSiC composite: the relative behavior can again be split into three phases. In the first part, at a volumic fraction below 0.35, the L-HS and WAV values were similar. The second phase corresponds to volume fractions ranging between 0.35 and 0.99. In this phase, the WAV value is midway between the two HS bounds. In the third phase, the WAV value tangents the U-HS bounds.

Comparison between Daubechies wavelets and B-Spline wavelets

We end this section with a discussion about the wavelets. 1 gives the bulk and shear moduli obtained with Al/SiC and resin/glass composites, using Daubechies and B-Splines wavelets. The volume fraction of the first component of the composite was equal to 20%, and the moduli were obtained by averaging on 100 draws.

In the case of the AlSiC composite (low contrast), the results given by the two discretization procedures (Daubechies and B-spline wavelets) were very similar (difference amounted to less than 3%); whereas with the resin/glass composite (large contrast), since the differences between the properties of the two materials were larger, the two discretization procedures did not give similar results. The differences increased by 7% in the bulk modulus and by 2% in the shear modulus. Moreover, the bulk modulus given by the Daubechies wavelets was outside the bounds, due to the lack of symmetry and regularity of Daubechies wavelets. Using the symmetric and C 1 B-spline wavelets, better results were obtained and the average bulk modulus computed was very near the lowest Hashin and Shtrikman bound (see Figs. [START_REF] Dumont | Meshfree methods and boundary conditions[END_REF][START_REF] Garboczi | An algorithm for computing the effective linear elastic properties of heterogeneous materials; threedimensional results for composites with equal phase Poisson ratio[END_REF]. We recall that the moduli have to be within the Hashin and Shtrikman bounds only in the asymptotic case of an infinite (sample length)/(particle length) ratio and purely isotropic macroscopic responses [START_REF] Hashin | A variational approach to the theory of the elastic behaviour of multiphase materials[END_REF]. This was not the case here where the ratio was equal to 32.

In addition, the results obtained with the B-spline wavelets were always within the Voigt and Reuss bounds, contrary to the results obtained with the Daubechies wavelets. This shows that the use of B-spline wavelets improves the numerical procedure.

Comment 5 Unlike the finite element method, the wavelet element method can be used with a single tool for the image processing, the image compression and the determination of the effective properties. In particular, this method does not require the use of Voronoi techniques. However, the solution of the linear systems can be compute less expensively using the finite element method [START_REF] Dumont | Wavelet-Galerkin method for periodic heterogeneous media[END_REF].

Concluding comments and perspectives

It is generally impossible to exactly determine the effective properties of random heterogeneous media. The most rigorous statements in the literature about these effective properties focus on the bounds. In this study, bounds of another kind were obtained in the form of an average and a standard deviation. One of the main results obtained here is the numerical convergence between the stochastic homogenization process with a mean solution, i.e., when the number of draws tends to infinity, the limit of the potential energy corresponding to a 'random lattice' is the potential energy of an isotropic material. The effective properties of this material depend on the volume fraction of each component. At lower volume fractions, the material has similar properties to those given by the lower Hashin and Shtrickman bounds and the self-consistent schema. A numerical theorem is obtained. This result is consistent with those obtained in [START_REF] Kanit | Determination of the size of the representative volume element for random composite: statistical and numerical approach[END_REF] using the finite element method and with an initial result obtained at low contrast in [START_REF] Lebon | Novel approaches in civil engineering[END_REF]. Note that this theorem is a generalization of rigorous mathematical stochastic homogenization results obtained in another context (scalar problems) [START_REF] Andrews | Stochastic homogenization of elliptic boundary-value problem with l p -data[END_REF][START_REF] Ay | Darcy equation for random porous media[END_REF].

The method developed here, based on a wavelet transform procedure is highly efficient and makes it possible to avoid meshing the microstructure after each draw. In addition, since the wavelet transform method was developed for image processing purposes, it can be used to determine the microstructure with a good resolution and a reasonable number of degrees of freedom.

This study could be extended in two ways: first, by determining the effects of the shape of the lattices. We propose to study the same problem investigated here with different morphologies: n-pixels in one direction, sets of n-pixels in the neighborhood of a point, etc. Another possible approach would consist in studying the effective properties in the case of more complex behavior such that involved in piezoelectric coupling, thermomechanical coupling and non linear behavior [START_REF] Lebon | Effective properties of nonlinear laminated composites with perfect adhesion[END_REF]29,32]. This method lends itself to dealing with problems of this kind because it has proved to be a particularly efficient tool.

Fig. 2

 2 Fig. 2 The scale function and the associated Daubechies wavelet of order M = 3

3

  The B-Spline wavelet of order M = 3 same regularity, one can take Ñ > 1.5N . Examples of functions ϕ, φ, ψ and ψ are given in[START_REF] Cohen | Biorthogonal bases of compactly supported wavelets[END_REF].

Fig. 4

 4 Fig.[START_REF] Beylkin | On the representation of operators in bases of compactly supported wavelets[END_REF] Particular form of the stiffness matrix discretized with a wavelet method

Fig. 5 Fig. 7 Fig. 8

 578 Fig.[START_REF] Christopoulos | Jpeg2000: the new still picture compression standard[END_REF] Homogenized plane elasticity tensor versus draws (volume fraction 70%). Each peak gives the result of one test

Fig. 9

 9 Fig. 9 Histogram in polar coordinates (volume fraction 70%, AlSiC)

Fig. 10 Fig. 11

 1011 Fig. 10 Bulk and shear moduli and Hashin-Shtrikman bounds (AlSiC)

Fig. 12 Fig. 13

 1213 Fig.[START_REF] Dumont | A numerical tool for periodic heterogeneous media. Application to interface in Al/SiC composites[END_REF] Histogram in polar coordinate system (volume fraction 20%, resin/glass composite)

Fig. 14

 14 Fig. 14 Bulk and shear moduli and Hashin-Shtrikman bounds (resin/glass composite)

Fig. 15

 15 Fig. 15 Bulk moduli and Hashin-Shtrikman bounds obtained with Daubechies and B-spline wavelets (resin/glass composite)

Table 1

 1 Comparison between the bulk and shear moduli (×10 3 Pa) obtained using Daubechies and B-Splines wavelets, in Al/SiC (low contrast) and B-Splines (large contrast), with a volume fraction of the first component equal to 20%

			Bulk mod.	Shear mod.
	Al/Sic	Daubechies	80.6	36.7
		B-Splines	83.0	36.9
	Resin/glass	Daubechies	5.83	1.74
		B-Splines	6.25	1.70
	Differences	Al/SiC	3%	0.5%
		Resin/glass	7%	2%
	Results are averages based on 100 draws
	Table

Acknowledgments

The authors gratefully acknowledge R. Peyroux for helpful discussions, enlightening comments and self-consistent schema software.