which both increase linearly with frequency. Recently, the static Darcy's law correctors have been established by [START_REF] Auriault | Filtration law in porous media with poor scale separation[END_REF] and the adiabatic acoustic regime was investigated in Boutin and Bazaille ͑2005͒.

The paper is organized as follows. Section II is devoted to the principle of the homogenization method and to the basic physical assumptions concerning the medium. Section III presents the macroscopic description up to third order, and focus on the physical meaning of the correctors. In Sec. IV, the perturbations of plane wave are examined. Finally, a simple analytical example is given.

II. HOMOGENIZATION APPLIED TO THE ACOUSTICS OF GAS SATURATED POROUS MEDIA A. Homogenization principle

The macroscopic representation of heterogeneous media makes sense only if there is a scale separation. This implies ͑Auriault, 1991͒ the following: ͑1͒ The material is regular enough to show a representative volume element. This is expressed by considering that the material is composed of repeated identical cells of characteristic size l. ͑2͒ The phenomenon must vary according to a size L larger than l. In acoustics, L is related to the wavelength by L = /2 ͑Boutin and Auriault, 1990͒.

In order to catch the variations at the well distinct lengths L and l, two space variables are introduced: x for the macrovariations and y for the microvariations, x and y being related by the scale ratio ⑀ = l / L Ӷ 1; y = ⑀ -1 x. It is worth mentioning that, for a given medium, the actual physical scale ratio ⑀ varies according to the wavelength, and therefore to the frequency. Moreover, the small parameter ⑀ suggests seeking variables ͑in the present case: the pressure p, the temperature , and the velocity v͒ in the form of asymptotic expansions in powers of ⑀:

p͑x,y͒ = ͚ 0 ϱ ⑀ i i p͑x,y͒, ͑x,y͒ = ͚ 0 ϱ ⑀ i i ͑x,y͒, ͑1͒ v͑x,y͒ = ͚ 0 ϱ ⑀ i i v͑x,y͒,
where ⑀ i means ⑀ to the power i, while i p͑x , y͒ is the ith term of the expansion of p͑x , y͒, etc. The scale separation and material periodicity induce the same periodicity for the physical quantities, thus the terms i p, i , and i v are ⍀-periodic according to the variable y. The homogenization proceeds in three steps:

͑1͒ Performing a physical analysis and rescaling the equations, using powers of ⑀ for expressing the order of magnitude of the dimensionless terms ͑see Sec. III B͒. ͑2͒ Two-scale expansions are introduced in the rescaled twoscale equations and the terms of same power in ⑀ are identified ͑developed in Sec. II C͒.

͑3͒ The resolution of the problems obtained in series ͑exposed Sec. III͒.

B. Harmonic gas motion porous media: Basic and rescaled equations

The gas saturating the motionless skeleton of porosity ͑ratio of the pore volume ⍀ f to the volume ⍀ of the periodic cell͒,

=

⍀ f ⍀ , is submitted to small harmonic perturbations ͑of frequency f = /2͒ from its equilibrium state ͑where the pressure, temperature, and density take the values P e , T e , and e ͒. The variables describing the perturbations are the variations of pressure, p, temperature, , density, , and the gas velocity v ͓D͑v͒ is the strain rate͔. The parameters governing the motion and the heat transfer are the gas viscosity, , the thermal conductivity, , the specific heat ratio, ␥, and the heat mass capacity c p . Considering air as a perfect gas, one has the relation c p ͑1-1/␥͒ = P e / T e e . The linearized equations governing harmonic oscillations are given by the following ͑here and in the following the term exp͑it͒ is omitted; ٌ stands for the gradient, ⌬ for the Laplacian, dot stands for contraction, double dots for double contraction, etc.͒: ͑1͒ In the pores ͑⍀ f ͒ of the periodic cell ͑Fig. 1͒:

Gas compressibility ͑mass balance͒:

div͑v͒ + i e = 0. ͑2͒

Navier-Stokes equation ͑momentum balance͒:

div͑2D͑v͒͒ -ٌpi e v = 0. ͑3͒ ͑2͒ On the gas-solid interface ⌫, the boundary conditions are ͑here and in the following the notation /⌫ means the value on ⌫͒: Adherence condition of the gas:

v /⌫ = 0. ͑6͒
Isothermal condition imposed by the weak thermal impedance of air, ͱ e c p , compared to that of the solid:

/⌫ = 0. ͑7͒

During wave propagation, the pressure and the density oscillate according to the wavelength whereas the velocity and temperature vary at the pore scale. Then, using the wavelength as reference length, the Navier-Stokes and Fourier equations ͑3͒ and ͑4͒ have to be rescaled, the mass balance ͑2͒ being unchanged. The richest case is obtained when the pore size is of the same order as the viscous and thermal skin depths, i.e., ͉␦ v ͉ = O͑l͒ where ␦ v = ͱ / i e and ͉␦ t ͉ = O͑l͒ where ␦ t = ͱ / i e c p = O͑l͒ ͑this is consistent with the fact that viscous and thermal layer thicknesses are of the same order for air͒. This situation leads one to replace Eqs. ͑3͒ and ͑4͒ by the two following rescaled equations ͓see, for instance, Boutin et al. ͑1998͔͒:

⑀ 2 div͑2D͑v͒͒ -ٌp -i e v = 0, ͑8͒ ⑀ 2 div͑ ٌ ͒ -i͑ e c p -p͒ = 0. ͑9͒
Note that the gas state equation enables the elimination of the density, and that the divergence of shear stresses, combined with the mass balance, may be expressed as follows:

div͑2D͑v͒͒ = ͓⌬͑v͒ + ٌdiv͑v͔͒ = ͫ ⌬͑v͒ -i ٌ ͩ p P e - T e ͪͬ.
Thus, the system is driven by the three differential operators acting on the variables p, v, and :

G͑p,v,͒ = 0, G͑p,v,͒ = div͑v͒ + iͩ p P e - T e ͪ, ͑10͒ N͑p,v,͒ = 0, N͑p,v,͒ = -ٌp -i e v + ⑀ 2 ͫ ⌬͑v͒ -i ٌ ͩ p P e - T e ͪͬ, ͑11͒ F͑p,͒ = 0, F͑p,͒ = i͑p -e c p ͒ + ⑀ 2 div͑ ٌ ͒. ͑12͒

C. Homogenization process

Using the two space variables x and y = ⑀ -1 x, the gradient ٌ is changed into ⑀ -1 ٌ y + ٌ x , and the Laplacian ⌬ becomes ⑀ -2 ⌬ y +2⑀ -1 ⌬ yx + ⌬ x , with ⌬ yx = ץ x i y i 2 . Consequently, the operators G, N, F expressed with the two space variables are changed into G xy , N xy , F xy , whose expansions are given by Eqs. ͑13͒-͑15͒:

Ά G xy ͑p,v,͒ = ⑀ -1 G -1 ͑v͒ + G 0 ͑p,v,͒, G -1 ͑v͒ = div y ͑v͒, G 0 ͑p,v,͒ = div x ͑v͒ + iͩ p P e - T e ͪ, ͑13͒ Ά N xy ͑p,v,͒ = ⑀ -1 N -1 ͑p͒ + N 0 ͑p,v͒ + ⑀N 1 ͑p,v,͒, + ⑀ 2 N 2 ͑p,v,͒ N -1 ͑p͒ = -ٌ y p, N 0 ͑p,v͒ = -ٌ x p -i e v + ⌬ y ͑v͒, N 1 ͑p,v,͒ = 2⌬ yx ͑v͒ -i P e ٌ y ͩp-P e T e ͪ , N 2 ͑p,v,͒ = ⌬ x ͑v͒ -i P e ٌ x ͩp-P e T e ͪ , ͑14͒ Ά F xy ͑p,͒ = F 0 ͑p,͒ + ⑀F 1 ͑p͒ + ⑀ 2 F 2 ͑p͒, F 0 ͑p,͒ = ip -i e c p + ⌬ y ͑͒, F 1 ͑͒ = 2⌬ yx ͑͒, F 2 ͑͒ = ⌬ x ͑͒.

͑15͒

Introduce now the expansions ͑1͒ of p, v, and in G xy , N xy , and F xy and identify terms of identical power of ⑀. These balance equations become series in power of ⑀ identically equal to zero whatever ⑀ Ӷ 1, so that each term must vanish. The so-derived equations at each order-combined with the adherence and isothermal conditions on ⌫ satisfied by each i v and i -lead to a series of problems to be solved recurrently. Section III presents the resolution up to the second order.

III. DESCRIPTION UP TO THE SECOND ORDER

A. Macroscopic mass balances

Before any resolution let us first mention that each mass balance of order i ജ 0, reads

͕G -1 ͑ i+1 v͒ + G 0 ͑ i p, i v, i ͒ = 0; i+1 v /⌫ = 0; i+1 v ⍀-periodic͖.
The local periodicity introduces a compatibility condition obtained by integrating the local mass balance on the pore cell volume:

͵ ⍀ f ͭ div y ͑ i+1 v͒ + div x ͑ i v͒ + iͩ i p P e - i T e ͪͮd⍀=0.
From the divergence theorem, the periodicity and the adherence condition on ⌫, the first term vanishes. Then, inverting y integration and x derivative,

div xͩ ͵ ⍀ f i vd⍀ ͪ + i ͩ͵ ⍀ f i p P e d⍀ -͵ ⍀ f i T e d⍀ ͪ = 0.
Introducing the physical ͑i.e., observable͒ macrovariables,

i V = 1 ⍀ ͵ ⍀ f ⑀ i i vd⍀ = ⍀ f ͵ ⍀ f ⑀ i i vd⍀, i P = 1 ⍀ f ͵ ⍀ f ⑀ i i pd⍀, i T = 1 ⍀ f ͵ ⍀ f ⑀ i i d⍀,
one obtains the macroscopic mass balance at order i,

div x ͑ i V͒ + iͩ i P P e - i T
T e ͪ=0.

͑16͒

B. Expressions of local fields

To go further we need to determine the local fields ͑ i p , i v , i ͒ by solving successively the harmonic viscothermal linear problems on the cell. At each order, the problem may be split into a leading set of equations related to the gas flow, and a set of equations related to the heat transfer induced by the gas pressure. The general scheme of the procedure and the first steps of the resolution are presented hereafter. The more technical steps are reported in Appendix A.

Pressure and temperature at the leading "zero… order

The momentum balance at order ⑀ -1 reduces to N -1 ͑ 0 p͒ =-ٌ y 0 p = 0, which gives 0 p͑x,y͒ = 0 P͑x͒.

This result introduced in the heat transfer at order ⑀ 0 , namely F 0 ͑ 0 P , 0 ͒ = 0, leads to the differential system ͕S t 0 ͖ set on ⍀ f , governing 0 :

S t 0 ͭ i 0 P͑x͒ -i e c p 0 + ⌬ y ͑ 0 ͒ = 0, 0 /⌫ = 0, 0 ⍀-periodic.
The solution of this linear problem, with 0 P͑x͒ as forcing term, may be expressed in the form: 0 ͑x,y͒

T e = 0 ͩ y ␦ t ͪ 0 P͑x͒ P e . ͑17͒
The temperature distribution 0 ͓solution for 0 P͑x͒ = P e / T e ͔ is complex and depends on the local variable and the frequency through the dimensionless variable y / ␦ t . In the isothermal regime, reached at low frequency ͑y / ␦ t → 0͒, the temperature vanishes and 0 ͑y / ␦ t ͒ → 0. In the adiabatic regime, reached at high frequency ͑y / ␦ t → ϱ ͒, the temperature tends to the uniform value ͑1-1/␥͒.

Local velocity, pressure, and temperature at the first order

The local velocity, 0 v, and the pressure of the first order, 1 p, are obtained from the following system of equations ͕S v 0 ͖ set on ⍀ f :

͕N -1 ͑ 1 p͒ + N 0 ͑ 0 P, 0 v͒ = 0, G -1 ͑ 0 v͒ = 0,
͕ 0 v /⌫ = 0, and 0 v, 1 p ⍀-periodic͖ which reads explicitly:

S v 0 Ά -ٌ y 1 p -ٌ x 0 P -i e 0 v + ⌬ y ͑ 0 v͒ = 0,
div y ͑ 0 v͒ = 0, 0 v /⌫ = 0, 0 v, and 1 p ⍀-periodic.

This set defines the linear dynamic permeability problem with ٌ x 0 P as forcing term, whose solution is as follows ͓see Auriault ͑1980͔͒:

0 v͑x,y͒ = 0 kͩ y ␦ v ٌͪ x 0 P, ͑18͒ 1 p͑x,y͒ = 1 ␣ͩ y ␦ v ٌͪ x 0 P + 1 p ˆ͑x͒.
The tensors 0 k and 1 ␣ are constituted by the three velocity and pressure distributions ͑-0 k i / , 1 ␣ i ͒ corresponding to unit pressure gradient in the three directions, ٌ x 0 P͑x͒ = e i . The solutions ͑ 0 k i , 1 ␣ i ͒ are complex and depend on the dimensionless variable y / ␦ v . In the quasistatic ͑low frequency͒ regime ͑l / ␦ v → 0͒, the inertial effect vanishes and 0 k tends to be real. At high frequency ͑l / ␦ v → ϱ ͒, the inertia dominates and 0 k tends to a purely complex tensor. Note that the pressure is defined up to a constant value 1 p ˆand that, by construction, 1 ␣ is of zero mean value on the cell. Now, the first corrector 1 of the temperature is derived from the heat transfer at order ⑀ governed by the differential system F 0 ͑ 1 p , 1 ͒ + F 1 ͑ 0 ͒ = 0. This leads to the system ͕S t 1 ͖ set on ⍀ f :

S t 1 ͭ i 1 p͑x,y͒ -i e c p 1 + ⌬ y ͑ 1 ͒ + 2⌬ yx ͑ 0 ͒ = 0, 1 /⌫ = 0, 1 ⍀-periodic.
The previously determined fields 0 and 1 p ͓Eqs. ͑17͒ and ͑18͔͒ introduce two independent forcing terms: ͑1͒ 1 p ˆ͑x͒, inducing an identical problem than ͕S t 0 ͖ except that 1 p ˆ͑x͒ replaces 0 P͑x͒ and ͑2͒ ٌ x 0 P, linked with 1 ␣ and ٌ y 0 ͓coming from the equality: ⌬ yx ͑ 0 ͒ = ͑T e / P e ٌ͒ y 0 . ٌ x 0 P͔.

Therefore, by linearity, the solution is in the form:

1 ͑x,y͒

T e = 0 ͩ y ␦ t ͪ 1 p ˆ͑x͒ P e + 1 ͩ y ␦ t , y ␦ v ͪ. ٌ x 0 P P e . ͑19͒
The tensor of temperature distribution 1 is constituted by the three solutions 1 i corresponding to pressure gradient in the three directions, ٌ x 0 P͑x͒ = ͑P e / T e ͒e i . These particular solutions are complex and depend on the local variable and frequency through the two dimensionless variables y / ␦ t and y / ␦ v . This corresponds to a nonlocal transient equilibrium with thermoviscous coupling.

Local fields at the two next orders

Following the process, the pressure and velocity, 2 p and 1 v, are determined from the set

͕S v 1 ͖ S v 1 Ά N -1 ͑ 2 p͒ + N 0 ͑ 1 p, 1 v͒ + N 1 ͑ 0 P, 0 v͒ = 0, G -1 ͑ 1 v͒ + G 0 ͑ 0 P, 0 v͒ = 0, 1 v /⌫ = 0, 1 v, and 2 p ⍀-periodic.
and the temperature 2 is solution of ͕S t 2 ͖:

S t 2 ͭ F 0 ͑ 2 p, 2 ͒ + F 1 ͑ 1 ͒ + F 0 ͑ 0 ͒ = 0, 2 /⌫ = 0, 2 ⍀-periodic.
Finally, 3 p and 2 v are derived from ͕S v 2 ͖:

S v 2 Ά N -1 ͑ 3 p͒ + N 0 ͑ 2 p, 2 v͒ + N 1 ͑ 1 p, 1 v͒ + N 2 ͑ 0 P, 0 v͒ = 0, G -1 ͑ 2 v͒ + G 0 ͑ 1 p, 1 v͒ = 0, 2 v /⌫ = 0, 2 v, and 1 p ⍀-periodic.
The resolution of ͕S v 1 ͖, ͕S t 2 ͖, and ͕S v 2 ͖ is reported in Appendix A. Let us underline that, at each order, the same differential set as the previous order is recovered, except for new additive terms, involving the gradient of the solutions of lower order. Thus, once ͑ j+1 p , j v , j ͒, 0ഛ j Ͻ i, have been determined, they become forcing terms in the problems related to ͑ i+1 p , i v , i ͒, and so on.

Local fields up to the second order

To sum up, the pressure, velocity, and temperature, read, up to the second order:

p͑x,y͒ = 0 P͑x͒ + ⑀͓ 1 p ˆ͑x͒ + 1 ␣ٌ x 0 P͑x͔͒ + ⑀ 2 ͫ 2 p ˆ͑x͒ + 1 ␣ . ٌ x 1 p ˆ͑x͒ + 2 ␣ . . ٌ x ٌ x 0 P͑x͒ + i P e 2 0 P͑x͒ ͬ + ¯, -v͑x,y͒ = 0 k . ٌ x 0 P͑x͒ + ⑀ ͫ 0 k . ٌ x 1 p ˆ͑x͒ + 1 k . . ٌ x ٌ x 0 P͑x͒ + i P e 1 n 0 P͑x͒ ͬ + ⑀ 2 ͫ 0 k . ٌ x 2 p ˆ͑x͒ + 1 k . . ٌ x ٌ x 1 p ˆ͑x͒ + 2 k ... ٌ x ٌ x ٌ x 0 P͑x͒ + i P e 2 n . ٌ x 0 P͑x͒ ͬ + ¯, ͑P e /T e ͒͑x,y͒ = 0 0 P͑x͒ + ⑀͓ 0 1 p ˆ͑x͒ + 1 . ٌ x 0 P͑x͔͒ + ⑀ 2 ͫ 0 2 p ˆ͑x͒ + 1 . ٌ x 1 p ˆ͑x͒ + 2 . . ٌ x ٌ x 0 P͑x͒ + i P e 2 0 P͑x͒ ͬ + ¯.
As previously, tensors i+1 ␣ and i k, respectively of rank i +1 and i + 2, are constituted by the 3 i+1 particular solutions of the purely viscoinertial transient problems ͕S v i ͖ under the 3 i+1 components of the forcing term ٌ͑ x ͒ i+1 0 P. All of them are complex and depend on the local variable and frequency through y / ␦ v only ͑cf. Appendix A͒. More precisely,

͑1͒ ͑-1 k mn / , 2 ␣ mn ͒ are the 9 solutions of ͕S v 1 ͖ with ٌ x ٌ x 0 P͑x͒ = e m e n . ͑2͒ ͑-2 k pqr / , 3 ␣ pqr ͒ are the 27 solutions of ͕S v 2 ͖ with ٌ x ٌ x ٌ x
0 P͑x͒ = e p e q e r .

Tensors i+1 and i n, of rank i -1 and i, are constituted by the 3 i-1 particular solutions of the coupled thermoviscoinertial transient problems ͕S v i ͖ under the 3 i-1 components of the forcing term ٌ͑ x ͒ i-1 0 P. Consequently, they are complex and depend on both y / ␦ t and y / ␦ v :

͑1͒ ͑-1 n / , 2 ͒ is the solution of ͕S v 1 ͖ with ͑i / P e ͒ 0 P͑x͒ =1. ͑2͒ ͑-2 n m / , 3 m ͒ are the 3 solutions of ͕S v 2 ͖ with ͑i / P e ٌ͒ x 0 P͑x͒ = e m .
Finally, tensors i of rank i are constituted by the 3 i solutions of ͕S t i ͖ under the 3 i unit components of the forcing term ͑T e / P e ٌ͒͑ x ͒ i 0 P. Except for 0 , solution of the purely thermoinertial transient problem ͕S t 0 ͖, i for i Ͼ 0, involves thermoviscoinertial coupling, so that they are complex and depend on both y / ␦ t and y / ␦ v . This is also the case of the scalar 2 , solution of ͕S t 2 ͖ under ͑i / P e ͒ 0 P͑x͒ =1. It is important to keep in mind that the fields observable in the reality are i p⑀ i , i v⑀ i , and i ⑀ i . Thus, the actual physical tensors-independent of the scale ratio !-are i k⑀ i , i ⑀ i , i ␣⑀ i , and 2 n⑀ 2 .

C. Continuum description and correctors

Knowledge of the local fields enables to express the macrovariables appearing in the macromass balances ͑16͒.

The averaging of the local fields ͑of nonzero mean value͒ on the cell enables definition of the macroscopic tensors whose dimension and magnitude are, respectively, i K ϳ O͑ 0 K͒l i , i ⌸ ϳ l i , 2 N ϳ l 2 , and 2 Z ϳ O͑1͒. They can be determined ͑nu-merically͒ as soon as the physical parameters of the gas, the frequency, and the pore geometry are given. Note that the tensors issued from viscothermal coupled problems depend on both dimensionless frequencies / c and / t ͓see Eqs. ͑21͒ and ͑22͒ for the definition of c and t ͔,

i Kͩ c ͪ= 1 ⍀ ͵ ⍀ f i k⑀ i d⍀, 0 ⌸ͩ t ͪ= 1 ⍀ f ͵ ⍀ f 0 ⑀ i d⍀ and for i Ͼ 1, i Nͩ c , t ͪ= 1 ⍀ ͵ ⍀ f i n⑀ 2 d⍀, i ⌸ͩ t , c ͪ= 1 ⍀ f ͵ ⍀ f i ⑀ i d⍀, 2 Zͩ t , c ͪ= 1 ⍀ f ͵ ⍀ f 2 ⑀ 2 d⍀.
Their main properties, especially symmetries, and relations enabling one to deduce effective tensors of the second order from local solutions of the first order, are established in Appendix B.

Biot-Allard modeling

The leading equations ͑i.e., zero order͒ are in agreement with the phenomenological approaches of Zwikker and Kosten ͑1949͒, Attenborough ͑1983͒, and Allard ͑1993͒. The equivalent continuum is described by the dynamic permeability tensor 0 K / and the effective compressibility ͓1 -0 ⌸͔ / P e ͑in the following, to save notations, the index x is omitted for the macroscopic derivatives͒:

Ά div͑ 0 V͒ + iͫ 0 P P e - 0 T T e ͬ=0, 0 V = -0 K . ٌ 0 P, 0 T T e = 0 ⌸ 0 P P e .

͑20͒

It is necessary for the following to recall the features of the dynamic permeability ͑Auriault et al., 1985͒. At low frequencies, viscous effects dominate and 0 K tends toward the real-valued intrinsic permeability, K. At high frequencies, inertia dominates and 0 K tends towards a pure imaginary value, / i e ␣ ϱ , where ␣ ϱ is the tortuosity. Low and high frequency domains are delimited by a critical frequency derived by equalizing viscous and inertial effects of the macroscopic flow:

c = K e ␣ ϱ . ͑21͒
As for the effective compressibility, at low frequencies the temperature variation vanishes, so that 0 ⌸ → 0 and the effective compressibility tends toward the isothermal compressibility / P e . At high frequencies, conduction effects are negligible, except in close proximity to the solid. The perturbations become quasiadiabatic so that 0 ⌸ → 1-1/␥, and the compressibility tends toward the adiabatic value / ␥P e . The thermal characteristic pulsation delimiting both regimes is related to the length ⌳ t , defined-up to a geometric factor of the order of one-as the ratio of the volume to the surface of pore ⍀ f / ⌫ ͑Champoux and Allard ͑1991͒͒:

t = ⌳ t 2 e c p . ͑22͒
Whereas viscous and thermal layer thicknesses are of the same order in air, the frequencies c and t can be rather different because the permeability essentially depends on the small ducts, while thermal transfer involves the whole pore volume. In consequence, c Ͼ t . Note also that the effective compressibility always contains a real part, and thus the change from isothermal to adiabatic conditions has consequences less drastic than the change from viscous to inertial regime.

Correctors to Biot-Allard modeling

The equations governing the two next orders provide the correctors to the continuum description.

First corrector:

Ά div͑ 1 V͒ + iͫ 1 P P e - 1 T T e ͬ=0, 1 V = -0 K . ٌ 1 P -1 K . . ٌٌ 0 P - i P e 1 N 0 P, 1 T T e = 0 ⌸ 1 P P e + 1 ⌸ . ٌ 0 P P e .

͑23͒

Second corrector:

Ά div͑ 2 V͒ + iͫ 2 P P e - 2 T T e ͬ=0, 2 V = -0 K . ٌ 2 P -1 K . . ٌٌ 1 P -2 K ... ٌٌٌ 0 P - i P e 2 N . ٌ 0 P, 2 T T e = 0 ⌸ 2 P P e + 1 ⌸ . ٌ 1 P P e + 2 ⌸ . . ٌٌ 0 P P e + i P e 2 Z 0 P P e . ͑24͒
These equations underline that a poor scale separation introduces deviations from the continuum description ͑strictly valid for infinite scale ratio͒ involving the successive gradients of pressure. The origin of those correctors lies in the terms neglected at preceding orders that become significant at the considered order. They bring unusual nonlocal terms and gas compressibility terms in both dynamic Darcy law and thermal behavior.

Isotropy or symmetric cell case: Cancellation of the first corrector of pressure, velocity, and temperature

In the case of macroscopic isotropy of the medium ͑up to the second order͒ or of cell presenting symmetry according to three orthogonal planes, tensors of odd rank, 1 K and 1 ⌸ cancel out. Therefore, Eqs. ͑23͒ driving the first correctors reduces to the same as that of the leading order ͑20͒. Consequently, without loss of generality, 1 P , 1 V , 1 T can be canceled out, and the effective correction of the three variables is reported to the second order, that is:

P = 0 P + 2 P + ... , V = 0 V + 2 V + ... , T = 0 T + 2 T + ... .
Furthermore, as any fourth rank isotropic tensor, 2 K takes the following form ͑a , b , c are three independent scalars and ␦ the Kronecker symbol͒:

2 K i jkl = a␦ i j ␦ k l + b␦ i k ␦ j l + c␦ i l ␦ k j giving: 2 K ... ٌٌٌ 0 P = ͑a + b + c͒ ٌ ͑⌬ 0 P͒.
Consequently, 2 K may be reduced to a scalar function 2 k = a + b + c and the velocity correctors read ͑ 0 k and 2 n are the scalar functions associated by the isotropic tensors 0 K and 2 N͒:

2 V = -0 k ٌ 2 P -2 k ٌ ͑⌬ 0 P͒ - i P e 2 n ٌ 0 P = -0 k ٌ 2 P + 2 k 0 k ⌬͑ 0 V͒ - i P e
2 n ٌ 0 P.

After adding the terms of zero order, we have

ͫ V -2 k 0 k ⌬͑V͒ ͬ = -0 k ͫ 1 + i P e 2 n 0 k ͬ ٌ P + O͑⑀ 3 ͒. ͑25͒
The effective dynamic Darcy's law is enriched by two terms: ͑a͒ One linked to 2 k is of viscoinertial nature-it corresponds to a generalization of the Brinkmans law in the dynamic range. ͑b͒ The other linked to 2 n involves the thermoviscous coupling induced by the gas compressibility-it introduces a correction of permeability.

As for the thermal corrector in isotropic case, one obtains ͑ 2 is the scalar function associated to the isotropic tensor 2 ⌸͒: This effective state equation for the gas includes nonlocal correction in the dynamic range with thermoviscous coupling effects linked with the compressibility.

General case: Cancellation of the first corrector of pressure

By combining the three equations of Eq. ͑20͒ governing the leading order, the velocity and the temperature may be eliminated, and the scalar wave equation expressed with the pressure 0 P only is derived:

-divͩ 0 K . ٌ 0 Pͪ + i ͓1 -0 ⌸͔ P e 0 P = 0. ͑26͒
Similarly, the first order set ͓Eq. ͑23͔͒ yields

-divͩ 0 K . ٌ 1 Pͪ + i ͓1 -0 ⌸͔ P e 1 P = S 1 ͑ 0 P͒, ͑27͒ S 1 ͑ 0 P͒ = divͩ 1 K . . ٌٌ 0 Pͪ + i P e 1 N . ٌ 0 P + i P e 1 ⌸ . ٌ 0 P.

͑28͒

It is worth mentioning that, if the same differential operator ͑left-hand side͒ applies for 0 P and 1 P, the difference between zero-and first-order equations ͑26͒ and ͑27͒ comes from the source term on the right-hand side. In fact, Eq. ͑27͒ means that the first corrector field 1 P results from radiation of a density of source S 1 ͑ 0 P͒ generated by the 0 P field. However, the first corrector presents generally a particular ener-getic property. In composite media, it was shown-in static or dynamic elasticity or for thermal transfer-that the work of the first corrective term under the zero-order field is null ͑Boutin and [START_REF] Boutin | Rayleigh scattering in elastic composite materials[END_REF][START_REF] Boutin | Microstructural influence on heat conduction[END_REF]. This is shown below in the frame of poroacoustics. The source term S 1 ͑ 0 P͒ reads

S 1 ͑ 0 P͒ = 1 1 K j pq0 P ,pqj + i P e ͓ 1 N p + 1 ⌸ p ͔ 0 P ,p .
However, the skew symmetry of 1 K and identity between 1 N and 1 ⌸ can be demonstrated, see Appendix B: 1 K j pq = -1 K p jq and 1 N p = -1 ⌸ p . ͑29͒

Combining this identity with 0 P ,pqj = 0 P ,jqp , one deduces that

S 1 ͑ 0 P͒ = 0 then 1 P = 0.
As a matter of fact, 1 P being governed by the same equations as 0 P, it can be canceled out without loss of generality. Nevertheless, contrary to the isotropic case, 1 V and 1 T do not necessarily disappear.

Effective correctors

To sum up, in any case, the governing equations ͓Eqs. ͑23͒ and ͑24͔͒ can be replaced by the following.

First corrector:

1 P = 0, 1 V = -1 K . ٌٌ 0 P - i P e 1 N 0 P, ͑30͒ 1 T T e = 1 ⌸ . ٌ 0 P P e .
Because of the identities ͑29͒, the macromass balance in Eq. ͑23͒ is necessarily satisfied.

Second corrector:

div͑ 2 V͒ + iͫ 2 P P e - 2 T T e ͬ=0, 2 V = -0 K . ٌ 2 P -2 K ... ٌٌٌ 0 P - i P e 2 N . ٌ 0 P, 2 T T e = 0 ⌸ 2 P P e + 2 ⌸ . . ٌٌ 0 P P e + i P e 2 Z 0 P P e . ͑31͒
The equation governing the actual pressure corrector 2 P derived from Eq. ͑31͒ takes the form:

-divͩ 0 K . ٌ 2 Pͪ + i ͓1 -0 ⌸͔ P e 2 P = S 2 ͑ 0 P͒ ͑ 32͒ with S 2 ͑ 0 P͒ = divͩ 2 K ... ٌٌٌ 0 Pͪ + i P e ͫ 2 N . . ٌٌ 0 P + ͩ 2 ⌸ . . ٌٌ 0 P + i P e 2 Z 0 Pͪ ͬ . ͑33͒
As shown previously, the term on the right-hand side acts as a density of source generated by the 0 P field. However, at this order the source S 2 ͑ 0 P͒ does not cancel out and, in turn, a corrective 2 P field is radiated ͑according to the zero-order macrobehavior͒: We are back to the Rayleigh scattering, in which the passing of a long wave through heterogeneity generates sources that radiate perturbations. The interest of the homogenization is to replace these sources and their precise radiation by mean sources and equivalent fields at the considered order. Furthermore, theoretical expressions are provided for determining the effective coefficients from the knowledge of the microstructure.

Simplification in adiabatic "or isothermal… regime

The above presented description ͓Eqs. ͑20͒, ͑30͒, and ͑31͔͒ may be simplified if the adiabatic regime is assumed in the gas. In that case, the conduction disappears and temperature, pressure, and density are related by

T e = ͑1 -1/␥͒ P P e , e =
P ␥P e , ͑34͒ so that 0 ⌸ =1-1/␥ and the tensors i ⌸, 2 Z, 1 N disappear. Tensors 0 K, 1 K, and 2 K are the same as previously noted, since they are independent of the heat transfer. However, the tensor 2 N is modified by the absence of thermal coupling. The isothermal case, which implies = 0, then i ⌸ = 0, and = P͑ e / P e ͒ leads to similar simplifications ͑with 2 N different than in the adiabatic case͒.

IV. RAYLEIGH SCATTERING OF PLANE WAVES

A. The continuum approximation

Consider a porous medium ͑isotropy is not needed͒ and investigate the scattering effect on an harmonic plane wave propagating in a given direction ͑of unit vector d͒. Denoting by h d = h r + ih i the complex wave number in this direction, the variations of the zero-order pressure take the form:

0 P = 0 P exp͑-ih d x . d͒exp͑it͒. ͑35͒
This plane wave is an eigenmode of the medium, and, from the zero-order balance Eq. ͑26͒,

͑-ih d ͒ 2 = i ͓1 -0 ⌸͔ P e K d
where

K d = d i 0 K i j d j .
The wave velocity C d , wavelength d , and attenuation d in direction d are related to the wave number by

h d = C d ͑1 -i d ͒, i.e., C d = h r , d = 2 h r , d = - h i h r
Because of the variations of the dynamic permeability with frequency recalled in Sec. III C 1, the P 2 wave shows a strong dispersion, evolving from ͑1͒ a diffusion wave at low frequency ͑ Ͻ c ͒:

C d Ӎ ͱ ͑ / c ͒ ͱ ͑P e / ␣ ϱ e ͒; d =1-O͑ ͱ ͑ / c ͒͒.
͑2͒ to a propagation wave at high frequency ͑ Ͼ c ͒:

C d Ӎ ͱ ͑␥P e / ␣ ϱ e ͒; d ϳ ͱ ͑ c / ͒.

B. Scattering effect

As 1 P = 0, the above-presented continuum description is valid up to a precision ⑀ 2 for the pressure and ⑀ for the flow and temperature. This is sufficient for long waves. For shorter waves the precision is improved by considering the higher order terms. With this aim, let us first calculate the source term S 2 ͑ 0 P͒. Replacing in Eq. ͑33͒ 0 P by its expression ͑35͒, gives

S 2 ͑ 0 P͒ = ͭ 1 2 K d ͑-ih d ͒ 4 + i P e ͭ ͓ 2 N d + 2 ⌸ d ͔͑-ih d ͒ 2 + i P e 2 Z ͮͮ 0 P exp ͑-ih d x . d͒ where 2 K d = d i 2 K i jkl d j d k d l , 2 N d = d i 2 N i j d j , 2 ⌸ d = d i 2 ⌸ i j d j .
This source is proportional ͑but not in phase͒ to 0 P. Thus, it acts as a forcing term loading the medium according to one of the eigenmodes, i.e., the plane wave in the direction d. This induces a self-resonance effect and subsequently the radiated field 2 P is linearly amplified as the wave progresses. Thus, the field 2 P satisfying Eq. ͑26͒ is in the form:

2 P = -Q d ͑-ih d ͒ 2 ͑-ih d x . d͒ 0 P exp͑-ih d x . d͒, ͑36͒
where the frequency-dependent complex coefficients, Q d and

Q d ͑-ih d ͒ 2 , read: 2Q d = 2 K d K d + 2 N d + 2 ⌸ d ͑1 -0 ⌸͒ + 2 Z K d ͓͑1 -0 ⌸͔͒ 2 ͑37͒ and 2Q d ͑-ih d ͒ 2 = i P e ͫ ͑1 -0 ⌸͒ 2 K d ͑K d ͒ 2 + 2 N d + 2 ⌸ d K d + 2 Z 1 -0 ⌸ ͬ .

͑38͒

Expression ͑36͒ means that the scattered field is coherent ͑consistently with the scale separation assumption͒ and radiates in the same direction as the zero-order field. Finally, up to the second order the macropressure reads:

0 P + 2 P = ͓1 -Q d ͑-ih d ͒ 2 ͑-ih d x . d͔͒ 0 P exp͑-ih d x . d͒ = P͑1 + O͑⑀ 3 ͒͒.
At the same level of approximation ͑apparently within a distance x such as O͉͑h d ͉ x . d͒ Ͻ 1/ ͱ ⑀, however this restriction can be removed by considering new boundary conditions beyond this distance͒, the macropressure can also be expressed as

P = 0 P ͑exp͑-ih d x . d͒ + O͑⑀ 3 ͒͒ with h d = h d ͓1 -Q d ͑-ih d ͒ 2 ͔.
Thus the diffraction modifies the wave number h d given by the continuum approximation. As expected, the correction is of the order of

Q d ͑-ih d ͒ 2 = O͑l 2 / ͑ /2͒ 2 ͒ = O͑⑀ 2 ͒.
This perturbation comes from the interference between the zero-order wave and the coherent and amplified wave generated by the induced source distribution. Writing Q d on the polar form: Q d = ͉Q d ͉ exp͑i͒, one derives the apparent velocity and attenuation factor by expressing explicitly that

h d = C d ͑1 -i d ͒ = C d ͑1 -i d ͒ ͫ 1 + ͉Q d ͉͑cos + i sin ͒ ϫ ͫ ͑1 -i d ͒ C d ͬ 2 ͬ .
Algebra leads to

C d = C dͫ 1 -͉Q d ͉ͩ C d ͪ 2 ͓͑1 -d 2 ͒͑cos + d sin ͒ + 2 d ͑sin -d cos ͔͒ ͬ , ͑39͒ d = d + ͉Q d ͉ͩ C d ͪ 2 ͓2 d ͑1 + d 2 ͒cos -͑1 -d 4 ͒sin ͔.

͑40͒

Those expressions point out the difference between the scattering effect for poroacoustic waves and for elastic or thermal waves. Let us recall that in the two latter cases the analogue of Q d is a purely real constant ͑i.e., =0͒ strictly related to the geometry of the microstructure ͓the reason is that-contrary to poroacoustics-in those cases, quasistatic conditions govern the local physics, see ͓Boutin and Auriault ͑1993͒, Boutin ͑1995͒, Chen and fish ͑2001͔͒. In elasticity, the zero-order velocity is a constant, and there is no attenuation ͑i.e., =0͒: The scattering leads to a velocity dispersion varying according to 2 and an apparent attenuation increasing according to 3 . The thermal waves being dispersive ͑C ϳ ͱ ͒ and damped ͑ =1͒, the scattering induces an addi- tional velocity dispersion and attenuation both varying linearly with .

In poroacoustics the rules for the physical consequences of scattering are not so easy to draw, first, because of the change in the P 2 wave from a diffusion to a propagation mode, and, second, because Q d is a frequency-dependent complex function. Consequently, the velocity and attenuation are both modified, but no general simple trends in the whole frequency range can be drawn for the frequency dependence.

In the viscous regime, the vanishing of the transient effects at the pore scale makes macroscopic tensors tend toward real value, then Q d → Q d0 and → 0 and one obtains a correction identical to that of a diffusion wave:

C d Ӎ C dͫ 1 + 2Q d0 ͩ C d ͪ 2 ͬ , d Ӎ d + 4Q d0 ͩ C d ͪ 2 ͑41͒ with C d Ӎ ͱ c ͱ P e ␣ ϱ e , d = 1 -Oͩͱ c ͪ.
In the inertial regime ͑below diffraction͒, i.e., c Ͻ Ͻ d , the wave propagates with damping and one may expect perturbation in both velocity and attenuation:

C d Ӎ C dͫ 1 -͉Q d ͉ͩ C d ͪ 2 ͑cos + 3 d sin ͒ ͬ , ͑42͒ d Ӎ d -͉Q d ͉ͩ C d ͪ 2 ͓sin -2 d cos ͔ with C d Ӎ ͱ ␥P e ␣ ϱ e , d ϳ ͱ c .
In the whole frequency range and for common pore morphology, the estimates of the tensors ͑Sec. III C͒ suggest the following assessment, where G is a dimensionless geometric factor, which accounts for the pore geometry and the ratio between pore and elementary representative volume sizes:

Q͑-ih͒ 2 ϳ G ͫ K 0 K ͬ i ␥P e , i.e., C C ˜1 -i 1 -i -1 ϳ G ͫ K 0 K ͬ i ␥P e .

C. Diffraction and critical frequency

Since the diffraction effects are O͑⑀ 2 ͒ they are actually significant when ⑀ → 1, i.e., for frequencies close to d . To identify the qualitative nature of the diffraction when it becomes significant, it is necessary to compare d to c .

According to the properties of P 2 waves, the wavelength may be assessed as ͑where C a = ͱ ␥P e / e is the sound veloc- ity in air͒

Ͻ c ⇒ /2 Ϸ C a ͱ c , Ͼ c ⇒ /2 Ϸ C a ͱ ␣ ϱ .
Applied at the diffraction frequency ͑ d is such that /2 = l͒, these relations give

d c Ͻ 1 ⇒ d c = ͩ C a l c ͪ 2 , d c Ͼ 1 ⇒ d c = C a l c ͱ ␣ ϱ = C a e K l .
If we now consider:

͑1͒ Air properties at the normal conditions, / C a e Ӎ 1 . 810 -5 Pa s ր 343 m / s1.2 kg/ m 3 Ӎ 4,3.10 -8 m and ͑2͒ The very rough permeability estimate K = O͑l 2 ͒, it appears that as soon as l Ͼ 4 ϫ 10 -8 m, then d Ͼ c . This means that the diffraction effect is significant in the dynamic regime, except for materials of pores smaller than about 50 nm.

The self-consistent approach for the intrinsic permeability of granular material-of characteristic size l =2R and porosity ͑Boutin, 2000͒-supplies the more realistic estimate:

K = R 2 ⌿͑͒, ⌿͑͒ = 1 3 ͫ 2 + 3͑1 -͒ 5/3 3 + 2͑1 -͒ 5/3 1 ͑1 -͒ 1/3 -1 ͬ .
Using this expression, Fig. 2 shows the ratio d / c versus the porosity for several realistic characteristic sizes. Clearly, the diffraction takes place in the viscous regime only for materials presenting rather fine pores ͑say less than 1 m͒ and/or small porosity. For most of the sound-absorbing materials, the diffraction frequency is to be expected within the dynamic regime.

V. A SIMPLE ACADEMIC EXAMPLE

To illustrate these results we consider a periodic array of parallel slits, Fig. 3. The porous medium consists in parallel impervious rigid plane plates of normal vector e 3 , of negligible thickness ͑i.e., Ӎ 1͒, separated by gas layers of constant thickness 2a. This medium is isotropic in the plane ͑e 1 , e 2 ͒ and impervious in direction e 3 . The one-dimensional geometry enables the exact determination ͑detailed in Appendix C͒ of the local fields and macroscopic tensors. For the sake of simplicity, the adiabatic regime is assumed. From the "cell," symmetry 1 K = 0 and the first corrector 1 V cancel out. According to Sec. III C 6 the description is, up to the second order, as follows.

Zero order:

div͑ 0 V͒ + i 0 P ␥P e = 0, 0 V = -0 K . ٌ 0 P ͑43͒
Second order ͑with adiabatic tensor 2 N͒:

Ά div͑ 2 V͒ + i 2 P ␥P e = 0, 2 V = -0 K . ٌ 2 P -2 K . ٌٌٌ 0 P - i P e 2 N . ٌ 0 P .

A. Macroscopic tensors

Using the notations,

͑recall that ␦ v = ͱ / i e ͒ u = a ␦ v , C = cosh͑u͒, S = sinh͑u͒, T = tanh͑u͒,
the zero-order dynamic permeability tensor takes the classical diagonal form ͑ 0 K i j = 0 for i j͒:

0 K 1 1 = 0 K 2 2 = ␦ v 2 ͫ 1 -T u ͬ , 0 K 3 3 = 0.
At the next order, as expected from the cell symmetry calculations show that 1 K vanishes, whereas some local fields 1 k are not zero:

∀i, j,k, 1 K i jk = 0.

After the algebra described in Appendix C, the nonvanishing components of 2 K are: ͑1͒ At low frequency, Q e 1 tends to the real value 2a 2 /5 =6K / 5 and, according to Eq. ͑41͒, the scattering perturbation is near to that of diffusion waves. ͑2͒ At high frequency, Q e 1 ϳ -i / : This leads to atypical effects, i.e., ͓cf. Eq. ͑42͔͒, a velocity dispersion varying according to ͱ , and an additional attenuation varying according to .

2 K 1 111 = 2 K 2 222 = ␦ v 4 ͫ 3 2 - 5T 2u -ͩ 1 6 - 1 u 2 ͪT 2 ͬ , 2 K 1 122 = 2 K 2 211 = ␦ v 4 ͫ 3 2 - 3T 2u - T 2 2 ͬ , 2 K 2 112 = 2 K 1 221 = ␦ v 4 ͫ -T u + ͩ 1 3 + 1 u 2 ͪT 2 ͬ ,

Oblique propagation

In oblique directions a significant anisotropic effect is induced in the whole frequency range by the magnification coefficient tan͑␤͒ 2 . Further, ͑1͒ At low frequency, Q d tends to a real constant ͑ϳK͒ and therefore diffraction will be similar to that of diffusion waves. ͑2͒ At high frequency, terms that vary according to i exactly compensate each other. Thus, Q d ϳ ͱ i, implying that Q d ͑ih d ͒ 2 varies as 2 ͱ i. Again, the consequences of diffraction differ from that of elastic waves: Both velocity dispersion and additional attenuation vary according to 2 ͱ i.

This simple case shows that the diffraction correction presents a rather complex frequency dependence. Figures 4-6 depict the strong quantitative and qualitative differences when ͑i͒ The diffraction frequency lies within the viscous, critical, or inertial regime ͑those cases are reached by varying the thickness of the gas layer͒ and ͑ii͒ the propagation is parallel to the slit or in oblique direction. It is worth noting that the atypical results at high frequency are strongly related to the one-dimensional-geometry and should not be generalized to other porous materials.

VI. CONCLUSION

To investigate the long wave scattering in air-saturated porous media, the periodic homogenization method was applied, extending the analysis up to the second-order terms. The actual enriched second-order description given by Eqs. ͑20͒, ͑30͒, and ͑31͒ accounts for the viscous, thermal, and inertia effects at the pore scale and is valid even for large concentration of solid.

In the Rayleigh domain, the improving of the continuum description implies a nonlocal behavior ͑associated with microstructural sources͒ and specific tensors depending on the microstructure and the frequency. It is shown that-up to the second order at least-the tensors associated the viscoinertial effects can be determined separately, while the thermal effects become coupled with viscous effects. A key result is that the first-order macropressure is demonstrated to vanish and that the correctors are reported at the second order. This points out the robustness of the usual continuum description, and explains in some manner the abrupt transition from the behavior in the Rayleigh domain to the behavior at higher frequencies.

The perturbation of plane wave propagation is derived from this description and interpreted as interferences between the zero-order field and the field-coherent and linearly amplified by a self-resonance mechanism-built up from the diffracted sources. The frequency dependence of the velocity's and attenuation's correction is rather complex due to the dispersion of the P 2 wave. Effects close to that observed for diffusion waves are recovered at low frequency. At high frequency an accurate discussion of these theoretical results would require numerical computation of high order tensors. An example is presented for the simple case of an array of slits. In that case, the atypical effect in high frequency range is strongly anisotropic and differs from what is observed with purely elastic waves.

Let us underline that those results are only valid in the scale separation frequency range. At higher frequencies, this approach becomes irrelevant since the scattered waves radiate in other directions than the initial wave. Another limitation is the periodicity assumption. Nevertheless, it may be inferred that this hypothesis does not modify the nature of the results for long wave in random media, provided that they can be characterized by a representative volume smaller than the wavelength. However, regarding shorter waves, strong divergences have to be expected, mostly because the diffracted field in random media loses its coherence.

APPENDIX A: HOMOGENIZATION PROCESS UP TO THE SECOND ORDER

This appendix gives the complement of the homogenization process not detailed in the main text. To save notations, the Einstein convention is used, the derivative ץf / ץy i is denoted by f ,i , the high indexes of the particular solutions are related to the forcing term, the low indexes indicate the components of the field and derivative, e.g., v = v i . e i ; div y ͑v͒

= v i,i ; ␦ i m is the Kronecker symbol ͑␦ i m =1 if i = m; =0
otherwise͒. The average on the pores domain is denoted by

͗-͘ = ͑1/⍀ f ͒͐ ⍀ f -d⍀.

The leading and first order

These two first steps are detailed in Sec. III B. The differential systems ͕S t 0 ͖ satisfied by the particular solutions 0 , is S t 0 ͑ 0 ͒ Ά i e c p 0 -⌬ y ͑ 0 ͒ = i P e T e , 0 /⌫ = 0, 0 ⍀-periodic.

Keeping in mind that P e / T e = e c p ͑1-1/␥͒, it clearly appears that 0 depends on the frequency throughout the single parameter ␦ t -2 = i e c p / . Thus, Ͻ 0 Ͼ = 0 ⌸, is a function of the thermal dimensionless frequency / t ͓ t being defined by Eq. ͑22͔͒.

As for the solutions ͑-0 k m / , 1 ␣ m ͒, they are governed by the set ͕S v 0 ͖:

S v 0 ͑ 0 k m , 1 ␣ m ͒ Ά i e 0 k i m -⌬ y ͑ 0 k i m ͒ -1 ␣ ,i m = ␦ i m , 0 k i,i m = 0, 0 k /⌫ m = 0, ͗ 1 ␣ m ͘ = 0, 0 k m , 1 ␣ m ⍀-periodic.
This form points out that ͑-0 k m , 1 ␣ m ͒ depends on the frequency throughout the single parameter ␦ v -2 = i e / . Therefore, ͗ 0 k͘ = 0 K is function of the viscous dimensionless frequency / c ͓ c being defined by Eq. ͑21͔͒.

The differential system ͕S t 1 ͖ satisfied by the particular solutions 1 m is given in the following:

S t 1 ͑ 1 m ͒ Ά i e c p 1 m -⌬ y ͑ 1 m ͒ = i P e T e 1 ␣ m + 2 0 ,m , 1 /⌫ = 0, 1 ⍀-periodic.
The presence of 1 ␣ m implies a viscothermal coupling. Then ͗ 1 ͘ = 1 ⌸ depends on both dimensionless frequencies / t and / c .

The second order

The gas flow problem ˆSv 1 ‰ The set of equations ͕S v 1 ͖ governing ͑ 1 v , 2 p͒ reads explicitly:

S v 1 Ά -i e1 v + ⌬ y ͑ 1 v͒ -ٌ y 2 p = ٌ x 1 p -2⌬ yx ͑ 0 v͒ + iٌ y 0 T e , div y ͑ 1 v͒ = -div x ͑ 0 v͒ + div x ͑͗ 0 v͒͘ + iͫ 0 -͗ 0 ͘ T e ͬ,
v /⌫ 1 = 0, 1 v, and 1 p⍀-periodic.

΅

Note the following:

͑1͒ Taking advantage of ͗div y ͑ 1 v͒͘ =0 ͑see Sec. III A͒ the average mass balance equation has been subtracted from the local mass balance. ͑2͒ The set ͕S v 1 ͖ is identical to the set ͕S v 0 ͖, except for the right-hand side of the equations acting as forcing terms induced by the solutions at previous order. When substituting Eqs. ͑17͒ and ͑18͒ of 0 , 0 v, and 1 p it appears that the forcing terms may be split into:

͑1͒ A term associated with ٌ x 1 p ˆ, which leads to an identical problem than at the preceding order. ͑2͒ Terms coming from 0 v, ͗ 0 v͘ and 1 p, associated with the second gradient of the zero-order pressure ٌ x ٌ x 0 P. ͑3͒ Terms coming from 0 and ͗ 0 ͘, associated to the zeroorder pressure 0 P.

As a consequence of the linearity, the solution is a linear combination of particular solutions associated with each of these forcing terms. Consequently the solution is in the form:

2 p͑x,y͒ = 2 p ˆ͑x͒ + 1 ␣ . ٌ x 1 p ˆ͑x͒ + 2 ␣ . . ٌ x ٌ x 0 P͑x͒ + i P e 2 0 P͑x͒, ͑A1͒ -1 v͑x,y͒ = 0 k . ٌ x 1 p ˆ͑x͒ + 1 k . . ٌ x ٌ x 0 P͑x͒ + i P e 1 n 0 P͑x͒. ͑A2͒
The solutions ͑-1 k mn / , 2 ␣ mn ͒ corresponding to unit second gradients of pressure, ٌ x ٌ x 0 P͑x͒ = e m e n , are governed by

S v 1 ͑ 1 k mn , 2 ␣ mn ͒ ϫ Ά i e 1 k i mn -⌬ y ͑ 1 k i mn ͒ -2 ␣ ,i mn = 1 ␣ m ␦ i n + 2 0 k i,n m , 1 k i,i mn = -0 k n m + ͗ 0 k n m ͘, 1 k /⌫ mn = 0;͗ 2 ␣ mn ͘ = 0,
1 k mn and 2 ␣ mn ⍀-periodic.

As previously, the single parameter is ␦ v -2 so that ͗ 1 k͘ = 1 K, only depends on / c .

The solution ͑-1 n / , 2 ͒ corresponding to a pressure, 0 P͑x͒ = P e / i, is driven by

S v 1 ͑ 1 n, 2 ͒ Ά i e 1 n i -⌬ y ͑ 1 n i ͒ -2 ,i = 0 ,i , 1 n i,i = 0 ⌸ -0 , 1 n /⌫ = 0,͗ 2 ͘ = 0,
1 n, and 2 ⍀-periodic.

The source terms brought a thermal coupling. Therefore ͗ 1 n͘ = 1 N depends on / t and / c .

The heat transfer problem ˆSt 2

‰

As for the temperature, 2 is solution of ͕S t 2 ͖:

S t 2 ͭ i e c p 2 -⌬ y ͑ 2 ͒ = 2 pi + 2⌬ yx ͑ 1 ͒ + ⌬ x ͑ 0 ͒, /⌫ 2 = 0, 2 ⍀-periodic.
Inserting Eqs. ͑17͒-͑19͒ and ͑A1͒ giving 0 , 1 , and 2 p, one notes that the forcing terms are constituted by the following:

͑1͒ Terms associated with 2 p ˆand ٌ x 1 p ˆ, which lead to identical problems than at the two preceding orders. ͑2͒ Terms associated with the second gradient of the zeroorder pressure ٌ x ٌ x 0 P. ͑3͒ Terms associated with the zero order pressure 0 P.

Again, from the linearity, the solution is a linear combination of particular solutions associated with each forcing term:

͑P e /T e ͒ 2 ͑x,y͒ = 0 2 p ˆ͑x͒ + 1 . ٌ x 1 p ˆ͑x͒ + 2 . . ٌ x ٌ x 0 P͑x͒ + i P e 2 0 P͑x͒.

͑A3͒

The temperature distributions 2 mn solutions of ͕S t 2 ͖ under second gradient of pressure ٌ x ٌ x 0 P͑x͒ = ͑P e / T e ͒e m e n are determined by

S t 2 ͑ 2 mn ͒ ϫ Ά i e c p 2 mn -⌬ y ͑ 2 mn ͒ = i P e T e 2 ␣ mn + 2 1 ,n m + 0 ␦ n m , 2 /⌫ mn = 0, 2 mn ⍀-periodic.
The temperature distribution 2 solution of ͕S t 2 ͖ under pressure 0 P͑x͒ = P e / i is determined by

S t 2 ͑ 2 ͒ Ά i e c p 2 -⌬ y ͑ 2 ͒ = i P e T e 2 ,
2 /⌫ = 0, 2 ⍀-periodic.

Note that the source terms introduce a viscous coupling in the first case, and a thermoviscous coupling in the second case. Then 2 ⌸ and 2 Z depend on / t and / c .

3. The third order: The gas flow problem ˆSv 2 ‰ Finally, the set ͕S v 2 ͖ governing 3 p and 2 v is

S v 2 Ά -i e2 v + ⌬ y ͑ 2 v͒ -ٌ y 3 p = ... =ٌ x 2 p -2⌬ yx ͑ 1 v͒ -iٌͩ y 1 p P e -ٌ y 1 T e ͪ-⌬ x ͑ 0 v͒ -iٌͩ x 0 P P e -ٌ x 0 T e ͪ div y ͑ 2 v͒ = -div x ͑ 1 v͒ + div x ͑͗ 1 v͒͘ -iͫ 1 p -͗ 1 p͘ P e - 1 -͗ 1 ͘ T e ͬ
v /⌫ 2 = 0, 2 v and 3 p ⍀-periodic.

΅

Following the same reasoning as mentioned earlier, one shows that the solution reads:

3 p͑x,y͒ = 3 p ˆ͑x͒ + 1 ␣ . ٌ x 2 p ˆ͑x͒ + 2 ␣ . . ٌ x ٌ x 1 p ˆ͑x͒ + 3 ␣ ... ٌ x ٌ x 0 P͑x͒ + i P e ͓ 2 1 p ˆ͑x͒ + 3 ٌ x 0 P͑x͔͒, -2 v͑x,y͒ = 0 k . ٌ x 2 p ˆ͑x͒ + 1 k . . ٌ x ٌ x 1 p ˆ͑x͒ + 2 k . . ٌ x ٌ x 0 P͑x͒ + i P e ͓ 1 n 1 p ˆ͑x͒ + 2 n . ٌ x 0 P͑x͔͒.
The solutions ͑-2 k mnp / , 3 ␣ mnp ͒ corresponding to unit third gradient of pressure, ٌ x ٌ x ٌ x 0 P͑x͒ = e m e n e p are governed by

S v 2 ͑ 2 k mnp , 3 ␣ mnp ͒ ϫ Ά i e 2 k i mnp -⌬ y ͑ 2 k i mnp ͒ -3 ␣ ,i mnp = 2 ␣ mn ␦ i p + 2 1 k i,p mn + 0 k i m ␦ np , 2 k i,i mnp = -1 k p mn + ͗ 1 k p mn ͘, 2 k /⌫ mnp = 0, ͗ 3 ␣ mnp ͘ = 0, 2 k mnp , 3 ␣ mnp ⍀-periodic.
Once more, the single physical parameter is ␦ v -2 so that ͗ 2 k͘ = 2 K, depends on / c only.

The solutions ͑-2 n l / , 3 l ͒ corresponding to gradients of pressure ٌ x 0 P͑x͒ = ͑P e / i͒e l are governed by

S v 2 ͑ 2 n l , 3 l ͒ ϫ Ά i e 2 n i l -⌬ y ͑ 2 n i l ͒ -3 ,i l = 1 ␣ ,i l -1 ,i l + ␦ i l -0 ␦ i l , 2 n i,i l = 1 ␣ l -1 l -͗ 1 l ͘, 2 n /⌫ l = 0; ͗ 3 l ͘ = 0,
2 n l and 3 l ⍀-periodic.

The source terms brought viscous and thermoviscous coupling, and then ͗ 2 n͘ = 2 N depends on / t and / c . The existence and uniqueness of the solutions of all these problems can be demonstrated from their variational formulation, not developed here.

APPENDIX B: TENSORS PROPERTIES

The leading idea to establish the properties of the macroscopic tensors is based on energy reciprocity. The principle is to calculate separately the virtual energies of a solution under another and reciprocally. Then, comparing the two expressions, and using the fact that the virtual kinetic and dissipated energy ͑or inertial and conduction energy for thermal problems͒ are identical in both cases, tensor identities are derived.

As a first example, take the scalar product of the momentum balance of the solution 0 k m ͓set S v 0 ͑ 0 k m , 1 ␣ m ͒, Appendix A͔ with the solution 0 k n and integrate on the pore volume. It gives

͗-1 ␣ ,i m0 k i n ͘ + i e ͗ 0 k m . 0 k n ͘ -͗⌬ y ͑ 0 k i m ͒ 0 k i n ͘ = ͗ 0 k m n ͘.
From the divergence theorem, the first term may be transformed as

͗-1 ␣ ,i m0 k i n ͘ = ͗ 1 ␣ m0 k i,i n ͘ - 1 ⍀ f ͵ ץ⍀ f 1 ␣ m0 k n . nds = 0
͓because of, ͑i͒ the zero divergence of 0 k n and ͑ii͒ the adherence and the periodicity conditions͔. Similarly, from the divergence theorem

͗⌬ y ͑ 0 k i m ͒ 0 k i n ͘ = ͗ 0 k i,jj m 0 k i n ͘ = ͗͑ 0 k i,j m 0 k i n ͒ ,j -0 k i,j m 0 k i,j n ͘ = 1 ⍀ f ͵ ץ⍀ f ٌ͑ y 0 k m . n͒ . 0 k n ds -ٌ͗ y 0 k m . . ٌ y 0 k n ͘
and the surface integral vanishes because of the adherence and periodicity. Thus finally:

i e ͗ 0 k m . 0 k n ͘ + ٌ͗ y ͑ 0 k m ͒ . . ٌ y ͑ 0 k n ͒͘ = ͗ 0 k m n ͘.
Inverting the role of 0 k m and 0 k n and comparing the two expressions lead directly to

͗ 0 k m n ͘ = ͗ 0 k n m ͘,
which demonstrate the symmetry of the dynamic permeability 0 K, ͓Auriault, 1980͔.

Skew symmetry of 1 K

Apply now this procedure to 0 k m and 1 k np . Consider first the virtual work of 0 k m under 1 k np :

͗-1 ␣ ,i m1 k i np ͘ + i e ͗ 0 k m . 1 k np ͘ -͗⌬ y ͑ 0 k i m ͒ 1 k i np ͘ = ͗ 1 k m np ͘.
Using the same argument as above, and taking into account ͑i͒ the expression of the divergence of 1 k np ͓set S v 1 ͑ 1 k mn , 1 ␣ mn ͒, Appendix A͔ and ͑ii͒ ͗ 1 ␣ m ͘ = 0, the first term becomes

͗-1 ␣ ,i m1 k i np ͘ = ͗ 1 ␣ m1 k i,i np ͘ + 0 = ͗ 1 ␣ m ͑-0 k p n + ͗ 0 k p n ͒͘͘ = -͗ 1 ␣ m0 k p n ͘
For the same reasons as previously:

-͗⌬ y ͑ 0 k i m ͒ 1 k i np ͘ = ٌ͗ y ͑ 0 k m ͒ . . ٌ y ͑ 1 k np ͒͘ and thus

-͗ 1 ␣ m0 k p n ͘ + i e ͗ 0 k m . 1 k np ͘ + ٌ͗ y ͑ 0 k m ͒ . . ٌ y ͑ 1 k np ͒͘ = ͗ 1 k m np ͘. ͑B1͒
In a second step, consider the virtual work of 1 k np under 0 k m . From the momentum balance ͓cf. set S v 1 ͑ 1 k mn , 1 ␣ mn ͒, Appendix A͔ we get

͗-2 ␣ ,i np0 k i m ͘ + i e ͗ 1 k np . 0 k m ͘ -͗⌬ y ͑ 1 k i np ͒ 0 k i m ͘ = ͗ 1 ␣ n0 k p m ͘ + 2͗ 0 k i,p n 0 k i m ͘. ͑B2͒
Here again, the first term is null because of the zero divergence of 0 k m , and -͗⌬ y ͑ 1 k i np ͒ 0 k i m ͘ may be changed into ٌ͗ y ͑ 0 k m ͒ . .ٌ y ͑ 1 k np ͒͘. Then, subtracting Eqs. ͑B1͒ and ͑B2͒ one obtains -͗ 1 ␣ m0 k p n ͘ + ͗ 1 ␣ n0 k p m ͘ + 2͗ 0 k i,p n 0 k i m ͘ = ͗ 1 k m np ͘. ͑B3͒

Note that this expression enables one to calculate ͗ 1 k m np ͘ from the solutions at previous order without calculating 1 k m np . Moreover, inverting m and n gives

͗ 1 k m np ͘ + ͗ 1 k n mp ͘ = 2͑͗ 0 k i,p n 0 k i m ͘ + ͗ 0 k i,p m 0 k i n ͒͘ = =2͗͑ 0 k n . 0 k m ͒ ,p ͘ = - 2 ⍀ f ͵ ץ⍀ f ͓ 0 k n . 0 k m ͔e p . nds = 0
͑because of periodicity and adherence conditions͒. This demonstrates the skew symmetry of 1 K relative to the two first indexes:

1 K m np = -1 K n mp .

Relation between 1 N and 1 ⌸

Similarly an expression of 1 N is derived by considering the following.

͑1͒

The virtual work of 0 k m under 1 n:

͗-1 ␣ ,i m1 n i ͘ + i e ͗ 0 k m . 1 n͘ -͗⌬ y ͑ 0 k i m ͒ 1 n i ͘ = ͗ 1 n m ͘,
which reads after transformation using the mass balance of S v 1 ͑ 1 n , 2 ͒, Appendix A:

͗ 1 ␣ m ͑ 0 ⌸ -0 ͒͘ + i e ͗ 0 k m . 1 n͘ + ٌ͗ y 0 k m . . ٌ y 1 n͘ = ͗ 1 n m ͘.
͑2͒ The virtual work of 1 n under 0 k m ͓cf. the momentum balance of S v 1 ͑ 1 n , 2 ͒, Appendix A͔:

i e ͗ 1 n . 0 k m ͘ -͗⌬ y ͑ 1 n i ͒ 0 k i m ͘ -͗ 2 ,i 0 k i m ͘ = ͗ 0 ,i 0 k i m ͘,
which reduces to i e ͗ 0 k m . 1 n͘ + ٌ͗ y 0 k m ٌ y 1 n͘ = 0.

Reporting this result in the first equality gives

͗ 1 n m ͘ = -͗ 1 ␣ m0 ͘.
Independently, an expression of 1 ⌸ is derived by considering:

͑1͒ The virtual work of 0 under 1 m , cf. set S t 0 ͑ 0 ͒, Appendix A: ͗i e c p 0 1 m ͘ -͗⌬ y ͑ 0 ͒ 1 m ͘ = i P e T e ͗ 1 m ͘.

͑2͒ The virtual work of 1 m under 0 , cf. set S t 1 ͑ 1 m ͒, Appendix A: ͗i e c p 1 m0 ͘ -͗⌬ y ͑ 1 m ͒ 0 ͘ = i P e T e ͗ 1 ␣ m0 ͘ + 2͗ 0 ,m 0 ͘.

Noting that 2͗ 0 ,m 0 ͘ = ͗͑ 0 ͒ ,m 2 ͘ = -2

⍀ f ͵ ץ⍀ f ͑ 0 ͒ 2 e p . nds = 0
͑because of the periodicity and isotherm boundary condi-tions͒, one deduces after usual transformations:

͗ 1 ␣ m0 ͘ = ͗ 1 m ͘.
Comparing with the above expression of ͗ 1 n m ͘, one obtains ͗ 1 n m ͘ = -͗ 1 m ͘, i.e., 1 N m = -1 ⌸ m .

Expression and symmetry of 2 K

Pursuing the analysis for 2 K, two identities can be established from:

͑1͒ The virtual work of 0 k m under 2 k npq and of 2 k npq under 0 k m :

͗ 2 k m npq ͘ = ͗ 0 k n . 0 k m ͘␦ pq -͗ 1 ␣ m1 k q np ͘ + ͗ 2 ␣ np0 k q m ͘ -2͗ 1 k i np0 k i,q m ͘.
͑2͒ The virtual work of 1 k np under 1 k mq :

-͗ 2 ␣ mq0 k np ͘ + ٌ͗ y ͑ 1 k mq ͒ . . ٌ y ͑ 1 k np ͒͘ + i e ͗ 1 k mq . 1 k np ͘ = ͗ 1 ␣ m1 k q np ͘ + 2͗ 1 k i np0 k i,q m ͘.

Combining these relations leads to ͗ 2 k m npq ͘ = ͗ 0 k n . 0 k m ͘␦ pq + ͗ 2 ␣ mq0 k p n ͘ + ͗ 2 ␣ np0 k q m ͘...

-ٌͩ͗ y ͑ 1 k mq ͒ . . ٌ y ͑ 1 k np ͒͘ + i e ͗ 1 k mq . 1 k np ͪ͘ , which provides the demonstration of the symmetry of 2 K relative to the two pairs of indexes:

2 K m npq = 2 K n mqp .

Expression of 2 N and 2 Z

As for 2 N, the comparison of the virtual works of 0 k m under 2 n n and reciprocally gives

2 N m n / = ͗ 2 n m n ͘ = ͗ 0 k n m ͑1 -0 ͒͘ + ͗ 1 ␣ m ͑ 1 ␣ n -1 n ͒͘.
Finally, the virtual works of 0 under 2 and inversely lead to

2 Z = ͗ 2 ͘ = ͗ 0 2 ͘.

FIG. 1 .

 1 FIG.1. Periodic cell ⍀ of porous media. ⍀ f is the pore volume filled of gas. ⌫ is the air-solid interface. = ⍀ f / ⍀ is the porosity.

  added to the zero-order terms: T T e = ͫ 0 ⌸ + i P e 2 Z ͬ P P e + 2 ⌬ͩ P P e ͪ+O͑⑀ 4 ͒.

FIG. 2 .

 2 FIG. 2. Logarithm of the ratio of diffraction and criticalfrequencies, Log͓ d / c ͔, vs porosity for media of different elementary representative volume size l. Lines from left to right corrrespond respectively to l = 1 cm, 1 mm, 100 m, 10 m, and 1 m.

FIG. 4 .

 4 FIG. 4. Diffraction effect within the viscous regime.The thickness 2a =5ϫ 10 -8 m gives a diffraction frequency d Ӎ 0.01 c , ͑f c Ӎ 100 GHz͒. The relative variation of velocity, ⌬C / C = ͑C ˜-C͒ / C, ͑plain line͒ and variation of attenuation, ⌬ = ˜-, ͑dashed line͒ are presented vs Log͓ / c ͔, for propagation parallel to the slits ͑tick line͒ and in oblique direction ␤ = /4 ͑normal line͒. Velocity and attenuation significantly increase, similarly in both directions.

Finally, tensor 2 N is diagonal ͑ 2 N i j = 0 for i j͒ and has the following components:

B. Plane wave diffraction

Consider now a plane wave propagating in a direction d inclined at an angle ␤ with the slits, i.e., d = cos͑␤͒e 1 + sin͑␤͒e 3 . In this direction, the wave number h d is given by

and the ͑adiabatic͒ coefficient Q d takes the form, see Eq. ͑38͒:

Thus the dimensionless corrective term Q d ͑ih d ͒ 2 reads:

Propagation parallel to the slits

When focusing on waves propagating parallel to the slits, the corrector reduces to

From expressions of 0 K 1 1 and 2 K 1 111 , it can be shown that, all over the whole frequency range, 2 K 1 111 ͑ 0 K 1 1 ͒ -2 is very close to a constant real value ͑varying from 7 / 5 at low frequency to 4 / 3 at high frequency͒. Consequently:

APPENDIX C: HIGH ORDER TENSOR FOR A PERIODIC ARRAY OF SLITS

Owing to the one-dimensional geometry, the fields depend locally on variable y 3 only, and the axes e 1 and e 2 play the same role. In this context, the problem ͕S v 0 ͖ governing the dynamic permeability is rewritten in the simple form:

The resolution is straightforward and gives, where

1 ␣ 3 =-y 3 , other 0 k j i and 1 ␣ i = 0. After integrating on the gas layer ͓-a + a͔ one obtains

These solutions enable one to tackle the next problem ͕S v 1 ͖, which takes the simplified form:

The resolution shows that the nonzero components read:

As the 1 k components are odd functions, they all are of zero mean value. Therefore,

The macroscopic tensors 2 K and 2 N are calculated with the expressions established in Appendix B that enable one to shun the explicit resolution of ͑S v 2 ͒:

and, within the adiabatic assumption:

From the previous results, the expressions of the nonzero components become ͓k is defined in ͑C1͔͒:

Integration on the gas layer provides the expressions given in Sec. V.