Laurent Guillot
email: guillot@ipgp.jussieu.fr

Yann Capdeville

Jean-Jacques Marigo

2D non periodic homogenization of the elastic wave equation -SH case

In the Earth, seismic waves propagate through 3-dimensional heterogeneities characterized by a large variety of scales, some of them much smaller than their minimum wavelength. Computing the wavefield in such media with the use of numerical methods, leads to high calculation costs. To lower these latter, but also to obtain a better geodynamical interpretation of tomographic images, we aim at calculating appropriate effective properties of heterogeneous and discontinuous media, by deriving convenient upscaling rules for the material properties and for the wave equation.

To progress towards this goal, we extend the successful work of Capdeville et al. (2009a), from 1-D to 2-D; basically, we first apply the so-called homogenization method -based on a two-scale asymptotic expansion of the field variables-, to model antiplane wave propagation in 2-D periodic media. These are characterized by short-scale variations of elastic properties, compared to the smallest wavelength of the wavefield. Seismograms are obtained using the 0th-order of this asymptotic expansion, plus a partial first-order correction, at a much lower computational cost. They are in excellent agreement with reference solutions calculated with spectral elements simulations, at least in the bulk of the medium. We then follow the suggestions of Capdeville et al. (2009a) to extend the homogenization of the wave equation, to 2-D nonperiodic, deterministic media.

INTRODUCTION

Are some parts of the Earth elastically anisotropic at the large scale? Even though some have performed sophisticated statistical analyses to prove that anisotropy is needed to explain data (e.g. the F-test in [START_REF] Suquet | Global anisotropic phase velocity maps for fundamental mode surface waves between 40 and 150 s[END_REF]), this is still an open question: the answer given by seismological data, often depends on the parameterization of the media that elastic wavefields have propagated through. Classically, synthetic data calculated in isotropic models with a fine spatial parameterization of elastic properties, fit real data as well as synthetic ones, calculated using anisotropic models with a coarser parameterization.

Nevertheless, seismological imaging often relies on relatively long-period data (e.g., [START_REF] Beucler | Computation of large anisotropic seismic heterogeneities (CLASH)[END_REF]), with respect to the size of the Earth heterogeneities, as can be found in the crust or even in the upper mantle [START_REF] Gudmundsson | Stochastic analysis of global traveltime data : mantle heterogeneity and random errors in the ISC data[END_REF]. But it is well-known that in that case, longwavelength wavefields naturally "uspscale" these media: what waves see, is a medium with effective properties (see Kennett (1983) or [START_REF] Chapman | Fundamentals of seismic wave propagation[END_REF] for instance). It is then of prime interest, to be able to understand the rules of this upscaling process, or, in other words, to understand how media containing small-scale heterogeneities, are "seen" by the the long-period part of an elastic wavefield. [START_REF] Backus | Long-wave elastic anisotropy produced by horizontal layering[END_REF] showed that finely-layered isotropic media are orthotropic for long-wave propagation. One of the important achievements of this article (and of the associated paper, see Capdeville et al., 2009b) will be to show that this kind of result can be generalized to more complex media, like randomly-generated ones. We shall therefore bring another interpretation to seismological data: perhaps the "real" Earth may be widely isotropic but highly heterogeneous, and the "effective" Earth, anisotropic.

As recalled by Capdeville et al. (2009b) in their introduction, numerous different kinds of methods have been used to determinate the effective properties of media with small-scale heterogeneities. In this article, we will tackle this issue using some tools of the so-called homogenization theory. Lots of studies have been devoted to this theory, either to its mathematical bases [START_REF] Murat | Calcul des variations et homogénéisation[END_REF][START_REF] Allaire | Homogenization and two-scale convergence[END_REF] or to some physical applications, often related to physical properties of composite materials [START_REF] Dumontet | Study of a boundary layer problem in elastic composite materials[END_REF][START_REF] Francfort | Homogenization and optimal bounds in linear elasticity[END_REF]. Unfortunately, most of these works were relative to periodic media, in the static case -remote from the preoccupations of a seismologist who study using wave properties, an Earth far from being a simple checkerboard. Some works were however devoted to time-depending equations (e.g Sanchez-Palencia 1980; [START_REF] Auriault | Dynamique des composites élastiques périodiques[END_REF][START_REF] Fish | Nonlocal dispersive model for wave propagation in heterogeneous media. part 1: One-dimensional case[END_REF][START_REF] Fish | Space-time multiscale model for wave propagation in heterogeneous media[END_REF][START_REF] Parnell | Dynamic homogenization in periodic fibre reinforced media. Quasi-static limit for sh waves[END_REF][START_REF] Lurie | On homogenization of activated laminates in 1D-space and time[END_REF][START_REF] Allaire | Diffractive behavior of the wave equation in periodic media: weak convergence analysis[END_REF] or to nonperiodic settings: [START_REF] Briane | Homogenization of a nonperiodic material[END_REF], who developed for elliptic equations, a theory based on a local, periodic approximation to nonperiodic materials; [START_REF] Shkoller | An approximate homogenization scheme for nonperiodic materials[END_REF], who applied Briane's theory to defects in fiber-reinforced composites; or more recently, [START_REF] Nguetseng | Homogenized structures and applications I[END_REF], who presented a quite complex mathematical treat- High-order homogenization in the dynamical, nonperiodic case, has been tackled by [START_REF] Capdeville | Second order homogenization of the elastic wave equation for non-periodic layered media[END_REF] and [START_REF] Capdeville | Shallow layer correction for spectral element like methods[END_REF] for wave propagation in stratified media. More recently, Capdeville et al. (2009a) proposed another method to understand upscaling of the wave equation, in 1D; in that paper, they suggested to apply a procedure that can be generalized to a higher dimension in space. In this article, we will apply this procedure to 2D-SH wave propagation (or antiplane elastodynamic motion). In a first part, we will recall general features of homogenization theory applied to 2D periodic media (as was done by [START_REF] Fish | Space-time multiscale model for wave propagation in heterogeneous media[END_REF] in the multidimensional case); and then we will generalize these results, to the nonperiodic case, using the heuristic procedure suggested in Capdeville et al. (2009a), or in the companion paper Capdeville et al. (2009b). We will show that order 0 homogenization theory gives results in very good agreement with solutions considered as reference ones (using the Spectral Element Method in very heterogeneous media; see e.g. [START_REF] Komatitsch | The spectral element method: an effective tool to simulate the seismic response of 2D and 3D geological structures[END_REF]).

PERIODIC CASE

Generalities

Let us consider a plane, infinite elastic surface (or a 2D finite domain surrounded by an absorbing boundary layer) whose elastic tensor and density are c 0 and ρ 0 , respectively. As in the 1D case (Capdeville et al, 2009), we consider the elastic plate to be infinite (or surrounded by an absorbing boundary layer, which is equivalent) in order to avoid the necessary treatment of boundary conditions.

This issue is postponed to a future work.

The position of a point in this space, is described by a vector x = t (x 1 , x 2 ), where t is the transposition operator.

The propagation of elastic waves in such a medium, results of the application of an external force f (x, t), that will be defined later (note that this quantity is a scalar in the case of antiplane wave propagation). The general, classical, linearized equation of motion and constitutive relation we have to solve for 2D wave propagation in an elastic medium are:

ρ 0 ∂ tt u -∇ • σ = f σ = c 0 • ∇ x u , ( 1
)
where u is the displacement, σ the incremental Cauchy stress tensor, c 0 the general stiffness tensor,

and

∇ x = t (∂ x1 , ∂ x2 ) the space gradient operator.
As we consider antiplane motion in 2D, the displacement generated by the application of a (scalar) source term, is a scalar quantity, which will be noted u. This displacement is perpendicular to the

ρ 0 ∂ tt u -∇ • σ = f σ = µ 0 • ∇ x u , (2)
where the elastic tensor c 0 reduces to

µ 0 =   µ 0 11 µ 0 12 µ 0 21 µ 0 22   (3)
defined by only 3 independent parameters, since µ 0 12 = µ 0 21 . Note that we do not employ in this article, the strain tensor in the definition of the constitutive law, as done by Capdeville et al. (2009b). This will imply a slightly different procedure to retrieve the so-called homogenized elastic coefficients, for we shall impose the symmetry of µ 0 .

The index notation will be extensively used throughout this article. Any tensor A of order p has components, in the index notation:

A i 1 i 2 ..ip .

Set up of the homogenization problem -Periodic case

Let us assume in this first section, that the physical properties of the surface, are -periodic, that is ρ 0 (x + ) = ρ 0 (x) and µ 0 (x + ) = µ 0 (x), for any x (see Fig. 1).

Let us also consider that the time dependence of the source term f (equation ( 2)) is characterized by a corner frequency f c in the spectral domain. This hypothesis underlies the existence of a minimum wavelength λ m for the propagating wavefield, away from the source.

In this section dedicated to antiplane wave propagation in periodic settings, will be assumed that

ε = λ m << 1 , (4)
where is the characteristic length (that can always be found) of the periodic pattern as seen by the propagating wavefield, at a given point.

Clearly, heterogeneities in the medium are considered to be much smaller than the minimum wavelength of the propagating wavefield -or, in other words, the scale of heterogeneities is small, in comparison with the scale of the wavefield.

Given λ m , the minimum wavelength of the wavefield in the direction of the wave-vector, we can define a sequence of wave propagation problems, by varying the initial periodicity , each problem being indexed by a specific ε. This approach is classical in the homogenization theory, the goal being then to find an asymptotic solution (u, σ) to the wave equation when ε tends towards 0 [START_REF] Sanchez-Palencia | Non homogeneous media and vibration theory[END_REF]. According to this point of view, we can rewrite equations (2) indexing each physical quantity by ε, at the exception of the source (which is in fact not exactly correct for a point source, as underlined in Capdeville et al. (2009a), for instance):

ρ ε ∂ tt u ε -∇ • σ ε = f σ ε = µ ε • ∇ x u ε .
(5)

Our purpose is to study large-scale wave propagation with respect to the scale of heterogeneities. A clever trick to explicitly take small-scale heterogeneities into account, is to re-scale them, in order for them to artificially have the same scale as the wavefield. A space variable which varies rapidly is then introduced:

y = t (y 1 , y 2 ) = x ε . ( 6
)
The variable y is usually called the microscopic variable and x, is the macroscopic one. It is easy to see that any change in y induces a very small change in x, when ε → 0. This observation leads to the idea that scales can be separated, and therefore, that x and y can be considered as being independent variables. A direct consequence of this, is the redefinition of the gradient operator:

∇ x → ∇ x + 1 ε ∇ y . ( 7
)
As underlined in the introduction, simple physical considerations (see for instance [START_REF] Backus | Long-wave elastic anisotropy produced by horizontal layering[END_REF] or [START_REF] Chapman | Fundamentals of seismic wave propagation[END_REF]) lead us to think that a wavefield "sees" small heterogeneities in an effective way at a large scale; but we can expect it should also be locally sensitive to rapid variations of elastic properties.

The homogenization theory is a tool that allows to catch both effects, by seeking solutions to the wave equations ( 5) under the form of asymptotic expansions in ε:

u ε (x, t) = i≥0 ε i u i (x, x/ε, t) = i≥0 ε i u i (x, y, t) , σ ε (x, t) = i≥-1 ε i σ i (x, x/ε, t) = i≥-1 ε i σ i (x, y, t) . (8)
In these equations, u i and σ i depend on both space variables x and y, and are chosen to be λ mperiodic in y. Note that the expansion for the stress starts at i = -1 because of the constitutive relation between u and σ, and the redefinition of the gradient operator in (7).

What remains to do before solving the wave equation, given the ansatz (8), is to define the y-dependence of material properties:

ρ(y) = ρ ε (εy) , µ(y) = µ ε (εy) . (9)
These quantities can be called, the unit cell density and stiffness tensor. They are of course, ε-independent and λ m -periodic.

The external source term f is assumed to be ε and y independent. This assumption is not obvious for a punctual source. The idea behind this assumption is to forget about the source during the asymptotic development and then to reintroduce it using energy considerations. This issue is discussed with more details in Capdeville et al. (2009b).

Finally, we define the so-called cell average (over the unit cell Y = [0, 1 /ε] × [0, 2 /ε], |Y| being its surface), for any (tensorial) function g(x, y) which is λ m -periodic in y:

g (x) = 1 |Y| Y g(x, y)dy . (10)
For any tensorial quantity g(x, y) of order q, λ m -periodic in y, it can easily be shown that:

∂ i g ii 1 ...i q-1 = 0 . ( 11
)
For any couple of tensorial functions g(x, y) and h(x, y), of order p and q, respectively, λ m -periodic in y, is verified the following property

∂ y i (g ii 1 ..i p-1 h j 1 ..jq ) = 0 , (12)
and then: In the following, the time dependence t is dropped to ease the notations.

(∂ y i (g ii 1 ..i p-1 )h j 1 ..jq = -g ii 1 ..i p-1 (∂ y i h j 1 ..jq ) . (13)
Introducing expansions (8) in equations (5), using (7) and identifying term by term in ε i , we readily obtain:

ρ∂ tt u i -∇ x • σ i -∇ y • σ i+1 = f δ i,0 , (14)
σ i = µ • (∇ x u i + ∇ y u i+1 ) . ( 15
)
These last equations have to be solved for each i.

• Equations ( 14) for i = -2 and (15) for i = -1 give

∇ y • σ -1 = 0 , σ -1 = µ • (∇ y • u 0 ) . (16)
These equations imply that

∇ y • (µ • ∇ y u 0 ) = ∂ y i (µ ij ∂ y j u 0 ) = 0 . (17)
Multiplying the last equation by u 0 , integrating over the unit cell, and then by part, taking account of the periodicity of u 0 and µ, leads to

Y µ ij (∂ y i u 0 )(∂ y j u 0 ) dy = 0 . ( 18
)
As µ is positive-definite, the unique solution to the above equation is u 0 (x, y) = u 0 (x). We therefore have

σ -1 = 0 , (19)
u 0 = u 0 (x) = u 0 . ( 20
)
This last equality simply implies that the order 0 solution in displacement is independent on the fast variable y, and does not oscillate with x/ε. This is a major result that confirms the well known fact that the displacement field is poorly sensitive to scales much smaller than its own scale (see [START_REF] Chapman | Fundamentals of seismic wave propagation[END_REF]). Anecdotically, this fact justifies the name of the theory: u 0 is an homogenized solution.

• Equations ( 14) for i = -1 and (15) for i = 0 give or rewritten using the summation convention:

∇ y • σ 0 = 0 , (21)
σ 0 = µ • (∇ y u 1 + ∇ x u 0 ) , (22)
∂ y i σ 0 i = 0 , (23)
σ 0 i = µ ij (∂ y j u 1 + ∂ x j u 0 ) . ( 24
)
Because the action of the divergence operator, equation ( 21) does not imply that σ 0 (x, y) = σ 0 (x) = σ 0 (x), as it it the case in 1D (Capdeville et al., 2009a), where the divergence becomes a simple gradient operator. Nevertheless we will next obtain a simple expression for σ 0 .

Both equations lead to this next equality:

∇ y • µ • ∇ y u 1 = -∇ y • µ • ∇ x u 0 , (25)
which should be verified, whatever the gradient of u 0 be. Because of this and of the linearity of the previous equation, we can look for a solution of the form

u 1 (x, y) = u 1 (x) + χ 1 (y) • ∇ x u 0 (x) . ( 26
)
The vectorial quantity χ 1 (y) (with components χ 1 k (y), k = 1, 2), is called the first order periodic corrector; it is λ m -periodic. Introducing (26) into (25), we obtain the equations of the so-called cell problems:

∇ y • µ • I + ∇ y χ 1 k = ∂ y i µ ij δ jk + ∂ j χ 1 k = 0 , (27)
that lead to find the first-order corrector components. To enforce the uniqueness of the solution for (27), we impose χ 1 k = 0. On the contrary to the 1D case where an analytical solution to the previous equation does exist, we can not solve it in 2D, but numerically. The special case of a stratified (and nonperiodic) model, which can be considered as a 1D medium, is treated in Appendix C.

In prevision to the next section, where homogenization in nonperiodic settings will be tackled, let us rewrite (27) under the equivalent form

∇ y • H = ∇ y • (µ • G) = 0 ( 28
)
where H and G = I + ∇ y χ 1 (I being the identity tensor) are second order tensors whose components have the dimension of a stress and of the gradient of a displacement, respectively. This observation will be useful in the section dedicated to nonperiodic homogenization.

Taking the cell average of (24) and using the ansatz (26), we obtain the order 0 constitutive relation: where µ * is the order 0 homogenized elastic tensor whose components are:

σ 0 = µ * • ∇ x u 0 , (29)
µ * ij = µ ik (δ jk + ∂ y j χ 1 k ) . (30)
Obviously, and because of ( 22), ( 29) and the definition of G, we obtain

µ * = H . ( 31
)
At the first sight, it is not obvious to physically interpret the expression for the effective elastic constants (30). As noticed by [START_REF] Papanicolaou | Boundary value problems with rapidly oscillating random coefficients[END_REF], they are the sum of their own average and of a correction term, which is the average of elementary stresses associated to displacements equal to the first-order corrector's components. This observation may lead to a practical method for the determination of effective elastic constants, as shown by [START_REF] Suquet | Global anisotropic phase velocity maps for fundamental mode surface waves between 40 and 150 s[END_REF] in the static case, and [START_REF] Grechka | Effective media: A forward modeling view[END_REF] in the dynamical one. We will use this kind of observation in the next section (dedicated to the nonperiodic case), in order to determinate conveniently, both the y-dependence of the stiffness tensor, and the effective, homogenized tensor associated to it.

• Equations ( 14) for i = 0 and ( 15) for i = 1 give

ρ∂ tt u 0 -∇ x • σ 0 -∇ y • σ 1 = f , (32)
σ 1 = µ • (∇ y u 2 + ∇ x u 1 ) . (33)
Applying the cell average on (32), using the property (13), the fact that u 0 does not depend on y, and the expression (29) for the average of σ 0 , we obtain the order 0 wave equation:

ρ * ∂ tt u 0 -∇ x • σ 0 = f σ 0 = µ * • ∇ x u 0 , (34)
where ρ * = ρ is the effective density, whose expression is the same as in 1D (Capdeville et al., 2009a). We can notice that this equation is simply the classical wave equation that can be solved using classical techniques. Because ρ * and µ * are constant, solving the order 0 wave equation is a much simpler task than in the original medium; no numerical difficulty arises, that is related to the rapid spatial variations of the physical properties in the plane.

Once u 0 is found, and because the cell problem ( 27) has been solved to determinate the effective elastic tensor, the first order correction χ 1 (x/ε) • ∇ x u 0 (x) can then be computed.

Obtaining the complete order 1 solution u 1 (x, y) using (26) (that means, finding the remaining term

u 1 (x)
) is much more difficult than in 1D, as we shall show now.

Subtracting ( 34) to (32) we obtain, which, together with (33), gives

∇ y • σ 1 = (ρ -ρ )∂ tt u 0 , (35)
∇ y • µ • ∇ y u 2 = -∇ y • µ • ∇ x u 1 + (ρ -ρ )∂ tt u 0 . (36)
Using ( 26), we easily obtain the following, expanded equation:

∂ y i µ ij ∂ y j u 2 = -∂ y i µ ij ∂ x j u 1 -∂ y i µ ij χ 1 k ∂ 2 x j x k u 0 + (ρ -ρ )∂ tt u 0 . ( 37
)
Exactly as was done in (26), using the linearity of the last equation, we can separate variables and look for a solution of the following form (in index notation)

u 2 (x, y) = u 2 (x) + χ 1 k (y)∂ x k u 1 (x) + χ 2 lm (y)∂ 2 x l xm u 0 (x) + χ ρ (y)∂ tt u 0 . ( 38
)
The components of the second-order tensor χ 2 , and of χ ρ , are solutions of the following partial differential equations:

∂ y i µ ij δ jl χ 1 k + ∂ y j χ 2 lk ) = 0 , (39)
∂ y i µ ij ∂ y j χ ρ = ρ -ρ , (40)
where χ 2 lk and χ ρ are λ m -periodic. To ensure the uniqueness of the solutions, we impose χ 2 lk = χ ρ = 0, whatever the couple (l, k) be. Note that χ 2 is obviously symmetric: 38) into (33) and averaging the expression over the unit cell, we obtain the contribution at order 1 to the constitutive relation:

χ 2 ij = χ 2 ji . Introducing (
σ 1 i = µ * ij ∂ x j u 1 + µ 1 * ilk ∂ 2 x k x l u 0 + µ ρ * i ∂ tt u 0 (41)
where

µ 1 * ilk = µ ij (δ jl χ 1 k + ∂ y j χ 2 lk ) and (42)
µ ρ * i = µ ij ∂ y j χ ρ (43)
Finally taking the average of equation ( 14) for i = 1 and using (41), gives the contribution of order 1 to the wave equation, in index notation:

ρ ∂ tt u 1 + ρχ 1 i ∂ x i ∂ tt u 0 -∂ x k σ 1 k = 0 σ 1 i = µ * ij ∂ x j u 1 + µ 1 * ilk ∂ 2 x k x l u 0 + µ ρ * i ∂ tt u 0 . ( 44
)
As was done by Capdeville et al. (2009a) for the 1D case, we can show that the previous equations can be recast, for µ ρ * i = ρχ 1 i (see appendix B):

ρ ∂ tt u 1 + µ ρ * i ∂ x i ∂ tt u 0 -∂ x k σ 1 k = 0 σ 1 i = µ * ij ∂ x j u 1 + µ 1 * ilk ∂ 2 x k x l u 0 + µ ρ * i ∂ tt u 0 . (45)
Unfortunately, the µ 1 * ilk coefficients are not null in general, therefore we can not further simplify the order 1 equations, to a form similar to that obtained for 1D wave propagation. However, as noted by [START_REF] Boutin | Microstructural effects in elastic composites[END_REF] -and this can easily verified numerically -, these coefficients are quite small in comparison to the order 0 homogenized coefficients µ * ij , and the associated term could be considered as negligible. Furthermore, they are identically null when homogenized elastic properties in the unit cell, are isotropic [START_REF] Boutin | Microstructural effects in elastic composites[END_REF]; this is the special case of a heterogeneous but macroscopically isotropic material.

Though we may solve up to a higher order (see for instance, [START_REF] Fish | Space-time multiscale model for wave propagation in heterogeneous media[END_REF] for a 2D periodic case), we stop here the expansion.

Practical resolution

Solving the homogenized wave equations derived above can be done in series, or by combining successive orders together, using one of different kinds of wave propagation solvers: normal-mode summation [START_REF] Capdeville | Second order homogenization of the elastic wave equation for non-periodic layered media[END_REF]; finite element methods [START_REF] Fish | Space-time multiscale model for wave propagation in heterogeneous media[END_REF]; or, as in Capdeville et al. (2009a) and this work, the Spectral Element Method (SEM).

In the most general case, we then want to solve

ûε,p (x) = u 0 (x) + ε u 1 (x) + ... + ε p u p (x) , (46) σε
,p (x) = σ 0 (x) + ε σ 1 (x) + ... + ε p σ p (x) , (47)
where the bracketed terms are the different homogenized fields, that we could calculate as in the previous section. Once these effective fields known, we could find the complete ones, ûε,p and σε,p , by applying a high-order corrector operator to it, as suggested by Capdeville et al. (2009a) in 1D -and it can be shown that the following approximations are then verified:

ûε,p = u ε + O(ε p ) , σε,p = σ ε + O(ε p ) . (48)
As noticed previously, the order 1 effective equation of motion and constitutive relation, can not be reduced in order to look like the same equations at order 0, unless µ 1 * is the null tensor. This implies that in the most general case, the obtention of u 1 (x) and σ 1 (x) involves not so light modifications in the SEM code. We could of course, proceed to homogenization to order 1, starting from a microscopically (at the scale of the unit cell) isotropic medium, and obtain a partial order 2 solution;

we have nevertheless not lead this work, the case of a isotropic and periodic medium, being irrelevant in geophysics.

As we shall see in the section 3.5 dedicated to some examples of wave propagation in a periodic setting, determining only a partial first-order homogenized solution seems enough to obtain quite accurate (2009b), we will only try to solve

ûε,1/2 (x) = ûε,0 (x) + χ 1 (x/ε) • ∇ x ûε,0 (x) , (49)
where the 1/2 superscript means "partial order 1". Of course, for it is only a partial order 1 solution, we do not have in general

u ε (x) = ûε,1/2 (x) + O(ε 2 ) , (50)
on the contrary of the 1D case (Capdeville et al., 2009a), unless u ε,1 is very small (or null).

External source term

In seismology, we generally consider that the source dimension is much smaller than the smallest wavelength of the (far) wavefield, and that a point source (located around a given x 0 ) is therefore a good approximation. This localization of the source leads to its mathematical formulation:

f (x, t) = δ(x -x 0 )g(t).
We claimed earlier in this article, that the external force should not depend on the ε-parameter. In general, as shown in Capdeville et al. (2009a), this is not true, because interactions of a (ideal) point source with its surrounding microscopic structure, must be accounted for.

In this article, as the asymptotic, homogenized solution we look for, is expanded at its sole, truncated order 1, and that only collocated forces located in homogeneous source regions, will be considered, the corrective source term of order 1 is not needed. Nevertheless the reader should keep in mind that this is an oversimplification, and that the procedure suggested in Capdeville et al. (2009a) or Capdeville et al. (2009b) -and based on energetic considerations-is necessary to treat source effects in a correct manner.

Homogenization in a periodic setting: an example

Let us turn now to a practical application of the homogenization theory in a periodic setting. A numerical experiment is led in a rectangular plane (see Fig. 2) of size 15 × 20 km 2 , surrounded by absorbing boundary layers (which is not important here, for the sole balistic arrival is of interest, the coda waves being non-existent). The microscopic structure is that of a stretched checkerboard, as shown in Fig. 1, with a horizontal periodicity of 120 m, and a vertical one of 200 m. Elastic properties and density in the elements of the periodic cell, are ±50% around a mean value of 60 GPa for µ 11 and µ 22 (µ 12 being null), and 2800 kg m -3 for ρ.

To obtain values for the first-order corrector components and their spatial derivatives, the cell problem (27) is solved using periodic boundary conditions. The differential equations are solved with a finite element method based on the same mesh and quadrature than the one that will be used to solve the wave equation. This then allows to compute the homogenized stiffness tensor µ * and density ρ * , and the 0th-order homogenized wave equation,

ρ ∂ tt u -∇ x • (µ * • ∇ x u) = f . ( 51
)
using the spectral element method (SEM, see for instance [START_REF] Komatitsch | The spectral element method: an effective tool to simulate the seismic response of 2D and 3D geological structures[END_REF]).

The collocated source is located at S = (x S 1 , x S 2 ) = (5 km, 10 km). Its time evolution is described by a Ricker with a central time shift of 0.4 s and a central frequency of 5 Hz (to which is associated a corner frequency of about 12.5 Hz). A receiver is located at R = (x R 1 , x R 2 ) = (10 km, 10 km), therefore, at the same vertical component as the source, which means, that the periodicity seen at this receiver point, is of 120 m. Considering the physical properties of the medium (with a mimimum wavelength of around 370 m), and the cut-off frequency of the source, the value of ε in the far field, is around 0.32. In order to properly compare the homogeneous solution with a reference one, calculated in the reference checkerboarding box, both solutions are computed with the same mesh (each element of this latter, corresponding to an element of the checkerboard), and the same time step (10 -3 s).

The seismogrames resulting of these numerical simulations are reported in Fig. 3, in which are shown the reference solution (grey line), the order 0 homogenized solution (black line) and a solution obtained when taking the arithmetic average of both density and elastic constants (dashed line). Clearly, this last solution (which could be seen as a "natural averaging" solution) is not in phase with the reference one, whereas there is an excellent agreement between this latter, and the order 0 homogenized seismogram. Snapshots along a line (L) (see Fig. 2) between the source and the receiver, are shown in Fig. 4.a for each of the simulations, at t =0.4 s. Once again, the agreement is excellent between reference and 0th-order homogenized solutions, and very poor with the "natural filtering" solution.

The residual between the reference solution and the order 0 homogenized solution u ε (x, t) -û0 (x, t), at t =0.4 s, is plotted in Fig. 4.b. The error amplitude is lower than the percent and contains fast variations.

Once done the 0th order homogenized simulation (and therefore, once known displacement gradients), we can compute the incomplete order 1 solution (49), again, at any time after the source excitation.

In Fig. 4.c is shown the partial order 1 residual u ε (x, t) -û1/2 (x, t). Comparing Fig. 4.b and Fig. 4.c, it appears that the partial order 1 periodic correction removes most of the fast variations present in the order 0 residual. The remaining fast variations are due to the neglected higher-order terms, in the expansion (48). The smooth remaining residual could be due to u 1 term (and also, to higher-order ones in 47) that is (are) not computed. In order to check the order of that this smooth remaining residual, we perform a test similar (at t =1.2 s) to the previous one, but with ε = 0.16 (which corresponds to a periodicity divided by a factor 2). Surprisingly, there is approximately a factor 4 in amplitude between both (see Fig. 4.d), which indicates a residual of order ε 2 , almost. That means that the u 1 term is quite close to being null; probably because the coefficients of µ 1 * are quite small in (45) -homogenized elastic properties in the unit cell, are close to being isotropic, which is in fact the case -and that (45) can be rearranged as in Capdeville et al. (2009a) ( u 1 being equal to 0 in 1D).

NON PERIODIC CASE

Let us now turn to the case of interest in geophysics: the Earth material properties being not spatially arranged periodically, we give up this drastic hypothesis and consider more complex patterns for heterogeneities -more precisely: we do not make any assumption on the spatial variability of density and elastic coefficients. After a theoretical treatment, we will apply it to the specific case of a plane in which properties are randomly generated around a constant mean value, as in Capdeville et al. (2009b), surrounded by a strip of progressively constant physical properties (to avoid any spurious reflection) and a PML (Perfectly Matched Layers, e.g. [START_REF] Festa | The newmark scheme as velocity-stress time-staggering: an efficient implementation for spectral element simulations of elastodynamics[END_REF]) layer (to avoid the treatment of boundary layers in our homogenization procedure), as shown in Fig. 5.a. Note that this problem we shall tackle, is not that of the homogenization of random structures as studied by [START_REF] Papanicolaou | Boundary value problems with rapidly oscillating random coefficients[END_REF]; in our example, the properties will be spatially known and unique.

Obviously we will still assume a minimum wavelength λ m (or a maximum wave-number k m = 1/λ m ) for the anti-plane wavefield u, far enough from the source. Therefore, to some sense, heterogeneities in the plane, whose the size is much smaller than λ m , should have a little, effective influence on the wavefield u and an homogenization procedure may be hopefully performed. As recalled by Aki & Richards (1980), there are 3 different wave propagation regimes (waves in a smoothly-varying body, coda waves, and the homogenized part of a wavefield) depending on the ratio of the wavefield characteristic scale, to the one of the heterogeneities.

But when properties are not periodic in space: what is the characteristic scale of the heterogeneities?

The whole difficulty we are faced to, is therefore to find a clear spatial scale delimitation (to apply an homogenization procedure), in order to catch wavefield properties in each of these regimes. In 1D, -b: order 0 residual, u ε (x, t) -û0 (x, t).

-c: partial order 1 residual, u ε (x, t) -û1/2,ε (x, t).

-d: grey line: partial order 1 residual for ε = 0.32. Black line: partial order 1 residual for ε = 0.16 with an amplitude multiplied by 4, at t =1.2 s . [START_REF] Capdeville | Second order homogenization of the elastic wave equation for non-periodic layered media[END_REF] suggested to apply a filter to the recast elastic operator in the frequency domain, or to some physical quantities known a priori (Capdeville et al., 2009a); this filtering operation allowed them to obtain very accurate "homogenized" solutions, when compared to reference ones. Unfortunately, in a most general case, we do not know to which physical quantities this filtering procedure has to be applied. We will suggest one way to proceed (as done in after having recalled some filtering notions, and having set up a heuristic homogenization procedure to wave propagation in 2D non-periodic domains.

Basic notions on spatial filtering

In the theoretical development that will follow, we shall desire to separate low from high wavenumbers k = 1/λ (the norm of the wavenumber vector k) of a spatial distribution of any given quantity g(x) (which can be a tensor of any order), around a given wavenumber k 0 . Note that the choice of k 0 will be determinant to accurately described waveform properties in each propagation regime.

To that purpose we shall introduce a low-pass space filter operator that takes the form, for any function g:

F k 0 (g) (x) = R g(x )w k 0 (x -x )dx , (52)
where w k 0 is a wavelet that would ideally be defined in the spectral domain as

wk 0 (k) =    1 for k ≤ k 0 ; 0 for k > k 0 , ( 53
)
where w is just the Fourier transform of w.

In practice, in order to have a wavelet w k 0 for which a compact support is a good approximation, we do not use a filter with such a sharp cutoff but one defined by a smoother transition from 1 to 0 around k 0 . There are many ways to design a filter with such a property. The one we shall use is characterized as:

w(k) =          1 for k ≤ k min ; 1 2 1 + cos π |k|-k min kmax-k min for |k| ∈]k min , k max [ ; 0 for |k| ≥ k max . , (54)
where k min and k max are two real numbers around 1, defining for the low pass filter the tapering zone from 1 to 0. Just notice that the following property is verified: R 2 w(x)dx = 1, and that we can easily define w k 0 (x) = k 0 w(xk 0 ), the same but contracted (if k 0 > 1) wavelet of corner spatial frequency k 0 , that also verifies: R 2 w k 0 (x)dx = 1. The choice of k min and k max is ad hoc and left to the user, this latter being aware of these 2 limitations: the perfectly sharp cutoff (k min = k max = 1) is characterized by an infinite support in the space domain, that can not be truncated with a convenient accuracy; and any other choice, for which this truncation can be safely applied, does not present the interesting property of a perfectly clear separation of scales for a given quantity g.

Set up of the homogenization problem in the nonperiodic case

In 1D, Capdeville et al. (2009a) applied the previous filtering operation to physical quantities that were explicitly known a priori. For wave propagation in higher dimensions, we do not have access to this kind of information, and do not know how to separate scales for density and elastic constants, in order to proceed to an homogenization of these quantities and of the wave equation. In other words, we do not know, for a given distribution of material properties, how to construct the x and y contributions of ρ and µ, from ρ 0 and µ 0 . To that purpose, we then present here, an original but heuristic procedure.

Classically a small parameter ε is introduced to solve the so-called two-scale homogenization problems:

ε = λ λ m , ( 55
)
where λ is a spatial wavelength, upon which will depend our asymptotic expansion (see below). For a periodic medium, λ would be the local length defining the periodicity of the model (called in the first section of this article). In the non periodic case, it is not as obvious (because of the absence of an isolated characteristic length), and the introduction of another parameter is required

ε 0 = λ 0 λ m . ( 56
)
λ 0 is the length below which a wavelength can be considered as belonging to the small scale (microscopic) domain, and reciprocally for the large (or macroscopic) scale. This λ 0 -parameter can be arbritrarily chosen; nevertheless it makes sense to assume that the wavefield does interact with heterogeneities whose size is smaller than λ m . Therefore, picking a ε 0 << 1, which means considering as microscopic, heterogeneities whose size is much smaller than the minimum wavelength, should be a good guess.

Exactly as done in the periodic case, we can now introduce the fast space variable y = x/ε, treat x and y as independent space variables when ε tends towards 0, and redefine the gradient operator as in (7).

We then introduce a wavelet w m (y) = w km (y) where w km is a low-pass filter as defined in 3.1 and

k m = 1/λ m .
We obviously assume that we can design it in such a way that its support in the space domain is contained in an interval [-αλ m , +αλ m ] 2 , where α is a positive real number.

Let Y 0 = [-βλ m , βλ m ] 2 be a square of R 2 , with β a positive real number much larger than α, and Y x the same square, but translated by a vector x/ε 0 . We define T = {g(x, y) : R 2 × Y 0 → R , Y 0 -periodic in y} the set of functions defined in y on Y 0 and extended to R 2 by periodicity. Associated to the wavelet w m , can be defined a filtering operator for any function g ∈ T : Finally let V be the set of functions g(x, y) such that, for a given x, the y part of g is periodic and contains only spacial frequencies higher than k m , plus a constant value in y:

F (g) (x, y) = R 2 g(x, y )w m (y -y )dy . (57)
V = {g ∈ T /F (g) (x, y) = g (x)} , (58)
where, similarly to the periodic case, the cell average of g over the newly defined periodic cell is

g (x) = 1 |Y 0 | Y 0 g(x, y)dy . ( 59
)
Of course, as in our procedure to treat wave propagation in nonperiodic media, the periodicity condition is kept, properties ( 11) and ( 13) are still valid. Furthermore, it is easy to show that, for any (tensorial) function g in V, the (spatial) partial derivatives of its components are also in V, and finally that:

∀h ∈ V with h = 0 and ∇g = h ⇒ g lies in V . (60)
In this section and the next one, we proceed in the same way as in Capdeville et al. (2009a). We first assume that we have been able to define (ρ ε 0 (x, y), µ ε 0 (x, y)) that verify

ρ ε 0 (x, x/ε 0 ) = ρ 0 (x) µ ε 0 (x, x/ε 0 ) = µ 0 (x) (61)
and that set up a sequence (as in the periodic case) of models indexed by ε

ρ ε 0 ,ε (x) ≡ ρ ε 0 (x, x ε ) , µ ε 0 ,ε (x) ≡ µ ε 0 (x, x ε ) . ( 62
)
We also guess that, with such a set of parameters, a solution to the homogenization problem described below exists. This is by far not obvious. The construction of such (ρ ε 0 (x, y) and µ ε 0 (x, y)) from (ρ 0 (x) and µ 0 (x)), which is therefore a critical issue, is left for section 3.4.

As in the periodic case (and after recasting of the wave equation, see appendix A), we then look for the solutions of the wave equation and constitutive relation

ρ ε 0 ,ε ∂ tt u ε 0 ,ε -∇ • σ ε 0 ,ε = f , σ ε 0 ,ε = µ ε 0 ,ε • ∇u ε 0 ,ε . (63)
A solution to the equations ( 63) is again sought as an asymptotic expansion in ε: with the additional hypothesis -in this specific, nonperiodic case -, that u ε 0 ,i and σ ε 0 ,i must belong to the space V. Introducing the expansions (64) in the wave equations ( 63) and using ( 7) we obtain the system of differential equations

u ε 0 ,ε (x, t) = ∞ i=0 ε i u ε 0 ,i (x, x/ε, t) = ∞ i=0 ε i u ε 0 ,i (x, y, t) σ ε 0 ,ε (x, t) = ∞ i=-1 ε i σ ε 0 ,i (x, x/ε, t) = ∞ i=-1 ε i σ ε 0 ,i (x, y, t) , (64)
ρ ε 0 ∂ tt u ε 0 ,i -∇ x • σ ε 0 ,i -∇ y • σ ε 0 ,i+1 = f δ i,0 (65)
σ ε 0 ,i = µ ε 0 • ∇u ε 0 ,i + ∇u ε 0 ,i+1 , (66)
which need to be solved for each i, up to a given i 0 . As in the periodic case (section 2.4), we will restrict this resolution to order 1/2, that means: i 0 = 0, plus the first-order periodic correction.

Resolution of the two scale homogenization problem

We follow the same procedure as for the periodic case. We work at ε 0 . Because the y periodicity is kept in V, the resolution of the homogenized equations is almost the same as in the periodic case. The major difference is that each physical quantity we look for, depends on ε 0 and on x.

• As for the periodic case, equations ( 65) for i = -2 and ( 66) for i = -1 give σ ε 0 ,-1 = 0 and u ε 0 ,0 = u ε 0 ,0 (x).

• Equations ( 65) for i = -1 and ( 66) for i = 0 give

∇ y • σ ε 0 ,0 = 0 , (67)
σ ε 0 ,0 = µ ε 0 • (∇ y u ε 0 ,1 + ∇ x u ε 0 ,0 ) . ( 68
)
Both previous equations lead to this next equality:

∇ y • µ ε 0 • ∇ y u ε 0 ,1 = -∇ y • µ ε 0 • ∇ x u ε 0 ,0 , (69)
Similarly to the periodic case, we can look for a solution of the form

u ε 0 ,1 (x, y) = u ε 0 ,1 (x) + χ ε 0 ,1 (x, y) • ∇ x u ε 0 ,0 (x) . ( 70
)
As u ε 0 ,1 must belong to V, and because u ε 0 ,0 is y-independent, the first order corrector χ ε 0 ,1 (x, y) (a vector of dimension 2) must also be in V. We will see in section 3.4, how to proceed to obtain a first-order corrector that has this property. Introducing ( 70) into (69), we obtain the equations of the cell problems:

∇ y • µ ε 0 • ∇ y y k + ∇χ ε 0 ,1 k = 0 . ( 71
)
To enforce the uniqueness of the solution, we again impose χ ε 0 ,1 k = 0.

As in the periodic case, we can rewrite (71) under the equivalent form with H ε 0 = µ ε 0 • G ε 0 , and G ε 0 = I + ∇χ ε 0 ,1 , where I is the unit tensor.

∇ y • H ε 0 = ∇ y • (µ ε 0 • G ε 0 ) = 0 , (72)
Taking the cell average of (68) and using the ansatz ( 70), together with the guess that the first-order corrector can be found in V, we obtain the order 0 homogenized constitutive relation:

σ ε 0 ,0 = µ ε 0 * • ∇ x u ε 0 ,0 , (73)
where µ * ε 0 is the order 0 homogenized elastic tensor whose components again are (as in the periodic case):

µ * ε 0 ij = µ ε 0 ik (δ jk + ∂ y j χ ε 0 ,1 k ) . ( 74
)
As in the periodic case, it is easy yo show the following equality

µ * ε 0 = H ε 0 , ( 75
)
and also that

σ ε 0 ,0 = H ε 0 • ∇u ε 0 ,0 , (76)
Obviously, because we look for σ ε 0 ,0 belonging to this space, H ε 0 mus lie inV. Furthermore, as it is needed -to obtain the desired asymptotic solution of our problem-that the first-order corrector also be in this space, G ε 0 , and because it is its gradient, must lie in V. These remarks will be fundamental, in our quest of the construction of a suitable µ ε 0 (x, y) in section 3.4.

• Finally, equation ( 65) for i = 0 gives

ρ ε 0 ∂ tt u ε 0 ,0 -∇ x • σ ε 0 ,0 -∇ y • σ ε 0 ,1 = f . ( 77
)
To be able to obtain σ ε 0 ,1 in V, the last equation implies that ρ ε 0 must lie in V. Taking the average of this equation, together with (68) leads to the order 0 homogenized wave equation

ρ ε 0 ∂ tt u ε 0 ,0 -∇ x • σ ε 0 ,0 = f (78) σ ε 0 ,0 = µ * ε 0 • ∇ x u ε 0 ,0 . ( 79
)
As we have seen for the periodic case, pushing the development further, leads to effective equations that are not those implemented in the spectral elements' codes. Furthermore, and more dramatically, in 1D (Capdeville et al., 2009a), a nonperiodic development as the one looked for in this article, is only valid for the order 0 and the first order corrector. It is valid for higher order only in certain specific cases. Therefore, we do not go further, and again looks for a partial first-order solution, as in the section of this article, dedicated to periodic settings. As for the periodic case, the different orders can be combined as shown in section 2.4. For our partial first-order solution, and using the same notations, we then have

ûε 0 ,ε,1/2 (x) = ûε 0 ,ε,0 (x) + χ ε 0 ,1 (x, x/ε) • ∇ x ûε 0 ,ε,0 (x) , (80)
where the 1/2 superscript again means "partial order 1". Of course, for it is only a partial order 1 solution, we do not have in general

u ε 0 ,ε (x) = ûε 0 ,ε,1/2 (x) + O(ε 2 ) , (81)
Let us notice finally, that practically, only the case where ε 0 = ε is of interest, as underlined by Capdeville et al. (2009a), because it is the only value of ε 0 that leads to the solution of the reference problem (see 61).

3.4 Construction of ρ ε 0 (x, y) and µ ε 0 (x, y)

Here we are at the heart of the problem: how to conveniently construct ρ ε 0 (x, y) and µ ε 0 (x, y), in order for u ε,0 , u ε,1 and σ ε,0 to be in V? As already said, µ ε 0 must be built such that the first order corrector χ ε 0 ,1 (equation 70) and the vector H ε 0 (equation 76) belong to V. Let us just recall that χ ε 0 ,1 being in V implies that G ε 0 also is (because it is its gradient). We will show that the reciprocity is verified.

We therefore look for ρ ε 0 (x, y) and µ ε 0 (x, y) such that (i) ρ ε 0 , H ε 0 and χ ε 0 ,1 are in V;

(ii) ρ ε 0 and µ ε 0 must be positive definite;

(iii) ρ ε 0 (x, x/ε 0 ) = ρ 0 (x) and µ ε 0 (x, x/ε 0 ) = µ 0 (x).
Let us introduce a given, starting ρ 0 s (x, y) = ρ 0 (ε 0 y) defined on R × Y x and then y-extended periodically to R 2 . Note that this starting density distribution, is x-dependent.

The construction of ρ ε 0 (x, y) is then trivial:

ρ ε 0 (x, y) = F ρ 0 s (x, x/ε 0 ) + (ρ 0 s -F ρ 0 s )(x, y) (82)
We indeed have ρ ε 0 (x, x/ε 0 ) = ρ 0 (x) and ρ ε 0 is in V and is a positive function with a well chosen wavelet w m .

It is by far not so easy to find a correct µ ε 0 . We shall follow the same procedure described in Capdeville (2009a) or Capdeville et al. (2009b). Let us recall that these fields are defined as:

G ε 0 = I + ∇ y χ ε 0 ,1 , (83)
H ε 0 = µ ε 0 • G ε 0 , (84)
∇ y • H ε 0 = 0 , (85)
G ε 0 = I . (86)
To find convenient G and H, we process in the following way:

(i)

Step 1: at a given x, let us define the area Y x , and build a starting µ ε 0 ,s (x, y) = µ 0 (x, ε 0 y) on it. Then solve for the (starting) periodic corrector solution of ( 85) with periodic boundary conditions on Y x , in order to find χ ε 0 ,1,s (x, y).

Then we can compute

G ε 0 ,s (x, y) = I + ∇ y χ ε 0 ,1,s , (87)
H ε 0 ,s (x, y) = µ ε 0 ,s (x, y) • G ε 0 ,s (x, y) . (88)
At this stage, we can already compute the effective elastic tensor for any x:

µ * ,ε 0 (x) = F (H ε 0 ,s ) • F (G ε 0 ,s ) -1 (x/ε 0 ) . (89) (ii)
Step 2: now let us compute for any y ∈ Y x ,

G ε 0 (x, y) = [(G ε 0 ,s -F (G ε 0 ,s )) (y)] • [F (G ε 0 ,s ) (x/ε 0 )] -1 + I , (90)
H ε 0 (x, y) = [(H ε 0 ,s -F (H ε 0 ,s )) (y) + F (H ε 0 ,s ) (x/ε 0 )] • [F (G ε 0 ,s ) (x/ε 0 )] -1 . (91)
The extension of G ε 0 and H ε 0 from Y x to R 2 in y in then done by periodic extension. At this stage, we can check that (H, G) ∈ V 2 , and that G = I.

This last property is interesting, because of the definition (83) of G ε 0 and the fact that this quantity belongs to V: it implies that ∇ y χ ε 0 ,1 also belongs to V; it has also the following property:

∇ y χ ε 0 ,1 = 0, and then, because of (60) , we finally derive the necessary condition: χ ε 0 ,1 ∈ V.

(iii) Step 3: From (84) we can now build

µ ε 0 (x, y) = H ε 0 • (G ε 0 ) -1 (x, y) , (92)
(iv) Step 4: once µ ε 0 (x, y) is known, the whole homogenization procedure can be pursued to find the different components of the 1st-order corrector vector.

This procedure may seems obscure and here is a tentative of interpretation (that is a bit different than the one given in Capdeville et al. (2009a) or Capdeville et al. (2009b)). The idea behind step one, is to first assume that all scales are fast. Then, the cell problem equations (71) (which are those of the static

y k + χ ε 0 ,1,s k
; therefore G ε 0 ,s only is the gradient of this displacement (as equation ( 83) indicates) and H ε 0 ,s the microscopic stress associated to this gradient. The idea behind any elastic homogenization technique, is to find the effective tensor that links an effective stress and an effective deformation: this is exactly what is done in equation ( 89) where F (H ε 0 ,s ) and F (G ε 0 ,s ) play the role of the effective stress and deformation. Equation ( 89) allows to directly determinate (after filtering) the homogenized elastic tensor, without solving for the "correct" other quantities, with the starting µ ε 0 ,s (y). Practically, in order to ensure the symmetry of this effective elastic tensor, we impose, in this procedure, the equality of µ * ,ε 0 12 and µ * ,ε 0 21 , adding then a condition in (89) and solving this latter using a singular value decomposition. Steps 2 and 3 allow to build µ ε 0 (x, y) by separating the scales of H ε 0 ,s and G ε 0 ,s using the filter F () and making sure that G ε 0 = I.

Following all these operations, we have by construction µ ε 0 (x, x/ε 0 ) = µ 0 (x) and µ ε 0 is positive definite for a well chosen wavelet w m . We also have checked that χ ε 0 ,1 also belongs to V (we have G ε 0 ∈ V and G ε 0 = I by construction (see step 2), which, using (60), indeed implies that χ ε 0 ,1 belongs to V ). This corrector is unique when we impose the classical averaging condition χ ε 0 ,1 = 0.

In order to determine χ ε 0 ,1 , we first find µ ε 0 with (92) and solve (71) once again.

As noticed by Capdeville et al. (2009b), we are not yet able to prove the symmetry of the homogeneous elastic tensor µ * ,ε 0 in the general case, using the procedure presented in this section. We prove it in the specific case of a nonperiodic but layered medium, in appendix C.

Homogenization in a non-periodic setting: an example

To validate it, we finally apply the homogenization procedure developed previously, for wave propagation throughout the random model shown in Fig. 5-a. This model is a square of size 20 × 20 km 2 , surrounded by a 1km thick strip of constant physical properties (µ 0 11 = µ 0 22 = 60 GPa, µ 0 12 = µ 0 21 = 0GPa, ρ = 2800 kg m -3 ), plus an additional PML to avoid reflections at the boundaries in long-time simulations. This square is divided in elements of size 100×100m 2 . Elastic values in each element are randomly generated within ±50% of the elastic values of the surrounding strip; note that it is decided to start with an isotropic model, so the reference values for µ 0 11 and µ 0 22 are equal. A reference solution can be computed in this model using SEM. A collocated source is located at S = (x S 1 , x S 2 ) = (100 m, 10 km). Its time evolution is described by a Ricker with a central time shift of 0.4 s and a central frequency of 5 Hz (to which is associated a corner frequency of about 12.5 Hz).

Seismograms are recorded at fourteen receivers as shown in Fig. 5-a in Fig. 8, in grey line with black dots.

As for the homogenization, and practically, we first solve the cell problem (71) on small domains Y X , with periodic boundary conditions, for a large number of x-values. We therefore obtain values for first-order correctors and for homogenized coefficients. The resolution of the cell problem is done using a finite element method, based on the same mesh and quadrature as the ones that is used in wave propagation simulations with SEM. In the following example, the polynomial order of the expansion of physical fields, is equal to 4.

In Fig. 5-b is shown the effective model corresponding to the previous random model Fig. 5-a; it is obtained by applying a homogenization procedure with a low-pass filter whose the cut-off corresponds to a ε 0 equal to 0.27. As can be seen, elastic quantities also show rapid spatial variations, but these variations are much smoother than in the original medium. This can be seen more precisely along a section, see Fig. 6. It is also possible to determinate the apparent anisotropy, as seen by the wavefield whose the smallest wavelength is much larger than the characteristic size of the (isotropic!) heterogeneities; here we decide to measure the anisotropy as the normalized difference between the homogenized values for µ 0 11 and µ 0 22 : (µ * ,ε 0 11 -µ * ,ε 0 22 )/µ * ,ε 0 11 (the other terms of the homogenized elastic tensor being very close to be null). As can be seen on Fig. 7, the anisotropy can be quite large, varying between -6% and 6%, with (absolute) average values lying between 1 and 2%. In Fig. 8 we compare the particle velocity (recorded at receiver 8, with the source located at S; see 5-a), calculated in the reference medium (grey line with black dots), to the one obtained with the sole order 0 homogenization ( uε 0 ,0 ), for different values of ε 0 , varying from 2.16 to 0.27 (black line). For large values of ε 0 , that means, when the homogenized model is too smooth with respect to the minimum wavelength of the wavefield, it clearly appears that the coda of direct waves is nonexistent or incorrectly calculated. The more ε 0 decreases, the more the coda correctly "built" is: for a ε 0 of 0.27, the agreement between the reference and the homogenized solutions, is excellent. Note that the arrival time of the ballistic wave is always rightly retrieved, whatever the value of ε 0 is; this is not the case when seismograms are recorded in effective media when a simple (and erroneous) filtering is applied to elastic quantities (see Fig. 9).

Let us now take a look at the convergence of the homogenized solution u ε 0 ,0 towards the reference one u ref .

As suggested in the companion paper (Capdeville et al., 2009b), let us introduce a measure of the error between a given field u and the reference one, at a given receiver k, as:

E k ( u) = 20 0 ( u -uref ) 2 (x k , t)dt 20 0 ( uref ) 2 (x k , t)dt , (93)
where the maximal bound for the time integral, indicates that seismic recordings of interest last 20 seconds (after that time the amplitude of the coda wave is completely negligible). Then a combined error E tot can be defined, taking into account the error at each individual receiver (see Fig. 5-a), as:

E tot ( u) = 1 14 14 k=1 E k ( u) . (94)
As in the P-SV case (Capdeville et al., 2009b), the error E tot ( uε 0 ,0 ) decreases quite slowly when ε 0 with large values, decreases (see Fig. 10). The reason is quite obvious when taking a look at Fig. 8:

the coda only begins to be fully constructed for values of ε 0 around 0.5. Quite surprinsingly, from this point, the convergence is then very fast: instead of a theoretical, expected convergence in ε 0 , the order 0 homogenized converges towards the reference one, in ε 2 0 . For this random example, this observation 94), as a function of ε 0 , between the order 0 homogenized solution (E(u ε0,0 ), straight line) or the partial first-order solution (E(u ε0,1/2 ), dashed line), and the reference one.

underlies the fact that high-order terms in the asymptotic expansion (64) are negligible in comparison with the leading term. This is confirmed when calculating the same kind of error, but this time between the partial first-order solution u ε 0 ,1/2 and the reference one (see Fig. 10): the effect of the correction term is null for large values of ε 0 , and quite weak for small ε 0 s -even if clearly observable. Note that Capdeville et al. (2009b) pursue the same kind of calculations for smallest values of ε 0 , for which the effect of the first-order correction, is more pronounced. Finally, we observe a quite good convergence even for high ε 0 values, because even if the source is located in a homogeneous zone, the area on which we calculate correctors and homogenized quantities, is larger than the homogeneous, surrounding strip of the box model (see Fig. 5) when ε 0 is equal to 2.16, therefore the first-order correction for the point source (see Capdeville et al. (2009a)) should be taken into account in this case (and the error associated with this large ε 0 should then be lower than the one calculated here).

The effect of the first-order correction can be observed in Fig. 11, where are reported the errors E( uε 0 ,0 ) and E( uε 0 ,1/2 ) as a function of x 1 , along a line joining all receivers plotted on Fig. 5-a.

The error almost always decreases with the supplemented corrective first-order term. We also observe, although it is less marked than in the P-SV case, that the error determined for the partial first-order homogenized solution, varies more slowly with x 1 , than the one for the order 0 solution. This effect can be explained by the fact that the first-order corrector term, depends on y (fast scale)-and then 94), for ε 0 = 0.27, between the order 0 homogenized solution (straight line) or the partial first-order solution (dashed line), and the reference one, measured at each receiver reported in Fig. 5, plotted as a fonction of their abscissa.

"corrects" for errors (between the order 0 homogenized and the reference solution) that vary rapidly in space. This "microscopic" dependence of the first-order correction term can be more precisely observed in Fig. 12, where is plotted this correction ( uε 0 ,1/2 -uε 0 ,0 ) against with the residual uref -uε 0 ,0 (normalized by the largest value of the velocity in each case), along a line connecting receivers shown in Fig. 5, with ε 0 =0.27, and at a time t = 12s. The fast oscillations are quite similar for both curves, and the differences only are the uncomputed higher-order terms in the asymptotic expansion (64).

They are quite small in amplitude with respect to the one of the reference wavefield (between 1 and 2% of the velocity amplitude), and quite close, which means that the partial first-order homogenized solution we suggested to use in this article, is an excellent approximation to the wavefield propagating throughout the initial, random model.

CONCLUSION

After having recalled classical results of the two-scale homogenization theory in heterogeneous but periodic media, a generalization of this theory has been introduced, as in the companion papers (Capdeville et al., 2009a;Capdeville et al., 2009b). It is valid for wave propagation in general, nonperiodic media -as ones typical of the Earth -, when the smallest wavelength of the wavefield, is much larger that the characteristic size of the heterogeneities. Excellent results, in comparison with that obtained in the reference medium (both using the Spectral Element Method as the wave propagation tool), are obtained using a partial first-order solution, including the order 0, leading term of the classical (in periodic media) or equivalent (in nonperiodic media) asymptotic two-scale expansion used in this theory, supplemented by a (though incomplete) first-order correction. Although it would be possible to add higher-order terms in the case of periodic media, at the cost of a change in the SEM algorithm, this is in general not possible in the case of general media, at less using the theory presented in this article.

Let us recall, as done extensively in Capdeville et al. (2009b), that one the major advantages of a homogenization (or upscaling) procedure, is to obtain effective media (and an associated effective wave equation) smoother than original ones (typically, the Earth'upper mantle and crust), implying for numerical simulations, much simpler and sparser meshes, resulting in an efficient reduction in calculation time; this smoothing effect in very heterogeneous media, has been shown in both parts of this article, although we did not insist on benefits relative to the meshing (see Capdeville et al. (2009b) and the application to P-SV propagation in the Marmousi model, which is a perfect illustration for that issue).

First-order corrections at the source (Capdeville et al., 2009a;Capdeville et al., 2009b), and boundary L. GUILLOT conditions' effects are issues that have not been tackled in this article; this latter, fundamental when studying surface-waves propagation for instance, will be the purpose of a future work, generalizing the results of [START_REF] Capdeville | Shallow layer correction for spectral element like methods[END_REF].

Let us underline two major results of this work. First, it appears that even far of the application range of the homogenization theory (that means, when effective media are constituted by large-scale heterogeneities, relative to the wavefield propagating through them -typically: the domain of application of ray theory), the homogenization procedure seems to allow to obtain a smooth medium in which time arrivals of ballistic phases, are very well retrieved (which is not the case when using a natural filtering of elastic properties). As it is well-known that this time arrival is sensitive to smooth effective models of the real ones, one may wonder if the homogenization theory would not be a convenient tool to obtain this effective model -as it seems the case in this article and its companion paper. Second, it is interesting to recall that (as shown by Backus in 1962 for layered media), an isotropic model (at the microscopic scale) is anisotropic when seen by a wavefield whose minimal wavelength is much larger than its characteristic size -and that this macroscopic anisotropy is far from being negligible (around 1 and 2% in average, with peaks at ±6%). This is a serious indication that models of the Earth, as looked for by seismologists when inverting for long-period data, should be considered as anisotropic (then involving, a specific parameterization).

Let us finally notice that the homogenization procedure presented in this paper (and specifically, the construction of the fast spatial dependence of the density and elastic tensor, see section 3.4), should be generalized without any difficulty to three-dimensional models, and that applications to the real Earth are then in sight. A patent [START_REF] Capdeville | Procédé de détermination d'un modèle élastique effectif[END_REF] has been filed on the non-periodic homogenization process (which is by no mean a restriction to any academic research on the subject).

-(∂ y j χ 1 k )µ ij (∂ y i χ ρ ) = ρχ 1 k , (B3)
which can be rewritten:

-

∂ y i (∂ y j χ 1 k )µ ij χ ρ + ∂ y i µ ij ∂ y j χ 1 k χ ρ = ρχ 1 k , (B4)
in which the first term is equal to 0.

Using successively the definition of the first-order corrector ( 27), and the property (13) the last equality can be transform as

-∂ y i (µ ij δ jk ) χ ρ = µ ij δ jk ∂ y i χ ρ = µ ik ∂ y i χ ρ = ρχ 1 k . (B5)
Knowing the definition of µ ρ * k (equation 43), the equality µ ρ * k = ρχ 1 k is then proved. One can first notice that, because of the symmetry of µ ε 0 ,s , the homogenized tensor also is symmetric, which is an important, theoretical result: our procedure, in the simplified case of a stratified medium, preserves the symmetry of the elastic tensor.

Second, it is interesting to note that each component of the homogenized tensor, is equal to some combinations of filtered quantities only. Furthermore, it appears that these combinations are exactly the same ones, as that obtained analytically for the periodic case (e.g., [START_REF] Cioranescu | An introduction to homogenization[END_REF]). We therefore have generalized this latter, well-known result, to stratified, nonperiodic media, the average over the unit cell being replaced by a filtering operation, but on the same quantities.

  equations in nonperiodic media.

Figure 1 .

 1 Figure 1. Sketching detail of one (very simple) periodic cell of an infinite plane. The periodicity vector =( 1 , 2 ) has a characteristic length, much smaller than the characteristic size of the propagating wavefield, Λ.

  seismograms, with respect to ones calculated in a reference medium. Therefore, as inCapdeville et al.

Figure 2 .

 2 Figure 2. Periodic model. The periodicity vector is = ( 1 , 2 ) = (120m, 200m). The point source is located at S, the receiver at R. The (L)-line is the one along which are calculated snapshots shown in Fig. 4.

Figure 3 .

 3 Figure 3. Displacement (in meters) recorded at the receiver R located at (10km,10km). The reference is plotted in grey, the order 0 homogenized solution in black, and the "natural averaging" (see text)solution is in dashed line.

Figure 4 .

 4 Figure 4. -a: grey line: snapshot of the displacement u ε (x, t) at t =0.4 s computed in the reference model. Black line: the order 0 homogenized solution û0 (x, t). Dashed line: solution computed in a model obtained by arithmetically averaging the elastic properties.

  equilibrium in mechanics) are solved. They describe the microscopic effect of the imposition of a unit gradient in displacement in the y k -direction: the associated (microscopic) displacement obviously is

Figure 5 Figure 6 .

 56 Figure 5. a-Square, random model. The point source is located at S, the 14 receivers are also indicated. Here are represented the values of the µ 0 11 coefficient. b-Homogenized (or effective) model for the square, random model. On this figure are reported the values for the homogenized coefficient µ * ,ε011 , with ε 0 equal to 0.27.

Figure 7 .

 7 Figure 7. Effective anisotropy for the random model, calculated as (µ * ,ε0 11 -µ * ,ε0 22 )/µ * ,ε0 11 , with ε 0 =0.27.

Figure 8 .

 8 Figure 8. Velocity traces computed for the collocated source S, at the receiver 8. The reference solution is in grey with black dots; the order 0 homogenized solution, in black. The ε 0 -parameter takes values 2.16 (a), 1.08 (b), 0.54 (c) and 0.27 (d).

Figure 9 .

 9 Figure 9. Velocity traces recorded at the receiver 8 (see Fig. 5-a). The reference solution is plotted in grey, and the "natural filtered" solution is in dashed line. Both are clearly out-of-phase.

Figure 11 .

 11 Figure 11. Total errors as defined in (94), for ε 0 = 0.27, between the order 0 homogenized solution (straight

Figure 12 .

 12 Figure12. Cut along a line joining receivers plotted in Fig.5, for the correction term uε0,1/2 -uε0,0 (straight line) and the residual between the reference solution and the order 0 homogenized one, uref -uε0,0 (dashed line), at a time t = 12s. The value of ε 0 is 0.27.

APPENDIX A:

Let us recall the classical equations in index notation, and neglecting the source term:

where the elastic constants satisfy the following symmetries:

In the antiplane case, the motion is considered as being perpendicular to, let us say, the (x 1 , x 2 ) plane (and parallel to the x 3 -axis), and physical properties only vary in this plane, so the previous equations can be rewritten:

Obviously the constant index can be erased, and we can redefine the elastic tensor, which now depends on 4 parameters only:

Note that the reduced elastic tensor µ is, as c, positive-definite, and that its components are symmetric.

Of course the dispacement is a scalar, and the stress, a simple vector. The equations we have to solve are therefore:

(A5)

APPENDIX B:

We want to show that the equality E ρ * = ρχ 1 , valid in 1D (see Capdeville et al. (2009a)), can be generalized:

This can be done in a similar manner, using the property (13) extensively. Let us start with the differential equation defining χ ρ :

Multiplying this last equation by χ 1 k , taking the cell average of the resulting expression, and using the fact that χ 1 k = 0, we obtain

APPENDIX C:

We show in this appendix, that in the specific case of a nonperiodic, layered medium (that means, one whose properties only change in one specific spatial direction), the effective elastic tensor µ * ,ε 0 as calculated in ( 89), is indeed a symmetric tensor.

Let us guess that material properties only evolve in the x 1 (y 1 )-direction, and let us omit the index from now on, as well as the x-dependence of physical quantities.

Because all derivatives relative to y 2 are null, equation (71) for the "starting" corrector as defined in the first step of the procedure described in section (3.4), leads to

We can solve both of these equations, for the derivative of the corrector's components:

Both constants c 1 and c 2 can be calculated by integrating the last equations, over the periodic cell Y , and by using the equalities χ ε 0 ,s 1 = χ ε 0 ,s 2 = 0. For instance, c 1 = 1/2β β -β (1/µ ε 0 ,s 11 (y))dy. Then we can easily obtain the expressions for G ε 0 ,s and H ε 0 ,s , according to (72): )


Sylvain Sorin
email: sorin@math.jussieu.fr

Xavier Venel
email: xavier.venel@sip.univ-tlse1.fr

Guillaume Vigeral
email: guillaumevigeral@gmail.com

Asymptotic Properties of Optimal Trajectories in Dynamic Programming

We show in a dynamic programming framework that uniform convergence of the finite horizon values implies that asymptotically the average accumulated payoff is constant on optimal trajectories. We analyze and discuss several possible extensions to two-person games.

Presentation

Consider a dynamic programming problem as described in Lehrer and Sorin [START_REF] Lehrer | A Uniform tauberian theorem in dynamic programming[END_REF]. Given a set of states S, a correspondence Φ from S to itself with non empty values and a payoff function f from S to [0, 1], a feasible play at s ∈ S is a sequence {s m } of states with s 1 = s and s m+1 ∈ Φ(s m ). It induces a sequence of payoffs {f m = f (s m )}, m = 1, ..., n, .... Recall that starting from a standard problem with random transitions and/or signals on the state, this presentation amounts to work on the set of probabilities on S and to consider expected payoffs.

Let v n (s) (resp. v λ (s)) be the value of the n stage program G n (s) (resp. λ discounted program G λ (s)) starting from state s. The asymptotic approach deals with asymptotic properties of the values v n and v λ as n goes to ∞ or λ goes to 0.

The uniform approach focuses on properties of the strategies that hold uniformly in long horizons. v ∞ is the uniform value if for each ε > 0 there exists N such that for each s ∈ S: 1) there is a feasible play {s m } at s with

1 n n m=1 f (s m ) ≥ v ∞ (s) -ε, ∀n ≥ N
2) for any feasible play {s ′ m } at s and any n ≥ N

1 n n m=1 f (s ′ m ) ≤ v ∞ (s) + ε.
Obviously the second approach is more powerful than the second (existence of a uniform value implies existence of an asymptotic value : the limit of v n exists) but it is also more demanding: there are problems without uniform value where the asymptotic value exists (see Section 2). Note that the condition for the existence of a uniform value implies that the average accumulated payoff on optimal trajectories remains close to the value.

We will prove that a similar phenomenon holds true under conditions that are stronger than the existence of an asymptotic value but weaker than the existence of a uniform value.

Say that the dynamic programming problem is regular if : i) lim v n (s) = v(s) exists for each s ∈ S. ii) the convergence is uniform. This condition was already introduced and studied in Lehrer and Sorin [START_REF] Lehrer | A Uniform tauberian theorem in dynamic programming[END_REF] (see Section 2).

We consider the following property P:

Date: October 2009.

For any ε > 0, there exists n 0 , such that for all n ≥ n 0 , for any state s and any feasible play {s m } ε-optimal for G n (s) and for any t ∈ [0, 1]:

(1) 3ε ≥ 1 n ( [tn] m=1 f m ) -tv(s) ≥ -3ε.
where [tn] stands for the integer part of tn.

This condition says that the average payoff remains close to the value on every almost-optimal trajectory with long duration (but the trajectory may depend on this duration). It also implies a similar property on every time interval.

Examples and comments

1) The existence of the asymptotic value is not enough to control the payoff as required in property P. An example is given in Lehrer and Sorin [START_REF] Lehrer | A Uniform tauberian theorem in dynamic programming[END_REF] (Section 2), where both lim v n and lim v λ exist on S but where the asymptotic average payoff is not constant on the unique optimal trajectory, nor on ε-optimal trajectories: in G 2n , an optimal play will induce n times 0 then n times 1 while v = 1/2. Note that this example is not regular: the convergence of v n to v is not uniform.

2) Recall that in the framework of dynamic programming, regularity is also equivalent to uniform convergence of v λ (and with the same limit), see Lehrer and Sorin [START_REF] Lehrer | A Uniform tauberian theorem in dynamic programming[END_REF] (Section 3). Note also that this regularity condition is not sufficient to obtain the existence of a uniform value, see Monderer and Sorin [START_REF] Monderer | Asymptotic properties in dynamic programming[END_REF] (Section 2).

3) General conditions for regularity can be found in Renault [START_REF]Uniform value in dynamic programming[END_REF].

Main result

Theorem 3.1. Assume that the program is regular, then P holds.

Proof

Let us start with the upper bound inequality in (1). The result is clear for t ≤ ε (recall that that the payoff is in [0, 1]). Otherwise let n 1 large enough so that n ≥ n 1 implies ||v n -v|| ≤ ε by uniform convergence. Then the required inequality holds for n ≥ n 2 with [εn 2 ] ≥ n 1 .

Consider now the lower bound inequality in (1). The result holds for t ≥ 1ε by the ε-optimal property of the play, for n ≥ n 1 . Otherwise we use the following lemma from Lehrer and Sorin [START_REF] Lehrer | A Uniform tauberian theorem in dynamic programming[END_REF] (Proposition 1). Lemma 3.1. Both lim sup v n and lim sup v λ decresase on feasible histories.

In particular, starting from s [tn] the value of the program for the last n -[tn] stages is at most v(s [tn] ) + ε for n ≥ n 2 , by uniform convergence, hence less than the initial v(s) + ε, using the previous Lemma. Since the play is ε-optimal in G n (s), this implies that (2)

[tn] m=1 f m + (n -[tn])(v(s) + ε) ≥ n(v n (s) -ε) ≥ n(v(s) -2ε)
hence the required inequality. For any ε > 0, there exists λ 0 , such that for all λ ≤ λ 0 , for any state s and any feasible play {s m } ε-optimal for G λ (s) and for any t ∈ [0, 1]:

(3) 3ε ≥ n(t;λ) m=1 λ(1 -λ) m-1 f m ) -tv(s) ≥ -3ε.
where n(t; λ) = inf{p ∈ IN; p m=1 λ(1λ) m-1 ≥ t}. Stage n(t; λ) corresponds to the fraction t of the total duration of the program. Theorem 4.1. Assume that the program is regular, then P' holds.

Proof

The proof follows the same lines than the proof of Theorem 3.1. Recall that by regularity both v n and v λ converge uniformly to v. Moreover the discounted sums

(1-λ) -N N m=1 λ(1-λ) m-1 f m belong to the convex hull of the averages 1 n n m=1 f m ; 1 ≤ n ≤ N . The counterpart of equation (2) is now (4) n(t;λ) m=1 λ(1 -λ) m-1 f m + (1 -t)(v(s) + ε) ≥ (v λ (s) -ε) ≥ v(s) -2ε 4.2. Continuous time.
Similar results holds in the following set-up: v T (x) is the value of the control problem Γ T with control set U where the state variable in X is governed by a differential equation (or more generally a differential inclusion) ẋt = f (x t , u t ) starting from x at time 0. The real payoff function is g(x, u) and the evaluation is given by:

1 T T 0 g(x t , u t )dt.
Regularity in this framework amounts to uniform convergence (on X) of V T to some V . (Sufficient conditions for regularity can be found in Quincampoix and Renault [START_REF] Quincampoix | On the existence of a limit value in some non expansive optimal control problems[END_REF]). The corresponding property is now P": For any ε > 0, there exists T 0 , such that for all T ≥ T 0 , for any state x and any feasible trajectory ε-optimal for Γ T (x) and for any θ ∈ [0, 1]:

(5) 3ε ≥ 1 T θT 0 g(x t , u t )dt -θV (x) ≥ -3ε.
Theorem 4.2. Assume that the optimal control problem is regular, then P" holds.

Proof

Follow exactly the same lines than the proof of Theorem (2).

Finally the same tools can be used for an evaluation of the form λ +∞ 0 e -λt g(x t , u t )dt.

Two-player zero-sum games

In trying to extend this result to a two-person zero-sum framework, several problems occurs.

Optimal strategies on both sides.

First it is necessary, to obtain good properties on the trajectory, to ask for optimality on both sides.

For example in the Big Match with no signals,

α β a 1 * 0 * b 0 1
where a * denotes an absorbing payoff, the optimal strategy of player 1 in the "asymptotic game" on [0, 1] is to play "a before time t" with probability t, see Sorin [START_REF] Sorin | A first course on zero-sum repeated games[END_REF] Section 5.3.2. Obviously, if there is no restrictions on player 2's moves the average payoff will not be constant. However, the optimal strategy of player 2 is "always (1/2, 1/2)" hence time independent on [0, 1]. It thus induces a constant payoff and it is easy to see that the property is robust to small perturbations in the evaluation of the payoff.

Player 1 controls the transition.

Consider a repeated game with finite characteristics (states, moves, signals, ...) and use the recursive formula corresponding to the canonical representation with entrance laws being consistent probabilities on the universal belief space, see Mertens, Sorin and Zamir [START_REF] Mertens | Repeated Games, CORE Discussion Papers[END_REF], Chapters III.1, IV.3. This representation preserves the values but in the auxiliary game, if player 1 controls the transition an optimal strategy of player 2 is to play a stage by stage best reply. Hence the model reduces to the dynamic programming framework and the results of the previous sections apply. A simple example corresponds to a game with incomplete information on one side where asymptotically an optimal strategy of the uniform player 1 is a splitting at time 0, while player 2 can obain u(p t ) at time t where u is the value of the non-revealing game and p t the martingale of posteriors at time t, see Sorin [START_REF] Sorin | A first course on zero-sum repeated games[END_REF], 3.7.2.

Example.

Back to the general framework of two person zero-sum repeated games, the following example shows that in addition one has to strengthen the conditions on the pair of ε-optimal strategies. We exhibit a game having a uniform value v but for some state s with v(s) = 0 one can construct, for each n, optimal strategies in Γ n (s) inducing essentially a constant payoff 1 during the first half of the game.

Starting from the initial state s, the tree representing the game Γ has countably many subgames Γ2n , the transition being controlled by player 1 (with payoff 0). In Γ2n there are at most n stages before reaching an absorbing state. At each of these stages of the form (2n, m), m = 1, ...n, the players plays a "jointly controlled" process leading either to a payoff 1 and the next stage (2n, m + 1) (if they agree) or an absorbing payoff x 2n,m with (m -1) + (2n -(m -1))x 2n,m = 0, otherwise. Hence every feasible path of length 2n in Γ2n gives a total payoff 0. Obviously the uniform value exists since each player can stop the game at each node, inducing the same absorbing payoff. The representation is as follows:

Notice that in the 2n + 1 stage game, after a move of player 1 to Γ2n , any play is compatible with optimal strategies, in particular those leading to the sequence of payoffs 2n times 0 or n times 1 then n times -1.

Conjectures.

A natural conjecture is that in any regular game (i.e. where v n converges uniformly to v): for any ε > 0, there exists n 0 , such that for all n ≥ n 0 , for any initial state s, there exists a couple (σ n , τ n ) of ε-optimal strategies in G n (s) such that for any t ∈ [0, 1]: A more elaborate conjecture would rely on the existence of an asymptotic game Γ * played in continuous time on [0, 1] with value v (as in Section 5.1). We use the representation of the repeated game as a stochastic game trough the recursive structure as above, see Mertens, Sorin, Zamir [START_REF] Mertens | Repeated Games, CORE Discussion Papers[END_REF], Chapter IV. The condition is now the existence of a couple of strategies (σ, τ ) in the asymptotic game that would depend only on the time t ∈ [0, 1] and on the current state s such that for any ε > 0, there exists η with the following property: in any repeated game where the (relative) weight of stage m is α m , with {α m } decreasing and less than η, thus defining a partition Π of [0, 1], the strategies (σ Π , τ Π ) induced in the repeated game by (σ, τ ) satisfies [START_REF] Sorin | A first course on zero-sum repeated games[END_REF].

(6) 3ε ≥ 1 n E s σn,τn ( [tn] m=1 f m ) -tv(s) ≥ -3ε. • s 0 • (2, 1) 0 • (4, 1) 0 0 • (2n, 1) 0 Γ 2 Γ 4 Γ 2n Figure 1. The game Γ starting from state s C A C A C A C A C A C 1 -→ 0* C 1 -→ x 2n,2 * C 1 -→ x 2n,m * C 1 -→ x 2n,n * C -1* -1* A 0* 0* A x 2n,2 * x 2n,2 * A x 2n,m * x 2n,m * A x 2n,n * x 2n,n * A -1* -1* (2n, 1) (2n, 2) • • • (2n, m) • • • (2n, n) -1 *

4 . Extensions 4 . 1 .

 441 Discounted case.A similar result holds for the program G λ corresponding to the evaluation ∞ m=1 λ(1λ) m-1 f m . Explicitly, one introduces the property P ′ :

Figure 2 .

 2 Figure 2. The subgame Γ 2n starting from state (2n, 1)

Acknowledgment:

This work was done while the three authors were members of the Equipe Combinatoire et Optimisation. Sorin's research was supported by grant ANR-08-BLAN-0294-01 (France).