M Louërat

M Ayouz

P Perré

Heat and moisture diffusion in spruce and wood panels computed from 3-D morphologies using the Lattice Boltzmann method

Keywords: anisotropy, effective value, fiberboard, parallel computing, wood, X-ray micro-tomography

Heat and moisture diffusion in spruce and wood panels computed from 3-D morphologies using the Lattice Boltzmann method M. Louërat, M. Ayouz, P.

Introduction

The use of bio-based materials in building construction is increasing because of their renewability and their insulation properties. As they are hygroscopic, understanding their hygrothermal behavior is necessary for appropriate building design. Over the last decade, it became obvious that the computational thermal codes used in buildings have to account for the coupling between heat and mass transfer [START_REF] Woloszyn | Tools for performance simulation of heat, air and moisture conditions of whole buildings[END_REF][START_REF] Delgado | A critical review of hygrothermal models used in porous building materials[END_REF]. These requirements gave rise to several experimental and theoretical works intended to deeply describe the effects of coupled heat and mass transfer in hygroscopic materials. These works are performed at the scales of wall assembly or building envelope [START_REF] Labat | Dynamic coupling between vapour and heat transfer in wall assemblies: Analysis of measurements achieved under real climate[END_REF][START_REF] Rafidiarison | Dataset for validating 1d heat and mass transfer models within building walls with hygroscopic materials[END_REF][START_REF] Zhang | Moisture buffering phenomenon and its impact on building energy consumption[END_REF][START_REF] Rouchier | Hygric characterization of wood fiber insulation under uncertainty with dynamic measurements and markov chain monte-carlo algorithm[END_REF] These new requirement of building simulation are likely to benefit from the mature domain of drying and processing of hygroscopic materials [START_REF] Perré | Multiscale modeling of drying as a powerful extension of the macroscopic approach: application to solid wood and biomass processing[END_REF][START_REF] Perré | A comprehensive dual-scale wood torrefaction model: Application to the analysis of thermal run-away in industrial heat treatment processes[END_REF]. With the increasing performances of the building envelops, accurate characterization of building materials in terms of heat and mass transfer is crucial. This can be done by clever characterization protocols. [START_REF] Lelievre | Hygrothermal behavior of bio-based building materials including hysteresis effects: Experimental and numerical analyses[END_REF] investigated heat and moisture transfer within hemp concrete by accounting for hysteresis and phase change effects. Moisture transport in plywood under isothermal and non-isothermal conditions has been measured by Glass [START_REF] Glass | Measurements of moisture transport in wood-based materials under isothermal and non-isothermal conditions[END_REF]. Effective thermal conductivity and mass diffusivity of bio-based materials are among the important properties engineers need to know to predict the energy consumption of buildings. An essential factor affecting such properties is the micro-structure of materials, which is particularly complex and difficult to characterize in the case of wood and derived materials. Therefore, a prediction approach that can compute the macroscopic properties from the medium morphology is of great interest.

Tomography allows a non-destructive description of the internal structure of a material. By its feature of non-destructive imaging and its ability to observe without the need to prepare a surface, the use of X-ray computed micro-tomography became a classical facility to observe various aspects of wood structure, for example to quantify anatomical features or determine density [START_REF] Steppe | Use of x-ray computed microtomography for non-invasive determination of wood anatomical characteristics[END_REF][START_REF] Van Den Bulcke | X-ray tomography as a tool for detailed anatomical analysis[END_REF][START_REF] Huang | Porosity estimation of phyllostachys edulis (moso bamboo) by computed tomography and backscattered electron imaging[END_REF], to observe musical instrument [START_REF] Osborne | Estimating the density of coast douglas-fir wood samples at different moisture contents using medical x-ray computed tomography[END_REF], to observe archaeological wood [START_REF] Stelzner | X-ray computed tomography for the anatomical and dendrochronological analysis of archaeological wood[END_REF] or to assess wood decay (Van den [START_REF] Van Den Bulcke | Three-dimensional imaging and analysis of infested coated wood with x-ray submicron ct[END_REF][START_REF] Hervé | Density mapping of decaying wood using x-ray computed tomography[END_REF]. [START_REF] Lux | 3D characterization of wood based fibrous materials: an application[END_REF] studied the morphology of wood-based fibrous materials using X-ray tomography.

In the last decades, the possibility to obtain sub-micrometric resolutions, at first using synchrotron facilities and more recently thanks to laboratory nano-tomographs, offered brand new possibilities. Indeed, with submicrometric resolution, the morphology of secondary cell wall can be clearly observed and quantified [START_REF] Trtik | 3D imaging of microstructure of spruce wood[END_REF][START_REF] Brodersen | Visualizing wood anatomy in three dimensions with high-resolution x-ray micro-tomography (µct)-a review[END_REF].

Computed tomography with very high resolution (sub-micrometric resolution) was used by [START_REF] Standfest | Determination of the Microstructure of an Adhesive-Bonded Medium Density Fiberboard (MDF) using 3-D Sub-micrometer Computer Tomography[END_REF] and Walther and Thoemen [START_REF] Walther | Synchrotron x-ray microtomography and 3d image analysis of medium density fiberboard (mdf)[END_REF] to investigate an industrially produced Medium Density Fiberboard (MDF). Beyond imaging, these tools are also used for in situ experiments, either with mechanical or hydric loading [START_REF] Derome | Hygromorphic behaviour of cellular material: hysteretic swelling and shrinkage of wood probed by phase contrast x-ray tomography[END_REF][START_REF] Zauner | Synchrotron-based tomographic microscopy (sbtm) of wood: development of a testing device and observation of plastic deformation of uniaxially compressed norway spruce samples[END_REF][START_REF] Baensch | Damage evolution in wood: synchrotron radiation micro-computed tomography (srµct) as a complementary tool for interpreting acoustic emission (ae) behavior[END_REF][START_REF] Zauner | Failure and failure mechanisms of wood during longitudinal compression monitored by synchrotron microcomputed tomography[END_REF] Therefore, it appears to be a great tool for micro-structural characterization of wood-based materials, and we used high-resolution X-ray microtomography in this paper. Using such 3-D description at high resolution involves a huge amount of information and complex geometrical shapes. Suitable computational methods are therefore needed to predict the equivalent macroscopic properties.

Based on the kinetic theory of gases, the Lattice Boltzmann method (LBM) is a computational method widely used in fluid dynamics [START_REF] Shan | Lattice boltzmann model for simulating flows with multiple phases and components[END_REF][START_REF] Chen | Lattice Boltzmann Method for Fluid Flows[END_REF][START_REF] Succi | The Lattice Boltzmann Equation: For Fluid Dynamics and Beyond[END_REF]. Unlike the traditional CFD methods, which solve the conservation equations of macroscopic properties after discretization, the LBM is an ascendant method where the macroscopic behavior of a fluid or a solid emerges from the discrete movement (consecutive propagation and collision processes) of particles [START_REF] Succi | The Lattice Boltzmann Equation: For Fluid Dynamics and Beyond[END_REF]. The LBM has become a general numerical method to solve partial differential equations. As no connectivity is needed in the grid used for the resolution, the LBM is usually classified in the family of so-called meshless methods [START_REF] Belytschko | Meshless methods: An overview and recent developments[END_REF][START_REF] Frank | The potential of meshless methods to address physical and mechanical phenomena involved during drying at the pore level[END_REF] with a local formulation. In the domain of heat or mass transfer, the scientific community remains very active in methodological developments able to open new configurations to LBM [START_REF] Hamila | Enthalpic lattice Boltzmann formulation for unsteady heat conduction in heterogeneous media[END_REF][START_REF] Chen | A lattice Boltzmann model for heat transfer in heterogeneous media[END_REF][START_REF] Demuth | Application and accuracy issues of TRT lattice Boltzmann method for solving elliptic PDEs commonly encountered in heat transfer and fluid flow problems[END_REF] Considering those characteristics, the LBM provides several advantages: simple numerical development, suitability for parallel computing and easy processing of complex morphologies. We would just emphasize here that the choose of LBM was guided by these advantages in using simply the 3-D morphologies of the materials. Other computational strategies to solve the solution such as Finite Elements or Control Volumes would have produced similar results.

For these properties, the LBM has been used to predict effective conductivity and diffusivity of heterogeneous media. Wang and Pan [START_REF] Wang | Predictions of effective physical properties of complex multiphase materials[END_REF] calculated the thermal conductivities of moist porous brick sands, in both frozen and unfrozen states. [START_REF] Wang | Thermal conductivity enhancement of carbon fiber composites[END_REF] focused on the effective thermal conductivity of carbon fiber composites. Xuan et al. [START_REF] Xuan | Investigation on mass diffusion process in porous media based on Lattice Boltzmann method[END_REF] investigated the effective mass diffusivity in porous media. [START_REF] Yablecki | Modeling the Effective Thermal Conductivity of an Anisotropic Gas Diffusion Layer in a Polymer Electrolyte Membrane Fuel Cell[END_REF] determined the anisotropic and heterogeneous effective thermal conductivity of the gas diffusion layer of a polymer electrolyte membrane. [START_REF] Walther | Cement paste morphologies and effective diffusivity: using the Lattice Boltzmann method[END_REF] predicted the effective diffusivity of cement paste during hydration.

However, in most studies, the effective properties of heterogeneous materials are calculated on model media or reconstructed morphologies and not on real morphologies (stochastic generation growth method (Wang et al., 2007b), uniformly distributed inclusions and randomly generated fibrous medium [START_REF] Jeong | Estimation of thermal and mass diffusivity in a porous medium of complex structure using a lattice Boltzmann method[END_REF], stochastic fibers placement algorithm [START_REF] Yablecki | Modeling the Effective Thermal Conductivity of an Anisotropic Gas Diffusion Layer in a Polymer Electrolyte Membrane Fuel Cell[END_REF] and correlated random field [START_REF] Walther | Cement paste morphologies and effective diffusivity: using the Lattice Boltzmann method[END_REF]). Wulf et al. proposed a full 3D approach applied to metallic foams [START_REF] Wulf | Experimental and numerical determination of effective thermal conductivity of open cell FeCrAl-alloy metal foams[END_REF]. [START_REF] Perré | New modelling approaches to predict wood properties from its cellular structure: imagebased representation and meshless methods[END_REF] studied both thermal and mass diffusivities in real morphologies of wood but only in two dimensions (2-D) whereas it is necessary to know the material properties in all directions because of their anisotropy.

In this paper, we investigate both heat and mass transfers in bio-based materials using the LBM to predict the effective thermal conductivity and mass diffusivity. The calculations are done on real 3-D morphologies obtained by micro-tomography, in each direction. The next section of this paper describes the materials studied, the imaging technique used, synchrotron X-ray micro-tomography and image processing. Then, we present the Lattice Boltzmann method in the case of diffusion and describe the computation of equivalent (or effective) properties. The methodology to reduce calculation time (choice of numerical parameters, determination of a convergence criterion and parallel computing) is explained in Section 4. Finally, the results for both heat and mass transfers are presented in Section 5.

Material morphologies

Several bio-based building materials are studied in this paper: spruce (Picea abies) earlywood (EW) and latewood (LW), Low-Density Fiberboard (LDF) and Medium-Density Fiberboard (MDF). These products are all woodbased materials but present a wide range of properties as a result of their internal structure (porosity, pore connectivity, directional orientation). Therefore, depending on their properties, these materials have different roles in building, such as mechanical structure, furniture, partition and insulation. These materials were chosen because they have the same solid phase, which will allow us to concentrate on the effect of morphology on macroscopic properties in the present work.

Sample preparation

To obtain 3-D images at high resolution, we had to prepare small cylindrical samples of our materials. Depending on the toughness of the product, we used two different methods. For spruce, we used a wood-turning machine to make cylinders and carefully reduce their diameter down to a few millimeters. This method does not work for LDF and MDF due to their low cohesion. Instead, we used a laser cutting device to prepare, without any mechanical stress, cylinders of 10 mm and 5 mm, respectively. However, we had to face another problem: for LDF, the density is so low, hence the air circulation so easy, that the samples burnt during laser cutting. To avoid this, we pre-cut larger samples with a hole saw and placed them in an open-top aluminum container free of oxygen thanks to an intense nitrogen flux during laser cutting.

Microtomography

X-Ray micro-tomography is a non-destructive and non-invasive 3-D imaging technique. It is used to characterize material micro-structure at a micron level spatial resolution. It has become a relatively common tool in materials science and is suitable for bio-based materials [START_REF] Perré | A Review of Modern Computational and Experimental Tools Relevant to the Field of Drying[END_REF].

Classical X-ray tomography is based on the attenuation of an X-ray beam by matter. A sample is scanned at different rotation angles by X-ray beams and the attenuation of the X-ray beam intensity is measured by a 2-D detector, which acquires a collection of projections. Then an algorithm reconstructs the 3-D object from the 2-D radio-density images.

Our samples have been scanned at ESRF (European Synchrotron Radiation Facility, line ID19) by Novitom, a company specializing in advanced analytical imaging powered by synchrotron technology. The tomography analysis mode used is phase contrast. In contrast to classical tomography, it is not the attenuation of the X-rays that is measured but the beam phase shift (transformed into variations in intensity) that is recorded by the detector. The exposure time is about 0.1 s per projection. With a total of 6000 projections, the time required per reconstructed volume is ca. 10 seconds. For wood samples, 2 or 3 volume were needed to scan the whole sample in the longitudinal direction. All scan were done with samples equilibrated with the surrounding air, i.e. at a moisture content of about 10%. We scanned at different resolutions: 3.05 µm for LDF,0.97 µm for MDF and 0.62 µm for spruce. The dimensions of the scanned volume are 12 × 12 × 13 mm 3 for LDF, 4.4 × 4.4 × 2.1 mm 3 for MDF and 1.4 × 1.4 × 2.6 mm 3 for spruce. The scanned spruce volume contains both early-and latewood.

Image processing and 3-D morphologies

Image processing of the reconstructed datasets is needed to define the 3-D morphologies of materials in a form suitable for the LBM code. The image processing is performed with Fiji, which is an open source image processing package based on ImageJ. For spruce (late-and earlywood), the sample is rotated to align the material directions of wood (longitudinal, radial and transverse) along the x-, y-and z-axes. The reconstructed datasets are cropped to have 200×200×200 voxels for fibers and around 100×150×400 voxels for spruce. To have enough tracheids with only 100 × 150 × 400 voxels, we used a 2×2×2 binning factor on earlywood. We apply a global threshold-based segmentation on each sample, the threshold level is determined by the moment method. Finally, the morphology is defined as a file with the coordinates and the corresponding grey level (0 for the solid phase or 1 for the gaseous phase). The four bio-based materials that have been scanned in this study exhibit contrasted morphologies, both in terms of porosity and in terms of pore connectivity (Fig. 1).

The solid fraction, s , is 0.29 for earlywood and 0.74 for latewood. This factor of ca. 2.5 is explained by two well-known cumulative effects: the larger cell wall thickness and the smaller radial extension of tracheids in latewood compared with earlywood. Considering a solid density of 1500 kg m -3 for the cell wall, the densities of the samples are 435 kg m -3 for earlywood and 1110 kg m -3 for latewood. These results are in close agreement with literature data [START_REF] Koubaa | Defining the transition from earlywood to latewood in black spruce based on intra-ring wood density profiles from X-ray densitometry[END_REF][START_REF] Trtik | 3D imaging of microstructure of spruce wood[END_REF].

For wood panels, we found a porosity value of 0.18 for LDF and 0.43 for MDF. With the same density of the cell wall, the calculated density values are respectively 270 kg m -3 and 645 kg m -3 . As the density of such panels depends mainly on the manufacturing conditions, it is difficult to compare with literature data. Instead, we directly measured their density by measuring the mass and volume of samples. We found 156 kg m -3 for LDF and 612 kg m -3 for MDF which is close to the computed densities. The density difference for LDF can be explained by the large heterogeneity inside the materials and the lower resolution of the tomography scan (3.05 µm for LDF instead of less than 1 µm for the other materials). This good result is an indirect validation of the whole procedure used here to define the 3-D morphology. Indeed, the excellent resolution, good contrast and low noise provided by synchrotron X-ray tomography ease the subsequent image processing and allow the solid phase to be accurately distinguished from the gaseous phase. We must keep in mind that such a segmentation is much difficult on classical 2-D images:

• for optical images, the resolution is too poor for the segmentation to be accurate. A calibration of the method is required in this case [START_REF] Gonçales | Mesure de la densité des bois sur coupe microscopique par analyse d'image : application quelques essences de la forêt amazonienne[END_REF],

• for scanning electron microscopy images, the resolution is much better, but several artifacts need to be addressed carefully: edge effects on the cell wall border due to electrical charges, sample deformation due to dehydration of the cell wall (except with ESEM [START_REF] Almeida | Free shrinkage of wood determined at the cellular level using an environmental scanning electron microscope[END_REF]).

The Lattice Boltzmann method

Presentation of the method

In this method, discrete and local equations are solved on a set of particles placed on a lattice. Particles are placed on the lattice nodes and can move in one time step towards the neighboring nodes according to discrete velocities. The macroscopic equation of interest depends on the rules applied on the lattice. At each time step, the particles exchange their energy by collision and stream along all possible directions to one of the nearest-neighboring nodes. The probability of a particle streaming in the i-direction is given by a distribution function f i that obeys the following equation:

f i ( x + c i δt, t + δt) -f i ( x, t) = Ω i ( x, t), (1)
where c i is the discrete particle velocity vector and Ω i represents the collision term. For this operator, we use the BGK (Bhatnagar-Gross-Krook) approximation [START_REF] Bhatnagar | A Model for Collision Processes in Gases. I. Small Amplitude Processes in Charged and Neutral One-Component Systems[END_REF]:

Ω i ( x, t) = -ω [f i ( x, t) -f eq i ( x, t)] , (2)
where ω = 1/τ is the relaxation factor, i.e. the inverse of the dimensionless relaxation time and f eq i ( x, t) is the equilibrium distribution function at site x and time t. Eq. ( 2) means that, after collision, the distribution function relaxes towards a discrete approximation of the equilibrium Maxwell-Boltzmann distribution in a single relaxation time.

LBM for the diffusion problem

The standard 3-D diffusion equation for the normalized mass or temperature field is:

∂Ψ ∂t = α ∂ 2 Ψ ∂x 2 + ∂ 2 Ψ ∂y 2 + ∂ 2 Ψ ∂z 2 , (3)
where α (assumed constant) is the local thermal or mass diffusion coefficient and Ψ is the normalized macroscopic variable (temperature field θ or mass X). For thermal diffusion, α = λ/(ρc p ) with λ the thermal conductivity, ρ the density and c p the specific heat capacity. We use the LBM to solve this equation in a 3-D parallelepiped porous medium

[x 1 , x m ] × [y 1 , y n ] × [z 1 , z p ].
The form of the equilibrium distribution function varies according to the problem. For a diffusion problem, it can be written as [START_REF] Mohamad | Applied lattice boltzmann method for transport phenomena, momentum, heat and mass transfer[END_REF]:

f eq i = w i Ψ, (4)
where w i are the constant weighting factors relative to the propagation direction, such that

i w i = 1.
In this paper, we use two velocity models: D 2 Q 5 for 2-D and D 3 Q 7 for 3-D. Indeed, it has been shown that five and seven discrete velocities, for 2-D and 3-D problems respectively, are sufficient to recover the diffusion equation with good accuracy [START_REF] Fang | Predictions of effective thermal conductivities for three-dimensional four-directional braided composites using the lattice Boltzmann method[END_REF] 2,3,4,5,6 c(0,0,0) c(±1,0,0),c(0,±1,0),c(0,0,±1) 1/4 1/8 4

. The D 2 Q 5 model is schematized Model i c i w i a D 2 Q 5 0 1, 2, 3, 4 c(0, 0) c(±1, 0), c(0, ±1) 1/3 1/6 3 D 3 Q 7 0 1,
Table 1: Characteristics of the velocity models used, c = δx δt = δy δt = δz δt with δx, δy and δz the space steps and δt the time step.

in Fig. 2. The characteristics of the velocity models used in this studyvelocities c i , weighting factors w i and the factor a of Eq. ( 5) -are given in Tab. 1.

In the literature, it has been proved that the diffusion equation can be derived from the LBM equation using the Chapman-Enskog expansion [START_REF] Higuera | Lattice gas method based on the Chapman-Enskog expansion[END_REF]. The relationship between macroscopic and mesoscopic scales is the following:

α = δx 2 aδt ( 1 ω - 1 2 ), (5)
where a is a factor depending on the velocities model; its values are given in Tab. 1. For stability reasons, ω must obey the following condition, 0 < ω < 2 [START_REF] Succi | The Lattice Boltzmann Equation: For Fluid Dynamics and Beyond[END_REF]. Moreover, the macroscopic normalized quantities, Ψ and related flux Φ, can be computed as follows:

Ψ = i f i , (6)
Φ = ( 1 ω - 1 2 )ω i e i f i , (7)
where the set of vectors e i are the unit velocity vectors.

Computation of equivalent properties using the LBM

Eqs. ( 1) and ( 2) represent the two main steps of the LBM, collision and streaming. They are illustrated in Fig. 2 in the case of a D 2 Q 5 velocity model. Then, boundary conditions are applied to complete one time step.

To compute the equivalent thermal conductivity of a heterogeneous medium, a constant macroscopic gradient between two opposite sides (x = x 1 and x = x m = L x ) is applied. Initially, a constant field is assumed: Ψ(x, y, z, 0) = 0. Then, we apply the Dirichlet boundary conditions in the x-direction and adiabatic boundary conditions in the other directions: Ψ(x 1 , y, z, t) = 1 and Ψ(x m , y, z, t) = 0 ; Φ(x, y, z, t). n = 0 at y = y 1 , y = y n , z = z 1 and z = z p .

In the LBM, boundary conditions are applied thanks to conditions imposed on f i . To implement the Dirichlet boundary conditions, the distribution functions on the boundaries are supposed to be at equilibrium and we compute the unknown f i according to Eq. ( 6). For adiabatic conditions, we use the free-slip boundary conditions [START_REF] Succi | The Lattice Boltzmann Equation: For Fluid Dynamics and Beyond[END_REF], where the unknown f i are determined by axial symmetry. For example, for a D 2 Q 5 model, we have on the side x = x 1 , f 1 = Ψ(x 1 , y, t) -(f 0 + f 2 + f 3 + f 4 ) and on the side y = y 1 , f 3 = f 4 . At the corners, we impose a combined Dirichlet/Neumann boundary condition: the degrees of freedom at each node (the distribution values f i ) allow us to impose both the temperature value and a zero flux. For example, at corner x = x 1 /y = y 1 , we first set f 3 = f 4 and then

f 1 = Ψ(x 1 , y, t) -(f 0 + f 2 + f 3 + f 4 ).
To apply the LBM on heterogeneous domains, two relaxation factors, ω s and ω g , are used for the collision step. They are related to the diffusivities of the solid phase α s and the gaseous phase α g through Eq. ( 5). According to this equation, the relaxation factor of the medium having the highest diffusivity is the lowest, hence named ω min . The relaxation factor of the other medium is therefore named ω max . In this work, we follow the conjugate condition [START_REF] Chen | A note on the solution of conjugate heat transfer problems using SIMPLE-like algorithms[END_REF] together with the "half lattice division" scheme proposed in Wang et al. (2007a). The conjugate condition implicitly assumes the product ρc p to be the same for the two media. This is not a problem in the present work because the equivalent properties are determined from steady-state fields.

To compute the equivalent macroscopic properties, the different stages of the LBM are iterated until the steady state is reached. The equivalent thermal conductivity and equivalent mass diffusivity are then computed as the average x-flux density, Φ x , divided by the macroscopic gradient imposed by the Dirichlet boundary conditions, ∆Ψ/L x :

Γ eq = Φ x L x ∆Ψ , (8)
where Γ eq represents λ eq or f eq for thermal conduction or mass diffusion, respectively.

Code validation

A Fortran90 LBM code was developed at LGPM (CentraleSupélec) to model the thermal and mass diffusion in a heterogeneous medium in 2-D and 3-D. To validate the code, we ran it for different simple morphologies (phases placed in parallel and series) and compared the computed equivalent conductivity with the analytical value. In the paper, we use the solid and gas thermal conductivities of the bio-based materials studied: λ s = 1 W m -1 K -1 and λ g = 0.023 W m -1 K -1 . For example, in the case of a D 2 Q 5 model, we get for parallel and series geometries λ = 0.51 W m -1 K -1 and λ ⊥ = 4.5 10 -2 W m -1 K -1 with a relative error lower than 10 -10 compared with analytical solutions.

We also ran the code on a more complex shape, with interfaces not parallel to the grid. The results are presented in Sec. 4.1. The error is lower than 1% compared with the reference value. Therefore, the code was successfully validated for both simple configurations and for morphological shapes representative of our materials.

Reducing the calculation time

Performing a 3-D calculation until steady state on quite large domains is CPU demanding. To reduce the total computational time, we combined different strategies.

Influence of the relaxation factors

To study the influence of the relaxation factors ω min and ω max , we ran our D 2 Q 5 code on an inclusion geometry. A disc having air properties is placed at the center of a matrix having wood properties, as illustrated in Fig. 3. The radius of the disc is equal to one-quarter of the length of the matrix. The equivalent conductivity was computed for different lattice node numbers m, such as δx = L x /(m-1), and different relaxation factors ω min . The reference solution of this problem was obtained using Comsol, a commercial finite element software. We used quadratic elements. With 'extra fine' (7 × 10 3 elements) and 'extremely fine' (3 × 10 4 elements) meshes, we get the same result with 6 digits after the decimal point. The reference value is λ eq,ref = 0.684071 W m -1 K -1 .

The LBM results are presented in Fig. 3. We can notice that the greater ω min is, the more accurate λ eq is. Moreover, when m increases, the influence of ω min decreases. However, the computational time increases drastically with ω min (and m of course). To limit the calculation time, we fix ω min to 1.4 as above this value, the error remains less than 1% whatever the lattice refinement.

Remark : with our strategy of image segmentation, the contours between solid and gaseous phase are stepped. Using a too small spatial resolution may alter the computed field locally. However, the results obtained with our reference solution tell us that this stepped description has very little effect on the effective property, which is rather a global effect of the morphology.

Steady-state criterion

As Γ eq is calculated at steady state, a convergence criterion must be specified to stop the time iterations. In our case, we computed at various times, on a spruce geometry, the calculation error on Γ eq , Γ , and the convergence parameter, ∆F/F :

Γ (t) = Γ eq,t -Γ eq,∞ Γ eq,∞ , (9)
∆F/F (t) = Φ x x 1 +δx,t -Φ x xm-δx,t Φ x t , (10)
where Γ eq,t is the equivalent property at time t, Γ eq,∞ is the reference property obtained after a large number of iterations (1 × 10 6 for heat, 5 × 10 6 for mass), Φ x t is the x-flux averaged over the whole domain at time t whereas Figure 3: Influence of the relaxation factor ω min and the lattice size m on the equivalent thermal conductivity λ eq and on the number of iterations. The disc with gaseous phase is depicted in red while the solid phase is in blue. Φ x x 1 +δx,t and Φ x xm-δx,t are the x-fluxes averaged along x-constant lines only (left and right sides of the domain, respectively). Indeed, according to the boundary conditions (adiabatic faces around all lateral faces), the mean flux along x is the same at every line parallel to the y-axis in the steady-state regime. Thus, when input Φ x (at the beginning of the disturbance) is close to output Φ x , we can stop the computation. Note that we choose x 1 + δx and x m -δx instead of x 1 and x m to avoid any boundary issues.

In this study, we ran the code twice: i) for thermal diffusion, by applying λ s on walls and λ g on pores; ii) for mass diffusion, by applying f s on walls and f g on pores (f is the dimensionless mass diffusivity, f s = 0.004 and f g = 1). The spruce morphology and the results are depicted in Fig. 4. When ∆F/F decreases, λ decreases also obviously. We notice the influence of averaging: the error on λ eq is lower than the difference between input and output fluxes. We fix the convergence criterion to ∆F/F = 10 -2 to have an accuracy better than 10 -4 on λ eq .

Parallel computing and successive relaxation factors

One of the advantages of the LBM is its suitability for parallel computing. To reduce the calculation time on 3-D morphologies with a 200 × 200 × 200 typical lattice size, we parallelized our code with OpenMP and run the code with successive relaxation factors.

∆F/F

ω min λ eq (10 -2 W m -1 K -1 ) Iteration no. 3,0.5,1,1.4 7.9562 120 × 10 3 1.08% 14

Table 2: Equivalent thermal conductivities, iteration numbers, error on λ eq and iteration gain factor obtained for two convergence criteria and for different ω min . The reference for the calculation of the error and of the gain factor is given by the case ω min = 1.4

Several parts of the code have been parallelized: collision step, streaming step and implementation of boundary conditions. According to tests realized with an LDF morphology of size 100 × 100 × 100, the calculation time is divided by five using eight cores.

Moreover, for each diffusion calculation, we run the code with successive and increasing relaxation factors using the following sequence: ω min = 0.1, 0.3, 0.5, 1, 1.4. When the convergence criterion defined in Sec. 4.2 is reached with one value of ω min , the code is re-run with the next ω min . The final result is obtained for ω min = 1.4 as defined in Sec. 4.1. Tab. 2 presents the different tests realized with the 100 × 100 × 100 LDF morphology in the case of heat diffusion for ∆F/F = 10 -1 and 10 -2 , with the reference ω min = 1.4 or five successive ω min . With ∆F/F = 10 -2 , the calculation time is divided by 14 thanks to this method and the loss in accuracy is only 1 % on λ eq compared with the results obtained with only one relaxation factor ω min = 1.4. Moreover, the choice of ∆F/F = 10 -2 is confirmed because the value ∆F/F = 10 -1 is not accurate enough.

Finally, thanks to the cumulative strategies, we reduced the calculation time by a factor of 70 (5 × 14).

Results and discussion

The optimized LBM code is used to compute the equivalent thermal conductivity and the mass diffusivity of the different bio-based materials presented in Section 2. The simulations were performed in the three main directions of the wood (radial R, tangential T and longitudinal L) and also in three directions of the fibre panels (p 1 and p 2 are the in-plane directions, t is the through-plane direction). Two selections of earlywood, one of latewood, four of LDF and five of MDF were studied. To check the suitable size of our Representative Elementary Volumes (REVs), the evolution of equivalent thermal conductivity and mass diffusivity in a given direction was computed for different REV sizes, using parts of the full volume of increasing length. For example, the results for earlywood in the radial direction are presented in Fig. 5. As a plateau is observed for the two properties when the radial length tends to the full dimension of the selected volume, the selected volume can be considered as representative.

A full LBM computational run on a grid of about 200 × 200 × 200 nodes requires one day using an Intel Xeon x5650 at 2.67 GHz. Note that the calculation time is much higher (a couple of weeks) when the connected phase has a low diffusivity (mass diffusion for wood, thermal diffusion for fibers).

Heat diffusion

For heat diffusion, the thermal conductivities of the solid and the gas are set to λ s = 0.7 W m -1 K -1 and λ g = 0.023 W m -1 K -1 , respectively. The value for the solid phase is a compromise between the two material directions, parallel and across the fibres [START_REF] Siau | Transport processes in wood[END_REF].

The computed equivalent values are compared with the lowest and highest possible values (heterogeneous phases placed in series and in parallel, respectively):

Γ ⊥ = 1 1-s Γg + s Γs , (11)
Γ = (1 -s )Γ g + s Γ s , (12)
where Γ represents λ for thermal conduction or f for mass diffusion and s is the solid fraction. Fig. 6 depicts examples of temperature fields and streamlines obtained in the case of thermal conduction in the radial and tangential directions in latewood. In the case of spruce, we can observe that the heat flux takes advantage of the connected and conductive phase (the solid phase), avoiding the inclusions of low conductivity (the lumens). Fig. 7 presents the equivalent thermal conductivity λ eq versus solid fraction s for the four types of morphologies. The results for all the selections of the bio-based materials studied are presented in Tab. 3. According to the morphology of spruce, it is not surprising to observe that the predicted macroscopic values for earlywood and latewood are close to the parallel model in the longitudinal direction. Moreover, spruce volume lengths are less than the fiber length, which amplifies the proximity with the parallel model.

As the cell walls are rather aligned along the radial direction in earlywood, a consequence of cambial cell divisions in trees, the thermal conductivity is higher in the radial direction than in the tangential direction. In contrast, the tracheids are flattened in latewood, which forms lumens having a large tangential extension. These inclusions, with low conductivity, block the thermal flux mainly in the radial direction, which explains the inverted anisotropic ratio. Those results have already been observed by [START_REF] Perré | New modelling approaches to predict wood properties from its cellular structure: imagebased representation and meshless methods[END_REF]. The value of λ s = 0.7 W m -1 K -1 chosen in this work, an averaged value of the fiber conductivity, can explained why the predicted transverse conductivities are slightly higher than published data [START_REF] Lagüela | Thermal conductivity measurements on wood materials with transient plane source technique[END_REF]. Note however that, due to the large variability of wood in terms of density and growth ring width, the experimental data available for wood also widely dispersed.

These new results obtained with 3D morphologies, allows us to update the correlations proposed in a previous work [START_REF] Perré | Determination of the Material Property Variations Across the Growth Ring of Softwood for Use in a Heterogeneous Drying Model. Part 2. Use of Homogenisation to Predict Bound Liquid Diffusivity and Thermal Conductivity[END_REF]. We kept the same expressions with optimized parameters to fit the LBM results (Fig. 8):

λ eqT = ( g λ g µ T + s λ s µ T ) 1/µ T , µ T = 0.47 λ eqR = α R λ + (1 -α R )λ ⊥ , α R = 0.62 (13) λ eqL = α L λ + (1 -α R )λ ⊥ , α R = 0.97
These expressions are likely to be applied to wood species with similar anatomy (most softwood species), including for the predicting of property profiles along the annual ring. They could also be used to refine the prediction of the directional conductivities of wood, as the fiber conductivity is often reported to be larger in the longitudinal direction than in the transverse direction [START_REF] Siau | Transport processes in wood[END_REF]. The value of the thermal conductivity for the solid phase could easily be dependent on the direction in expressions (13).

In agreement with published data [START_REF] Czajkowski | Thermal properties of wood-based panels: thermal conductivity identification with inverse modeling[END_REF], the fibreboard panels have always a larger conductivity parallel to the plane than perpendicular to the plane (a factor close to 2 for LDF and close to 1.5 for MDF). This is a consequence of the panel manufacturing with tends to align the fibres along the panel plane. In the case of LDF, a good insulating material, the through-plane conductivity λ eq is very close to the series model. This result is consistent with the LDF morphology : due to the fibres alignment, the connections between them is rather poor along the perpendicular direction. Our results for fiberboard panels are in good agreement with measured data [START_REF] Troppová | Influence of temperature and moisture content on the thermal conductivity of wood-based fibreboards[END_REF][START_REF] Lagüela | Thermal conductivity measurements on wood materials with transient plane source technique[END_REF].

Mass diffusion

For mass diffusion, we use the dimensionless mass diffusivity f , which accounts for the diffusion resistance relative to binary diffusion in air [START_REF] Perré | Determination of the Material Property Variations Across the Growth Ring of Softwood for Use in a Heterogeneous Drying Model. Part 2. Use of Homogenisation to Predict Bound Liquid Diffusivity and Thermal Conductivity[END_REF]. Therefore, this value ranges from 0 to 1. The values of f used for the solid and the gaseous phases are f s = 0.004 and f g = 1, respectively. These values are representative of bound water diffusion in wood at moisture content value of ca. 12% [START_REF] Siau | Transport processes in wood[END_REF].

As for heat diffusion, we compare the computed values with the lowest and highest possible values, f ⊥ and f . Fig. 9 depicts an example of streamlines obtained in the case of mass diffusion in LDF. Only one-eighth of the computed volume is shown for better visibility. As the gaseous phase diffuses more than the solid phase, i.e. fibers, the streamlines go around the fibers.

Fig. 10 presents the dimensionless equivalent diffusivity f eq versus solid fraction s for the four types of morphologies. The results for all the selections of the bio-based materials studied are presented in Tab. 5. In the case of mass diffusivity, the situation is opposite for spruce because the connected phase, the solid phase, has a low diffusivity. Therefore, the macroscopic values in the radial and tangential directions are close to the series model. In contrast to what has been observed in 2D [START_REF] Perré | New modelling approaches to predict wood properties from its cellular structure: imagebased representation and meshless methods[END_REF], the anisotropic ratio is not inverted between earlywood and latewood. This is due to the presence of rays in the radial direction that were considered as plain solid in 2-D. Indeed, they present a large porosity, aligned along the radial direction, and act significantly to vapor diffusion in that direction. This contribution is fully accounted for thanks to the high resolution 3-D morphology. In the longitudinal direction, results are still very close to the parallel model.

The fiberboard panels, whatever their density, have a larger equivalent diffusivity f eq than wood. These facts result from the pore morphology of these materials: the macro-porosity, the conductive phase for mass diffusion, is fully connected in low-density fiberboard panels. Regarding the low diffusivity of the solid phase, the contribution of the solid phase can be neglected. Consequently, the equivalent mass diffusivity is simply obtained by the ratio of porosity over the electrical tortuosity τ e :

f eq = g τ e (14)
Note that the concept of electrical tortuosity must be distinguished from the geometrical or hydraulic tortuosity values [START_REF] Ghanbarian | Tortuosity in porous media: A critical review[END_REF]. The values of electrical tortuosity found for LDF are quite large, leading to low dimensionless diffusivities with respect to porosity (Table 4). The tortuosity values, as determined from the LBM results, were compared to the empirical law proposed by [START_REF] Archie | The electrical resistivity log as an aid in determining some reservoir characteristics[END_REF] :

τ e = g (1-m) (15)
The best fit for the parameter m was found to be equal to 3.6 in the direction parallel to the panel plane and 2.1 through the plane. Once again, the effect of the diffusion direction is due to panel manufacturing which tends to align the fibers parallel to the panel plane. These values of m are quite large compared to other porous media. This observation is certainly due to the existence of a closed porosity (within fibers/bundle of fibers) that should not be considered in expression ( 15).

Yet, the diffusivity values computed in the present work remain slightly larger than experimental values reported in the literature [START_REF] Rémond | Mass diffusivity of low-density fibreboard determined under steady-and unsteady-state conditions: Evidence of dualscale mechanisms in the diffusion[END_REF]. We believe that the most plausible explanation for this difference lies in the difficulty of experimental determination of moisture diffusivity, especially when determined in transient state [START_REF] Kumaran | Moisture diffusivity of building materials from water absorption measurements[END_REF]. Indeed, due to the latent heat of evaporation, an intricate coupling between heat and mass transfer takes place as soon as moisture moves. This coupling always slows down the mass transfer and underestimates the mass diffusivity if the experiment is not analyzed using a coupled heat and mass transfer model [START_REF] Perré | The proper use of mass diffusion equations in drying modeling: introducing the drying intensity number[END_REF]. This issue is particularly acute for LDF panels, that cumulate a large moisture diffusivity and a low thermal conductivity. Works currently in progress in our team seem to confirm this assumption : a rigorous analysis of experimental results gives dimensionless diffusivities in good agreement with our LBM results.

LDF

Morphology families of Fiberboard panels

To distinguish the effect of the solid fraction from the effect of the phase morphology, we further analyzed the results obtained for LDF and MDF. A mixing law was chosen to distinguish these effects in a very simple way:

Γ µ = ( g Γ g µ + s Γ s µ ) 1/µ . (16)
For each selection and each physical phenomenon (thermal or mass diffusion), two values are studied: the transverse property Γ t and the mean property in the two tangential directions Γ p = (Γ p 1 + Γ p 2 )/2. The factor µ of Eq. ( 16) is determined by minimizing the difference between the computed values Γ i and the values Γ µ determined by Eq. ( 16).

The results for heat and mass diffusion are presented in Fig. 11, where the equivalent properties of the different selections of LDF and MDF are plotted versus the solid fraction. In the case of mass diffusion, all values can be nicely predicted using one single family whatever the panel density:

• µ t = 0.18 for mass transfer in the transverse direction;

• µ p = 0.43 for mass transfer in the in-plane direction.

In the case of heat diffusion, the situation is not as simple:

• a unique value µ p = 0.53 in the in-plane direction works well whatever the panel density;

• two values are required for a good prediction in the transverse direction µ t = 0.11 for LDF and µ t = 0.22 for MDF.

This is due to a percolation effect: in LDF, the contacts between fibers (the conductive phase for heat transfer) are inefficient because the low density imposes a very high tortuosity on this conductive phase. In contrast, in the case of MDF the higher density imposes numerous contacts between fibres in the transverse direction, which drastically reduces the tortuosity of the conductive phase. Indeed, for heat conduction in the transverse direction, the morphology cannot be assumed to be identical for LDF and MDF.

However, this analysis proves that it is possible to define a synthetic variable characterizing the morphology of a heterogeneous medium. This is particularly useful for simulation codes because, once the parameter µ is known, the material property in a given direction is simply predicted by its density.

As final discussion, one has to keep in mind that the physical properties of hygroscopic materials depend on both temperature and moisture content. The relationships proposed in this work allow the effect of morphology to be distinguished from the property of the solid phase. This is quite powerful as any reported effect of temperature or MC on the solid phase property could be used in our expressions to predict the effective property. For example, Vololonirina et al. [START_REF] Vololonirina | Characterization of hygrothermal properties of wood-based products-impact of moisture content and temperature[END_REF] proposed an expression for the thermal conductivity of wood fibres (the solid phase of the present work) as a function of temperature and moisture content. In the case of mass diffusivity, the huge diffusivity contrast between the gaseous and solid phase renders the task even simpler : we have just to account for the temperature level on the binary diffusion of water vapour in air. For both transfer mechanisms, the effect of shrinkage/swelling on the phase partition has to be taken into account:

• in the case of wood, previous works reported that the lumen size is almost constant in spite of shrinkage/ swelling [START_REF] Mariaux | La section transversale de fibre observée avant et après séchage sur bois massif[END_REF]. This allowed a simple expression to be derived to address this question of phase allocation [START_REF] Perré | Determination of the material property variations across the growth ring of softwood for use in a heterogeneous drying model. part 2. use of homogenisation to predict bound liquid diffusivity and thermal conductivity[END_REF],

• in the case of fiberboards panels, shrinkage/swelling of the solid phase (fibres) is likely to induce dimensional changes in the panel thickness, as the fibre length is very stable. In the case of a building assembly, the constrained panel thickness would turn into an increase of the solid phase fraction and an equivalent decrease of the gaseous phase. Yet, in both cases, the proposed expressions could be used.

Regarding the effect of temperature, our study is certainly limited to low/medium temperature level, before the contribution of evaporation/condensation become important in the global heat transfer [START_REF] Azizi | Approche expérimentale et théorique de la conductivité thermique des milieux poreux humides-i. expérimentation[END_REF]. According to the results published in [START_REF] Troppová | Influence of temperature and moisture content on the thermal conductivity of wood-based fibreboards[END_REF], in the worse case of high moisture content M C ≈ 15%, this effect is already important at 30 o C.

Conclusion

This paper presents a method to predict diffusion properties of different wood-based materials from anatomical images. The digital description of the 3-D morphologies of spruce (latewood and earlywood), LDF and MDF were derived from micro-tomography scans. The Lattice Boltzmann method was used to simulate heat and mass transfer in those materials and thus to compute their macroscopic properties: equivalent thermal conductivity and mass diffusivity. We chose this numerical method because it can account quite easily for complex geometries and is suitable for parallel computing. We developed a strategy to limit computation time and assure good accuracy with the use of successive and increasing relaxation factors. Finally, morphology families were identified and simple expressions are proposed to predict the equivalent properties as a function of phase properties and solid fraction.

Thanks to the impressive possibilities of nano-tomography, the chain of tools proposed in the present work is now used to compute macroscopic properties of several heterogeneous materials : wood-based materials, light concrete, plaster, brick and wood-cement composites. Further studies are also in progress to couple heat and mass transfers, which is an important phenomenon occurring in bio-based materials.

  Figure 1: 3-D reconstruction of selected regions of the bio-based materials scanned by X-ray microtomography. The size, voxel number and solid fraction are given for each 3-D structure.

Figure 2 :

 2 Figure 2: The main steps of the LBM approach in the case of a D 2 Q 5 model (the thickness of the arrows represents the values of the distribution function f i ).

Figure 4 :

 4 Figure4: Calculation of the error on Γ eq and of the convergence parameter for different time steps on a spruce morphology. For heat diffusion, a point is plotted every 2 × 10 4 iterations, for mass diffusion, a point is plotted every 12 × 10 4 iterations.

Figure 5 :

 5 Figure5: Equivalent thermal conductivity λ eq and mass diffusivity f eq of earlywood in the radial direction versus radial length of the computed volume. The length is expressed as a percentage of the radial dimension L R of the selected volume.

Figure 6 :

 6 Figure 6: Temperature field and streamlines for heat diffusion in the radial (left) and tangential (right) directions in spruce latewood. Only half of sample is shown in the tangential direction.

Figure 7 :

 7 Figure7: Equivalent conductivity λ eq versus solid fraction s for three directions of spruce early-and latewood, LDF and MDF. Series (dashed line) and parallel (solid line) theoretical models are represented.

Figure 8 :

 8 Figure8: Prediction of the equivalent conductivity λ eq of spruce as a function of the solid fraction s using equation 13.

Figure 9 :

 9 Figure 9: Streamlines for mass diffusion in the x-direction in LDF. Only one-eighth of the computed volume is shown for better visibility.

  diffusion: µ(LDF t ) = 0.11 ; µ(MDF t ) = 0.22 ; µ p = 0diffusion: µ t = 0.18 ; µ p = 0.43

Figure 11 :

 11 Figure 11: Identification of morphology families for fiberboard according to the direction considered and the type of diffusion.

Table 3 :

 3 Equivalent thermal conductivities λ eq ( W m -1 K -1 ) of the bio-based materials studied in this work.

Table 4 :

 4 Tortuosity calculated for LDF panels from LBM results and with equation (15).Figure10: Equivalent dimensionless mass diffusivity f eq versus solid fraction s for three directions of spruce early-and latewood, LDF and MDF. Series (dashed line) and parallel (solid line) theoretical models are represented. The inset depicts the radial and tangential diffusivity for the latewood on an enlarged scale.

		1
								parallel
		0.9						series
								EW1
		0.8						LW LDF1
		0.7	p 1	p 2	L			MDF1
		0.6
	f eq	0.5		t
		0.4				p 1
		0.3
		0 0.2 0.1			R T	p 2 t	R	L T		0.03 0.02 0.01	R	T
		0		0.2	0.4	0.6		0.8	1	0.7 0	0.75	0.8

Table 5 :

 5 Equivalent dimensionless mass diffusivities f eq of the different bio-based materials studied.

Acknowledgement

The authors wish to thank Erika Jean-Bart (EM2C, CentraleSupélec) for her assistance in laser cutting. Financial supports from Grand Reims, the Marne department and the Grand Est region are gratefully acknowledged.


Luc Abergel

A simple theory to obtain asymptotics of such integrals

Keywords: Asymptotics, Laplace method. MSC class. 40D05, 44A10

Contents 1 Definitions. 2 2 An adjonction between I π and J. 2 3 Comparison relationships on J, then on I π . 2 4 Example of applications. 4

Objectives.

The objective is to provide a method of asymptotic evaluation of integrals such as

I(x, f ) = b a f x (t) dt et I π (x, f ) = b a f x (t)π(t) dt when x tends to +∞.
This is a classic and important subject in mathematics : the evaluation of a fine and often delicate phenomenon. The proof we find are, in general, very technical and difficult to follow. These difficulties complicate the approach considerably in unconventional or even finer contexts.

We will therefore try here to give a simple approach and give a tool allowing access to this problem. This is the central result established in 3, i. e. it is sufficient to know an equivalent of f -1 in the vicinity of the point where f reaches its maximum to obtain an equivalent of these integrals.

We will then give classical or not classical examples of such equivalents, ending with 3 delicate examples, with an infinite number of terms asymptotic for the last example.

In the following, for x quite large, f x π is assumed to be integrable on the integration interval, an interval that may not be a segment.

1 Definitions.

Let's start with some definitions :

Definition :

Under the condition of existence, one poses

I π (x, f ) = b a f x (t)π(t) dt and J(x, f ) = 1 0 f (t x ) dt.
We will abregate I π (x, f ) en I(x, f ) si π = 1.

We will study I π (x, f ) for x → ∞ and J(x, f ) for x → 0.

It should also be noted that the integration interval is not specified for I π (x, f ), when it's supposed to be [0.1] for J(x, f ).

Definition :

It is said that a function f : [a, b] → R + satsifies the (H)-hypothesis if : f defines an homeomorphism on a neighborhood of a such as [a, a + α], f defines a C 1 -diffeomorphism on a neighborhood of a such as [a, a + α], f (a) = 1 and ∀c ∈]a, b], ∃ρ ∈]0, 1[, f ≤ ρ sur [c, b].
If futhermore, f is a decreasing diffeomorphism from ]a, b] to [0, 1[, we'll say that f satisfies the (H )hypothesis.

2 An adjonction between I π and J.

Proposition :

Let Π the primitive of π which vanishes in a. If f satisifies (H ) then I π (x, f ) = J( 1x , Π • f -1 ).

Proof :

We use the variable change u = f (t) :

I π (x, f ) = - 1 0 u x .(f -1 ) (u).π • f -1 (u) du = -[u x Π • f -1 (u)] 1 0 =0 + x 1 0 u x-1 π • f -1 (u) du.
We then put u x = v :

I π (x, f ) = 1 0 Π • f -1 (v 1
x ) dv as wished.

This writing is effective for obtaining comparison relationships for the study of I π (x, f ) in +∞ as shown in the following proposal :

3 Comparison relationships on J, then on I π .

Proposition :

If f = 1 o(g) (resp. f = 1 O(g)) with J(x, g) = 0 e o( 1 x ) , then J(x, f ) = 0 o(J(x, g)) (resp. J(x, f ) = 0 O(J(x, g)).
Proof :

J(x, f ) = a 0 f (v x ) dv + 1 a f (v x ) dv. But f (v) ≤ εg(v) if 1 -α ≤ v ≤ 1.
We then choose a = (1 -α)

1

x and we have :

|J(x, f )| ≤ a 0 f (v x ) dv ≤aN∞(f ) + ε 1 a g(v x ) dv ≤J(x,g) . But (1 -α) 1 x = o(J(x, g)) because J(x, g) = e o( 1 x ) and so aN ∞ (f ) ≤ εJ(x, g) if x is close to 0. So, for such x, |J(x, f )| ≤ 2εJ(x, g). The case f = 1 O(g) is treated in a similar way.
Let's start by establishing the relationship J( 1x , Π

• f -1 ) = 0 e o( 1 x ) if f satisifies (H) :
Proprosition : x) as soon as f satisfies the (H)-hypothesis and that π is never identically zero on any a-neighbourhood.

[I π (x, f )] 1/x -→ x→+∞ 1, or else I π (x, f ) = ∞ e o(
Proof :

-First of all, [I π (x, f )]

1 x ≤ N 1 (π) 1 x → x→∞ 1.
-To minorize, even if it means reducing the interval [a, b], we can assume that f ≥ 1 -ε and π ≥ α for a certain α > 0. We thus obtain I π (x, f )

1 x ≥ [α(b -a)] 1 x (1 -ε) therefore lim inf x→∞ [I π (x, f )] 1 x ≥ 1.
So we've established that :

I π (x, f ) = ∞ exp(o(x)).
We also know, by domination, that :

I π (x, f ) -→ x→+∞ 0.
Proprosition :

If f satisfies (H) then, for all a < c we have

Consequence :

Simplification of the problem.

By changing f out of a a neighbourhood into a new f 1 function, we get

I π (x, f ) = I π (x, f 1 ) + O(ρ x ) with ρ < 1.
We can therefore assume that f verifies (H ) to obtain an equivalent of x) ).

I π (x, f ) (since I π (x, f 1 ) = ∞ e o(
The integration interval can also be reduced. We can therefore treat the case b = ∞ and truncate the integral to give an equivalent.

For example :

∞ 0 f x (t)π(t) dt ∼ ∞ 1 0 f x (t)π(t) dt if f only reaches its sup in 0.
As an application, we obtain the following fundamental result :

Theorem :

If Π • f -1 ∼ 1 Π • g -1
with π not identically zero around a, and if f and g satisfy (H), then

I π (x, f ) ∼ ∞ I π (x, g).
It is an effective result for obtaining research equivalents. It is a consequence of the previous statements.

Remarks :

(1) If π(t) ∼ 0 at α , then :

f -1 ∼ 1 g -1 ⇒ Π • f -1 ∼ 1 Π • g -1 Indeed, lim 1 f -1 = 0 et Π(t) ∼ 0 a α + 1 t α+1 pour tout α = -1, Π(t) ∼ 0 a ln(t) for α = -1,
so an equivalent of f -1 is enough to get an equivalent of

[c -α, c + α] ⊂ [a, b], b a π(t)f x (t) dt ∼ c+α c-α π(t)f x (t) dt because c-α a π(t)f x (t) dt and b c+α π(t)f x (t) dt are dominated by ρ x for a real ρ ∈]0, 1[. If moreover the behavior of Π • f -1 is symmetrical around h = 0, that means that if Π • f -1 (c -h) ∼ h→0 Π • f -1 (c + h) then : c c-α π(t)f x (t) dt ∼ 1 2 b a π(t)f x (t) dt and c+α c π(t)f x (t) dt ∼ 1 2 b a π(t)f x (t) dt.
4 Example of applications.

Wallis integrals and generalizations.

We assume that around 0,

f (t) = 1 -pt d + o(t d ) avec d > 0.
We put g(t) = e -pt d sur[0, ∞[, so :

g -1 (y) ∼ 1 f -1 (y), thus I(x, f ) ∼ +∞ ∞ 0 e -xpt d dt. then : I(x, f ) ∼ +∞ Γ(1 + 1 d ) (px) 1 d
The case d = 2 is the classical example of Wallis integrals In addition, we have the asymptotic

I n = π 2 0 sin n (t) dt.

Equivalent of I(x) =

f (1 + h) = 1 - h 2 2e + o(h 2 ).
Using the symmetry of the behavior of f as well as its asymptotic around h = 0, (cf the simplification of the problem in paragraph 3 and the result of paragraph 4.1), on déduit :

1 0 f x (s) ds ∼ +∞ ∞ 1 f x (s) ds ∼ +∞ eπ 2x , et donc +∞ 0 x t t t dt ∼ +∞ 2π e √ x exp x e
4.3 Formule de Stirling.

It's still a classic problem, at least in the case x ∈ N. It is a question of giving an equivalent of

x! = Γ(x + 1) = ∞ 0 t x e -t dt.
The function f (t) = t x e -t reaches its sup in t = x.

We therefore put t = x(1 + h) in the integral, which gives

x! = x ∞ -1 exp [xln(x) + xln(1 + h) -x -xh] dh, so x! = x x e -x x ∞ -1 f x (h) dh where f (h) = e ln(1+h)-h = h→0 1 -h 2 2 + o(h 2 ).
Here again, by symmetry of the f behavior, as well as its asymptotic around h = 0 :

∞ -1 f x (h) dh ∼ ∞ 2 Γ( 3 2 )
x 2

, and so :

x! ∼ ∞ x x e -x √ 2πx
4.4 Other examples.

We are doing to study an asymptotic of I(x, f ) =

1 0 f x (t) dt with f (t) = 1 + t ln(t) then f (t) = 1 + t ln(t) .
1 Case f (t) = 1 + t ln(t).

Here

f -1 (y) ∼ 1 1-y -ln(1-y) . I(x, f ) = J( 1 x , f -1 ) ∼ ∞ 1 0 1 -t 1 x -ln(1 -t 1 x ) dt.
We then pose u = 1-t 1

x thereafter u = v x and we find

I(x, f ) ∼ ∞ xK(x-1) where K(x) = 1 x 2 x 0 v(1 -v x ) x ln(x) -ln(v)
dv.

We then rewrite K(x)

∼ ∞ 1 x 2 ln(x) ∞ 0 ve -v dv.
Hence :

1 0 (1 + t ln(t)) x dt ∼ ∞ 1 x ln(x) 2 Case f (t) = 1 + t ln(t) .
Here

f -1 (y) ∼ 1 -(1 -y) ln(1 -y).
Then in the same way

I(x, f ) ∼ ∞ xK(x -1) where K(x) = 1 x 2 x 0 (ln(x) -ln(v))v(1 - v x ) x dv, K(x) ∼ ∞ ln(x) x 2 ∞ 0 ve -v dv,
and finally :

1 0 1 + t ln(t) x dt ∼ ∞ ln(x)
x

3 Asymptotic to any order of

1 0 t x (1 -t) α-1 ln β 1 1 -t dt = I α,β (x) pour x → ∞.
By changing variables,

I α,β (x) = 1 0 (1 -t) x t α-1 ln β (1/t) dt ∼ 1 0 e -xt t α-1 ln β (1/t) dt = 1 x α x 0 e -u u α-1 ln β (x/u) dt,
by the theorem of paragraph 3.

An equivalent could then easily be deduced, namely Γ(α) ln β (x)

x α , but we'll give a more precise asymptotic here by analyzing the error.

Technical preliminaries :

1 ∞ x e -u u α-1 ln k (u) du ≤ e -(1-ε)x ∞ 0 e -εu u α+k-1 du = O(e -x(1-ε) ) at k fixed and x ≥ 1 (by ln(u) ≤ u over [1, ∞[)
2 By Taylor's formula with integral remainder, let's show that

|(1 + h) β - p-1 k=0 β k h k | ≤ C p h 0 h -u 1 + u p (1 + u) β (h -u) du ≤ a|h| p+c + a |h| p+c
where a, a , c, c design constants : xt 2 e -xt t α-1 ln β (1/t) dt for the integration interval ]0, t 0 ], and in O(ρ x ) with ρ < 1 for the interval [t 0 , x[. We majorize on the integration interval ln β (1/t) by M t ε , whitch gives

Let's write A = (h -u) p+1 (1 + u) β-p = h-u 1+u p-a (1 + u) β+a (h -u) 1+a , and choose a ≥ 0 such as β + a ≥ 0. By h-u 1+u ≤ |h| if u ∈ [0, h], if p ≥ a, we have |A| ≤ |h| p-a+1+a (1 + |h|) β+a , but (1 + |h|) β+a ≤ C(1 + |h| β+a ) (the quotient being bounded on R), hence the desired result. 3 The error between I α,β (x) = 1 0 (1 -t) x t α-1 ln β (1/t) dt and 1 0 e -xt t α-1 ln β (1/t) dt is in O 1 x α+1-ε : |e -xt -(1 -t) x | = e -xt |1 -e x(ln(1-t)+t) | ≤ Cxt
E ≤ C x t0 0 e -xt t α+1-ε dt + O(ρ x ) ≤ C x ∞ 0 e -xt t α+1-ε dt + O(ρ x ) = O 1 x α+1-ε
, which will be enough here.

Application :

Let

I 1 = x 0 e -u u α-1 ln β (x/u) du = ln β (x) x 0 e -u u α-1 1 + ln(1/u) ln(x) β du.
By point 2 :

Let ε = I 1 -ln β (x) x 0 e -u u α-1 p-1 k=0 β k ln k (1/u) ln k (x) du .
We thus have ε ≤ C ln β (x) Finally by point 3 :

I α,β -I 1 = O 1 x α+1-ε .
So that

1 0 t x (1 -t) α-1 ln β 1 1 -t dt = 1 x α p-1 k=0 β k ln β-k (x)(-1) k Γ (k) (α) + O ln β-p (x) x α
Remarks :

-Γ (k) (α) ∼ (-1) k k! α k+1 and therefore the series

∞ k=0 α k ln β-k (x)(-1) k Γ (k) (α) is divergent
(otherwise absolute convergence would have immediately given the asymptotic by truncating the series).

-By the same method, we can establish or even, give a more precise asymptotic.

  t)π(t) dt = O(ρ x ) for a ρ < 1. The proof is immediate by writing that f x0 π is integrable on [c, b[ for a certain x 0 , and on the other hand that f x π ≤ ρ x-x0 f x0 π over [c, b[.

  t)π(t) dt at +∞. (2) The case f reaches its sup on ]a, b[ : If : f reaches its maximum in c ∈]a, b[ f defines an homeomorphism on neighbourhoods of c such as [c -α, c] and [c, c + α] f defines a C 1 -diffeomorphism on neighbourhoods of c such as ]c -α, c[ and ]c, c + α[ ∀d < c( resp. d > c), ∃ρ ∈]0, 1[, f ≤ ρf (c) over [a, d] (resp. [d, b]) π is not identically 0 on any neighborhood of c then for all α > 0 such as

  We write I(x) = I(x) = +∞ 0 e tln(x)-tln(t) dt. By puting t = xs e et en notant f (s) = e s-1-sln(s) e , we obtain I(x) = x e e x/e ∞ 0 f x (s) ds.

  2 e -xt si t ≤ t 0 by a Taylor formula. Note the error between the two integrals. So she's bounded by C t0 0

x 0 e

 0 -u u α-1 a | ln p+c (u)| ln p+c (x) + a | ln p+c (u)| ln p+c (x) du, then ε ≤ C a ln -p-c+β (x) ∞ 0 e -u | ln(u)| p+c du + a ln -p-c +β (x) ∞ 0 e -u | ln(u)| p+c du = O ln -p-c"+β (x) .Futhermore, by point 1 : ρ x ) for one ρ < 1.

  t a-1 ln b (1/t) ln c (ln(1/t)) dt ∼ ∞ Γ(a) ln b (x) ln c (ln(x))x a

Professor, Lycée Janson de Sailly, PARIS, email : lucabergel@cegetel.net