Didier Lasseux

Francisco J Valdés Parada

Mark L Porter

An improved macroscale model for gas slip-flow in porous media

We report on a refined macroscopic model for slightly compressible gas slip-flow in porous media developed by upscaling the pore-scale boundary value problem. The macroscopic model is validated by comparisons with an analytic solution on a 2D ordered model structure and with direct numerical simulations on random microscale structures. The symmetry properties of the apparent slip-corrected permeability tensor in the macroscale momentum equation are analyzed. Slip-correction at the macroscopic scale is more accurately described if an expansion in the Knudsen number, beyond the first order considered so far, is employed at the closure level. Corrective terms beyond the first order are a signature of the curvature of solid-fluid interfaces at the pore-scale that is incompletely captured by the classical first order correction at the macro-scale. With this expansion, the apparent slip-corrected permeability is shown to be the sum of the classical intrinsic permeability tensor and tensorial slip-corrections at the successive orders of the Knudsen number. All the tensorial effective coefficients can be determined from intrinsic coupled, but easy-to-solve, closure problems. It is further shown that the complete form of the slip boundary condition at the microscale must be considered and an important general feature of this slip-condition at the different orders in the Knudsen number is highlighted. It justifies the importance of slip-flow correction terms beyond the first order in the Knudsen number in the macroscopic model and sheds more light on the physics of slip-flow in the general case, especially for large porosity values. Nevertheless, this new non-linear dependence of the apparent permeability with the Knudsen number should be further verified experimentally.

Introduction

The interest for studying slightly compressible gas slip-flow in channels of characteristic dimension comparable to the gas mean free path at the pressure and temperature under consideration is tremendous for many applications that encompass gas flow in microchannels and nanofluidic systems [START_REF] Porodnov | Experimental investigation of rarefied gas flow in different channels[END_REF][START_REF] Harley | Gas flow in micro-channels[END_REF][START_REF] Karniadakis | Microflows and nanoflows. Fundamentals and simulation[END_REF][START_REF] Cai | Gas flows in microchannels and microtubes[END_REF], characterization of low permeable porous materials [START_REF] Lasseux | Permeability measurement of graphite compression packings[END_REF]Jannot & Lasseux 2012;[START_REF] Profice | Permeability, porosity and Klinkenberg coefficient determination on crushed porous media[END_REF] involved in processes ranging from gas production [START_REF] Darabi | Gas flow in ultra-tight shale strata[END_REF], gas or nuclear waste storage, filtration and separation [START_REF] Chmielewski | Aerosol filtration with slip flow[END_REF], composite manufacturing [START_REF] Zhang | Cvd in weakly rarefied rotating disk flows[END_REF], among many others. Consequently, a great deal of interest may also be focused on the prediction and estimation of the coefficients governing the physics of gas transport at the macroscale. The emergence of powerful imaging methods together with the increase of computational capabilities make this prediction from direct numerical simulation, based upon images performed on a Representative Elementary Volume (REV) of the porous structure at the microscale, a realistic and promising approach to avoid sophisticated experiments [START_REF] Ghaddar | On the permeability of unidirectional fibrous media: A parallel computational approach[END_REF][START_REF] Nakashima | Estimate of transport properties of porous media by microfocus X-ray computed tomography and random walk simulation[END_REF][START_REF] Marino | Gas flow in a permeable medium[END_REF][START_REF] Prodanović | 3D image-based characterization of fluid displacement in a Berea core[END_REF].

To progress towards this goal, our purpose is to carefully derive an accurate macroscopic model for gas flow in homogeneous porous media in the slip regime. The existing model obtained so far in this context [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF][START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF]) is inaccurate. In fact, we show that the existing macroscopic slip correction filters out some important pore-scale topological information related to the curvature of solid-fluid interfaces, which can be of significant impact in some circumstances at the macroscopic level. The model derived in this work avoids such a simplification and yields a macro-scale model together with associated closures that are easy-to-solve boundary value problems defined at the microscale on a REV, which provide the required macroscopic coefficients.

The development proposed hereafter is organized as follows. The problem statement is presented first in section 2 together with the main steps of the upscaling leading to the non-closed macroscopic form. The procedure to close the macroscopic momentum equation is reported in sections 3.1 and 3.2 providing the closure problem that defines the apparent slip-corrected permeability, K s , the non-intrinsic effective coefficient appearing in the macroscopic model. Validation of the upscaled model is presented in section 3.3. Solutions obtained from Direct Numerical Simulations (DNS) of incompressible slip-flow performed on random 2D structures are used in the predictive macroscopic averaged model to identify K s that is further compared to the solution of the closure problem. In addition, the analytic estimation of K s is compared to the values obtained from numerical solutions of the closure over a periodic unit cell of the structure. In both cases, K s depends nonlinearly on the Knudsen number at the macroscale, Kn, based on the characteristic length-scale within the fluid-phase, a behavior that can not be predicted from porous media slip-flow models reported so far [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF][START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF].

In order to elucidate the nonlinearity and to obtain a macroscale model involving coefficients that are intrinsic to the micro-structure, an expansion in Kn is provided in section 4. Effective coefficients in this model correspond to the intrinsic permeability tensor at the zeroeth order and to slip-correction tensors at the successive higher orders in Kn, all of them being obtained from the solution of coupled closure problems. This is detailed in section 4.1. The existing Darcy-Klinkenberg model extensively used for gas slip flow in porous media, corresponds to the expansion up to the first order, together with an ideal gas law. In section 4.2, it is shown that a more accurate description of the macroscopic behavior, capturing the nonlinearity in Kn previously observed, requires higher order terms of the expansion. These terms are an explicit signature of the solidfluid interface curvature that is not present in the first-order term. This is assessed in section 4.3, where important features of the slip boundary condition with respect to the successive orders of the expansion are reported, shedding more light on the physics of slip-flow far beyond the application to porous media. Symmetry properties of the tensor K s are investigated in Appendix A. Conclusions are drawn in section 5.

Problem statement and upscaling

The starting point of the analysis is the initial boundary value problem at the pore scale with characteristic length β (see Fig. 1) describing the isothermal, Newtonian and slightly compressible flow of a barotropic fluid β occupying the sub-domain V β inside an averaging domain V of radius r 0 .

We recall the problem statement as

∂ρ β ∂t + ∇ • (ρ β v β ) = 0 in V β (2.1a) 0 = -∇p β + µ β ∇ 2 v β in V β (2.1b) ρ β = F (p β ) in V β (2.1c) B.C.1 v β = -ξλ β (I -nn) • n • ∇v β + ∇v T β at A βσ (2.1d ) B.C. 2 v β = f (t) at A βe (2.1e) I.C. ρ β = ρ β,0 when t = 0 (2.1f )
In these equations, ρ β , v β and p β are the density (of initial value ρ β,0 ), velocity and pressure in the β-phase obeying the state equation (2.1c), µ β is the dynamic viscosity, taken as a constant, λ β is the mean free path, ξ is a parameter related to the tangential momentum accommodation coefficient, σ v , given by ξ = 2-σv σv [START_REF] Agrawal | Survey on measurement of tangential momentum accommodation coefficient[END_REF][START_REF] Selden | Analysis of accommodation coefficients of noble gases on aluminum surface with an experimental/computational method[END_REF]. In practice, ξ ranges from roughly 1.3 to 1.7 [START_REF] Suetin | Poiseuille flow at arbitrary Knudsen numbers and tangential momentum accomodation[END_REF][START_REF] Arkilic | Mass flow and tangential momentum accommodation in silicon micromachined channels[END_REF][START_REF] Perrier | Mass flow rate measurements in microtubes: From hydrodynamic to near free molecular regime[END_REF]). The unit normal vector at the solid-fluid interface A βσ , directed from the β-to the solid σ-phase, is denoted by n while A βe represents the entrances and exits surfaces of the β-phase at the boundaries of the averaging domain on which the fluid velocity is f (t). The first order slip boundary condition expressed in Eq. (2.1d) and originally derived by Navier and Maxwell [START_REF] Maxwell | On stresses in rarefied gases arising from inequalities of temperature[END_REF], includes the complete strain rate ∇v β + ∇v T β . The ∇v T β term plays an important role in general [START_REF] Einzel | Boundary condition for fluid flow: curved or rough surfaces[END_REF][START_REF] Panzer | The effect of boundary curvature on hydrodynamic fluid flow: calculation of slip length[END_REF][START_REF] Lockerby | Velocity boundary condition at solid walls in rarefied gas calculations[END_REF][START_REF] Barber | Isothermal slip flow over curved surfaces[END_REF][START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF], although it has been often omitted for porous media flow [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF][START_REF] Pavan | Gas/solid heat transfer in gas flows under Klinkenberg conditions: comparison between the homogenization and the kinetic approaches[END_REF]. More physical insight will be provided into this boundary condition in section 4.2 below.

While writing the initial boundary value problem given by Eqs. (2.1), a set of constraints were assumed that can be listed as

ρ β ∆ρ β v β t * ρ β (1 + ξKn) L ρ β 1
(2.2)

ρ β ∆ρ β 1 1 + ξKn L ρ β 2 1
(2.3)

ρ β (1 + ξKn) 2 β µ β t * v 1 (2.4)
The first two constraints are required to neglect viscous terms related to compressible effects in the momentum equation whereas the last one is necessary to neglect the temporal acceleration; in both cases the comparisons are made with respect to the viscous diffusion term. In addition, to neglect the convective acceleration term, we require

Re = ρ β β v β µ β 1 (2.5)
In the above constraints, which are easily met in practice and were thoroughly discussed in a recent article [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF], v β represents the order of magnitude of v β whereas L ρ and t * ρ are the characteristic length-and time-scales over which ρ β experiences significant variations ∆ρ β , Kn = λ β β is the Knudsen number, which is assumed to remain smaller than about 0.1 to 0.4, characterizing the slip-flow regime in straight channels [START_REF] Harley | Gas flow in micro-channels[END_REF][START_REF] Maurer | Second-order slip laws in microchannels for helium and nitrogen[END_REF]. The left hand side in the relationship (2.4) is nothing else than the frequency parameter of the flow with t * v being the characteristic time over which v β experiences significant variations and Re in (2.5) is the Reynolds number of the flow.

Upscaling the problem given by Eqs. (2.1) can be performed with the aid of the volume averaging method [START_REF] Whitaker | The method of volume averaging[END_REF] for which the intrinsic and superficial averages of any physical variable ψ β are defined as

ψ β β = 1 V β V β ψ β dV (2.6) ψ β = ε ψ β β = 1 V V β ψ β dV (2.7)
ε denoting the porosity of the medium. The development of the averaged equations makes use of the averaging theorem (or Leibnitz rule) [START_REF] Howes | The spatial averaging theorem revisited[END_REF] given by

∇ψ β = ∇ ψ β + 1 V A βσ nψ β dA (2.8)
along with the Reynolds theorem [START_REF] Truesdell | The classical field theories[END_REF] and the classical decomposition [START_REF] Gray | A derivation of the equations for multiphase transport[END_REF])

ψ β = ψ β β + ψ β (2.9)
With this decomposition, it must be noted that ψ β β has a typical length-scale of variation L (see Fig. 1), while ψ β varies at the scale β . The upscaling is subject to the scale hierarchy defined by β r 0 L, which can be used to prove that ψ β β β = ψ β β and consequently that ψ β β = 0 [START_REF] Whitaker | The method of volume averaging[END_REF]. When the procedure is applied to the initial boundary value problem in Eqs. (2.1), the following averaged mass, momentum and state equations are obtained [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF])

∂ ρ β β ∂t + ∇ • ρ β β v β β = 0 (2.10) 0 = -∇ p β β + µ β ∇ 2 v β β + 1 V β A βσ n• (-I p β + µ β ∇ṽ β ) dA (2.11) ρ β β = F p β β (2.12)
To arrive at this form, the assumption of a slightly compressible flow, that is expressed through the constraint

ρ β ρ β β (2.13)
was employed, along with the homogeneity of the medium that allows treating ε as a constant.

The macroscale momentum equation (2.11) is not closed as it involves pressure and velocity deviations that are essentially varying at the pore-scale. Therefore, one needs to derive and solve the associated closure problem that is written in terms of p β and ṽβ over a periodic unit cell representative of the macroscopic region as [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF])

∇ • ṽβ = 0 in V β (2.14a) 0 = -∇ p β + µ β ∇ 2 ṽβ - 1 V β A βσ n• (-I p β + µ β ∇ṽ β ) dA in V β (2.14b) ṽβ + ξλ β (I -nn) • n • ∇ṽ β + ∇ v T β = -v β β at A βσ (2.14c) pβ β = 0 (2.14d ) Periodicity ṽβ (r + l i ) = ṽβ (r) , p β (r + l i ) = p β (r) , i = 1, 2, 3 (2.14e)
This problem simply results from the subtraction of the averaged equations from their initial counterparts expressed in Eqs (2.1). In Eq. (2.14c), λ β represents the mean free path at the average density, which defines the macro-scale Knudsen number, Kn = λ β β that will be used throughout this work. This definition of the Knudsen number is the most physically relevant as it reflects the wall effect occurring in the Knudsen layer (of thickness roughly given by the mean free path) relative to the size of the channel where the flow is taking place. Upon assuming molecular collisions between hard spheres and considering the slightly compressible hypothesis, λ β can be expressed as [START_REF] Cowling | Molecules in motion[END_REF])

λ β = M √ 2πN A ρ β β δ 2 (2.15)
Here, N A is the Avogadro's number, M the molar mass of the gas and δ the gas particle diameter, πδ 2 being the effective collision section. While Eq. (2.14d) is necessary to define a well-posed problem, the analog for the velocity deviation ṽβ β = 0 is not, although this last relationship will be employed later. Equations (2.14) are subject to the slightly compressible flow hypothesis expressed in Eq. (2.13) and to the constraint of scale hierarchy β L similar to that employed in homogenization [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF].

At this point, a closure procedure relating deviations to average quantities in the boundary value problem (2.14) is required in order to derive a closed macroscopic model. This is carried out in the next section, where a Darcy-like average momentum equation is obtained in which the apparent slip-corrected permeability, K s , is determined from the associated closure problem.

Closure and macroscopic model

A procedure to close the averaged momentum equation (2.11) from the deviations equations (2.14) was proposed and validated earlier, yielding a macroscopic model at O ξKn , which involves the intrinsic permeability tensor K and a slip-flow correction tensor S [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF]. Both tensors are obtained from the solution of intrinsic closure problems that are coupled through a 'slip-like' boundary condition. Whereas the closure problem for K can be written in a 'Stokes-like' form, S derived from an integro-differential closure. In terms of solution strategy and physical interpretation, this represents a real difficulty that is removed in the development presented below.

Closure problem

A convenient closure is now derived by noticing that the boundary value problem for the deviations (2.14) at the pore scale is made non-homogeneous due to the macroscopic source term v β β in the right hand side (rhs) of the boundary condition (2.14c). Since this problem is incompressible in nature, λ β can be regarded as a parameter and does not represent a source for p β and ṽβ despite its dependence on ρ β β . Due to the linearity of the problem, the solution can be sought as a linear combination of the source, namely

ṽβ = C • v β β (3.1) p β = µ β c • v β β (3.2)
Any additive constant in the representations of ṽβ and p β can be shown to be unimportant in the final macroscopic model, just as for the upscaling of the classical incompressible Stokes flow with no slip [START_REF] Whitaker | The method of volume averaging[END_REF]. Using these representations in the closure equations (2.14), while treating v β β as a constant due to the contrast of scale between β and L, yields the following closure problem in terms of the closure variables C and c

∇ • C = 0 in V β (3.3a) 0 = -∇c + ∇ 2 C - 1 V β A βσ n• (-Ic + ∇C) dA in V β (3.3b) C + ξλ β (I -nn) • n • ∇C + (∇C) T 1 = -I at A βσ (3.3c) c β = 0 (3.3d ) Periodicity C (r + l i ) = C (r) , c (r + l i ) = c (r) , i = 1, 2, 3 (3.3e)
In Eq. (3.3c), the superscript T 1 designates the transpose of a third order tensor that permutes the two first indices †. Note that Eq. (3.3d) follows from Eq. (2.14d) along with the pressure representation in Eq. (3.2) and is required for Eqs. (3.3) to form a well-posed problem.

At this point, a simple transformation, similar to that employed while upscaling the incompressible creeping flow with no slip leading to Darcy's law, can be used to obtain a pure differential form of this closure [START_REF] Barrère | On the closure problem for Darcy's law[END_REF][START_REF] Whitaker | The method of volume averaging[END_REF]. Letting

K -1 s = -ε -1 1 V β A βσ n• (-Ic + ∇C) dA (3.4)
defining D and d as

D = ε -1 (C + I) • K s (3.5) d = ε -1 c • K s (3.6)
and returning to Eqs. (3.3) yields the following local closure problem

∇ • D = 0 in V β (3.7a) 0 = -∇d + ∇ 2 D + I in V β (3.7b) D = -ξλ β (I -nn) • n • ∇D + (∇D) T 1 at A βσ (3.7c) d β = 0 (3.7d ) D = K s (3.7e) Periodicity D (r + l i ) = D (r) , d (r + l i ) = d (r) , i = 1, 2, 3 (3.7f )
From this tractable form having the structure of an incompressible Stokes problem with a slip boundary condition, K s can be computed from a simple average on D as indicated in Eq. (3.7e). This last equation is not required in the solution procedure but is used as a direct consequence of the definition of D in Eq. (3.5) together with the intrinsic average of the velocity decomposition in Eq. (3.1), which yields C β = 0. Although local, the closure (3.7) is non-intrinsic due to the presence of ξλ β in the boundary condition of Eq. (3.7c) and we shall come back to this in section 4.

Closed macroscopic model

When the representation of the deviations given in Eqs. (3.1) and (3.2) are reported in the unclosed from of the averaged momentum equation (2.11), one readily obtains a Darcy-like macroscopic law

0 = -∇ p β β -µ β K -1 s • v β + µ β ε -1 ∇ 2 v β (3.8)
which, along with Eqs. (2.10) and (2.12) forms the closed averaged model for the slightly compressible slip-flow considered in this work. In equation (3.8), K s , defined in Eq. (3.7e), is identified as the apparent slip-corrected permeability.

The last term in the rhs of Eq. (3.8) is referred to as the Brinkman correction. An order of magnitude analysis can be used to show that this macro-scale viscous diffusion term has a negligible contribution in the bulk of the porous medium. Indeed, as can be inferred from the boundary condition in Eq. (2.14c)

ṽβ = O v β β 1 + ξKn (3.9)
where v β β represents the leading order of v β β . This shows that

C = O 1 + ξKn -1 (3.10)
and consequently, from the definition of K s in Eq. (3.4), the order of magnitude of the Darcy term in Eq. (3.8) can be estimated as

µ β K -1 s • v β = O µ β ε -2 -1 β 1 + ξKn -1 a v v β β (3.11) In this estimate, a v = A βσ
V represents the interfacial area per unit volume of the medium that is expected to scale as -1 β . Similarly, the order of magnitude of the Brinkman term can be estimated to be

µ β ε -1 ∇ 2 v β = O µ β L -2 v β β (3.12)
which, when compared to Eq. (3.11), clearly indicates that the Brinkman term is completely negligible in the context of slip-flow, i.e. when ξKn 0.1. The macroscopic model can hence be written, to within an approximation O ( β /L), as

∂ ρ β β ∂t + ∇ • ρ β β v β β = 0 (3.13a) d -phase -phase
Figure 2. Unit cell for the square lattice of cylinders of circular cross section.

0 = -∇ p β β -µ β K -1 s • v β (3.13b) ρ β β = F p β β (3.13c)
It must be noted that the apparent slip-corrected permeability, K s , in the Darcy-like form of the macroscopic momentum equation (3.13b) is non-intrinsic as it lumps together topological properties, intrinsic to the microstructure of the medium, and slip effects that depend on the slip-length ξλ β . Moreover, a detailed analysis, provided in Appendix A, indicates that, unlike the intrinsic permeability that is the fundamental coefficient in the classical Darcy's law when no slip occurs, K s is not symmetric in the general case. Quasi-symmetry is observed under the constraint

ξKn O (ε) (3.14)
that is also detailed and illustrated in Appendix A, completing a previous analysis of this property [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF]. An important consequence is that the full characterization of the tensor K s from a measurement point of view is significantly simplified when the constraint in (3.14) is satisfied. In section 3.3 below, the validity of the macroscopic model and of the associated closure problem (3.7) yielding K s is verified by comparison with analytic solutions in a simple configuration and with direct numerical simulations on a more complex medium.

Validation

The finite element solver Comsol Multiphysics 4.4 was used to carry out the numerical solutions of both the pore-scale equations (2.1) and the closure problem (3.7). Meshing tests were carried out to reach convergence and led to use nearly 13,000 domain mesh elements with more than 500 boundary elements. This mesh was built by adding a boundary layer at the solid-fluid interface with a stretching factor of 1.2. As a first validation, a simple structure made of a periodic square pattern of parallel circular cylinders was considered, as represented in Fig. 2.

The slip-flow orthogonal to the cylinders axes is such that K s = k s I, where k s can be computed from the solution of the projection of the closure problem either on e x or e y . Due to the unit cell symmetry, the projection leads to a boundary value problem equivalent to the initial incompressible pore-scale problem when symmetry on the velocity and Numerical results on k * s = k s / 2 are reported in Figure 4a versus ξλ * β = ξλ β / , for ε = 0.8. As discussed elsewhere [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF], a Chang's unit cell, composed of the solid cylinder immersed in a circular fluid shell, may be used, along with a zero vorticity boundary condition at the outer edge, to derive an approximate analytic solution. For this relatively large value of the porosity, it is given by

k * s = 1 8π 1 + 2ξλ * β π φ -ln φ - 3 2 + 2φ - φ 2 2 +2ξλ * β π φ -ln φ - 1 2 + φ 2 2 (3.15) where φ = 1 -ε. 10 -4 10 -3 10 -2 10 -1 2 2.5 3 •10 -2 ξλ * β k * s Computed -ε = 0.8 Prediction 0 5 • 10 -2 0.1 0.15 0.2 0.25 1 1.2 1.4 1.6 ξKn k * s k * Computed -ε = 0.8 Prediction a) b)
1

k * s k * versus ξKn = ξλ * β * β
, * β being the dimensionless distance between parallel plates that would exhibit the same intrinsic permeability,

i.e. * β = 2 3k * ε .
The agreement of the prediction from Eq. (3.15) with our numerical results is excellent as shown in Figure 4. The relative error on k * s , taking the computed value as the reference, is less than 0.4% over the interval 10 -4 ξλ * β 0.13. Within this range of ξλ * β , ξKn remains smaller than ∼ 0.23 when * β = β is estimated from the slit aperture of the equivalent bundle of regularly spaced parallel plates given by * β = 2 3k * ε , k * = k 2 being the dimensionless intrinsic permeability, i.e. the permeability without slip. As clearly indicated in Figure 4b, representing

k * s k * versus ξKn, k *
s exhibits a non linear dependence on the Knudsen number even though the flow is expected to remain in the slip regime at the microscale. This important remark will be analyzed in more detail in section 4.2.

An additional validation is now investigated in the case of a more complex 2D structure featuring no special symmetry properties and for which no analytic solution is available. In that case, validation of both the upscaling procedure and macroscopic model can be achieved through the comparison of DNS of the pore-scale flow problem (Eqs. (2.1)) with the numerical solution of the closure problem (Eqs. (3.7)).

The structure at scale L for DNS is obtained by duplicating n times, in the direction of the applied macroscopic pressure gradient, a periodic unit cell of size made of randomly placed cylinders of circular cross sections with randomly chosen radii within a prescribed interval [r min , r max ] according to a log-normal distribution. Placement of cylinders centers is constrained to ensure periodicity on the x and y directions and is repeated until a target porosity is achieved, while controlling cylinders overlapping.

Two random unit cells, of porosity ε = 0.375 and ε = 0.804, were generated using respectively r min = 0.023 , r max = 0.062 and r min = 0.016 , r max = 0.031 . They were used to compute the xx dimensionless component, k * sxx = ksxx 2 , of the apparent slip-corrected permeability tensor from the closure problem solution. Additionally, DNS were performed on structures at scale L obtained from these generic unit cells with n, the number of unit cell replicas, ranging from 3 to 9. Constant pressures were prescribed at the entrance (x = 0) and exit (x = L) of the structure while periodic boundary conditions were applied in the y direction. Superficial and intrinsic averages of v β • e x and p β were performed on the central unit cell (n + 1)/2 to compute k * sxx . Computations were carried out for 10 -4 ξλ * β 10 -2 . The DNS procedure is illustrated in Figure 5 where we have represented an example of the two L-scale structures, along with the two generic unit cells on which maps of v β • e x were superimposed. In this figure, we have also reported the pressure profiles computed on the left vertical edge of the central unit cell. The evident non-uniformity of p β explains why DNS can not be carried out on a single unit cell with predefined entrance/exit constant Dirichlet boundary conditions on the pressure.

As shown in table 1 and in Figure 6a, results on k * sxx obtained from DNS and the closure problem solution are in excellent agreement since the relative error is less than 1.3% and decreases when n increases, confirming the validity of the macroscopic model and associated closure.

In figure 6b, a clear non-linear dependence of k * sxx on ξKn, mostly pronounced for the largest value of the porosity, is again highlighted for the two unit cells of Figure 5b, although ξKn remains smaller than about 0.23, a range typical of the slip-flow regime as justified by some measurements for even larger Knudsen numbers in straight channels [START_REF] Harley | Gas flow in micro-channels[END_REF]. This confirms the behavior reported above for the regular square pattern of cylinders.

At this point it is worth mentioning that, in the seminal work of [START_REF] Klinkenberg | The permeability of porous media to liquids and gases[END_REF], the linear dependence of the apparent permeability with the Knudsen number was regarded as an approximation resulting from a representation of the porous medium as a system of straight capillaries. As a matter of fact, [START_REF] Klinkenberg | The permeability of porous media to liquids and gases[END_REF] noticed some non-linearities in his experimental results with the same tendency as the one observed here. However, the porous media that he considered were relatively tight with intrinsic permeabilities that could go as low as 2.4 mDarcy. As a consequence, the porosity is expected to be quite small and this may explain why the observed non-linearity remained not too drastic in his experimental results. As shown in figure 6b), our predictions indicate that the non linearity decreases as porosity decreases and becomes almost insignificant for a porosity of 0.25, which is in a qualitative agreement with Klinkenberg's results.

n ξλ * β = 10 -4 ξλ * β = 10 -2 ε = 0.375 ε = 0.804 ε = 0.375 ε = 0.804 k * sxx × 10 -6 Error % k * sxx × 10 -4 Error % k * sxx × 10 -5 Error % k * sxx × 10 -
Before considering in more detail the above mentioned non-linearity, the impact of an incomplete version of the slip boundary condition, in which the ∇v T β term is omitted, shall be illustrated. When such an incomplete boundary condition is employed, the analytic form of the apparent slip-corrected permeability of the periodic square pattern of parallel circular cylinders, for a sufficiently large porosity (ε 0.55), is given by

k * s = 1 8π 1 + 2ξλ * β π φ -ln φ - 3 2 + 2φ - φ 2 2 +ξλ * β π φ -ln φ + 1 2 -2φ + 3φ 2 2 (3.
* = 1 8π -ln φ -3 2 + 2φ -φ 2 2
is again the dimensionless intrinsic permeability and the relative error between the two, taking the former as the reference, is given by -2 ln φ-3+4φ-φ 2 4(1-2φ+φ 2 )

. This relative % error, independent of ξKn, is represented versus the porosity in Figure 7a showing an increasingly significant error with ε, as it can reach about 62% for ε = 0.9.

When the ∇v T β term is omitted in the pore-scale slip boundary condition, the closure problem (Eqs. (3.7)) is unchanged, except the term ∇D T that is no longer present in Eq. (3.7c). This form of the closure was solved on the unit cell of the random pack of parallel cylinders (see Fig. 5b) for ε = 0.804, yielding k * sxx that turns out to be larger than that obtained with the complete boundary condition. The slip correction can again be estimated from k * sxx k * -1 computed with and without the complete shear rate in the boundary condition and the relative error can be formed taking the former as the reference. This % error, represented versus ξKn in Figure 7b, remains small at exceedingly small values of the Knudsen number (about 7% for ξKn 0.04), but strongly increases with ξKn as it reaches roughly 26% for ξKn 0.2. This clearly evidences that the complete strain rate must be kept in the expression of the slip at the solid boundary. More insight on the impact of the form of the boundary condition and on the non-linear behavior mentioned above will be provided in section 4, following a reformulation of the closure problem.

Reformulation and effective intrinsic coefficients

In the macroscopic model developed above, both viscous and slip effects are lumped together in the apparent slip-corrected permeability K s that is determined by a nonintrinsic closure problem. To isolate the Knudsen contribution to the flow at the macroscopic scale and to identify intrinsic macroscopic coefficients, the closure problem in Eqs.

(3.7) needs to be further developed.

Expansion in ξKn

The development is carried out on the dimensionless form of the closure given by the starred dimensionless quantities

d * = d/ β , D * = D/ 2 β , ∇ * = β ∇ as ∇ * • D * = 0 in V β (4.1a) 0 = -∇ * d * + ∇ * 2 D * + I in V β (4.1b) D * = -ξKn (I -nn) • n • ∇ * D * + (∇ * D * ) T 1 at A βσ (4.1c) d * β = 0 (4.1d ) D * = K s -2 β (4.1e) Periodicity D * (r * + l * i ) = D * (r * ) , d * (r * + l * i ) = d * (r * ) , i = 1, 2, 3 (4.1f )
Since ξKn remains smaller than unity in the context of slip-flow, d * and D * can be developed up to the m th order in Maclaurin series expansions that we shall write as

d * = d * 0 + m j=1 ξKn j e * j + R dm (4.2)
and

D * = D * 0 + m j=1 ξKn j E * j + R Dm (4.3)
where R dm and R Dm are the residuals at the m th order. When these expansions are introduced back into Eqs. (4.1), the closure problem can be split in a series of closures that ensue from the identification at the successive orders in ξKn. Returning to the dimensional form, the closure problem at the 0 th order takes the form 0 th order :

∇ • D 0 = 0 in V β (4.4a) 0 = -∇d 0 + ∇ 2 D 0 + I in V β (4.4b) D 0 = 0 at A βσ (4.4c) d 0 β = 0 (4.4d ) Periodicity D 0 (r + l i ) = D 0 (r) , d 0 (r + l i ) = d 0 (r) , i = 1, 2, 3 (4.4e)
This closure problem exactly corresponds to that defining the intrinsic permeability, K, while carrying out the upscaling of the incompressible Stokes problem without slip [START_REF] Whitaker | The method of volume averaging[END_REF][START_REF] Barrère | On the closure problem for Darcy's law[END_REF], and this means that

D 0 = K (4.5)
Taking this relationship into account, along with the expansion of Eq. ( 4.3) and the definition of K s in Eq. (3.7e), allows writing, at the m th order

K s K •   I + m j=1 ξλ β j S j   (4.6)
Here S j (j = 1, • • • , m) is the macroscopic coefficient obtained from the closure problem at the j th order given, in its dimensional form, by j th order (j = 1, • • • , m):

∇ • D j = 0 in V β (4.7a) 0 = -∇d j + ∇ 2 D j in V β (4.7b) D j = -(I -nn) • n • ∇D j-1 + (∇D j-1 )
T 1 at A βσ (4.7c)

d j β = 0 (4.7d ) D j = K • S j (4.7e) Periodicity D j (r + l i ) = D j (r) , d j (r + l i ) = d j (r) , i = 1, 2, 3 (4.7f )
In this j th order problem, d j and D j represent the closure variables defined from their analog, e j = β e * j and

E j = 2 β E * j (j = 1, • • • , m)
, after applying the rescaling

d j = e j / ( β ) j (4.8) D j = E j / ( β ) j (4.9)
As a consequence, the macroscopic model at the m th order in ξKn takes the form

ε ∂ ρ β β ∂t + ∇ • ρ β β v β = 0 (4.10a) v β - 1 µ β K •   I + m j=1 ξλ β j S j   • ∇ p β β (4.10b) ρ β β = F p β β (4.10c)
in which ξλ β j S j is the macroscopic j th order slip-flow correction. The tensors K and S j are respectively given by Eqs. (4.5) and (4.7e) once the closure problems of Eqs. (4.4) and (4.7) are solved. All these macroscopic coefficients are intrinsic as they derive from boundary value problems that only depend on the microstructure of the medium. At the first order, the macroscopic model is formally identical to that previously derived [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF], and generalizes the classical Darcy-Klinkenberg form [START_REF] Klinkenberg | The permeability of porous media to liquids and gases[END_REF]) when the gas is ideal. However, it must be noted that the form of the first order correction can be more complex than the dependence on 1/ p β β predicted by the classical Klinkenberg correction if a different state law is to be considered.

One must be clear that the expansion carried out in ξKn at the closure level is an alternative representation of the original closure problem given by Eqs. (3.7) that is basically local, and this is quite a different approach from the one followed by [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF] that requires an additional constraint on the Knudsen number. The approach, in this last reference, potentially allows to take into account non local effects that could play a role while studying gas slip momentum transport near macroscopic boundaries.

Under the form of Eq. (4.10b), the macroscopic momentum equation clearly indicates that the apparent slip-corrected permeability, up to the first order, remains linear in ξKn. As a consequence, the non linear behavior observed in Figures 4b and6b can only be captured by higher order terms in the macroscopic model, and this is a major result of the present development. In this perspective, a more thorough analysis of the macroscopic slip-flow correction, beyond the first order, is proposed in section 4.2.

Effective coefficients

All the closure problems providing the tensors K and S j have an incompressible Stokes structure, the latter being coupled to each other and to the 0 th order problem through the slip-like boundary condition, which makes them non-homogeneous. Therefore, a unique numerical procedure can be used to determine all the macroscopic coefficients in a simple manner.

The same computational tool as the one employed for the determination of the apparent slip-corrected permeability, K s , (see sections 3.1 and 3.3), was used to solve the closure problems (4.4) and (4.7) up to m = 3 on the model structure represented by the unit cell of Figure 2 and ε ranging from 0.25 to 0.8. Due to the symmetries in this particular case, all the macroscopic coefficients are spherical tensors (K = kI, S j = s j I), thus requiring the solution of the projection of the corresponding closure problems on e x (or e y ) only. Fields of the xx component of the tensors D 0 , D 1 , D 2 and of the x components of the vectors d 0 , d 1 , d 2 are reported in their dimensionless forms in Figure 8.

In Figure 9a, we have represented k * = k 2 versus ε, which shows an excellent agreement with our previous results [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF] as well as with predictions available in the literature in the range of small [START_REF] Bruschke | Flow of generalized Newtonian fluids across a periodic array of cylinders[END_REF] and large [START_REF] Kuwabara | The forces experienced by randomly distributed parallel circular cylinders or spheres in a viscous flow at small Reynolds numbers[END_REF] porosities.

Dimensionless slip-flow corrective terms, k * s * j = k * ( ) j s j are represented versus ε in Figure 9b for j = 1 to 3. In this figure, we have also reported the O ξKn slip-flow corrective term k * s * appearing in the classical macroscopic model as the one reported earlier [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF]. A perfect agreement with the first order term of the present model can be observed, showing the equivalence of the two approaches restricted to the first order in ξKn. Higher order terms in the present model are such that s j is positive when j is odd and negative otherwise.

As can be inferred from the non-linear dependence of k s on the Knudsen number reported above, the first order macroscopic slip-correction might become relatively inaccurate in some circumstances. This is highlighted in Figure 10 where we have reported the apparent slip-corrected permeability, k * s , versus ξKn for ε = 0.25, 0.4, 0.6 and 0.8., keeping ξKn smaller than ∼ 0.19 so that slip-flow is expected to remain physically relevant. The apparent slip-corrected permeability was computed from the closure problem (3.7), on the one hand, and from Eq. (4.6) (i.e. from the solutions of problems (4.4) and (4.7)) for m up to 3, on the other hand. Figure 10 shows that the inaccuracy of the classical first order approximation increases with ξKn and becomes much more significant while increasing ε.

For a more quantitative picture of the above, the relative error on the slip correction in k s at the j th order, given by

ks k -1- m j=1 (ξλβ) j sj ks k -1 = k * s -k * 1+ m j=1 (ξλ * β ) j s * j k * s -k *
, is reported, up to the third order, versus ξKn in Figure 11, for the four values of ε. When k s is estimated at the first order, the relative error on the slip correction is ∼ 3.7% for ε = 0.25 and ξKn ∼ 0.19 but reaches ∼ 20.3% for ε = 0.6 and ξKn ∼ 0.1. The estimation of the slip correction is significantly improved while taking into account the second order term, or further the third order term. In fact, for m = 3 and the structure under consideration, the relative error on the slip correction falls down to approximately 0.6% for ε = 0.25 and ξKn ∼ 0.19 or 1.7% for ε = 0.6 and ξKn ∼ 0.1. This is clearly illustrated in Figure 11. If a different Knudsen number, say Kn σ , based on the characteristic size, σ , of the solid-phase is used instead of Kn, one would observe that the non-linearity appears for ξKn σ = β σ ξKn 0.02, a value smaller than, or, when ε is large, almost equal to ξKn, as β σ takes the values 1.98 × 10 -2 , 1.49 × 10 -1 , 4.21 × 10 -1 and 1.07 for ε = 0.25, 0.4, 0.6 and 0.8 respectively.

One must however be clear about the physical meaning of the macroscopic slipcorrection terms for j 2 obtained from the development leading to Eq. (4.10b). In fact, it must be kept in mind that the pore-scale physical model from which the macroscale balance equations derive, involves a slip-boundary condition and a momentum equation that are both first-order accurate in ξKn. As a consequence, the macroscopic slip correction terms beyond the first order in ξKn account for the contribution of microstructural effects, i.e. the signature of pore topology on the averaged slip effect that can significantly depend upon local constrictions, enlargements and curvature of pore walls through the strain rate at A βσ . One can note, in particular, that the macroscopic model and associ-ated closures provide an exact macroscopic solution in the case of slip-flow in a bundle of straight capillary tubes or slits made of plane parallel plates. Indeed, this solution exactly corresponds to that obtained by averaging the pore-scale flow solution. In this situation with no curvature at the pore-scale in the direction of the flow, which is 1D, the macroscopic slip-correction remains at the first order that is exactly obtained from the closure problem (4.7) for j = 1, all the higher-order closure problems yielding D j = 0 (j 2). To be more precise, specific general properties of the slip boundary condition and the consequences on the closure problems are highlighted in the following section.

Before moving on, it is important to mention that some results obtained from a numerical solution of an approximate version of the Boltzmann equations (i.e. Bhatnagar Gross and Krook (BGK) model) †, under a linearized approach, predict a linear dependence of the apparent permeability with the Knudsen number over the whole range of ξKn investigated here for ε = 0.8. The contrast with the present results would suggest that the range of application of the slip flow model should be drastically reduced as the porosity increases. This result is however difficult to justify physically. Moreover, the linear dependence of k * s with the Knudsen number resulting from this linearized BGK approach seems to hold for Knudsen number values up to about 0.5 fitting the predictions of Chang's unit cell solution that is however expanded at the first order in the limit of ξKn 1; such a surprising behavior remains unclear. Some possible sources for explaining this linear behavior are the use of a linearized BGK model, along with the numerical method used to solve it.

The above remarks should not be considered as a definitive argument to accept a nonlinear dependence of the apparent permeability with the Knudsen number. In this regard, the finding of this behavior should motivate the performance of slightly-compressible rarefied gas flow experiments in highly permeable porous media.

Some important features of the slip boundary condition

An incomplete slip-boundary condition, in which the ∇v T β term has been dropped, has often been used for porous media flow [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF][START_REF] Pavan | Gas/solid heat transfer in gas flows under Klinkenberg conditions: comparison between the homogenization and the kinetic approaches[END_REF], and we have shown, in accordance with some previous references [START_REF] Panzer | The effect of boundary curvature on hydrodynamic fluid flow: calculation of slip length[END_REF][START_REF] Lockerby | Velocity boundary condition at solid walls in rarefied gas calculations[END_REF][START_REF] Barber | Isothermal slip flow over curved surfaces[END_REF], that this physically unjustified form can lead to a significant error in some circumstances (see section 3.3). Together with the non-linear dependence of the slip correction on ξKn, this suggests a close attention to the slip-like boundary condition of Eq. (4.7c) in the j th -order closure problem.

The analysis starts with an alternative expression of the second term in the rhs of Eq. (4.7c) given by

(I -nn) • n • (∇D j-1 ) T 1 = (I -nn) • ∇ (n • D j-1 ) -(I -nn) • (∇n • D j-1 ) at A βσ (4.11)
For all j, n • D j-1 = 0 at A βσ since D j-1 is either 0 (j = 1) or purely tangential to A βσ at this interface (j 2). Consequently, ∇ (n • D j-1 ) is purely normal to A βσ and this allows rewriting Eq. (4.11) as

(I -nn) • n • (∇D j-1 ) T 1 = -∇ s n • D j-1 at A βσ (4.12)
where ∇ s denotes the surface differentiation operator defined by ∇ s ≡ (Inn) • ∇.

For j = 1, this result can be further simplified due to the no-slip boundary condition on D 0 (see Eq. (4.4c)) leading to

(I -nn) • n • (∇D 0 ) T 1 = 0 at A βσ (4.13)
and this proves that the ∇v T β term in the slip boundary condition does not play any role up to the first order of the expansion in ξKn at the closure level. At higher orders, Eq. (4.12) can not be further simplified and indicates that the contribution of (Inn)

• n • (∇D j-1 )
T 1 depends on ∇ s n that is related to the curvature of A βσ . In 2D, for instance, it can be easily demonstrated (see Appendix B) that ∇ s n = -κ (Inn) where κ = -∇ s • n is exactly the curvature of A βσ . This clearly shows that the (∇D j-1 )

T 1 term does play a role in the slip-like boundary condition (4.7c) of the closure problem at the j th order, j 2, in the general case. Its contribution vanishes only in some special geometrical situations where the flow is 1D and the curvature of A βσ is zero in the flow direction. This last remark applies to the boundary condition of Eq. (2.1d) at the porescale and justifies why the incomplete form of the slip boundary condition can be used in the case of straight tubes, slits or channels that have been extensively used for the study of Knudsen effects in the slip regime [START_REF] Lauga | A note on the stability of slip channel flows[END_REF][START_REF] Fishman | Viscosity and slip velocity in gas flow in microchannels[END_REF][START_REF] Shen | A kinetic theory based first-order slip boundary condition for gas flow[END_REF]).

The analysis is pursued with an alternative expression of the first term in the rhs of Eq. (4.7c), which can be written as

(I -nn) • (n • ∇D j-1 ) = ∇ s • (nD j-1 ) + κD j-1 at A βσ (4.14)
This finally allows to express the slip-like boundary condition in Eq. (4.7c) for the j thorder closure problem as

D j = -∇ s • (nD j-1 ) + (∇ s n -κI) • D j-1 at A βσ (4.15)
In 2D, this last expression takes the form D j = -∇ s • (nD j-1 ) -2κD j-1 . For j = 1, it simplifies to

D 1 = -(I -nn) • (n • ∇D 0 ) = -∇ s • (nD 0 ) at A βσ (4.16)
Equations (4.15) and (4.16) show that the slip-like boundary condition in the j thorder closure problem contains an explicit dependence upon the curvature of A βσ that is however filtered out at the first order in ξKn. This explains why the classical firstorder macroscopic slip-flow model remains inaccurate for systems presenting a significant curvature of the solid-fluid interface. For the periodic square pattern of circular cylinders, it clearly justifies why the error on the slip-correction, when restricted to the first-order, increases with the porosity, as curvature scales as (1ε)

-1/2 in that case. This also justifies the necessity of taking into account higher order terms in the macroscopic model as in Eq. (4.10b).

Finally, it must be noted that the analysis carried out in this section may be straightforwardly applied to an expansion in the Knudsen number performed on the original initial boundary value problem given in Eqs. (2.1). In that way, the slip boundary condition properties at the successive orders in Kn may be generalized far beyond the strict context of porous media flow, i.e. for any situation allowing such an expansion in Kn.

Conclusion

Macroscopic modeling of slightly compressible gas slip-flow in homogeneous porous media was thoroughly revisited and several important results were highlighted.

Under a specified set of constraints, the macroscopic model was shown to have a Darcylike form in which the apparent non-intrinsic slip-corrected permeability tensor is given by a non-intrinsic closure problem having the structure of an incompressible Stokes flow problem with a slip-like boundary condition. The validity of the macroscopic model and associated closure was verified through comparisons with analytic solutions and DNS on model porous structures. The apparent slip-corrected permeability was shown to be a non-symmetric tensor in the general case, whereas a sufficient condition for quasisymmetry to occur was derived. A non-linear dependence of the apparent slip-corrected permeability upon the averaged Knudsen number was evidenced that can become increasingly significant while increasing the porosity and the Knudsen number on the model structure under concern.

A reformulation of the macroscopic model and closure was derived using an expansion of the closure variables in the Knudsen number. This reformulation identifies the viscous and slip contributions that are, otherwise, lumped together in the apparent slip-corrected permeability. To the 0 th order, the macroscopic model corresponds to the classical flow problem without slip, the associated macroscopic coefficient being the intrinsic permeability tensor. The successive higher order terms in the macroscopic model are characterized by intrinsic slip-flow correction tensors that are all determined by easy-to-solve intrinsic closure problems that have basically the same Stokes structure and are coupled to each others and to the 0 th order closure problem through their slip-like boundary condition. The first-order correction was found to be in excellent agreement with a previous existing and validated model [START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF]). Higher-order corrections account for the nonlinear behavior of the apparent slip-corrected permeability versus the Knudsen number observed previously. Correction beyond the first order must be understood as a signature of pure geometrical effects of the microstructure, through the strain rate at the solid-fluid interface involved in the microscale slip boundary condition which is sensitive to local interface curvatures. Within a simple model structure, at large values of the porosity, errors up to 45% on the slip-correction were found when restricted to an estimation at the first order.

The closure problem providing the first order slip-correction was shown to be insensitive to an incomplete slip-flow boundary condition where the ∇v T β term is omitted and that has been often abusively used in the statement of the pore-scale physical problem as reported in many references [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF][START_REF] Pavan | Gas/solid heat transfer in gas flows under Klinkenberg conditions: comparison between the homogenization and the kinetic approaches[END_REF]. This result does not hold for higher-order slip-correction terms. Moreover, it was shown that the explicit dependence of the slip-like boundary condition on the curvature of the solidfluid interfaces is filtered out at the closure level when the model is restricted to a first-order approximation, explaining the inaccuracy of the existing model reported in the literature and justifying the importance of higher-order slip-correction terms. This is an important general feature of the slip boundary condition, that applies beyond the context of porous media flow, providing a new physical insight into slip-flow.

Supplementary work is needed for a further comparison of the present investigation with experimental data in particular for highly porous structures.

∇ • ∇ D T • D = 1 V A βσ n • ∇ D T • D dA (A 16)
We can now make use of the general identity

u • ∇A • B = (∇A) T 1 : (uB) T 1 (A 17)
valid for any two second order tensors A and B and any vector u, and when A is identified to D T , B to D and u to n, this provides an alternative form of the interfacial integral term in the rhs of Eq. (A 16). When replaced in Eq. (A 13), the apparent slip-corrected permeability, K s , can finally be recast into the following expression

K s = (∇D) T 3 : ∇D - 1 V A βσ (∇D) T 3 : (nD) dA (A 18)
While the tensor (∇D) T 3 : ∇D in this expression can be easily shown to be symmetric, the tensor represented in the area average term is not symmetric in the general case, except if the interface A βσ has specific symmetry properties. This is the case, for example, of an ordered porous structure for which the solid σ-phase exhibits symmetries about the three planes parallel to the edges of the periodic unit cell and passing through its centroid. In such circumstances, the off-diagonal terms of K s are all zero. The expression of K s in Eq. (A 18) shows that this tensor is generally not symmetric [START_REF] Skjetne | Homogenization of wall-slip gas flow through porous media[END_REF].

The investigation can be further pursued by examining the condition under which K s is close to a symmetric tensor. A sufficient condition is when the non-symmetric part in Eq. (A 18) remains small compared to the symmetric part and this can be expressed in terms of the orders of magnitude of these two respective terms as The order of magnitude of the area average can be estimated to be These two estimates allow us to formulate the constraint (A 20) as

O    1 V A βσ
ξλ β O ε a v (A 23)
If a v is thought to vary as -1 β as is often accepted [START_REF] Whitaker | The method of volume averaging[END_REF], one would be left with

ξKn O (ε) (A 24)
as a sufficient condition for K s to remain quasi symmetric. It should be noted that this is consistent with the fact that, when no slip occurs, corresponding to Kn 0, K s = K, which is a symmetric tensor as predicted by Eq. (A 18) when ξλ β = 0, i.e. when D = 0 at A βσ .

A short illustration of the symmetry properties can be provided from direct computation of K s over a unit cell that does not possess any particular geometrical symmetry as the one represented in Figure 12 for which ε 0.515 and a * v = a v 11.935. The apparent slip-corrected permeability was determined from the solution of the closure problem of Eqs. (A 1) above and for a cell Knudsen number, ξλ * β = ξλ β / , ranging from 10 -4 to 1.

Results on the components of K * s = K s / 2 , reported in table 2, clearly show that the contrast between the off-diagonal terms of K * s remains small when the constraint expressed in (A 23) 1.00 10 -4 1.26 10 -4 5.87 10 -6 5.80 10 -6 9.28 10 -5 1.15 0.01 2.13 10 -4 6.87 10 -6 9.09 10 -6 1.59 10 -4 27.80 0.10 4.95 10 -4 6.91 10 -6 1.57 10 -5 3.40 10 -4 78.00 1.00 7.24 10 -4 5.82 10 -6 1.87 10 -5 4.62 10 -4 105.00

Table 2. Components of the dimensionless apparent slip-corrected permeability tensor for the unit cell of Figure 12 Coming back to Eq. (B 4), it follows that

∇ s • n = ∇ • n = -κ (B 8)
which leads to the expected relationship expressed in Eq. (B 1), valid in 2D, when the result of Eq. (B 3) is taken into account.

Figure 1 :

 1 Figure 1: Averaging domain and characteristic lengths of the system.

Figure 1 .

 1 Figure 1. Averaging domain and characteristic lengths of the system.

Figure 3 .

 3 Figure 3. Dimensionless closure variable fields computed on the unit cell of Figure 2. ε = 0.8. a) Dxx/ 2 , ξλ * β = 10 -4 . b) Dxx/ 2 , ξλ * β = 0.1. c) dx/ , ξλ * β = 10 -4 . d) dx/ , ξλ * β = 0.1.

Figure 4 .

 4 Figure 4. Dimensionless apparent slip-corrected permeability k * s = ks 2 for the unit cell of Figure 2, ε = 0.8, k * = k 2

Figure 5 .

 5 Figure 5. Slip-flow DNS on 2D packs of cylinders. a) Large scale structures obtained from periodic replicas along ex of random unit cells. Velocity maps in the x direction are superimposed. b) Random unit cells of respective porosity ε = 0.375 and ε = 0.804 used to generate large scale structures along with x-velocity maps. c) Pressure profiles at the left edge of the central unit cells.

Figure 6

 6 Figure 6. a) Comparison of computed values of k * sxx obtained from DNS and volume averaging for the two unit cells of Figure 5b). b) Dependence of k * sxx k * xx on ξKn = ξλ * β β = ξλ * β * β for ξKn 0.2. ξKn is estimated from the dimensionless distance between plane parallel plates, * β = 2 3k * ε , that would lead to the same intrinsic permeability k * = k * xx . k * xx 2.44 × 10 -6 (ε = 0.375); k * xx 1.60 × 10 -4 (ε = 0.804).

Figure 7 .

 7 Figure 7. Relative error on the slip-correction while using an incomplete slip boundary condition (i.e. omitting the ∇v T β term). a) Unit cell of Figure 2. The error is computed from the analytic expressions (3.15) and (3.16). b) Random pack of parallel cylinders (see Fig.5b)), ε = 0.804. The error is computed from the closure problem solutions using complete and incomplete boundary conditions at A βσ .

  16) instead of Eq. (3.15) above, with again φ = 1ε. This last expression underestimates k * s obtained from Eq. (3.15). The two slip flow corrections can be extracted from Eqs. (3

Figure 8 .Figure 9 .

 89 Figure 8. Dimensionless closure variable fields up to the second order computed on the unit cell of Figure 2. ε = 0.8. a) D0xx/ 2 . b) D1xx/ 2 . c) D2xx/ 2 . d) d0x/ . e) d1x/ . f ) d2x/ .

Figure 10 .

 10 Figure 10. Comparison of the apparent slip-corrected permeability, k * s = ks/ 2 , computed from the solution of the closure problem (3.7), on the one hand, and estimated at the 1 st , 2 nd and 3 rd orders with Eq. (4.6) and the solutions of expanded closure problems (4.4) and (4.7), on the other hand. Periodic model structure of Figure 2. a) ε = 0.25. b) ε = 0.4. c) ε = 0.6. d) ε = 0.8. Kn is the macro-scale Knudsen number based on the fluid-phase characteristic length-scale: ξKn = ξ λ β β .

Figure 11 .

 11 Figure 11. Relative error on the slip flow correction in ks estimated at the 1 st , 2 nd and 3 rd orders taking the slip correction extracted from ks computed from the solution of the closure problem (3.7) as the reference (see complete expression in the text). Periodic model structure of Figure 2. a) ε = 0.25. b) ε = 0.4. c) ε = 0.6. d) ε = 0.8.

  of the boundary condition in Eq. (A 1c), the above constraint can be expressed as

  Figure 12. Unit cell for the computation of the full apparent slip-corrected permeability tensor. ε 0.515, a * v 11.935.

  and four values of ξλ * β . The contrast between k * syx and k * sxy is estimated in the last column.0, which leads to the straightforward relationship∇n • n = 0 (B 2)When this result is employed in the expression of the product (Inn) • ∇n, along with the fact that n is a unit vector, one ends up with∇ s n = (Inn) ∇ • n (B 3)in the specific 2D case.Let us now consider the expression of the curvature in the general case (2D or 3D) and write ∇ s • n from its definition as∇ s • n = ((Inn) • ∇) • n = ∇ • n -(nn • ∇) • n (B 4)The last term in the rhs of Eq. (B 4) can be explicitly re-written as((nn • ∇) • n) ij = n i n j ∂n i ∂x j = n j (∇n) ji n i (B 5)where we have implicitly used Einstein notation, and this means that(nn • ∇) • n = n • ∇n • n (B 6)However, due to Eq. (B 2), the last result simplifies to (nn • ∇) • n = 0 (B 7)

Table 1 .

 1 Comparison of the predictions of k * sxx from DNS with those from volume averaging. n is the number of unit cell replicas. k * sxx values are those obtained from DNS taking the average at the (n + 1)/2 unit cell. The error percent is computed relative to the DNS predictions.

	4 Error %

† Email address for correspondence: didier.lasseux@u-bordeaux.fr

† In the reference[START_REF] Lasseux | A macroscopic model for slightly compressible gas slip-flow in homogeneous porous media[END_REF]) (see Eq. (83c)), we simply used the superscript T for this transpose and this was ambiguous.

† These results were provided by an anonymous reviewer and are available from the authors.

Acknowledgements

FVP expresses his gratitude to Fondo Sectorial de Investigación para la educación from CONACyT (Project number: 12511908; Arrangement number: 112087) for the financial aid provided. Authors are thankful to L. Mieussens, B. Dubroca, S. Brull and P. Charrier from University of Bordeaux for fruitful discussions about the BGK method. Special thanks are due to an anonymous reviewer who brought to our attention the results obtained from this numerical approach.

Appendix A.

In this appendix, the symmetry properties of the apparent slip-corrected permeability tensor, K s , in the Darcy-like form of the macroscopic momentum equation (3.13b), are analyzed. This apparent permeability is explicitly given by the closure problem (3.7)

The analysis starts with a pre-multiplication of Eq. (A 1b) by D T and when the superficial average of the result is formed, one obtains

The first term in the rhs of the above equation can be equivalently expressed as

and, since D is a divergence-free tensor as stated in Eq. (A 1a), this is equivalent to

Making use of the averaging theorem (see Eq. (2.8)) in its divergence form, one obtains

Because D and d are periodic

and, since D has no normal component at A βσ as indicated by the boundary condition in Eq. (A 1c),

From Eq. (A 5), it follows that

Returning to equation (A 2), we hence have

The first term in the rhs of this last relationship can be rewritten as

in which we have used the trivial property ∇ 2 D T = ∇ 2 D T . Equation (A 10) can be arranged according to the fact that, for any two second order tensors A and B,

T 1 and, when A is identified to D T and B to D, this yields

In the above expressions, we have used the superscript T1 to designate the transpose of a third order tensor that permutes the two first indices and we adopted the classical nested convention for double dot product. The last equation can be equivalently written as

where the superscript T 3 denotes the transpose of a third order tensor that permutes the first and third indices. When this result is inserted back into Eq. (A 9), we have

The attention is now focused on the first term in the rhs of the last relationship, which, upon use of the averaging theorem, can be equivalently written as

Since D is periodic, ∇ D T is also periodic and hence

The objective of this appendix is to demonstrate that, in the 2D case and with the notations used throughout the paper,

where κ = -∇ s • n is the curvature of A βσ and ∇ s the surface differentiation operator defined by ∇ s ≡ (Inn) • ∇.

We shall start by noticing that, since n is a unit vector, n • n = 1 and hence ∇ (n • n) =


Keywords: Cellular technology (LTE, 5G, etc.) Bluetooth Radar Wi-Fi GPS Wireless Battery Charging

Autonomous (Self-driving) Vehicles (2/3) Autonomous (Self-driving) Vehicles (3/3) DSRC (760 MHz, 5.9 GHz)

Impact of vehicles density on the achieved channel busy ratio and collision rate:

P&A-A vs. SuRPA and ETSI schemes (Kirchberg map)

Ø Traffic congestion and accidents information sharing Ø Support for faster travel and emergency service delivery • Explosion of content retrieval and data dissemination in the Internet

• IP video will be 82% of all IP global traffic by 2022

[Cis17] Cisco, "Cisco Visual Networking Index: Forecast and methodology, 2017-2022", White paper, 2017.

• This growth impacts the infrastructure of ISPs and CPs/SPs Typical CDN architecture: Hierarchical caches structure

Modern ecosystem (2/3)

• Internet usage has evolved to become mainly used for content consumption

• There's a need for an upgraded infrastructure better suited to today's requirements

• Caching should go deeper into the network, to be closer to the end-users.

• Issues with multiple operators (wholesale CDNs)

• Alternative approach: Content-Centric Networking (CCN) [Jac09]

• Interesting idea but it still doesn't take ... The question of caches' placement becomes central • Limited size of the video library

• Popularity Zipf (α = 1.2) [Fricker12]
• Rapid decrease of the popularity Relevance of caching for UGC

• ANR ViPeer [Li12]
• Partners: INRIA, TB, Orange, …

• Objective: Relevance of caching Is caching really useful?

• Clearly beneficial for all parties:

• Customers: we improve their quality of experience • CPs/SPs: their customers will be satisfied and their origin servers less solicited • ISPs: less congested • However,

• Mandatory cooperation between CPs/SPs and ISPs

• Already done (Open connect appliance, Netflix) How a cache works?

Focus on a single cache

Important performance metric will be the hit probability

• Goal: find a policy to maximize the hit probability.

Hit probability = number of items served from cache total number of items served

Which item to replace?

• Classical cache replacement policies:

• RAND • FIFO • LRU • LFU • LRU-K Content store

Application

Cache R e q u e s t ( A )

Miss R e q u e s t ( A )

Replace one item

Even when the popularity is constant, LFU is not optimal.

• LFU is optimal under IRM (it maximizes the steady-state hit probability).

• LFU statically places in the cache the C most popular objects (on the basis of the average request rate of contents)

• Higher complexity than LRU

• LFU is no longer optimal under a general distribution

• The content of the cache is never adapted to instantaneous traffic conditions, resulting into suboptimal performance.

The online problem is not easy …

LRU is motivated by the heuristic that the element that has not been used for the longest time is likely to also not be needed for the longest time.

heuristic works well even in a worst-case analysis

Content store Application Cache R e q u e s t ( A ) Miss R e q u e s t ( A )

Replace one item

Which probabilities to use for the contents' request?

• Zipf distribution has been widely tested against the distributions collected from the real traces • The probability to request the content of rank 𝑟 is

𝑝 # = 𝑟 $% ∑ &'( ) 𝑖 $%
• where 𝛼 is the skew of the distribution; It depends on the type of the accessible objects

What is the impact of caching within a wireless and mobile edge (e.g. CAVs)?

• When considering only the control plane:

• Mobility pattern could be considered as it impacts the traffic pattern • Future requests' prediction (prefetching) could improve the performance (cooperation between CAVs, and CAVs-Infrastructure) [Ghada16]

• When considering the data plane:

• Link state variation can cause the user to choose a better or worse quality, which directly impacts the caching strategies (and the control plane) Cache placement in 5G and B5G

Network of caches

Cache placement: challenges

• Assume we have a fixed amount of cache to distribute over a multi-cache network

• Placing most of the cache at the core:

• Allows to reduce the solicitation of origin servers (use less peering links) • Does not solve the problem of congestion within the ISP network

• Distributing the cache at the network edge:

• Minimizes the delay for the most popular content • Less effective caching due to the content redundancy (higher demands on the origin servers)

• How to allocate efficiently the resources of cache? (320,1696,2080) (1952,1040,1104) (224,208,592) (1264,1184,2672) Perf. evaluation (%) Cache allocation 1216,832,2048) (1728,960,1408) (320,384,320) (64,3200,1856) Perf. evaluation (%)

Cache allocation

• The choice of a specific cache allocation depends on specific desired objectives

ØØè

  Localization: where am I? Ø Sensing the surroundings: what's happening around me? Ø Perception (fusion of sensor data) Ø Reasoning and decision making Ø Motion control Autonomous (Self-driving) Vehicles (1/3) Ø Promise to solve many of today's road traffic problems by • Improving road users' safety • Shortening their trip times • Enhancing their travel experience • Reducing air pollution Vehicular Networks Ø Wireless communication between • Vehicles (V2V) • Vehicles and the road infrastructure (V2I) • Vehicles and X? Ø Safety information sharing and dissemination 1 Ø Cooperative manoeuvres• Safe and efficient lane change 2 Hidden obstacles notification • … 1 -Safer driving 1 S. Zemouri, S. Djahel, and J. Murphy. A Fast, Reliable and Lightweight Distributed Dissemination Protocol for Safety Messages in Urban Vehicular Networks. Ad Hoc Networks Elsevier, Vol. 27, pp. 26-43, 2015. 2 A. Lissac, S.Djahel and J. Hodgkiss. Infrastructure Assisted Automation of Lane Change Manoeuvre for Connected and Autonomous Vehicles. IEEE ISC2, Oct. 14-17, 2019.Safer driving -Related challenges Ø How to ensure fast & reliable safety information dissemination? Ø How to avoid periodic beacons congestion? • At microscopic level • At macroscopic level Beacons congestion control Ø Safety and critical information (periodic broadcast) • Inter-vehicular communication (infrastructure-free) • Urban scenarios Ø Two types of messages • Event driven warning messages • Control messages (beacons, , CAMs, etc.) 19 Wireless congestion problem Ph.D. thesis "Enhanced Wireless Congestion Control in Urban Vehicular Environments", Sofiane Zemouri, University College Dublin, Ireland , July 2016 Safety messages dissemination Ø Critical information, multi hop transmission, high periodicity, limited timeframe, short to medium length messages Ø Goal, limit the amount of data inserted in the network • Get rid of unnecessary retransmissions after failures • Channel condition (shadowing, Doppler Effect,…), packet collisions, integrity of the message, sender receiver distance 20 è Who's the most suitable next forwarder? Safety messages dissemination Road-Casting Protocol (RCP) Ø Distance factor (D) • A function of the distance between the sender & the potential next forwarder, the next junction and the potential next forwarder, and the next junction to the transmission range of the sender Ø Link Quality factor (LQ) • A function of the signal quality (SQ), the channel quality (CQ) and the collision probability (CP) 21 S. Zemouri, S. Djahel, and J. Murphy. A Fast, Reliable and Lightweight Distributed Dissemination Protocol for Safety Messages in Urban Vehicular Networks. Ad Hoc Networks Elsevier, Vol. 27, pp. 26-43, 2015 of vehicles density on Packet Delivery Ratio and End-to-End Delay: RCP vs. distance based forwarding, link based forwarding and flooding Safety messages dissemination Ø Microscopic congestion control ✓ • Amount of data inserted in the network on a point-to-point level Ø Macroscopic congestion control • The extent of data inserted in the network (end-to-end) Ø Determine the Region of Interest (RoI) for multi-hop safety messages Ø Keep the amount of data inserted in the network minimal The set of critical junctions (red ticks) delimitating the zone where vehicles are already blocked due to an incident Road map connectivity check performed by each receiver of an emergency message: using simplified A* algorithm ratio, under different vehicle densities, of vehicles located inside and outside the RoI which did not receive the emergency message (Algiers map) of vehicles density on the accuracy of the RoI definition: RCP scheme RoI vs.

FasterØØ

  How to design efficient re-routing strategies to efficiently deal with unexpected events on the road Ø i.e., re-routing protocols that minimize the extra delay caused by the unexpected event Ø Answer: Next Road Re-routing (NRR)Ø Truck platooning requires new algorithms to determine:• Optimal frequency of platooning and length of a platoon • Smooth and safe integration of platoons in the surrounding traffic • Optimal speed • Gaps • Distance to surrounding traffic Ø Evaluating its impact on the infrastructure CACC -Related challenges 4 -Traffic congestion reduction • CR-TMS Impact of the percentage of buses on the achieved ATT 7*7 grid map with 25% of bus lanes 5 How to ensure efficient resources provision to ITS applications? • Rapid and dynamic change of resources availability due to vehicles mobility • Modeling and accurately predicting users' request patterns • Handover between clusters of CAVs service providers Ø Solution: AI enabled techniques such as Federated learning Ø Bandwidth scarcity of the DSRC spectrum • TV White Space? LTE-A/5G/6G? Ø Congestion and awareness control • Microscopic & Macroscopic -Machine learning Ø Electric vehicles (EVs) • Route planning & re-routing, charging planning Ø Security, privacy and liability • Cyber-attacks, location privacy, crash, etc. Emerging challenges for CAVs • Flexible and low cost solution with modest infrastructure needs Ø Ground to Air and Air to Air communication Ø Key challenges for cooperative UAVs • Optimal routing and recharging scheduling • Lightweight dynamic re-planning • Optimal tasks redistribution strategies Ø CAVs assisted cooperative UAVs? Ø UAVs assisted CAVs? Emerging challenges for CAVs on the Sky Ø CAVs is an emerging technology facing many challenges & offering many opportunities Ø Successfully overcoming these challenges promises a significant societal and economic impact Ø Most experts foresee that CAVs will be a key enabler of Smart Mobility, thus addressing its major technical challenges is essential Conclusion Ø How CAVs could support Smart Mobility • What is CAVs?

[•

  Jac09] V. Jacobson et al., "Networking named content", CoNEXT, 2009. [Ravindran17] R. Ravindran, et al., "5G-ICN: Delivering ICN Services over 5G Using Network Slicing," in IEEE Com. Mag., vol. 55, no. 5, pp. 101-107, May 2017 Virtualization of the network services and the infrastructures • Enabled by NFV and SDN • Possible expansion of virtualization to the field of content distribution • 5G/B5G network slicing • will offer new opportunities for ISPs to dynamically manage their shared storage assets • Mobile edge could be included • ICN-as-a-slice within 5G [Rahman18] [Rahman18] A. Rahman, et al., "Deployment Considerations for Information-Centric Networking (ICN)", IETF Draft draft-irtf-icnrg-deployment-guidelines-01 (May 2018) [Ravindran17] R. Ravindran, et al., "5G-ICN: Delivering ICN Services over 5G Using Network Slicing," in IEEE Com. Mag., vol. 55, no. 5, pp. 101-107, May 2017 5G ICN Architecture [Ravindran17]

•

  Number of active users > Number of actives videos [Dario13]

  T. Bottger et al. , « Open Connect Everywhere: A Glimpse at the Internet Ecosystem through the Lens of the Netflix CDN ». ACM SIGCOMM Computer Communication Review, (April 2018) 64

•

  Solutions produced by GRASP are close to the set of dominant pointsTotal cache = 25% / α = 1.0 / LCE Cache allocation solutions using GRASP

•

  The links congestion impacts the QoE [Rozhnova15][Ghada16] Ghada Moualla, Pantelis A. Frangoudis, Yassine Hadjadj Aoul, Soraya Ait Chellouche: « A Bloom-Filter-based socially aware scheme for content replication in mobile ad hoc networks ».

CCNC 2016 [Rozhnova15] Natalya Rozhnova, « Congestion control for Content-Centric Networking ». Pierre and Marie Curie University, Paris, France, 2015 77

-Clean or Green mobility Ø Support the adoption of Electric Vehicles (EVs) Ø Cooperative Adaptive Cruise Control (CACC) Ø Cars share info on speed changes in real-time • More efficient Adaptive cruise control system • Reduced carbon emissions

[Jac09] V. Jacobson et al., "Networking named content", CoNEXT, 2009. [Ravindran17] R. Ravindran, et al., "5G-ICN: Delivering ICN Services over 5G Using Network Slicing," in IEEE Com. Mag., vol. 55, no. 5, pp. 101-107, May 2017

• this will make caching possible (today we can only cache very few objects with the widespread use of https)

• Potentially a mobile edge (wireless and/or ad hoc)

• CAVs could play an important role

• sufficient processing power is available (e.g. availability of GPUs, …) • potential extension of the coverage to mobile users • potential caching with why not a sharing of caching capabilities • Everything remains to be invented ...

Mobile Edge Architecture

Before placing caches, let's ask proper questions ...

• Is caching really useful?

• How to measure the performance of a cache or a network of caches?

• How to place caches? Miss R e q u e s t ( A )

Replace one item

The offline problem is easy … If you know the sequence of requests (not realistic):

• Min policy

• At time t, if X t is not in the cache, evict an item in the cache whose next request occurs furthest in the future.

• Theorem (Maston et al. 1970) • MIN is optimal Content store Application Cache R e q u e s t ( A ) Miss R e q u e s t ( A )

Replace one item

How to compare different policies?

• Case 1: We can use trace-based simulations/experimentations Average hit ratio of the network using MACS • Good accuracy in estimating the overall cache hit performance (even when considering a multi-root tree structure)

• The cache becomes less efficient as its capacity increases How to compensate the lack of "traces"?

• We can model requests arrivals as stochastic processes (Started with [King 1971, Gelenbe 1973])

• Independent reference model (IRM) • ISPs are developing more and more caching capabilities:

• To relieve the load on their networks and improve the QoE of users

• Key questions:

• Where to strategically place caches for a network operator?

• For a cache budget requested by the content provider, where should the resources be reserved?

• Contradictory results in the state of the art [Sah17]:

• Caching at the edge is more effective • Caching at the core of the network is better

• Major drawbacks of current solutions dealing with the cache allocation problem [Sah17]:

• Only one single metric or objective is used to place the cache Setting up an efficient caching system (1/2)

• Cache allocation configuration of a total cache budget 𝑻 𝒄 defined as 𝑋 = (𝑥 # , … , 𝑥 $ ):

• 𝑥 ! represents the amount of storage capacity placed in node 𝑣 ! (𝑖 ∈ {1, . . 𝑀})

• Two performance metrics to minimize:

• f 1 (X): Content provider load à Lower cost • f 2 (X): Average distance ratio à Better Quality of Experience • f 1 (X) and f 2 (X) can be obtained analytically using our model MACS Setting up an efficient caching system (2/2)

Core

• Comes down to computing the weak composition of an integer 𝑇 % into 𝑀 parts

• Problem formulated as an integer nonlinear program Solving the cache allocation problem using GRASP

• Using Greedy Randomized Adaptive Search Procedure (GRASP) [Feo95]

• An iterative metaheuristic algorithm commonly applied to combinatorial optimization problems

𝑹𝑪𝑳 (Restricted Candidate

List): will contain some of the best elements from 𝑪𝑺 𝑋 = 𝑥 ! , 𝑥 " , 𝑥 # 𝑇 $ = 100 𝑆 = 0,0,0 𝑅𝐶𝐿 ← 𝑋 ∈ 𝐶𝑆 𝑓 𝑋 ≤ 𝑓 %&' + 𝜆 𝑓 %() -𝑓 %&' 𝜆𝜖 0,1 (𝜆 = 0: purely greedy / 𝜆 = 1: purely random) 𝑅𝐶𝐿 = { 10,0,0 , (0,10,0)} 𝑆 = 0,10,0 𝐶𝑆 = { 10,10,0 , 0,20,0 , (0,10,10)} Partial cache: 𝑃 $ = 10 𝐶𝑆 = { 10,0,0 , 0,10,0 , (0,0,10)} 𝑓 10,0,0 = 𝑣𝑎𝑙 ! 𝑓 0,10,0 = 𝑣𝑎𝑙 " 𝑓 0,0,10 = 𝑣𝑎𝑙 #

GRASP: local search phase

Create neighborhood

Select firstimproving/bestimproving solution

Contructed Solution

Repeat until no better solution is found

Update current solution

Improved Solution 𝑆 = 20,40,40 𝑁(20,40,40) = { 10,50,40 , 10,40,50 , 30,30,40 , 20,30,50 , 30,40,30 , (20,50,30)} 𝑆 = 30,30,40 𝑁(30,30,40) = { 20,40,40 , 20,30,50 , 40,20,40 , 30,20,50 , 40,30,30 , (30,40,30

Challenges

• Improving cache placement in 5G and B5G using real-world infrastructures

• Deep learning could be efficient to solve such a problem

• Placing caches in the context of multi-tenant networks • Placing in a multi-domain context

• To solve the issue of cooperation between several operators

• Placement caches in the mobile edge

Conclusions

• Caching is effective in improving performance given the incessant increase in demand • Understanding the traffic pattern helps improving the effectiveness of caching strategies • Cache placement is not a simple problem where you have to put the caches at the edge or in the core of the network

• The distribution of caches along the infrastructure is the key to getting closer to an optimum

• Including a moving edge makes the task even more complex • There are still many challenges that need to be addressed to make cache slicing a reality