M Ben Youssef

F Lavergne
email: francis.lavergne@cerema.fr

K Sab

K Miled

J Neji

Upscaling the elastic stiffness of foam concrete as a three-phase composite material

Keywords: foam concrete, homogenization, Young modulus, microstructure, numerical simulations

The stiffness of foam concrete depends primarily on the added porosity. Nevertheless, by performing 3D elastic numerical simulations on artificial unit cells in the frame of periodic homogenization, it is shown that describing foam concrete as a porous material is not sufficient to explain the experimental measurements of the Young modulus for added porosity higher than 40%. Indeed, introducing sand as a third phase enables to recover accurate estimates of the Young Modulus. Furthermore, for highly porous concrete foams, it is shown that the stress concentrates in thin members deprived of stiff sand particles, thus leading to a softer overall stiffness.

Introduction

Foam concrete, also named cellular concrete, offers many advantages as a building material. Indeed, it features a low density, a low thermal conductivity and a valuable fire and water resistance [START_REF] Cong | Properties of a foamed concrete with soil as filler[END_REF][START_REF] Sayadi | Effects of expanded polystyrene (eps) particles on fire resistance, thermal conductivity and compressive strength of foamed concrete[END_REF]. Lastly, it can easily be manipulated, pumped and cast into various shapes [START_REF] Gambhir | Lightweight foamed or aerated concrete[END_REF]. Regarding the formulation of foam concrete, the volume fraction of void, hereby called the porosity, can be tuned to trade cost, density and insulation properties for durability, strength and stiffness [START_REF] Miled | Particle size effect on {EPS} lightweight concrete compressive strength: Experimental investigation and modelling[END_REF][START_REF] Tang | Mechanical and drying shrinkage properties of structural-graded polystyrene aggregate concrete[END_REF][START_REF] Mydin | Mechanical properties of foamed concrete exposed to high temperatures[END_REF][START_REF] Gambhir | Lightweight foamed or aerated concrete[END_REF][START_REF] Cong | Properties of a foamed concrete with soil as filler[END_REF][START_REF] Hajek | Elasticity modulus and flexural strength assessment of foam concrete layer of poroflow[END_REF]. Typical uses include floors, trench fills, roof insulation and masonry units [START_REF] Neville | Properties of Concrete[END_REF].

For lightweight concrete containing soft EPS beads formulated in [START_REF] Miled | Particle size effect on {EPS} lightweight concrete compressive strength: Experimental investigation and modelling[END_REF], extending [START_REF] Roy | Taking into account the inclusions' size in lightweight concrete compressive strength prediction[END_REF], it has been shown that the Young modulus does not depend on the size of the beads in the range 1mm-6.3mm, for porosities in the range 0%-50%. As a result, the Young modulus mainly depends on the porosity, exhibiting a one-to-one relationship in this domain. In addition, different mean-field homogenization schemes have been defined with the view to explain the dependence of the Young Modulus to the added porosity, thus aiming at a scientific modelling of the stiffness. For instance, the fact that the experimental measurements of the Young Modulus are compliant with the Hashin-Shtrikman upper bound [START_REF] Hashin | A variational approach to the theory of the elastic behaviour of multiphase materials[END_REF] has been checked and the trisphere model [START_REF] Roy | Déformation instantanées et différées des bétons à hautes performances[END_REF][START_REF] De Larrard | Structures granulaires et formulation des b étons -Concrete Mixture -Proportionning -A scientific approach[END_REF] has been successfully applied to provide a better estimate of the stiffness [START_REF] Roy | Taking into account the inclusions' size in lightweight concrete compressive strength prediction[END_REF]. Furthermore, the differential scheme [START_REF] Mclaughlin | A study of the differential scheme for composite materials[END_REF] ensures a good agreement on a wider range of porosity (0≤ p ≤56%) as the Young Modulus varies directly with (1 -p) 2 [START_REF] Miled | Effective elastic properties of porous materials: Homogenization schemes vs experimental data[END_REF]. This exponent of 2 is consistent with expected and measured trends for foams where the bending is mainly attributed to the bending of the cell ribs [START_REF] Gibson | The mechanics of foams: basic results[END_REF][START_REF] Choi | Analysis of elastic modulus of conventional foams and of re-entrant foam materials with a negative poisson's ratio[END_REF]. Nevertheless, for foam concretes featuring high porosities (up to 75%), power laws featuring exponents higher than 2 (E = 24GPa(1 -p) 2.5 [START_REF] Youssef | Mechanical properties of non-autoclaved foam concrete: analytical models vs. experimental data[END_REF] or E = 32.9GPa(1 -p) 2.8 [START_REF] Mydin | Mechanical properties of foamed concrete exposed to high temperatures[END_REF]) have been fitted to experimental results, thus questioning the use of the differential scheme for porous materials.

As an alternative to Eshelby-based mean-field schemes, full-field numerical elastic simulations can be performed on an artificial unit cell of periodic composite materials to estimate the effective overall stiffness [START_REF] Povirk | Incorporation of microstructural information into models of two-phase materials[END_REF][START_REF] Gusev | Representative volume element size for elastic composites: A numerical study[END_REF][START_REF] Roberts | Elastic moduli of model random three-dimensional closed-cell cellular solids[END_REF][START_REF] Kanit | Determination of the size of the representative volume element for random composites: statistical and numerical approach[END_REF]. Two steps are involved in the process. First, an artificial periodic unit cell complying with the observed microstructure of the composite material must be built. Parameters such as volume fractions, pore size distributions and minimum spacing may be accounted for in this step. Then, elastic computations help estimating the effective stiffness of the composite material by relying on the periodic homogenization theory.

Postprint of Cement & Concrete Reserch doi.org/10.1016/j.cemconres.2018.04.021

The description of foam concrete as a biphasic material will be questioned by performing numerical simulations. First, concrete foam is considered as a biphasic material and the Young Modulus will be estimated for different volume fractions, pore size distributions and parameters driving the geometry of the microstructure. These estimates will be compared to experimental measurements to show that foam concrete can hardly be described as a biphasic composite material at high porosity. Then, sand is introduced as the third phase and the estimated Young Modulus is compared to the experimental measurements.

3D numerical modelling

Generation of artificial unit cells

So as to perform 3D numerical simulations, unit cells complying with the morphology of the considered material must be generated. By mean of X-ray tomography, foam concrete can visibly be described as porous matrix-inclusion material, where the pores are rather spherical (Fig. 1.1). In addition, the sieve curve of the sand and the pore size distribution are such that the sand particles are smaller than the pores (Fig. 1.2). Nevertheless, there is no separation of scale between the sand and the pores. In the present section, a procedure to obtain unit cells satisfying this description is proposed. Figure 1.1: A 6×6mm 2 slice extracted from a X-ray tomographic image of a foam concrete of density 0.6 and porosity 70% [START_REF] Youssef | Mechanical properties of non-autoclaved foam concrete: analytical models vs. experimental data[END_REF], the pixel size being 6µm. Dark areas are pores, grey uniform areas are sand particles, white dots are likely related to unhydrated cement particles.

On the one hand, for matrix-inclusion composite materials, the random sequential adsorption algorithm [START_REF] Feder | Random sequential adsorption[END_REF] or the Lubachevsky-Stillinger algorithm [START_REF] Lubachevsky | Geometric properties of random disk packings[END_REF] are often applied to pack spherical inclusions and reach the targeted volume fractions. On the other hand, the microstructure of closed-cell foams is often modelled by using Voronoi or Laguerre diagrams. In the present work, a continuum of microstructures between these geometries is defined to model closed-cell foams of any volume fractions. First, the Lubachevsky-Stillinger algorithm [START_REF] Lubachevsky | Geometric properties of random disk packings[END_REF] is adopted to produce periodic unit cells featuring thousands of spherical inclusions and volume fractions up to about 60% [START_REF] Lavergne | Effects of grain size distribution and stress heterogeneity on yield stress of polycrystals: A numerical approach[END_REF]. This algorithm is briefly recalled in the next section.

The Lubachevsky-Stillinger algorithm for periodic unit cells

The Lubachevsky-Stillinger algorithm [START_REF] Lubachevsky | Geometric properties of random disk packings[END_REF][START_REF] Kansal | Computer generation of dense polydisperse sphere packings[END_REF] is an event-driven algorithm designed to pack hard spherical particles in a given volume. As input, it requires the size of the unit cell and both the targeted volume fractions φ i and size distributions for each kind of particles i ∈(sand,pore).

The size distribution of sand particles or pores, named the sieve curve, is defined by the volume fraction of particles P i (v) of volume lower than the volume v (Fig. 1.2). For the shake of simplicity, it is modelled as a piecewise linear function on intervals [v j , v j+1 ], j ∈ 0..N s i , the volume fraction of particles in these range being p i,j . N s i is the number of sieves. A list of targeted volumes v i,k , k ∈ 0..N u i is randomly picked according to the sieve curve so as to fill a volume k v i,k = V φ i , where V is the volume of the rectangular unit cell. N u i is the number of spheres of kind i to be placed in the unit cell. The Lubachevsky-Stillinger algorithm is seeded by randomly placing dots in the unit cell at time t = 0. Each dot becomes a sphere as time goes by, the growth of its radius being linear. Consequently, the growth rate is set to a i,k = (3/(4π)v i,k ) 1/3 so that all spheres reach their targeted volume v i,k at time t = 1. The Lubachesky-Stillinger [START_REF] Lubachevsky | Geometric properties of random disk packings[END_REF][START_REF] Kansal | Computer generation of dense polydisperse sphere packings[END_REF] enforces that the growing spheres do not overlap. Indeed, the spheres are allowed to move at uniform velocity and collisions are handled by changing these velocities at the time of collisions. The time of the potential collision can be analytically computed since the velocities are uniform and the growth rate is linear. Consequently, the algorithm can jump from one collision to the next without time-stepping. Other potential events are encounters with the boundaries of the unit cell, where instances of the sphere must be introduced to ensure periodicity [START_REF] Ghossein | A fully automated numerical tool for a comprehensive validation of homogenization models and its application to spherical particles reinforced composites[END_REF][START_REF] Lavergne | Effects of grain size distribution and stress heterogeneity on yield stress of polycrystals: A numerical approach[END_REF]. Indeed, either the sphere hits a boundary and a new connected instance of the sphere must be created at the opposite boundary or an instance leaves the unit cell and it can be destroyed to save memory. An additional improvement consists in splitting the domain into rectangular sectors to reduce the number of potential events. Consequently, if N instances of different spheres overlap with a sector, N (N -1)/2 + 6N potential events can occur in that sector.

The output of the Lubachevsky-Stillinger algorithm specifies the center, the volume and the kind of each sphere in the unit cell. The algorithm is very efficient at packing spherical particles as it can reach volume fractions of 74% for an unimodal pore size distribution. Nevertheless, the degree of order increases monotonically with the jammed packing fraction [START_REF] Torquato | Is random close packing of spheres well defined?[END_REF] and reaching high volume fractions (≥ 64.5%) requires more time as the growth rates need to be reduced [START_REF] Kansal | Computer generation of dense polydisperse sphere packings[END_REF]. To avoid these shortcomings, an additional step allowing to reach any volume fractions of pores is introduced in the next section, though the pores become non-spherical.

A continuum of microstructures between spherical inclusions and Laguerre diagrams

Foams are often modelled by considering Laguerre diagrams, also called Power diagrams. These partitions of space into convex polyhedral cells are built starting from n seeds characterized by their position x k and a scalar, named power, r k . The cell i is defined as the set of points x such that:

||x -

x i || 2 -r 2 i ≤ ||x -x k || 2 -r 2 k ∀k (1)
If all powers r k are equal, the resulting tessellation is a Voronoi diagram. Seeding a Laguerre diagram with the centers and radii of non-overlapping spheres is a practical way to gain some control on the volume and shape of each cell of the foam (or polycristal) [START_REF] Fan | Simulation of polycrystalline structure with voronoi diagram in laguerre geometry based on random closed packing of spheres[END_REF][START_REF] Benabbou | Geometrical modeling of granular structures in two and three dimensions. application to nanostructures[END_REF][START_REF] Jafari | Microstructure generation of severely deformed materials using voronoi diagram in laguerre geometry: Full algorithm[END_REF][START_REF] Hitti | Precise generation of complex statistical representative volume elements (rves) in a finite element context[END_REF][START_REF] Lavergne | Effects of grain size distribution and stress heterogeneity on yield stress of polycrystals: A numerical approach[END_REF]. Indeed, by doing so, the Laguerre cells obviously contain their corresponding spheres. 3

Finally, to complete the modelling of the closed cell foam, a thickness b ik is associated to each solid planar member of the Laguerre diagram [START_REF] Benabbou | Geometrical modeling of granular structures in two and three dimensions. application to nanostructures[END_REF], providing that the total solid volume is consistent with the targeted porosity φ. Nevertheless, the generated microstructures features convex polyhedral pores and it clearly conflicts with the round pores of foam concrete due to the equilibrium between pore pressure and surface tensions. This issue is now to be tackled.

Starting from a sphere packing at an intermediate volume fraction φ pore ≤ φ, each spherical pore can be dilated until the minimum thickness b ik becomes the limiting factor. A dilation parameter d is introduced and a point x is in a pore i if and only if three conditions are met:

1. The point is in the corresponding Laguerre cell i 2. The point is not in the closed-cell foam defined by the minimum thickness b ik 3. The point is in the dilated sphere:

||x -x i || 2 -(r i + d) 2 ≤ 0 (2)
The value of the dilation parameter d must be set according to the targeted porosity by mean of a dichotomy algorithm. Lastly, a unique minimum thickness b = b ik ∀i, k is set for all planar members so as to limit the number of parameters. Consequently, the complete procedure is driven by a couple of parameters:

The size of the unit cell l

The targeted volume fractions of sand φ sand

The grain size distribution of sand v sand,j , p sand,j

The targeted porosity φ

The pore size distribution v pore,j , p pore,j

The intermediate volume fraction of pores φ pore

The minimum thickness of planar members b

The only requirement on l is that it must be big enough for the estimated effective stiffness to be representative of the composite material. Volume fractions and size distributions can be extracted from the formulation or tomography images. Then, parameters φ pore and b are yet to be set. If the difference φ -φ pore is significant, the microstructure deviates from well-dispersed spherical inclusions. Indeed, large planar members are created between pores which are close to one another. Then, if the minimum thickness b were null, the pores would percolate and the solid skeleton could be fragmented. Moreover, the minimum thickness b divides the matrix into two parts: the matrix is either mainly in the planar member (large b) or concentrated at the foam nodes.

It can be noticed that the modelling of highly porous foams by such a procedure is not perfect, due to the fact that the Plateau's laws are almost always infringed: angles between the planes defining edges are likely not equal to 120°and angles between edges can be very different from arccos(-1/3) ≈ 109.47. In addition, a given non-spherical pore can feature both a flat face and curved boundaries. Since the curvature is related to the pressure drop across the interface, and since the pressure within a pore is uniform, the solid skeleton is assumed to be able to withstand some deviatoric stress during setting, through viscous effects for instance.

The resulting unit cell defines the microstructure of a periodic foam concrete (Fig. 1.3). The next section addresses the numerical computation of the effective elastic properties of this material.

3D numerical simulations by the Fast Fourier Transform (FFT) algorithm

Periodic boundary conditions are applied, a mean stress is applied to the microstructure and the local strain is computed. Finally, a component of the overall stiffness is obtained by integrating the local elastic energy. By generating a unit-cell respecting periodic constrains, boundary effects are avoided and the estimated overall stiffness corresponds to that of the material defined by periodizing the unit cell. Nevertheless, the unit cell must be large enough to ensure an estimated Young modulus representative of the composite material. Different algorithms can be applied to solve the elastic problem, including the finite element method and the FFT method. The latter is preferred in the present study due to its versatility and ease of use.

The FFT algorithm

The FFT algorithm is recognized as one of the most powerful tool for computing the overall response of composites featuring complex periodic microstructures [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF][START_REF] Eyre | A fast numerical scheme for computing the response of composites using grid refinement[END_REF][START_REF] Moulinec | Comparison of three accelerated fft-based schemes for computing the mechanical response of composite materials[END_REF]. As input, the local elastic stiffness C(x) must be provided at each point of a regular grid representing the unit cell. Then, either a macroscopic strain E or a macroscopic stress Σ is applied to the unit cell. As output, the algorithm provides a value of the local strain ε(x) and stress σ(x) for each voxel x of the grid as well as averages of local stresses σ and strains ε. Finally, the elastic energy can be computed according to Hill's lemma. Since the accelerated scheme of Eyre and Milton [START_REF] Eyre | A fast numerical scheme for computing the response of composites using grid refinement[END_REF] has been pointed as the most efficient regarding the convergence rate at high contrast [START_REF] Moulinec | Comparison of three accelerated fft-based schemes for computing the mechanical response of composite materials[END_REF], it is applied without modification in the present study. If a macroscopic strain E is applied to the unit cell, its iteration writes:

ε i+1 = ε i + 2(C + C 0 ) -1 : C 0 : Γ 0 * (C 0 -C) : ε i -ε i + E (3)
where * denote a convolution and Γ 0 is the Green operator describing the strain triggered by an eigenstress in an infinite homogeneous material of an isotropic reference stiffness C 0 . This convolution is turned into a product in the Fourier space. Indeed, the Green operator Γ 0 is explicitly known in the Fourier space for any frequency ξ [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF]:

Γ0 ijkh (ξ) = 1 4µ 0 ||ξ|| 2 (δ ki ξ h ξ j + δ hi ξ k ξ j + δ kj ξ h ξ i + δ hj ξ k ξ i ) - λ 0 + µ 0 µ 0 (λ 0 + 2µ 0 ) ξ i ξ j ξ k ξ h ||ξ|| 4 (4)
where λ 0 and µ 0 are the Lamé coefficients of the isotropic reference material C 0 . The fastest convergence rate is obtained for a reference material defined as [START_REF] Eyre | A fast numerical scheme for computing the response of composites using grid refinement[END_REF][START_REF] Moulinec | Comparison of three accelerated fft-based schemes for computing the mechanical response of composite materials[END_REF]:

k 0 = max x k(x) min x k(x) µ 0 = max x µ(x) min x µ(x) (5)
where k(x) and µ(x) are respectively the bulk and shear modulus of the isotropic stiffness C(x). The complete algorithm [START_REF] Moulinec | Comparison of three accelerated fft-based schemes for computing the mechanical response of composite materials[END_REF] writes :

ε 0 (x) = E ∀x while not converged τ (x) = (C 0 -C(x)) : ε i (x) ∀x τ = FFT(τ ) ê(ξ) = Γ0 (ξ) : τ (ξ) ∀ξ = 0 ê(0) = E e = FFT -1 (ê) e(x) = e(x) -ε i (x) ∀x ε(x) i+1 = ε i (x) + 2(C(x) + C 0 ) -1 : C 0 : e(x) ∀x σ(x) i+1 = C(x) : ε(x) i+1 ∀x convergence tests (6)
If a macroscopic stress Σ is applied, the overall strain is to be updated at each iteration as [START_REF] Moulinec | A numerical method for computing the overall response of nonlinear composites with complex microstructure[END_REF]:

E = C -1 0 : (Σ -σi ) + εi (7)
The previous algorithm requires little memory since only two fields need to be stored if the Fast Fourier Transform (FFT) is computed in place. The complexity of the iteration is that of FFT, that is n log(n) where n is the size of the grid. The number of iteration needed to ensure convergence varies directly with the square root of the maximum contrast (max k/ min k, max µ/ min µ). Convergence criteria are related to the balance of the stress, compatibility of the strain and respect of the applied overall loading as defined in [START_REF] Moulinec | Comparison of three accelerated fft-based schemes for computing the mechanical response of composite materials[END_REF]. Among improvements of the FFT algorithms that are not considered in the present article, the use of a finite difference based Green operator [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] or a Green operator consistent with voxel-wise constant polarization fields [START_REF] Brisard | Fft-based methods for the mechanics of composites: A general variational framework[END_REF][START_REF] Brisard | Combining galerkin approximation techniques with the principle of hashin and shtrikman to derive a new fft-based numerical method for the homogenization of composites[END_REF] may lead to more accurate local fields [START_REF] Willot | Fourier-based schemes for computing the mechanical response of composites with accurate local fields[END_REF] or provide bounds on the macroscopic properties based on variational principles [START_REF] Brisard | Fft-based methods for the mechanics of composites: A general variational framework[END_REF][START_REF] Brisard | Combining galerkin approximation techniques with the principle of hashin and shtrikman to derive a new fft-based numerical method for the homogenization of composites[END_REF].

Handling of pores and grey voxels

A porous material is considered: the Young modulus of the solid skeleton is set to E m = 24GPa, its Poisson's ratio being 0.2. Therefore, the bulk and shear modulus of the matrix respectively are K m = 13.33GPa and G m = 10.00GPa. The accelerated scheme of Eyre and Milton [START_REF] Eyre | A fast numerical scheme for computing the response of composites using grid refinement[END_REF] is not able to handle pores as the contrast becomes infinite. Consequently, a non-null stiffness must be associated to the pores so as to recover a finite contrast. The larger the contrast is, the softer the pore are and the more accurate the estimated overall response becomes. Nevertheless, the computation time increases as the square root of the contrast. As a consequence, a sensitivity study to the contrast is to be performed: contrasts of 100, 1000 and 10000 are tested for a given porous unit cell featuring a porosity of 63.61% (Fig. 1.4). While changing the contrast from 100 to 1000 leads to a significant decrease of the estimated moduli, increasing the contrast from 1000 to 10000 induces a smaller decrease of the estimated moduli. In addition, the sensitivity of the fineness of the discretization is to be checked. Lastly, some of the voxels overlap on both the pores and the solid skeleton: a rule must define the stiffness of such voxels. To this end, an estimate of the porosity in the voxel is produced by using 64 sensing points in each voxel. The stiffness is then set according to an averaging rule, by applying either the Voigt bound (arithmetic average of stiffness) or the Reuss bound (harmonic average of stiffness). As the Reuss bound is applied, it is observed on figure 1.4 that the estimated moduli vary directly with the size of the pixel, that is the volume fraction of grey voxels. On the contrary, it is also shown on figure 1.4 that choosing the Voigt bound leads to an estimate of the overall stiffness that is less sensitive to the discretization. It must be noticed that the Reuss bound could have been chosen for the very same reason if the inclusions were stiffer than the matrix [START_REF] Lavergne | Investigation of the effect of aggregates' morphology on concrete creep properties by numerical simulations[END_REF]. Other mixing rules [START_REF] Kabel | Use of composite voxels in fft-based homogenization[END_REF][START_REF] Gélébart | Filtering material properties to improve fft-based methods for numerical homogenization[END_REF] have also been proposed: a macroscopically consistent rules defining the stiffness of grey voxels would likely be non-local, as introduced in [START_REF] Bignonnet | Macroscopically consistent non-local modelling of heterogeneous media[END_REF]. Scaling the Young modulus of the solid skeleton would not affect the outcome of the sensitivity study due to the fact that the estimated overall bulk and shear moduli vary directly with it. In the sequel, all computations are performed using the Voigt bound for grey voxels and a contrast of 10000 to account for the pores.

Representativity of the unit cell

It must be checked that the generated unit cell is a representative elementary volume (REV). Generating the unit-cell of a periodic matrix-inclusion microstructure and applying periodic boundary conditions ensures that the estimated overall response is representative of the effective stiffness of that particular microstructure. Indeed, the overall response is not plagued by boundary effects which occur as displacement or force boundary conditions are considered, as investigated in [START_REF] Sab | On the homogenization and the simulation of random materials[END_REF]. Nevertheless, if the unit cell is too small, the estimate of the stiffness can be erroneous since the unit cell does not account for the size distribution and other statistics of the material. For a maximum error of 5% on the estimated overall modulus of matrix-inclusion materials, the unit cell must be at least twice as large as the reinforcement diameter [START_REF] Drugan | A micromechanics-based nonlocal constitutive equation and estimates of representative volume element size for elastic composites[END_REF].

Furthermore, the statistical variability of the estimated stiffness can be plotted as a function of the size of the unit cell to justify the choice of this size [START_REF] Povirk | Incorporation of microstructural information into models of two-phase materials[END_REF][START_REF] Gusev | Representative volume element size for elastic composites: A numerical study[END_REF][START_REF] Kanit | Determination of the size of the representative volume element for random composites: statistical and numerical approach[END_REF][START_REF] Zeman | Numerical evaluation of effective elastic properties of graphite fiber tow impregnated by polymer matrix[END_REF]. Unit cells of periodic microstructures featuring a porosity of 61% and a pore size distribution given in figure 1.2 (maximum diameter D = 2.21mm) are generated by using the Lubachevsky-Stillinger algorithm. Unit cells of sizes 4 × 4 × 4 mm 3 , 6 × 6 × 6 mm 3 and 8 × 8 × 8 mm 3 are discretized on grids of respective sizes 80 3 , 120 3 and 160 3 so as to discard the effect due to the pixel size. Ten unit cells of each size are created so as to compute the average and the standard deviation of the estimated bulk and shear moduli. It is shown on figure 1.5 that the standard deviation is reduced as the size of the unit cell increases. Moreover, all the estimated bulk moduli and shear moduli respectively fall into ranges (2.63GPa, 2.69GPa) and (1.72GPa, 1.85GPa). Consequently, it can reasonably be assumed that the relative error on the estimated bulk and shear moduli due to representativity is less than 5% for a given unit cell of size 6 × 6 × 6 mm 3 or larger.

Comparison between numerical estimates and experimental measurements

Foam concrete as a porous material

Various unit cells featuring different pore size distributions, different intermediate porosities and minimum thicknesses are considered.

Influence of the pore size distribution

For each porosity, three pore size distributions are considered:

An unimodal pore size distributions, where all pores feature an equivalent diameter of 1.5mm.

A bimodal pore size distribution, where half the porosity corresponds to pores of equivalent diameter 2mm and the other half to pores of diameter 0.5mm.

A multimodal pore size distribution, given in figure 1.2 (maximum diameter D = 2.21mm). For generating the unit cells, the intermediate volume fraction is set to 61% if the porosity is larger than 61% and the minimum thickness is set to 0. The size of the unit cells is 6 × 6 × 6 mm 3 and these unit cells are discretized on 120 3 grids for the FFT algorithm. For different porosities, the overall moduli estimated by the 3D numerical computations are compared to those stemming from mean field schemes in figure 2.1. It is shown that these numerical estimates comply with the Hashin-Shtrikman upper bound and could be approximated by using the differential scheme. In addition, the different pore size distributions cannot be distinguished starting from the estimated moduli, keeping in mind their limited accuracy. As a result, the effect of the pore size distribution can be described as small and non-significant. Lastly, the experimental measurements of the Young moduli of foam concrete cannot be properly approximated by any of the computed estimates. For instance, the measured Young modulus at a porosity of 67% is about 5% of that of the solid matrix, while the ratio estimated by 3D numerical simulations is about 12.7%. That ratio becomes 11.0% for the differential scheme and 19.7% for the Hashin-Shtrikman upper bound. In the next section, the intermediate porosity and the minimum width will be tuned to try to match the measured moduli.

Influence of the intermediate porosity and minimum thickness

The pore size distribution being set to that given in figure 1.2 (maximum diameter D = 2.21mm), the effect of the intermediate porosity and minimum thickness are to be investigated. As the intermediate volume fraction becomes significantly lower than the targeted porosity, the pores tend to bundle, thus forming larger pores of complex shapes. Consequently, this clustering induces a significant decrease of the estimated stiffness (Fig. 2.2), thus getting closer to the experimental measurements of the Young modulus of foam concrete. Nevertheless, explaining the experimental measurements by introducing a single intermediate volume fraction is hardly possible. Indeed, even if the intermediate volume fraction is null, the estimated stiffness is still stiffer than that measured at high porosities. In addition the generated unit cell does not look like the tomographic image of a real sample of foam concrete. Indeed, as the intermediate porosity becomes close to zero, the resulting microstructure becomes similar to that obtained by randomly placing overlapping spheres and the pore network is almost entirely percolated. Lastly, introducing a minimum thickness raises the estimated stiffness and does not solve the mismatch between numerical estimates and measurements.

Depending on the porosity, the following conclusions arise. On the one hand, for porosities lower than 40%, the upper Hashin-Strickman bound accurately estimates the overall stiffness. In addition, there is a little dependence to the pore size distribution, to the intermediate volume fraction or the minimum thickness. Indeed, the clustering of pores only results in bigger pores without affecting the stiffness. On the other hand, for porosities higher than 40%, the upper Hashin-Strickman bound slightly overestimates the effective stiffness of the considered unit cells. Hence, the differential scheme becomes a valuable alternative. Nevertheless, the estimated stiffness proves sensitive to the intermediate porosity and the minimum thickness and the mean-field schemes should be used with care.

However, the estimated Young Modulus is significantly higher than that experimentally measured on foam concrete. Thus, introducing the sand as a third phase could improve the numerical estimates and resolve the mismatch.

Foam concrete as a three-phase composite material 2.2.1. Questioning the scale separation

Scale separation of sand and pores is not ensured. Indeed, pores and sand particles features comparable radii, as shown on tomographic images and size distributions. Even if the sand particles were much smaller than the pores, scale separation would be questionable due to the small spacing δ between the pores. The order of magnitude of the spacing between the pores can be estimated by considering the mortar as a uniform coating on the pores. The thickness δ/2 of that coating must be consistent with the volume fraction of the solid skeleton, thus leading to the following equation:

N s pore -1 j=0 p j v pore,j = φ N s pore -1 j=0 p j v pore,j 1 + 3 4π 3v j δ 2 (8)
If a single diameter d pore is provided, this equation boils down to δ = d pore (φ -1/3 -1). If 61% of pores of equivalent diameter 1mm are considered, the spacing δ is about 179µm. If 67% of pores of equivalent diameter 0.5mm are considered, this thickness decreases to about 71µm. Therefore, the characteristic spacing between the pores is comparable to the size of the smaller sand particles (≥ 65µm). Consequently, it is likely that there are little or no sand particles in the thin members of the solid skeleton, where the stress concentrates.

As the cement paste is softer than the sand particles, neglecting this phenomenon is likely responsible for the overestimation of the overall stiffness encountered as foam concrete is considered as a two phase material. Therefore, sand is to be introduced as a third phase in the homogenization process. On the contrary, for the EPS concrete formulated in [START_REF] Miled | Particle size effect on {EPS} lightweight concrete compressive strength: Experimental investigation and modelling[END_REF], featuring porosities lower than 50% and EPS beads larger than 1mm, the thickness δ becomes larger than 260µm, which is comparable to the maximum diameter of the considered rounded quartz fine sand (300µm). As a consequence, the mortar can be considered as an homogeneous matrix and the differential scheme provides a valuable estimate of the Young modulus.

The present reasoning can be applied to other foam composite materials. For instance, the solid skeleton of the composite metal foam considered in [START_REF] Elbir | Sic-particulate aluminum composite foams produced by powder compacts: Foaming and compression behavior[END_REF] is made of an aluminum matrix reinforced by 22µm silicon carbide particles. The porosity is about 83% and the size of the pores is about 3mm. As a result, the spacing δ = 280µm is much larger than the particles and the particles are clearly visible on BSE and SE micrographs of the cell walls [START_REF] Elbir | Sic-particulate aluminum composite foams produced by powder compacts: Foaming and compression behavior[END_REF]. Another example is the case of reinforced polyurethane foams dedicated to acoustic insulation [START_REF] Sung | Influence of filler surface characteristics on morphological, physical, acoustic properties of polyurethane composite foams filled with inorganic fillers[END_REF]. The typical porosities and cell size are 95% and 350µm. The corresponding thickness is about 6µm. It is to be compared to the size of the fillers : 13µm for talc, 10µm for Zinc Borate and 33µm for Aluminum Hydroxide. The fact that Zinc Borate particles lead to a better improvement of the storage modulus than Aluminum Hydroxide particles [START_REF] Sung | Influence of filler surface characteristics on morphological, physical, acoustic properties of polyurethane composite foams filled with inorganic fillers[END_REF] is likely related to their smaller size: more of these particles are located in the thin cell walls and ribs, where the strain concentrates during dynamic tests. Nevertheless, while talk particles are slightly bigger than the Zinc Borate particles and despite the fact that they induce smaller pores, the storage modulus is improved. As explained in [START_REF] Sung | Influence of filler surface characteristics on morphological, physical, acoustic properties of polyurethane composite foams filled with inorganic fillers[END_REF], these particles are likely trapped in the thinner cell walls during the synthesis of the foam due to their particular hydrophobic property. Consequently, the size and microstructural locations of the particles are clearly influencing the elastic properties of the reinforced polyurethane foam, which is modelled as a three-phase material in [START_REF] Siegmann | Mechanical behavior of reinforced polyurethane foams[END_REF].

Regarding foam concrete, the segregation of sand particles at foam nodes could be attributed to a capillary effect. Indeed, in a fresh mortar, the occurrence of a big particle in a thin wall would lead to an important curvature of the mortar-pore interface. The resulting capillary pressure would push the particle toward the nodes of the foam so as to recover smooth spherical pores or flat interfaces. This phenomenon, investigated in [START_REF] Louvet | Capture of particles in soft porous media[END_REF][START_REF] Guignot | Liquid and particles retention in foamed suspensions[END_REF][START_REF] Khidas | Capture-induced transition in foamy suspensions[END_REF] ( see [START_REF] Wang | A critical review of the growth, drainage and collapse of foams[END_REF] for a review), led to the definition of a confinement parameter λ which is close to the ratio d p /δ, where d p is the diameter of the solid particles. The confinement parameter for a porosity of 61%, a pore diameter of 1mm and a particle diameter of 0.5mm is about λ ≈ 1.8. Therefore, the sand particles are likely trapped at the nodes of the considered foam concrete.

The segregation of sand particles at foam nodes is made obvious by analysing the 3D tomographic images (Fig. 1.1). Nevertheless, distinguishing the siliceous sand particles from the cement paste by processing the image is hardly achievable because the X-ray linear attenuations of both materials are similar within the considered X-ray energy range [START_REF] Lavergne | Contribution to the study of the time-dependent strains of viscoelastic composite materials[END_REF][START_REF] Stamati | Phase segmentation of concrete x-ray tomographic images at meso-scale: Validation with neutron tomography[END_REF]. As an alternative, the image is treated so as to segment the volume of the solid skeleton that could be occupied by spherical sand particles larger than a given volume (Fig. 2.3). To this end, a median filter is first applied to reduce noise (Fig. 2.3a) followed by an Ostu threshold to segment the solid skeleton (Fig. 2.3b). The segmented volume fraction for the foam concrete of porosity 70% is 42%. Therefore, the segmented solid skeleton likely overestimates the extend of the actual solid skeleton. Then an opening morphological filter is applied so as to define the volume that could be occupied by sand grains larger than a given diameter. Hence, the structuring element of the morphological filter is a ball featuring the considered diameter. As the diameter of the ball increases from 60µm to 144µm, the remaining volume accessible to sand grains significantly decreases from 37% to 24%. Hence, a significant volume of the solid skeleton is made of thin members of cement paste deprived of sand grains, since the considered sieve curve starts at an equivalent diameter of 250µm (Fig. 1.2). Lastly, the different connected components of the remaining volume are labelled (Fig. 2.3d). All these treatments are performed on the 3D image using the Insight Segmentation and Registration Toolkit (ITK) [START_REF] Johnson | The ITK Software Guide[END_REF][START_REF]The insight segmentation and registration toolkit[END_REF]. Unsurprisingly, right after thresholding, the solid skeleton is found to be made of a single continuous component, apart from some small spurious components near the edges of the image. On the contrary, once the thin members deprived of sand grains are removed, the solid skeleton becomes disconnected. Indeed, at diameter 144µm, in an abstract of size 1.8mm×3mm×3mm (Fig. 2.3), the three largest components respectively occupy volume fractions of 7.8%, 3.5% and 3.4% and components disconnected from the edges of the domain reach a volume fraction of 1.1%. It proves that these thin members are essential to the stiffness of the foam. First a median filter is applied, then an Ostu threshold filter, an opening by a ball and a labelling of the components. As the diameter of the ball increases, the colored components depicting the volume potentially occupied by sand grains of diameter bigger than that of the ball become disconnected.

3D numerical estimates of the stiffness

Both sand particles and pores are now accounted for. The volume fraction of sand in the mortar is f s/m = 35% and the size distribution of sand is plotted in figure 1.2 (maximum diameter D = 2mm), along with that of pores (maximum diameter D = 2.21mm). Hence, φ being the final porosity, the volume fraction of sand in foam concrete is f s/m (1 -φ). The total volume fraction of particles and pores is therefore f s/m (1 -φ) + φ. An intermediate volume fraction is to be introduced for the pores: it is set to 61% -f s/m (1 -φ) and sand particles are not modified during the dilation process, thus preserving their spherical shape. Now that the sand and the cement paste are distinguished, it is required to set a stiffness for each of these phases. The Young modulus and Poisson's ratio of the quartz sand particles are respectively set to E s = 95.5GPa and ν s = 0.078 [START_REF] Bass | Elasticity of Minerals, Glasses, and Melts[END_REF] (see also [START_REF] Heyliger | Elastic constants of natural quartz[END_REF]). The Poisson's ratio of the cement paste is set to 0.2 and its Young modulus to 14.1GPa so as to recover a Young modulus of E m = 24GPa for the mortar by using a Mori-Tanaka scheme. The bulk modulus and shear modulus estimated by the Mori-Tanaka scheme respectively are K m = 12.5GPa and G m = 10.2GPa. Such a Young modulus of 14.1GPa for a cement paste of water-to-cement ratio of 0.38 is very low compared to that measured on pure cement paste of similar age and water-to-cement ratios, around 20-25GPa [START_REF] Haecker | Modeling the linear elastic properties of portland cement paste[END_REF]. The difference can partly be attributed to the mismatch between static moduli and dynamic moduli, since dynamic moduli are measured in [START_REF] Haecker | Modeling the linear elastic properties of portland cement paste[END_REF] while the moduli considered in the present study are static moduli measured by performing three loading cycles between 0% and 30% of the compressive strength. Furthermore, it may be attributed to a slow hydration and curing conditions: the cement (CEM I 42.5N) is rather coarse as it features a Blaine fineness of 3200 cm 2 /g and the samples are unmolded at 6-8 hours and stored in open air. However, the occurrence of an ITZ or microcracking of the cement paste in the mortar are not visible on tomographic images of foam concrete. Lastly, though the Young modulus of the mortar estimated by using the Mori-Tanaka scheme is set so as to be consistent with its experimental measurement, the Poisson's ratio is estimated at 0.18 without experimental calibration.

The stiffness of grey voxels containing three phases must be set according to a mixing rule for elastic linear simulations using the FFT algorithm. First, the cement paste and sand grains are merged using a Reuss bound, as performed in [START_REF] Lavergne | Investigation of the effect of aggregates' morphology on concrete creep properties by numerical simulations[END_REF]. Then, this mix is merged with the pores using a Voigt bound, as justified in section 1.2.2.

First, the minimum thickness is set to zero. The porosity being set, the overall estimated moduli are softer that those estimated by assuming the mortar to be homogeneous (Fig. 2.4). As a result, the estimates of the Young modulus accounting for sand grains are more consistent with the experimental measurements of the Young moduli on foam concrete samples. In addition, the artificial unit cells seem closer to the tomographic images and the thin members squeezed between pores are actually deprived of sand grains (Fig. 1.1,2.5). Displaying components of the local stress unveils stress concentrations patterns (Fig. 2.5). On the one hand, at low porosity (29%), the stress concentrates in the sand particles, thus leading to an efficient stiffening of the concrete foam. On the other hand, at high porosity (67%), the stress concentrations are driven by the geometry of the solid skeleton. Indeed, the stress concentrates in thin members made of soft cement paste, thus explaining the lower overall estimated moduli: these thin members acting as soft microscopic pillars. The occurrence of large stress concentrations at high porosities could also explain the reduced compressive strength of concrete foam [START_REF] Youssef | Mechanical properties of non-autoclaved foam concrete: analytical models vs. experimental data[END_REF].

Finally, the estimated moduli are found to be sensitive to the minimum thickness of the members for high porosities. Indeed, introducing a minimum thickness of 24µm induces a 16% increase of the estimated Young moduli at 67% porosity. Nevertheless, the tomographic images do not feature such a measurable minimum thickness. In addition, keeping a null minimum thickness leads to the most accurate estimates of the overall moduli Young moduli.

Discussion

The proposed model of foam concrete does not account for the potential occurrence of an interfacial transition zone (ITZ) around the sand particles that might affect the elastic stiffness [START_REF] Neubauer | A three-phase model of the elastic and shrinkage properties of mortar[END_REF][START_REF] Scrivener | The interfacial transition zone (itz) between cement paste and aggregate in concrete[END_REF]. Indeed, the typical width of the ITZ is similar to the size of the cement grains (50µm), that is comparable to the voxel size (6mm/120) in the elastic simulations: representing the ITZ in elastic simulations at the scale of mortars would likely result in reducing the stiffness of grey voxels around the sand particles [START_REF] Grondin | How to consider the interfacial transition zones in the finite element modelling of concrete?[END_REF], replacing each sand particle by a bigger softer particle incorporating the ITZ [START_REF] Gal | Meso-scale analysis of FRC using a two-step homogenization approach[END_REF] or handle the ITZ as an imperfect interfaces [START_REF] Bary | Analytical and 3d numerical analysis of the thermoviscoelastic behavior of concrete-like materials including interfaces[END_REF] since the discretization is too coarse to represent the ITZ as an effective layer [START_REF] Garboczi | Elastic moduli of a material containing composite inclusions: effective medium theory and finite element computations[END_REF][START_REF] Lavergne | Investigation of the effect of aggregates' morphology on concrete creep properties by numerical simulations[END_REF]. The occurrence of an ITZ is often attributed to the wall effect: based on geometrical considerations, the volume fraction of big cement particles near sand grains must be very low. Hence, the water to cement ratio is these areas is supposed to be higher, thus triggering a lower water to cement ratio far from the sand grains. In the case of foam concrete featuring a high porosity, if the ITZ were limited to sand particles, the thin members of the foam deprived of sand particles would become stiffer. Since the stress concentrates in these thin members, the overall stiffness of the foam concrete could be improved.

On the contrary, the ITZ may not be limited to sand particles. Indeed, if it is admitted that the capillary forces are strong enough to exclude small sand particles from the thin members of the foam, this phenomenon might also affect big cement particles. It would result to a capillary-induced ITZ around the pores: the water to cement ratio in the thin member of the foam would be higher, thus lowering the stiffness where the stress concentrates. As a consequence, the overall stiffness would be lower.

Conclusion

The present study is dedicated to foam concrete, showing that the scales of sand, pores and cement paste can hardly be uncoupled to upscale the elastic stiffness. While 3D numerical simulations and mean field schemes provide consistent estimates of the overall stiffness of a biphasic composite material made of spherical pores embedded in a matrix, the comparison to some moduli measured on foam concrete samples is puzzling, notably for high porosities. To resolve the mismatch, foam concrete can be modelled as a three-phase material by considering the cement paste, the sand and the air voids. Indeed, at high porosity, the width of the mortar members becomes so small that the big sand particles are expected to be segregated at the foam nodes. Since the stiffness of the foam is mainly driven by the ribs and walls where the strain concentrates, the effective stiffness is lower than the one expected if the mortar was homogeneous. The thickness of the foam members can be estimated according to the equation 8 and compared to the maximum diameter of the sand to decide whether considering the mortar as homogeneous is acceptable. This procedure could be extended to various porous materials. Regarding the formulation of foam concrete, the present study unveils that a fine sand is often required to ensure an homogeneous mortar and an optimized stiffness.

The elastic stiffness is not the only physical property of interest as concrete foam is considered. Whether the sand is to be considered or not as other physical linear phenomena (electrical or thermal conductivity) are investigated remain unanswered. Indeed, the pores and the sand may both provide insulation and contrasts of electric or thermal conductivities can be different from the contrast of elastic stiffnesses. The numerical simulations can be adapted to provide estimates of the overall electrical of thermal conductivities if corresponding conductivities are provided for the air pores, the cement paste and the sand particles. Finally, the increase of local stress concentrations as porosity increases could be related to a lower overall strength of foam concrete.

Figure 1 . 2 :

 12 Figure 1.2: Pore size distributions are extracted from tomographic images at the Navier Laboratory and the size distribution of sand grains is obtained by using sieves. The pore size distribution slightly depends on the amount of foam added to the mix, i. e. the porosity.

Figure 1 . 3 :

 13 Figure 1.3: Top left: spheres are packed at the intermediate volume fraction of 61% according to the considered size distribution. The corresponding Laguerre diagram is colored in blue. Top right: the distance to the nearest pore is displayed. Bottom: the resulting periodic unit cells of size 6mm×6mm×6mm at porosities of 67%, the minimum thickness being b = 0 (left) and b = 24µm (right). A binning has be applied to provide an estimate of the local porosity.

Figure 1 . 4 :

 14 Figure 1.4: Effect of the discretization of the unit cell on the estimated stiffness. For a given cubic unit cell featuring 44 spherical inclusions and a volume fraction of 63.61%, the effective bulk modulus (left) and shear modulus (right) are estimated by using the FFT algorithm. Parameters to be tested are the discretization (64 3 , 96 3 , 128 3 , 144 3 , 160 3 ), the contrast (100,1000,10000) and the rule defining the stiffness of the grey voxels (Voigt or Reuss).

Figure 1 . 5 :

 15 Figure 1.5: Effect of the size of the unit cell on the estimated stiffness. For cubic unit cells of sizes 4 × 4 × 4 mm 3 , 6 × 6 × 6 mm 3 and 8 × 8 × 8 mm 3 and a porosity of 61%, the scaled effective bulk modulus K/K m (left) and shear modulus G/G m (right) are estimated by using the FFT algorithm. Error bars represent the standard deviation of each series of 10 unit cells.

Figure 2 . 1 :

 21 Figure 2.1: The scaled effective bulk modulus K/K m (left), shear modulus G/G m (middle) and Young modulus E/E m (right) estimated by numerical simulations are compared to those stemming from the Hashin-Shtrikman upper bound and the differential scheme. Different pore size distributions are considered for the 3D numerical simulations.

Figure 2 . 2 :

 22 Figure 2.2: The scaled effective bulk modulus K/K m (left), shear modulus G/G m (middle) and Young modulus E/E m (right) estimated by numerical simulations are compared to those stemming from the Hashin-Shtrikman upper bound and the differential scheme. Different intermediate volume fractions are considered for the 3D numerical simulations.

  (a) Median filter (b) Ostu filter (c) Opening, 60µm diameter (d) Opening, 144µm diameter

Figure 2 . 3 :

 23 Figure 2.3: 3mm×3mm slices of a treated tomographic abstract of size 1.8mm×3mm×3mm.First a median filter is applied, then an Ostu threshold filter, an opening by a ball and a labelling of the components. As the diameter of the ball increases, the colored components depicting the volume potentially occupied by sand grains of diameter bigger than that of the ball become disconnected.

Figure 2 . 4 :

 24 Figure 2.4: The scaled effective bulk modulus K/K m (left), shear modulus G/G m (middle)and Young modulus E/E m (right) estimated by numerical simulations are compared to those stemming from the Hashin-Shtrikman upper bound and the differential scheme for a porous material. The cement paste, the sand particles and the pores are distinguished for the 3D numerical simulations. The moduli estimated by the Mori-Tanaka scheme, using the cement paste as the matrix, are also displayed.

Figure 2 . 5 :

 25 Figure 2.5: Top: Generated three-phase artificial microstructures are discretized on regular grids for porosities of 29% (left) and 67% (right). Pores are black, sand grains are white and the cement paste is grey. Bottom: the component σ xx (x) of the local stress is displayed for a hydrostatic load Σ xx = Σ yy = Σ zz = 1..


Thomas Krajewski
email: thomas.krajewski@cpt.univ-mrs.fr

A renormalisation group approach to the universality of Wigner's semicircle law for random matrices with dependent entries

In this talk, we show that if the non Gaußian part of the cumulants of a random matrix model obey some scaling bounds in the size of the matrix, then Wigner's semicircle law holds. This result is derived using the replica technique and an analogue of the renormalisation group equation for the replica effective action.

Introduction

Random matrix theory (see the classical text [START_REF] Mehta | Random matrices[END_REF]) first appeared in physics in Wigner's work on the level spacing in large nuclei. Since then, it has proved to have multiple applications to physics and other branches of science, see for instance [START_REF]Oxford Handbook of Random Matrix Theory[END_REF]. Most of these applications rely on the universal behaviour of some of the observables for matrices of large size. A simple example is Wigner's semicircle law for the eigenvalue density that holds in the large N limit for matrices whose entries are independent and identically distributed.

Understanding the universal behaviour of eigenvalue distributions and correlations ranks among the major problems in random matrix theory. In this respect, the renormalisation group turns out to be a powerful technique. Introduced in the context of critical phenomena in statistical mechanics by K. Wilson to account for the universality of critical exponents, the latter has also proved to be useful in understanding probability theory. For instance it leads to an insightful proof of the central limit theorem, see the review by G. Jona-Lasinio [START_REF] Jona-Lasinio | renormalisation group and probability theory[END_REF] and references therein.

The renormalisation group has been used to derive the semicircle law for random matrices in the pioneering work of E. Brézin and A. Zee [START_REF] Brézin | Universal correlations between eigenvalues of large random matrices[END_REF]. In the latter approach, the renormalisation group transformation consists in integrating over the last line and column of a matrix of size N + 1 to reduce it to a size N matrix. This leads to a differential equation for the resolvent G(z) = 1/N Tr (z -M ) -1 in the large N limit whose solution yields the semicircle law.

In this talk, we follow a different route: We first express the resolvent as an integral over replicas and introduce a differential equation for the replica effective action. This differential equation is a very simple analogue of Polchinski's exact renormalisation group equation [START_REF] Polchinski | Renormalization and Effective Lagrangians[END_REF]. It is used to derive inductive bounds on the various terms, ensuring that the semicircle law is obeyed provided the cumulants of the original matrix model fulfil some simple scaling bounds in the large N limit.

This talk is based on some work in collaboration with A. Tanasa and D.L. Vu in which we extend Wigner's law to random matrices whose entries fail to be independent [START_REF] Krajewski | Wigner law for with dependent entries: A perturbative approach[END_REF], to which we refer for further details. There have been other works on such an extension, see [START_REF] Schencker | Semicircle law and freeness for random matrices with symmetries or correlations[END_REF], [START_REF] Hofmann-Credner | Wigner theorems for random matrices with dependent entries: Ensembles associated to symmetric spaces and sample covariance matrices[END_REF] and [START_REF] Friesen | The Semicircle Law for Matrices with Dependent Entries[END_REF].

What are random matrices ?

A random matrix is a probability law on a space of matrices, usually given by the joint probability density on its entries,

ρ(M ) = ρ(M11, M12, . . . ) (1)
Thus a random matrix of size N is defined as a collection of N 2 random variables. However, there is a much richer structure than this, relying notably on the spectral properties of the matrices.

Here we restrict our attention to a single random matrix. Note that it is also possible to consider several random matrices, in which case the non commutative nature of matrix multiplication plays a fundamental role, leading to the theory of non commutative probabilities.

There are two important classes of probability laws on matrices.

• Wigner ensemble: The entries are all independent variables,

ρ(M ) = i,j ρij(Mij), (2)
up to the Hermitian condition M ij = Mji.

• Unitary ensemble: The probability law is invariant under unitary transformations

ρ(U M U † ) = ρ(M ), (3)
for any unitary matrix U ∈ U(N ).

The only probability laws that belong to both classes are the Gaußian ones

ρ(M ) ∝ exp -1 2σ 2 Tr(M ) 2 , (4)
up to a shift of M by a fixed scalar matrix.

The main objects of interest are the expectation values of observables, defined as

O = dM ρ(M )O(M ).
(5)

Among the observables, the spectral observables defined as symmetric functions of the eigenvalues of M , play a crucial role in many applications. This is essentially due to their universal behaviour: In the large N limit, for some matrix ensembles and in particular regimes, the expectation values of specific spectral observables do not depend on the details of the probability law ρ(M ).

Universality is at the root of the numerous applications to physics and other sciences, since the results we obtain are largely model independent. Among the applications to physics, let us quote the statistics of energy levels in heavy nuclei, disordered mesocopic systems, quantum chaos, chiral Dirac operators, ...

Wigner's semicircle law

In this talk, we focus on the eigenvalue density, defined as

ρ(λ) = 1 N 1≤i≤N δ λ -λi M √ N . (6)
In particular, a universal behaviour is expected in the large N limit for some ensembles. For a Gaußian random Hermitian matrix ρ(M ) ∝ exp -1 2σ 2 Tr(M 2 ), the eigenvalue density obeys Wigner semicircle law,

lim N →∞ R dλ λ k ρ(λ) = 1 2πσ k+2 2σ -2σ dλ λ k √ 4σ 2 -λ 2 if k is even, 0 if k is odd. (7)
Empirically, ρ(λ) may be determined by plotting the histogram of eigenvalue of a matrix taken at random with a given probability law, see figure 1.

The derivation of Wigner's semicircle in the large N limit is based on the resolvent (also known as the Green function)

G(z) = lim N →+∞ 1 N Tr z - M √ N -1 = z 2σ 2 1 -1 - 4σ 2 z 2 . ( 8
)
Then, the density of eigenvalues is recovered as

ρ(λ) = G(λ -i0 + ) -G(λ + i0 + ) 2iπ , (9)
1 M = 1 Nz Tr 1 1 M z = k=0 1 Nz k+1 hTr(M k )i.
(2.5) our computation of the average value of the trace and its function very helpful. e average density of state for the Gaussian Unitary Ensemble GUE mit an that in the large N limit, the

(2k)! k!(k + 1)!z 2k+1 1 r 1 4 z 2 ◆ (2.6) e obtain z 2 if z 2 [ 2, 2] otherwise (2.7)
ar law holds in fact for a much e GUE. Chapter 3 is devoted to ge N limit, the distribution of eigenvalues might become continuous. ntinuous probability measure ⇢ on the real line, instead of Green's f Stieltjes transform which is defined by

S ⇢ (z) = Z R d ⇢( ) z . (2.8)
fined and analytic for z outside the support of ⇢, in particular in f the complex plane. If the support of the probability measure ⇢ is ments then ⇢ can be obtained from its Stieltjes transform through

(z) = i 2⇡ [S ⇢ (z + i✏) S ⇢ (z i✏)]. (2.9)
es transform of a probability measure can also be regarded as the ts

⇢( ) = Z R d 1 z (1 + z + 2 z 2 + ...)⇢( ) = 1 X k=0 1 z k+1 E k .
(2.10) 15 where we have used the relation

1 x ± i0 + = p.v. 1 x ∓ 2iπδ(x). (10)
In the large N limit, for the Gaußian model, the resolvent obeys the self consistency equation (also known as the Schwinger-Dyson equation), see for instance [START_REF] Zee | Quantum field theory in a nutshell[END_REF], section VII.4,

G(z) = ∞ k=0 σ k G k (z) z k+1 = 1 z -σG(z) . (11)
Its solution that behave as 1/z for large z is

G(z) = z 2σ 2 1 -1 - 4σ 2 z 2 . ( 12
)
Taking the cut of the square root on the negative real axis, we obtain the Wigner semicircle law (7) in the large N limit. The semicircle law is not limited to the Gaußian case, it also holds for Wigner matrices in the large N limit. A random Hermitian N × N matrix is a Wigner matrix if

• real and imaginary parts of upper diagonal elements are independent and identically distributed (i.i.d.) with mean 0 and variance σ;

• diagonal elements are i.i.d. with finite mean and variance and independent of the off diagonal ones.

Then, in the limit N → +∞, the eigenvalue distribution of M √ N is the semicircle law [START_REF] Schencker | Semicircle law and freeness for random matrices with symmetries or correlations[END_REF].

The original proof is of combinatorial nature and involves the expectation of the moments lim

N →+∞ 1 N k/2+1 Tr(M k ) = (2l)! (l!) 2 (l+1) for k = 2l even 0 for k odd. ( 13
)
To derive this result, the idea is to first factorise ρ for a Wigner ensemble as

ρ(M ) = i ρ (Mii) i<j ρ (Re Mij)ρ (Im Mij), ( 14
)
where ρ is the common probability density of the real diagonal terms and ρ the common probability density of the real and imaginary parts of the off diagonal terms. Then, we expand the trace and integrate over the independent real variables Mii, Re Mij and Im Mij. The power of N in the expectation of a given moment arises from the denominator 1 N k/2+1 and from the number of independent indices in the summations. In the large N limit, the only configurations that survive are counted by Catalan numbers, C l = (2l)! (l!) 2 (l+1) . Since the latter also appear in the following Taylor expansion

z 2σ 2 1 -1 - 4σ 2 z 2 = l≥0 (2l)! (l!) 2 (l + 1) σ 2l z 2l+1 , (15)
we conclude that

G(z) = lim N →+∞ 1 N Tr z - M √ N -1 (16) = lim N →+∞ ∞ k=0 1 z k+1 1 N k/2+1 Tr(M k ) (17) = z 2σ 2 1 -1 - 4σ 2 z 2 . ( 18
)
This is the form of the resolvent that leads to Wigner's semicircle law.

Here, we see universality at work: In the large N limit, the eigenvalue density is given by the semicircle law, whatever the probability densities ρ and ρ are. However, this result relies on the independence of the matrix elements. In the next section, we will extend it to matrices whose entries are not necessary independent.

Wigner's law beyond Wigner ensembles

Let us introduce the cumulants, defined through their generating function

Mi 1 j 1 • • • Mi k j k c = ∂ ∂Jj 1 i 1 . . . ∂ ∂Jj k i k log exp Tr(M J) J=0 . (19)
In the physics terminology, these are the connected correlation functions.

In particular, the Gaußian cumulants vanish beyond the quadratic term

ρ(M ) ∝ exp -1 2σ 2 Tr(M 2 ) ⇒ MijM kl c = σ 2 δ il δ jk vanish otherwise . ( 20
)
Therefore, cumulants of degree higher than 2 are a measure of the deviation from the Gaußian case.

Turning back to the general case, for each cumulant we construct an oriented graph as follows (see figure 2 for some examples):

• vertices are distinct matrix indices in the cumulant, • there is an edge from i to j for every Mij.

i j k l i j Since non quadratic cumulants measure deviations from the Gaußian case, if the perturbation is small it is reasonable to expect that the semicircle law is still obeyed.

(M ij ) 2 M jk M ll c M ij M ji c
To state this result, recall that an oriented graph is Eulerian if every vertex has an equal number of incoming and outgoing edges. Equivalently, it means that every connected component admits an Eulerian cycle, i.e. an oriented cycle that passes through all edges, respecting the orientation. Furthermore, let us denote v(G), e(G), c(G) the number of vertices, edges and connected components of G. Theorem 1 (Wigner's law for matrices with dependent entries). Let ρN be a probability law on the space of Hermitian N ×N matrices M such that its cumulants can be decomposed as CG = C G + C G , with C G a Gaussian cumulant and C G a perturbation such that, uniformly in the vertex indices i1, ..., i v(G) (i.e. all constants involved should not depend on these indices),

• lim

N →∞ N v(G)-c(G)-e(G)/2 C G (i1, ..., i v(G) ) = 0 if G is Eulerian, • N v(G)-c(G)-e(G)/2 C G (i1, ..., i v(G) ) bounded if G is not Eulerian.
Then, the moments of the eigenvalue distribution of the matrix M √ N converge towards the moments of the semicircle law, with σ given by the Gaussian cumulant MijM kl c = σ 2 δ il δ jk , lim

N →∞ R dλ λ k ρN (λ) = 1 2πσ k+2 2σ -2σ dλ λ k √ 4σ 2 -λ 2 if k is even, 0 if k is odd. (21)
For instance, for the graph i j k l which is not Eulerian, with v = 3, e = 4 and c = 2, the cumulant should obey

1 N (Mij) 2 M jk M ll c ≤ K, (22)
with K a constant that does not depend on the indices i, j, k and l. On the other hand, for the graph i j , which is Eulerian, with v = 2, e = 4 and c = 1, we impose

lim N →+∞ 1 N (Mij) 2 (Mji) 2 c = 0, (23)
uniformly in i and j.

As an illustration, we recover the case of Wigner matrices (with finite moments). Indeed,

• there is no graph with v ≥ 3 (independence of off diagonal matrix elements);

• for v = 1 and v = 2, e ≥ 3, bounds are satisfied because of 1/N e/2 and all moments are assumed to be finite;

• C i j (i, j) = Mij c = Mij = 0 (off diagonal elements have mean value 0); • C i j (i, j) = MiiMjj c = MiiMjj -Mii Mjj = 0 (indepen- dence of diagonal elements); • C i j (i, j) = MijMjj c = MijMjj -Mij Mjj = 0 (
independence of diagonal and off diagonal elements);

• C i j (i, j) = MijMij c = MijMij -Mij Mij = (ReMij) 2 -(ImMij) 2 + 2i ReMijImMij = 0
(independence of real and imaginary parts and equality of their distributions with mean value 0);

• C i j
= σ 2 is the Gaußian cumulant leading to the semicircle law.

The case of unitarily invariant matrices is critical since the bounds are saturated, see [START_REF] Krajewski | Wigner law for with dependent entries: A perturbative approach[END_REF]. This is consistent since we known the semicircle law is not obeyed by unitary non Gaußian ensembles [START_REF] Brezin | Planar Diagrams[END_REF].

is possible to give a combinatorial proof of this result based on the relation between moments and cumulants,

Mi 1 j 1 • • • Mi k j k = I 1 ,...,Ip partition of {(i 1 ,j 1 ),••• ,(i k ,j k )} ij∈I 1 Mij c • • • ij∈Ip Mij c (24)
In the moment method, we have to estimate

1 N k/2+1 Tr(M k ) = 1 N k/2+1 1≤i 1 ,...,i k ≤N Mi 1 i 2 • • • Mi k i 1 . (25)
Then, we express the moments in (25) in terms of cumulants using (24) and represent each cumulant as a graph. Because of the trace, one has to draw Eulerian cycles on the graphs after some vertex identifications.

Then, the scaling bounds on the cumulants can be used to show that only Gaussian terms survive.

Proof based on the replica effective action

Let us give a renormalisation group proof of this result based on the replica effective action. The use of replicas in random matrix theory is a classical subject, see for instance [START_REF] Kanzieper | Replica Approach in Random Matrix Theory[END_REF] or [START_REF] Zee | Law of addition in random matrix theory[END_REF]. To begin with, let us note that

Tr z - M √ N -1 = ∂ ∂z log det z - M √ N . ( 26
)
It is convenient to express the logarithm using the replica method. First observe that

log(A) = lim n→0 A n -1 n . (27)
Then, we express the n th power of the determinant as a Gaußian integral over n replicas of a complex vector of size N (with a factor of π nN included in the measure),

1 det n (z -M ) = dX exp -Tr X † (z -M )X , (28)
which fit into a N × n complex matrix X = (Xi,a) 1≤i≤N 1≤a≤n

.

The limit n → 0 may be worrisome, its meaning is as follows. Because of U(n) invariance, any perturbative result in powers of 1/z is a polynomial in n, from which we retain only the linear term. Of course, this may not hold beyond perturbation theory, where replica symmetry breaking can occur.

After averaging over M with the random matrix density ρ(M ), we obtain the following expression for the resolvent,

G(z) = - 1 N ∂ ∂z dX exp -Tr(X † X) + V0(X) order 1 in n , (29)
where the replica potential is

V0(X) = log exp Tr X † M √ N X . (30)
Because of the logarithm, the potential involves the cumulants and can be expanded over graphs as

V0(X) = G oriented graph 1 |Aut(G)|N e(G)/2 1≤i 1 ,...,i v(G) ≤N all different CG(i1, . . . , i v(G) ) e edge (XX † )i s(e) i t(e) , (31)
where s(e) is the source of edge e and t(e) its target.

Let us introduce a replica effective action, obtained by a partial integration

V (t, X) = log dY exp - Tr(Y † Y ) t + V0(X + Y ) -N n log t. (32)
The parameter t ranges between 0 (where we have no integration, V (t = 0, X) = V0(X)) and t = 1/z.

The effective potential obeys a semi-group property that follows from Gaußian convolution (see for instance [START_REF] Zinn-Justin | Quantum field theory and critical phenomena[END_REF], section A10.1),

V (t + s, X) = log dY exp - Tr(Y † Y ) s + V (t, X + Y ) -N n log s.
(33) For small s = dt, it translates into the following renormalisation group equation, which is a simple version of Polchinski's exact renormalisation group equation [START_REF] Polchinski | Renormalization and Effective Lagrangians[END_REF],

∂V (t, X) ∂t = i,a ∂ 2 V (t, X) ∂Xi,a∂Xi,a + ∂V (t, X) ∂Xi,a ∂V (t, X) ∂Xi,a . (34)
The first term on the RHS is referred to a the loop term, since it creates a new loop in the Feynman graph expansion of the effective action while the second inserts a one particle reducible line and is referred as the tree term, see figure 3. Taking into account the boundary condition V (t = 0, X) = V0(X), it is convenient to write (34) in integral form

∂ ∂t = +
V (t, X) = V0(X)+ t 0 ds i,a ∂ 2 V (s, X) ∂Xi,a∂Xi,a + ∂V (s, X) ∂Xi,a ∂V (s, X) ∂Xi,a . (35)
This allows us to derive inductive bounds in powers of t = 1/z. From a physical point of view, we evaluate the effective potential by a large succession of small partial integrations, with a total weight given by t. Let us stress that in our context this differential equation is merely a tool to control the t dependence of the effective action after integrating with a t dependent propagator.

The effective potential also admits an expansion over graphs,

V (t, X) = G oriented graph 1 |Aut(G)|N e(G)/2 1≤i 1 ,...,i v(G) ≤N all different CG(t; i1, . . . , i v(G) ) e edge (XX † )i s(e) i t(e) . (36)
This leads to a graphical interpretation of the action of the two differential operators in the renormalisation group equation, see figure 4. Indeed, in the expansion (36), an edge joining a vertex carrying label i to a vertex carrying j is equipped with a factor a Xi,aXj,a, with a a replica index. Then, the differential operator ∂ ∂X i,a (resp.

∂ ∂X j,a

) removes the outgoing (resp. incoming) half edge. Finally, the remaining half edges are reattached and the vertices identified to yield a new graph on the RHS of (35), with one less edge. These operations are performed on the same graph for the loop term and on distinct ones for the tree term.

∂ 2 V ∂X i,a ∂X i,a (loop term) : v G → v G ∂V ∂X i,a ∂V ∂X i,a (tree term) : v 1 G 1 v 2 G 2 → v 1 = v 2 G

Figure 4: Action of differential operators on the vertices of the effective action

Let us decompose the effective cumulants appearing in (36) into Gaußian ones and perturbations, and expand both in a power series in t = 1/z,

CG(t) = ∞ k=0 t k C (k) G Gaußian + C (k) G perturbation . (37)
The Gaußian terms are those that are constructed using only the Gaußian term in the initial potential V0(X). Even is V0(X) is quartic in X, this does not hold for the Gaußian part of Vt(X), that contains terms of all orders. The perturbation collects all the remaining terms, they contain at least one non Gaußian perturbation from V0(X).

The renormalisation group equation (35) allows us to prove inductively on k that the perturbations C (k) G obey the same scaling bound imposed on C (0) G = C G (0) and that the purely Gaußian terms do not grow to fast,

• lim N →∞ N v(G)-c(G)-e(G)/2 C (k) G order 0 in n = 0 for G Eulerian, • N v(G)-c(G)-e(G)/2 C (k) G order 0 in n bounded for G not Eulerian, • N v(G)-c(G)-e(G)/2 C (k) G order 0 in n bounded for any G.
This involves a combinatorial discussion based on the graphical interpretation of figure 4 that can be found in [START_REF] Krajewski | Wigner law for with dependent entries: A perturbative approach[END_REF]. Let us simply mention that the terms that may violate the bounds are of higher order in n. Thus, they are harmless when taking the limit n → 0 before the limit N → +∞.

Finally, using (29) and the renormalisation group equation (35), the resolvent can be expressed as

G(z) = 1 z + 1 N 3/2 z 2 1≤i≤N C i (1/z; i) order 0 in n . (38)
The scaling bounds for the non Gaußian cumulants impose, perturbatively in 1/z, lim

N →∞ 1 √ N C i (1/z; i) order 0 in n = 0. (39)
Therefore, only the Gaußian cumulants contribute and we recover Wigner's semicircle law.

Conclusion and outlook

this talk, we have argued that Wigner's semicircle law remains valid for matrices with dependent entries. The deviation from the independent case is measured by the joint cumulants of the entries, which are assumed to fulfil some scaling bound for large N . To establish this result, we have introduced an effective action for the replicas. This effective action obeys a renormalisation group equation that allowed us to prove perturbative bounds on the effective cumulants. As a consequence of these bounds, only the Gaußian terms contribute in the large N limit, thus establishing the validity of Wigner's semicircle law.

It may also be of interest to investigate the case of the sum of a random matrix M and a deterministic one A, see for instance [START_REF] Zee | Law of addition in random matrix theory[END_REF] where such a model is discussed. In this case, the resolvent is expressed as

G(z) = - 1
In our context, the deterministic matrix A induces a non trivial kinetic for the replicas. In particular, if A is a discrete Laplacian, it yields a non trivial renormalisation group flow that bears some similarities with the QFT renormalisation group. In this case, we expect to exploit the true power of the renormalisation group equation, with a discussion of fixed points and scaling dimensions.

Figure 2 . 1 :

 21 Figure 2.1: Eigenvalue histogram of a 1000⇥1000 normalized Gaussian matrix

Figure 1 :

 1 Figure 1: Spectrum of a large (N = 1000) Hermitian matrix with a Gaußian distribution (σ = 1)

Figure 2 :

 2 Figure 2: Examples of graph associated to cumulants

Figure 3 :

 3 Figure 3: Graphical interpretation of the renormalisation group equation

  exp -Tr(X † (A + z)X) + V0(X)

Acknoledgements I wish to thank the organisers of the workshop, Mario Gattobigio, Thierry Grandou and Ralf Hoffmann for their kind invitation as well as Stan Brodsky for insightful remarks. I am also grateful to my collaborators, Dinh Long Vu and Adrian Tanasa.