Keywords: propagation d'ondes, condensation de masse, schéma numérique, réseau, arbre, Kirchhoff

Résolution numérique de l'équation des ondes sur un réseau de fentes Résumé : Dans ce rapport technique, nous présentons un modèle théorique et numérique pour simuler la propagation des ondes dans des réseaux finis de fentes avec des conditions de Kirchhoff classiques et améliorées aux noeuds des réseaux. Nous commençons par décrire le cadre continu, puis nous discrétisons le problème en utilisant la technique des éléments finis avec condensation de masse, introduite par G. Cohen et P. Joly (voir [3]). Finalement, nous montrons une implémentation du schéma numérique dans un code écrit en C++ en collaboration avec K. Boxberger, quelques résultats et des estimations d'erreurs.

Introduction

In this technical report, we extend the results about propagation of acousic wave in a general network of thin slots. Study of waves in networks is not so recent: one can check for example the works of B. Dekoninck and S. Nicaise [5], and in a more recent case the works of B. Maury, D. Salort and C. Vannier [START_REF] Maury | Trace theorems for trees, application to the human lungs[END_REF]. It's also quite easy to do some numerical simulations in some particular case of networks (one see for example the work of Y. Achdou, C. Sabot and N. Tchou [1]). But actually, it is not possible to find an efficient numeric code which solves wave propagation problem on a general network. To avoid this, a subject was proposed for an internship during second trimester 2009. During this internship which was supervised by P. Joly and A. Semin, K. Boxberger studied the propagation of acoustic wave in "semi-homogeneous" networks, i.e. networks whose geometrical and physical parameters are constant functions on each edge of those networks, but are not constant on the whole network. The approach used in this stage was with a finite differences method. The approach allowed us to get some good numeric scheme, but we can see that we lost one order of error estimate when dealing with boundary conditions. The idea is then to use a finite element approach. One could say that numerically, finite element method is slower than finite differences method, that's why we will use the mass-lumping method.

The aim of this technical report is the following: in section 1, we present the continuous problem we want to solve. In section 2, we present the discrete problem associated to the continuous problem. We first use space discretisation, then time discretisation. In section 3, we present an implementation of the numeric code netwaves (whose writing has during the previously cited internship -the code can be soonly found at http://gforge.inria.fr/projects/netwaves/). In section 4, we give some error estimates for both problems introduced in sections 1 and 2. Finally, in section 5, 1 Setting of the continuous problem

Geometry of the problem

In this part, we recall some notations about the geometry we consider, in the spirit of [START_REF] Maury | Trace theorems for trees, application to the human lungs[END_REF], but taking the geometry in ℝ , with = 2 or = 3. Reader can find more information in [2] or in [START_REF] Joly | Construction and Analysis of Improved Kirchoff Conditions for Acoustic Wave Propagation in a Junction of Thin Slots[END_REF] about the choice of the conductances. Let = ( , ℰ, ) denote a graph: is a finite set of vertices in ℝ , ℰ subset of × such that, for (( , ), ( , )) ∈ ℰ × ℰ with ( , ) ∕ = ( , ), the open segments ( ) and (

) do not intersect in ℝ and the the intersection of the closed segments [ ] and [ ] is a subset (possibly empty) of , and ∈ (0, +∞) ℰ a conductance field. Note that the condition we put on ℰ let us allow say that edges are counted only once in ℰ, that is to say ( , ) ∈ ℰ ⇒ ( , ) ∕ ∈ ℰ. Such example of sets are given by figure 1.1.

We now have to introduce some definitions to describe more precisely the particularities of graph : Here we introduce some definitions to describe more precisely the structure of our network.

• Edges connected to a vertex. Given ∈ , the set of edges ( ) connected to is the subset of ℰ given by elements whose one of the extremities is . In other words,

( ) := {( , ) ∈ ℰ / = or = }
• Inner and outer vertices. A vertex ∈ will be called an outer vertex if the number of edges connected to is equal to 1. Otherwise, it will be called an inner vertex. We call the set of outer vertices, and the set of inner vertices. In other words,

= { ∈ / # ( ) ⩾ 2} = { ∈ / # ( ) = 1}
Remark. In this definition, one assume that there's not any non-connected vertex (i.e. a vertex satisfying # ( ) = 0). Now one wants to define the wave equation on this kind of structure: we have to consider parametrization on each edge of the graph, and to define precisely the equation we consider on each edge ∈ ℰ and conditions we put on each vertex ∈ . : → + -Definition 1.3 (Sobolev spaces on a edge). Given ∈ ℕ and ∈ ℰ, a function defined on will be in the Sobolev space H ( ) if and only if the function ∘ -1 belongs to H (]0, [).

Definition 1.4 (Discontinuous Sobolev spaces on a graph). Let = ( , ℰ, ) be a graph and ∈ ℕ given. A function defined on will be in the Discontinuous Sobolev space H d ( ) if and only if the restriction of to each edge ∈ ℰ belongs to H ( ).

Continuous problem

We want to solve the following problem: find

∈ C 0 (ℝ + , H 1 ( )) ∩ C 1 (ℝ + , H 0 ( )) such that: ⎧          ⎨          ⎩ ∂ 2 ∂ 2 -Δ = 0 in , ∀ ∈ ℰ and ∀ ∈ ℝ * + ∑ ∈ ( ) ( ) ∂ ∂ , = 0 for each ∈ and ∀ ∈ ℝ * + = for = 0 and ∀ ∈ ℰ ∂ ∂ = for = 0 and ∀ ∈ ℰ (1.1)
where:

• Δ means the second derivative of with respect to the curvilinear abscissa, i.e.

Δ = ∂ 2 ( ∘ ) ∂ 2 ∘ -1 (1.2)
• ∂ ∂ , means the inner derivative of with respect to the curvilinear abscissa on an edge ∈ ( ), for ∈ . To be more precise, let ∈ and ∈ ( ). Then we define ∂ ∂ , by

∂ ∂ , ( ) = ⎧  ⎨  ⎩ ∂ ( ∘ ) ∂ ∘ -1 ( ) if -1 ( ) = 0, - ∂ ( ∘ ) ∂ ∘ -1 ( ) otherwise.
Remark. One can see that the sign of this quantity does not depend of the direction of the parametrization .

RT n°0369

A. Semin

Of course, one can define the outer derivative of on an edge ∈ ( ), for ∈ , as the opposite of the inner derivative. More precisely,

∂ ∂ , ( ) = - ∂ ∂ , ( )
• ∈ H 1 (ℰ) and ∈ H 1 (ℰ) vanish near any vertex.

The first line of this system is the "classical" wave equation (one can see that, thanks to (1.2), this line becomes the classical 1D wave equaion in ℝ * + × (0, )). The second line, associated to the fact that is continuous at the nodes, is the so-called Kirchhoff conditions, named in this way by analogy with the conditions Kirchhoff wrote in 1906 for electrical networks. We will see in section 5 what changes if we use other type of conditions.

To complete the system (1.1), one has to add boundary conditions on the outer vertices:

• (In)homogeneous Dirichlet condition: on a given 0 ∈ , one has

( , 0 ) = 0 ( ) where 0 ∈ C 1 (ℝ * + ) • (In)homogeneous Neumann condition: on a given 0 ∈ , one has ∂ ∂ ( , 0 ) = 0 ( ) where 0 ∈ C 0 (ℝ * + ) Remark. One denotes ∂ ∂ instead of ∂ ∂ , .
This is due to the fact that on outer vertices, by definition, there's only one edge connected (there's no ambiguity about ). We omit also for convenience (the convention for normal derivatives is to consider outer normal derivatives).

• Outgoing condition -this is the particular case of the Sommerfield radiation condition: on a given 0 ∈ , one has

∂ ∂ ( , 0 ) + ∂ ∂ ( , 0 ) = 0
We call respectively , , , and , the set of outer vertices with Dirichlet conditions, Neumann conditions and Outgoing conditions. With these conditions, multiplying the first line of (1.1) with the product of ( ) and a test function ∈ H 1 ( ) which vanishes on each vertex ∈ , , and making some integration over the term (-Δ ) , by noting

∇ = ∂ ( ∘ ) ∂ ∘ -1 ,
one gets, by summing the integrals over all edges:

∑ ∈ℰ ∫ ( ) ∂ 2 ∂ 2 ( , ⋅) + ( )∇ ( , ⋅)∇ + ∑ ∈ , ( ) ( ) ( ) + ∑ ∈ , ( ) ∂ ∂ ( , ) ( ) = 0 (1.3)
Remark. The notation ( ) is an abusive notation. We should rather write ( ( )), where ( ) is the sole edge connected to .

INRIA 2 Discretization of the continuous problem

Discretization in space

Once the continuous problem under its variationnal form (1.3) is defined, one can discretize this problem by using ℙ finite elements with mass lumping. Moreover the fact that mass lumping allows lesser storage and improved computation time, one will see that this method also allows us to uncouple computation on open edges and computation on vertices. In the following we will take = 1.

On a given edge ∈ ℰ, we introduce a local space step Δ such that /Δ is a (great) positive integer, denoted . We introduce then, for 0 ⩽ ⩽ , the point , given by , = ( Δ )

Note that there are some couples ( , ) which give the same point ∈ ℝ : these points correspond to inner vertices of the graph.

Once we introduced discretization points, we can introduce on each edge ∈ ℰ the basis functions of ℙ 1 called (Φ , ) 0⩽ ⩽ such that:

• for any ∈ {1, . . . , }, Φ , ∘ is affine on the open segment

( -1 ( , -1 ), -1 ( , )
) ,

• for any ∈ {0, . . . , }, Φ , ( , ) = where is the Kronecker symbol.

We approximate our functions and by writing, on each edge ∈ ℰ,

| ( , ) = ∑ =0 , ( )Φ , ( ) and | ( ) = ∑ =0 , Φ , ( )
which can be rewritten under the condensed form, by denoting = ( ,0 . . . . , , ) for ∈ { , , Φ}:

( , ) = ∑ ∈ℰ Φ ( ) ( ) and
( ) = ∑ ∈ℰ Φ ( ) = ∑ ∈ℰ Φ ( ) (2.1)
We put up the vectorial approximations (2.1) in the variationnal formulation (1.

3), and one can see that two natural set of matrices appear:

• rigidity matrices ( ) ∈ℰ are given by

= ∫ ∇Φ ∇Φ = 1 Δ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 -1 0 . . . 0 -1 2 . . . . . . . . . 0 . . . . . . . . . 0 . . . . . . . . . 2 -1 0 . . . 0 -1 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (2.2) RT n°0369
A. Semin

• mass matrices ( ) ∈ℰ are given by

= ∫ Φ Φ = Δ 6 ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 2
⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (2.3)
Instead of using the mass matrices (2.3), we use mass lumping technic (we give here only for the order 1, the reader may see [3] for mass lumped matrices for any ℙ ). We set that

= Δ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 2 0 . . . . . . 0 0 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 0 0 . . . . . . 0 1 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (2.4)
Moreover, for any point ∈ , we introduce:

• ( ) the sole element belonging to ( ),

• the ( ) + 1 square matrix whose is identically equal to zero, except for: -the coefficient on the first line and first column if -1 ( ) = 0, -the coefficient on the last line and last column otherwise.

And finally, for any point ∈ (yes, even for inner vertices), for any ∈ ( ), we introduce the two following functions:

• 0 defined as 0( , ) = 0 if -1 ( ) = 0 and 0( , ) = otherwise,

• 1 defined as 1( , ) = 1 if -1 ( ) = Δ and 1( , ) = -1 otherwise.

Remark. For these two functions, we will omit the index when ∈ .

Discretization in time

For the time discretization, we introduce here a global time step Δ , and one denotes the approximation of ( ) at time = Δ . We discretize first order time derivative by using a centered scheme:

∂ ∂ ( Δ ) ⇒ +1 - -1 2Δ (2.5)
and second order time derivative by using classical scheme:

∂ 2 ∂ 2 ( Δ ) ⇒ +1 -2 + -1 Δ 2 (2.6) INRIA 2.

Writing the discete problem and resolution

Finally, using notations and matrices of section 2.1 with time discretizations (2.5) and (2.6) leads to the numerical scheme:

∑ ∈ℰ ( ) ( +1 -2 + -1 Δ 2 + ) + ∑ ∈ , ( ) ( Δ )+ ∑ ∈ , ( ) +1 ( ) - -1 ( ) 2Δ = 0
(2.7) under the two following constraints:

• continuity of the solution at the inner vertices -for any ∈ , for any , ′ ∈ ( ), one has

+1 ,0( , ) = +1 ′ ,0( , ′ ) (2.8)
• Dirichlet condition on some outer vertices: for any ∈ , , one has +1 ( ),0( ) = (( + 1)Δ )

(2.9)

Computation at initial time

One has to compute the solution at iterations = 0 and = 1 (corresponding respectively to times = 0 and = Δ ). Te be more precise about the "near node vanishing" Cauchy datas of the continuous problem (1.1), we suppose that:

Hypothesis 2.1. For any ∈ , for any ∈ ( ), for any ∈ { , }, one has

,0( , ) = ,1( , ) = 0 (2.10)
On the following, we treat separately the contribution from the Cauchy data and the contribution from the boundary conditions.

Contribution from the Cauchy data For 0 , we simply write that 0 = for each edge ∈ ℰ. For 1 , one starts for the D'Alembert formula for the solution of the 1D time-wave equation:

( , ) = 1 2 ( ( -) + ( + )) + 1 2 ∫ + - ( ) (2.11)
We can use this formula for any inner discretization point on each edge -one can see easily that 1 ( ) = 0 for any ∈ . By using our interpolation formula in (2.11), one gets (assuming that Δ ⩽ Δ for any edge ):

∀ ∈ ℰ, ∀1 ⩽ ⩽ -1, 1 , = Δ -Δ Δ , + Δ 2Δ ( , -1 + , +1 ) + Δ 2Δ ( Δ 2 , -1 + ( 2 - Δ Δ ) , + Δ 2 , +1
)

(2.12)

RT n°0369

A. Semin

Contribution from the boundary conditions For 0 , we simply take into account the Dirichlet part: ∀ ∈ , , 0 ( ),0( ) = (0) For 1 , we have a little more work to do. We have three parts to consider:

• Dirichlet part: ∀ ∈ , , 1
( ),0( ) = (Δ )
• propagation of the Dirichlet part: using (2.7) for = 0 (assuming that -1 is identically equal to 0 near vertices) gives that

∀ ∈ , , 1
( ),1( ) = Δ 2 Δ 2 (0)
Remark. Saying that -1 is identically equal to 0 is available if and only if the functions and are support in ℝ + and if the Cauchy data vanish near vertices.

• Neumann part: still using (2.7) for = 0, one gets that

∀ ∈ , , 1
( ),0( ) = - 2Δ 2 Δ (0)

Computation over time

For computation at other times, one simply starts from the numeric scheme (2.7). Knowing -1 and , it will be very simple, to compute +1 . The fact that each matric is a diagonal matrix let us allow to compute each discretization point independantly of the others.

Computation of inner edges discretization points One simply has

∀ ∈ ℰ, ∀1 ⩽ ⩽ -1, +1 , -2 , + -1 , Δ 2 - , +1 -2 , + , -1 Δ 2 = 0
Computation of inner vertices Given 0 ∈ . Because of condition (2.8), we denote +1 0 the common value of +1 at vertex 0 . One has, using the numeric scheme (2.7), and extracting parts corresponding to this vertex,

∑ ∈ ( ) ( )Δ 2Δ 2 ( +1 0 -2 0 + -1 0 ) + ( ) Δ ( 0 -,1( 0) ) = 0 (2.13) which gives +1 0 = 2 0 --1 0 + 2Δ 2 ∑ ∈ ( 0) Δ ( ,1( 0) -0 ) ∑ ∈ ( 0) Δ (2.14) INRIA

Computation of outer vertices

We have to treat Neumann and Outgoing conditions. Dirichlet condition is given explicitely by (2.9).

• Neumann condition: one gets, for any ∈ , :

Δ ( ) 2Δ 2 ( +1 ( ),0( ) -2 ( ),0( ) + -1 ( ),0( ) ) + ( ),0( ) -( ),1( ) Δ + ( Δ ) = 0 (2.15)
• Outgoing condition: one gets, for any ∈ , :

Δ ( ) 2Δ 2 ( +1 ( ),0( ) -2 ( ),0( ) + -1 ( ),0( ) ) + ( ),0( ) -( ),1( ) Δ + +1 ( ),0( ) - -1 ( ),0( ) 2Δ = 0
(2.16)

Implementation of the discrete problem

Here, we show the implementation of the numeric scheme dscribed in the numerical code netwaves. This program uses mainly two types of structures: a structure of geometrical graph, and a structure of problem. Actually, this code considers that Dirichlet and Neumann conditions are homogeneous conditions.

Remark. In this section and associated sub-sections, we'll introduce some classes and some methods (we won't introduce all the classes and methods used, one can check them by retrieving the sources and run Doxygen, available at http://www.stack.nl/~dimitri/doxygen/). We decided to prefix all members of our classes by an underscore, to differientiate with local variable in our class methods (in addition of using the object * this).

Geometrical graph structure

ª ¦ ¥

The two associatives maps _vertices and _edges respectively contains and ℰ. More precisely:

RT n°0369

• index of map _vertices is the label of vertices of the graph, and value correspond to coordinates of those vertices,

• index of map _edges is the label of edges of the graph, and value correspond to the two vertices linked with each edge.

To describe more quickly the following functions, we introduce the two following typedef: § t y p ed ef s t d : : map<i n t , s t d : : v e c t o r <r e a l > > t a b u l a r _ v e r t i c e s ; t y p ed ef s t d : : map<i n t , s t d : : v e c t o r <i n t > > t a b u l a r _ e d g e s ;

ª ¦ ¥

ª ¦ ¥

Each associative map (except the map _tabular_boudary_conditions which index is linked to index of the map _vertices of member _geom) are indexes linked to index of the map _edges of member _geom. We stock

• in map _hspaces the couples ( , Δ ), for ∈ ℰ,

• in map _relative_width the couples ( , ( )) for ∈ ℰ (we called this map as it because the conductance is linked to some geometrical parameters of our network),

• in maps _Unminus and _Un the solution computed on discretization points at iteration -1 and .

INRIA

Set up the problem

Compute initlal data

Assume that we only know the structure of the graph and we know the desired value of space step on each edge, named Δ . For each edge :

1. we retreive the coordinates of the two vertices associated with this edge, 2. we compute the length , 3. we compute = ⌊ / Δ ⌋ and then Δ = / . This ensures that Δ ⩾ Δ , and if the CFL condition is satisfied for the choice Δ , it will be satisfied for the choice Δ , 4. we create a pointer on a vector whose size is equal to + 1 and whose elements are equal to 0, ª ¦ ¥ 6. we compute the initial data by using formula (2.12) by using this method (thanks to hypothesis 2.

Computation over time

Assume that _Unminus and _Un contains the vectors -1 and . Then we compute +1 as described in section 2.3.2. The computation is the following (for the sequential form -one can see it would be easy to parallelize the code by using shared memory): §

ª ¦ ¥

Some remarks about this listing:

• on line 51, we get in the map tabular_edges_connected the set ( ), for each ∈ (we do not differientiate if ∈ or not),

• on lines 62-71, since the values -1 and does not depend on the edge, we look for these values on the first edge (in sense of the std :: less order associated to int) , and we write that

-1 = -1
,0( ) and = ,0( ) lines 68-71 are for the case 0( ) = 0, lines 63-66 are for the other case.

• on lines 124-133, we compute, for ∈ 0 , the values of 0( ) and 1( ), then we use this values in various boundary conditions implemented in lines 134-163.

Error estimates

In this section we give some error estimate associated to both continuous framework (section 1) and discrete framework (section 2). One has to define some energy, and then to check error estimates both theoretically and numerically.

Continuous framework

One can check in (1.3) that taking = ∂ ∂ leads to the relation

1 2 ∂ ∂ [ ∑ ∈ℰ ( ) ∫ ∂ ∂ ( , ⋅) 2 + |∇ ( , ⋅)| 2 ] + ∑ ∈ , ( ) ( ) ∂ ∂ ( , )+ ∑ ∈ , ( ) ∂ ∂ ( , ) 2 = 0
(4.1) From (4.1), one can define some continuous energy by the following:

INRIA Definition 4.1.
The energy of solution of (1.3) at time is defined by

( ) = 1 2 [ ∑ ∈ℰ ( ) ∫ ∂ ∂ ( , ⋅) 2 + |∇ ( , ⋅)| 2 ] (4.2)
The idea is to get some a priori error estimates. For this case, we suppose that there's only Neumann or Outgoing conditions (for the Dirichlet conditions, one has to introduce some function locally supported the graph satisfying the Dirichlet conditions, and to check up which problem is satisfied by -). We will suppose also that the Neumann functions belong to the set H 1 (ℝ * + ). One starts from (4.1), integrate over ∈ [0, ], and by using the fact that some terms are positive, one gets:

( ) -(0) ⩽ - ∑ ∈ , ( ) ∫ 0 ( ) ∂ ∂ ( , ) (4.3)
Using some intergration by parts in (4.3) and using that (0, ) = 0 for each ∈ , one gets

( ) -(0) ⩽ - ∑ ∈ , ( ) ( ) ( , ) + ∑ ∈ , ( )
∫ 0 ′ ( ) ( , ) (4.4)
To be able to conclude, one has to use the two following lemmas Lemma 4.1 (Poincare-Wirtinger inequality). For any ∈ , , by taking = ( ), there exists a constant which depends only of such that, for any

∈ ℝ + | ( , )| ⩽ ∥ ∥ H 1 ( ) ⩽ ( ∥ ( , ⋅)∥ L 2 ( ) + ∥∇ ( , ⋅)∥ L 2 ( ) )
Lemma 4.2. For any ∈ ℝ + , for any ∈ ℰ, one has

∥ ( , ⋅)∥ L 2 ( ) ⩽ 2 ( ∥ (0, ⋅)∥ L 2 ( ) + √ ∫ 0 ∂ ∂ ( , ⋅) L 2 ( ) )
Proof of lemma 4.2. One starts from

∥ ( , ⋅)∥ 2 L 2 ( ) = ∫ ( , ) 2 = ∫ ( , ) ( (0, ) + ∫ 0 ∂ ∂ ( ,
) ) ⩽ ∫ | ( , )| ⎛ ⎝ | (0, )| + √ ( ∫ 0 ∂ ∂ ( , ) 2
) 1/2 ⎞ ⎠ Then using Cauchy-Schwartz inequality gives the relation of the lemma.

RT n°0369

By using both lemmas 4.1 and 4.2, there exists two constant , ′ and " depending on time such that

( ) ⩽ (0) + ( ) + ∫ 0 ′ ( ) √ ( ) + "( ) √ ( ) (4.5)
Remark. One can easily see that all these constants are equal to zero when all Neumann conditions are homogeneous Neumann conditions. Now using the fact that

"( ) √ ( ) ⩽ ( "( )) 2 2 + 1 2 ( )
in (4.5) gives that there exist two constants ˜ ( ) and ′ ( ) (which is the same constant as before) such that

ℰ( ) ⩽ ˜ ( ) + ∫ 0 ′ ( ) √ ( ) (4.6)
Finally, using Gronwall inequality leads to the following a priori following majoration:

Theorem 4.3. There exist functions ˜ and ′ which depends on time and on boundary conditions of problem (1.1) such that, for all ∈ ℝ + ,

( ) ⩽ ( √ ˜ ( ) + 1 2 ∫ 0 ′ ( ) ) 2
where ˜ and ′ are the previously defined constants.

Remark. When all Neumann conditions are homogeneous Neumann conditiosn, theorem 4.3 gives that energy is bounded over time (see back the remark associated to equation (4.5)).

Discrete framework

Here we will suppose for convenience that we have only homogeneous Neumann conditions. We recall the numeric scheme (2.7) with our hypothesis:

∑ ∈ℰ ( ) ( +1 -2 + -1 Δ 2 + ) = 0
One has to mulltiply this line by +1 --1 , and one gets, by definying the following discrete energy

= ∑ ∈ℰ ( ) ( +1 - Δ ) ( +1 - Δ ) + ∑ ∈ℰ ( )( +1 ) (4.7)
INRIA the fact that this energy is constant over time. However, one can also check that not all the terms of energy defined by (4.7) are positive ones. To avoid that, one remarks that, thanks to the symmetry of each ,

( +1 ) = 1 2
( ( +1 ) +1 + ( ) ) - Δ 2 2 ( +1 - Δ ) ( +1 - Δ ) (4.8)
Using (4.8) in (4.7) let us allow to rewrite the discrete energy in a different way: one has

= ∑ ∈ℰ ( ) ( +1 - Δ ) ( - Δ 2 2 ) ( +1 - Δ ) + ∑ ∈ℰ ( ) 2
( ( +1 ) +1 + ( ) ) (4.9)
Let us define, for each edge ∈ ℰ, the two follwing matrices One has then, thanks to (4.10) and (4.11),

˜ = 2 Δ = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1
- Δ 2 2 = Δ 2 ( ˜ - Δ 2 Δ 2 ˜ ) (4.
12)

The advantage to introduce the matrices ˜ and ˜ is that their coefficients do not depend on the discretisation. Let us now resolve the generalised eigenvalue problem: find ( ,

) ∈ ℝ × ℝ such that ˜ = ˜ (4.13) RT n°0369
Then one can easily say that, for all

∈ ℝ , ( ˜ - Δ 2 Δ 2 ˜ ) ⩾ (1 - Δ 2 Δ 2 ) ˜ (4.14)
where is the largest eigenvalue of problem (4.13). Rouch calculus show that ⩽ 1. Then one can say that Theorem 4.4. The energy contains positive terms (and then numeric scheme is stable) if and only if

∀ ∈ ℰ, Δ Δ ⩽ 1
Remark. One can also check by computation that the numeric scheme is not stable (with exponential increasing of the computed solution) if the condition of theorem 4.4 is not satisfied on at least one edge.

Dealing with Improved Kirchhoff conditions

Instead of writing problem (1.1) with Kirchhoff conditions, we want to solve the following problem: find

∈ C 0 (ℝ + , H 1 d ( )) ∩ C 1 (ℝ + , H 0 ( )) such that: ⎧     ⎨     ⎩ ∂ 2 ∂ 2 -Δ = 0 in , ∀ ∈ ℰ and ∀ ∈ ℝ * + = for = 0 and ∀ ∈ ℰ ∂ ∂ = for = 0 and ∀ ∈ ℰ (5.1)
with some improved Kirchhoff conditions, introduced by P. Joly and A. Semin in [START_REF] Joly | Propagation of an acoustic wave in a junction of two thin slots[END_REF] for the particular case of two slots and in [START_REF] Joly | Construction and Analysis of Improved Kirchoff Conditions for Acoustic Wave Propagation in a Junction of Thin Slots[END_REF] for the general case (which interest us). However, we have to define some notations. In particular, the notation here will differs from the notation introduced in section 2.3.2.

Improved Kirchhoff conditions

We first have to define an order function, which says us how we locally number our edges.

Definition 5.1. An order function on the graph is a function defined on × ( ) = { , ( )} ∈ with values in ℕ * such that for any ∈ , ( , ⋅) is a bijection between ( ) and {1, . . . , # ( )}. We will then denote -1 ( , ⋅) : {1, . . . , # ( )} ∈ ℰ( ) the reciprocal bijection.

Once we defined our order function , for any ∈ , one can denote: this vector is not proportionnal to the vector (1, . . . , 1) , since the function is not continous on each vertex ∈ .

• ∂ ( ) the vector in ℝ # ( ) ( ) defined by

∂ ( ) = ⎛ ⎜ ⎜ ⎝ ∂ ∂ -1 ( ,1), ( , ) . . . ∂ ∂ -1 ( ,# ( )), ( , ) ⎞ ⎟ ⎟ ⎠ (5.3)
Improved conditions will be given as follow: for any ∈ , there exist two matrices and in ℳ # ( ) (ℝ) such that

∂ ( ) = ( + ∂ 2 ∂ 2 ) ( ) (5.4)
The problem with Improved Kirchhoff conditions is (5.1) with (5.4), i.e.: find

∈ C 0 (ℝ + , H 1 d ( ))∩ C 1 (ℝ + , H 0 ( )) such that: ⎧         ⎨         ⎩ ∂ 2 ∂ 2 -Δ = 0 in , ∀ ∈ ℰ and ∀ ∈ ℝ * + ∂ ( ) = ( + ∂ 2 ∂ 2 ) ( ),
∈ C 0 (ℝ + , H 1 d ( )) ∩ C 1 (ℝ + , H 0
( )) satisfying inhomogeneous Dirichlet condition such that, for any ∈ H 1 d ( ) satisfying associated homogeneous Dirichlet, one has (denoting the vector for defined as (5.2) for ):

∑ ∈ℰ ∫ ( ) ∂ 2 ∂ 2 ( , ⋅) + ( )∇ ( , ⋅)∇ + ∑ ∈ , ( ) ( ) ( ) + ∑ ∈ ( ) + ∂ 2 ∂ 2 ( ) + ∑ ∈ , ( ) ∂ ∂ ( , ) ( ) = 0 (5.6) RT n°0369

Discretization

The discrete problem associated to the new variationnal formulation (5.6) is the same as the discrete problem described in section 2.3, except for inner vertices which are dealt as below: for any ∈ , one introduces

• vector is defined as the vector containing the values of on each discretization point associated to vertex on each edge in ( ), i.e. = ⎛ ⎜ ⎝ -1 ( ,1),0( ) . . . Then, for inner vertex, instead of computing (2.14), one has to solve this system:

1 2Δ 2 Δ ( +1 -2 + -1 ) + Δ -1 ( -,Δ )
+ +1 + 2 + -1 4
+ +1 -2 + -1 Δ 2 = 0
(5.7)

INRIA

Remark. One can check in (5.7) that the term associated to has been treated implicitely. This is not classical, and this is due to the fact that if we do not implicit this term, there might be some resonance frequencies (linked to the amplitude of matrix ) that will corrupt CFL condition of theorem 4.4 -the case for two slots has been proved in [START_REF] Semin | Propagation d'ondes dans des jonctions de fentes minces[END_REF].

Implementation

In our class struct_problem, we add some associative maps, whose indexes will be indexed on indexes of _geom->_vertices. § template<c l a s s g eo metry , c l a s s r e a l > c l a s s s t

ª ¦ ¥

The class Matrix is the class of matrices of the external library newmat, written by R.B. Davies (see [4] for instance).

For initialization of the problem, the matrices _J and _A are read in a data file, as the vector _order, and one checks that the dimensions of those objects match with the number of edges connected to associated vertex. Moreover, one checks (when building vector _order), that there's no element associated with two differents keys, with the following code: § s t d : : map<i n t , i n t > o r d e r 2 ; s t d : : v e c t o r <i n t > c o n s t _ i t e r a t o r i t = t h i s ->_order->f i n d ( v e r t e x _ l a b e l )->b e g i n ( ) , i t e = t h i s ->_order->f i n d ( v e r t e x _ l a b e l )->b e g i n ( ) f o r ( ; i t != i t e ; ++i t ) { o r d e r 2 [ i t ->s e c o n d ] = i t -> f i r s t ; } i f ( _ t h i s->o r d e r ->s i z e ( ) != o r d e r 2 . s i z e ( ) ) { // I n t h i s c a s e , t h e r e was two d i f f e r e n t s i t e r a t o r s w i t h // t h e same v a l u e . We e x i t w i t h f a i l e x i t ( EXIT_FAILURE ) ;

RT n°0369

A. Semin

Conclusion

We showed the existence of a numeric scheme solving propagation of acoustic wave in networks. Moreover, we can say that this numeric scheme is stable under some classical conditions, and implementation of this scheme is easy -one does not need to use matrices to implement this scheme. However we can see that the time step used is a global time stepone expects to use a local time step, in the spirit of the works of J. Rodriguez.

Figure 1 . 1 :

 11 Figure 1.1: On the left: an example of set in ℝ 2 . On the right: an example of set ℰ associated to .

INRIA Definition 1 . 2 (

 12 Parametrization). Given = ( , ) ∈ ℰ, we introduce the parametrisation :]0, [→ , where is the length (in the sense of the Euclidian norm) of , by

3. 2

 2 Problem structure For the problem structure, one defines also a template class struct_problem based on a float point class and on a geometrical graph class. The complete syntax is: § template<c l a s s g eo metry , c l a s s r e a l > c l a s s s t r u c t _ p r o b l e m { p r o t e c t e d : g e o m e t r y * _geom ; s t d : : map<i n t , r e a l > * _ hspaces s t d : : map<i n t , r e a l > * _ r e l a t i v e _ w i d t h s t d : : map<i n t , s t d : : v e c t o r <r e a l > > * _Unminus ; s t d : : map<i n t , s t d : : v e c t o r <r e a l > > * _Un ; s t d : : map<i n t , i n t > * _ t a b u l a r _ b o u d a r y _ c o n d i t i o n s ; }

5 .

 5 we set up in our problem structure the space step Δ and we initialize vectors _Unminus and _Un by the following instructions: § 1 ( * t h i s ->_ hspaces ) [ l a b e l _ e d g e ] = delta_x_e ; 2 ( * t h i s ->_Unminus ) [ l a b e l _ e d g e ] = * U ; 3 ( * t h i s ->_Un) [ l a b e l _ e d g e ] = * U ;

1 ): § 1

 11 template<typename struct_ g eo m , c l a s s r e a l > 2 v o i d s t r u c t _ p r o b l e m<struct_ g eo m , r e a l > : : i n i t _ d a t a ( 3 c o n st r e a l dt , 4 r e a l ( * p o i n t e r _ f ) ( r e a l , r e a l , r e a l ) , 5 r e a l ( * p o i n t e r _ g ) ( r e a l , r e a l , r e a l ) ) 6 { 7 s t d : : v e c t o r <r e a l > c o o r d i n a t e s ( 3 ) ; 8 s t d : : v e c t o r <r e a l > c o o r d i n a t e s _ p r e v i o u s ( 3 ) ; 9 s t d : : v e c t o r <r e a l > c o o r d i n a t e s _ n e x t ( 3 ) ; 10 11 // Loop o v e r e d g e s 12 13 typename t a b u l a r _ e d g e s : : c o n s t _ i t e r a t o r 14 e d g e _ i t e r = t h i s ->_geom->g e t _ e d g e s ( )->b e g i n ( ) , 15 edge_end = t h i s ->_geom->g e t _ e d g e s ( )->end ( ) ; 16 A. Semin f o r ( ; e d g e _ i t e r != edge_end ; ++e d g e _ i t e r ) { // Get i n f o r m a t i o n s a b o u t edg e i n t l a b e l _ e d g e = e d g e _ i t e r -> f i r s t ; i n t l a b e l _ p o i n t 1 = ( e d g e _ i t e r ->s e c o n d ) [ 0 ] ; i n t l a b e l _ p o i n t 2 = ( e d g e _ i t e r ->s e c o n d ) [ 1 ] ; // R e t r i e v e v e r t i c e s typename t a b u l a r _ v e r t i c e s : : c o n s t _ i t e r a t o r v e r t e x 1 = t h i s ->_geom->g e t _ v e r t i c e s ( )->f i n d ( l a b e l _ p o i n t 1 ) , v e r t e x 2 = t h i s ->_geom->g e t _ v e r t i c e s ( )->f i n d ( l a b e l _ p o i n t 2 ) ; // R e t r i e v e number o f unknowns and a s s o c i a t e D e l t a x_e i n t N = t h i s ->_Unminus->f i n d ( l a b e l _ e d g e )->s e c o n d . s i z e ( ) ; r e a l dx = t h i s ->_hspaces->f i n d ( l a b e l _ e d g e )->s e c o n d ; // Compute l e n g t h o f edg e r e a l l e n g t h = 0 . ; f o r ( i n t j =0; j <3; j ++) { l e n g t h +=(( v e r t e x 1 ->s e c o n d ) [ j ] -( v e r t e x 2 ->s e c o n d ) [ j ] ) * ( ( v e r t e x 1 ->s e c o n d ) [ j ] -( v e r t e x 2 ->s e c o n d ) [ j ] ) ; } l e n g t h = s t d : : s q r t ( l e n g t h ) ; // I n i t i a l i z e v e c t o r s s t d : : v e c t o r <r e a l > * U1 = new s t d : : v e c t o r <r e a l >(N , 0 . ) , } r e a l Fn = ( * p o i n t e r _ f ) ( c o o r d i n a t e s [ 0 ] , c o o r d i n a t e s [ 1 ] , c o o r d i n a t e s [ 2 ] ) ; r e a l Fnm = ( * p o i n t e r _ f ) ( c o o r d i n a t e s _ p r e v i o u s [ 0 ] , c o o r d i n a t e s _ p r e v i o u s [ 1 ] , c o o r d i n a t e s _ p r e v i o u s [ 2 ] ) ; r e a l Fnp = ( * p o i n t e r _ f ) ( c o o r d i n a t e s _ n e x t [ 0 ] , c o o r d i n a t e s _ n e x t [ 1 ] , c o o r d i n a t e s _ n e x t [ 2 ] ) ; r e a l Gn = ( * p o i n t e r _ g ) ( c o o r d i n a t e s [ 0 ] , c o o r d i n a t e s [ 1 ] , c o o r d i n a t e s [ 2 ] ) ; r e a l Gnm = ( * p o i n t e r _ g ) ( c o o r d i n a t e s _ p r e v i o u s [ 0 ] , c o o r d i n a t e s _ p r e v i o u s [ 1 ] , c o o r d i n a t e s _ p r e v i o u s [ 2 ] ) ; r e a l Gnp = ( * p o i n t e r _ g ) ( c o o r d i n a t e s _ n e x t [ 0 ] , c o o r d i n a t e s _ n e x t [ 1 ] , c o o r d i n a t e s _ n e x t [ 2 ] ) ; ( * U1 ) [ e v a l p o i n t ] = Fn ; ( * U2 ) [ e v a l p o i n t ] = ( dx-d t ) / dx * Fn+d t / ( 2 . * dx ) * (Fnm+Fnp ) + 0 . 5 * d t * ( 2 . -d t / dx ) * Gn + 0 . 5 * d t * d t / (2 * dx ) * (Gnm+Gnp ) ; } ( * t h i s ->_Unminus ) [ l a b e l _ e d g e ] = * U1 ; ( * t h i s ->_Un) [ l a b e l _ e d g e ] = * U2 ; U1->c l e a r ( ) ; d e l e t e U1 ; U2->c l e a r ( ) ; d e l e t e U2 ; } } ª ¦ ¥

1

  template<typename struct_ g eo m , c l a s s r e a l > 2 v o i d s t r u c t _ p r o b l e m<struct_ g eo m , r e a l > : : r e s o l v e ( c o n st r e a l d t ) e t up a new map ( t h a t w i l l be U^(n+1) ) s t d : : map<i n t , s t d : : v e c t o r <r e a l > > * Up = new s t d : : map<i n t , s t d : : v e c t o r <r e a l > >( * ( t h i s ->_Un) ) ; // S e t up map o f s p a c e s t e p s s t d : : map<i n t , i n t > * Ne = new s t d : : map<i n t , i n t >; // F i r s t p a r tc o m p u t a t i o n o v e r e d g e s typename t a b u l a r _ e d g e s : : c o n s t _ i t e r a t o r e d g e _ i t e r = t h i s ->_geom->g e t _ e d g e s ( )->b e g i n ( ) , edge_end = t h i s ->_geom->g e t _ e d g e s ( )->end ( ) ; f o r ( ; e d g e _ i t e r != edge_end ; ++e d g e _ i t e r ) { e d g e _ l a b e l = e d g e _ i t e r -> f i r s t ; r e a l dx = t h i s ->_hspaces->f i n d ( e d g e _ l a b e l )->s e c o n d ; r e a l c f l = d t / dx ; assume ( c f l <= 1 . ) ; ( * Ne ) [ e d g e _ l a b e l]= i n t ( ( * ( t h i s ->_Un) ) [ e d g e _ l a b e l ] . s i z e ( ) ) -1; f o r ( i n t m=1; m < ( * Ne ) [ e d g e _ l a b e l ] ; m++) { ( * Up ) [ e d g e _ l a b e l ] [ m]=2 * (1 -c f l * c f l ) * ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ m] + c f l * c f l * ( ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ m+1] +( * t h i s ->_Un) [ e d g e _ l a b e l ] [ m-1] ) -( * t h i s ->_Unminus ) [ e d g e _ l a b e l ] [ m ] ; } } // Second p a r tc o m p u t a t i o n o v e r v e r t i c e s // S e t up t a b u l a r o f e d g e s c o n n e c t e d t o a g i v e n v e r t e x // C o n v e n t i o n o f b o o l e a n : t r u e i f v e r t e x i s t h e r i g h t p o i n t , // f a l s e i f v e r t e x i s t h e l e f t p o i n t s t d : : map<i n t , bool> t a b u l a r _ e d g e s _ c o n n e c t e d ; typename t a b u l a r _ v e r t i c e s : : c o n s t _ i t e r a t o r v e r t e x _ i t e r = t h i s ->_geom->g e t _ v e r t i c e s ( )->b e g i n ( ) , v e r t e x _ e n d = t h i s ->_geom->g e t _ v e r t i c e s ( )->end ( ) ; INRIA RT n°0369 { // Outer v e r t e x e d g e _ l a b e l = t a b u l a r _ e d g e s _ c o n n e c t e d . b e g i n ( )-> f i r s t ; i n t BC = t h i s ->_ t a b u l a r _ b o u n d a r y _ c o n d i t i o n s->f i n d ( v e r t e x _ l a b e l )->seco nd , i n t z e r o , one ; i f ( t a b u l a r _ e d g e s _ c o n n e c t e d . b e g i n ( )->s e c o n d ) { z e r o = ( * n b _ i t e r ) [ e d g e _ l a b e l ] ; one = ( * n b _ i t e r ) [ e d g e _ l a b e l ] -1 ; } e l s e { z e r o = 0 ; one = 1 ; } s w i t c h (BC) { c a se b c _ d i r i c h l e t : ( * Up ) [ e d g e _ l a b e l ] [ z e r o ] = 0 ; break ; c a se bc_neumann : ( * Up ) [ e d g e _ l a b e l ] [ z e r o ] = 2 * ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ z e r o ] + ( 2 * d t * d t / dx / dx ) * ( ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ one ] -( * t h i s ->_Un) [ e d g e _ l a b e l ] [ z e r o ] ) -( * t h i s ->_Unminus ) [ e d g e _ l a b e l ] [ z e r o ] ; break ; c a se bc_ o utg o ing : r e a l c1 , c2 , c3 ; c1 = dx / ( 2 . * d t * d t ) + 1 . / d t ; c2 = dx / ( d t * d t ) ; c3 = -dx / ( 2 . * d t * d t ) + 1 . / d t ; ( * Up ) [ e d g e _ l a b e l ] [ z e r o ] = ( c2 * ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ z e r o ] + c3 * ( * t h i s ->_Unminus ) [ e d g e _ l a b e l ] [ z e r o ] + ( ( * t h i s ->_Un) [ e d g e _ l a b e l ] [ one ] -( * t h i s ->_Un) [ e d g e _ l a b e l ] [ z e r o ] ) ) / c1 ; break ; RT n°0369 A. Semin d e f a u l t : s t d : : c e r r << "ERROR␣ : : ␣ r e s o l v e ( ) ␣ : : ␣ C o n d i t i o n ␣ " << BC << " n o t ␣ y e t ␣ i m p l e m e n t e d " << s t d : : e n d l ; s t d : : c e r r << " A b o r t i n g . . . " << s t d : : e n d l ; e x i t ( EXIT_FAILURE ) ;

  0 . . . . . .

(

  , ) defined on the edge -1 ( , 1) . . . ( , ) defined on the edge -1 ( , # ( ))

  classical conditions (Dirichlet, Neumann, Outgoing). The variationnal formulation is: find

- 1 (

 1 ,# ( )),0( ) ⎞ ⎟ ⎠• vector ,Δ is defined as the vector containing the values of on each discretization point next to discretization points associated to vertex on each edge in ( ), i.e.

  r u c t _ p r o b l e m { p r o t e c t e d : // . . . p r e v o u s l y d e f i n e d members s t d : : map<i n t , s t d : : v e c t o r <i n t > > _ o rder ; s t d : : map<i n t , M a t r i x * > _J ; s t d : : map<i n t , M a t r i x * > _A; }

Acknolowedgment I would like to thanks P. Joly for his knowledge which helped me in some technical points, and K. Boxberger with helped me to provide the code netwaves. INRIA


local str="${1}" if [[ -z "$str" ]]; then return 1 # error else # return by var ret="$str" # return by stdout printf "$str" return 0 # OK fi } printf $(g "First") if g "Second"; then printf "$ret" fi #!/bin/bash cmd x=$(cmd) x=$((1+2))

Shell calls:

run this script with bash call cmd capture cmd standard output integer arithmetics

Flavors of UNIX:

-Linux -FreeBSD -Darwin (MacOSX) -Cygwin -WSL (Windows) -…

UNIX Shell Cheatsheet

Flavors of Shell:

-Bash -Zsh -Fossil -…

UNIX wisdom: everything is a file.

Shell wisdom: everything is a text stream.

Script wisdom: everything is a string

-cmd [options] ARG [OPTARG] -Short: -h -v [or] -al -Long: --help --verbose -Arguments: --opt[= ]<value>

Options convention:

Examples: $ ls -a -l .

string is empty if [[ -w "$fname" ]]; then # file is writable fi -r -f -d -x File conditions: is readable is a file is a directory is executable $0 $1 … $* $# $IFS

Special variables:

program name 1 st argument … all arguments number of arguments Internal File Separartor for x in "$X"; do # iterate over split X done for ((i=0; i<10; i++)); do # numeric counter i done

  detachable shell(s) detach the shell(s) get back the detached shell(s) background=depends on the session detached =does not depends on session alias f="cmd -opt" function f

  ' read -r line || [[ -n "$line" ]]; do printf "$line\n" done<"$filename"


Nicolas Verzelen
email: nicolas.verzelen@math.u-psud.fr

Technical appendix to "Adaptive estimation of stationary Gaussian fields"

Keywords: AMS 2000 subject classifications: Primary 62H11; secondary 62M40 Gaussian field, Gaussian Markov random field, model selection, pseudolikelihood, oracle inequalities, Minimax rate of estimation

.

Proof of Proposition 8.1

Proof of Proposition 8.1. First, we recall the notations introduced in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF]. Let N be a positive integer. Then, I N stands for the family of subsets of {1, . . . , N } of size less than 2. Let T be a set of vectors indexed by I N . In the sequel, T is assumed to be a compact subset of R (N (N +1)/2)+1 . The following lemma states a slightly modified version of the upper bound in remark 7 in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF].

U i U j t {i,j} + N i=1 t {i} + t ∅ ,
where U 1 , . . . , U N are independent Rademacher random variables. Then for any x > 0,

P {T ≥ E[T ] + x} ≤ 4 exp - x 2 L 1 E[D] 2 ∧ x L 2 E , (1)
where D and E are defined by:

D := sup t∈T sup α: α 2 ≤1 N i=1 U i j =i
α j t {i,j} , E := sup t∈T sup α (1) ,α (2) , α (1) 2≤1 α (2) ≤1 N i=1 j =i t {i,j} α

i α

(2) j

.

Contrary to the original result of [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF], the chaos are not assumed to be homogeneous. Besides, the t {i} are redundant with t ∅ . In fact, we introduced this family in order to emphasize the connection with Gaussian chaos in the next result.

A suitable application of the central limit theorem enables to obtain a corresponding bound for Gaussian chaos of order 2.

Lemma 1.2. Let T be a supremum of Gaussian chaos of order 2.

T := sup t∈T {i,j} t {i,j} Y i Y j + i t i Y 2 i + t ∅ , (2)
where Y 1 , . . . , Y N are independent standard Gaussian random variable. Then, for any x > 0,

P {T ≥ E[T ] + x} ≤ exp - x 2 E[D] 2 L 1 ∧ x EL 2 , (3)
where

D := sup t∈T sup α∈R N α 2 ≤1 i,j Y i (1 + δ i,j )α j t {i,j} , E := sup t∈T sup α1, α1 2 ≤1
sup α2, α2 2 ≤1 i,j α 1,i α 2,j t {i,j} (1 + δ i,j ) .

The proof of this Lemma is postponed to the end of this section. To conclude, we derive the result of Proposition 8.1 from this last lemma. For any matrix R ∈ F , we define the vector t R ∈ R nr(nr+1)/2+1 indexed by I nr as follows t R {(i,k),(j,l)} := δ k,l (2δ i,j )

R[i,j] n , t R {(i,
Let us now turn the constant

E E = sup t R ∈T sup α 1 , α 2 ∈ R nr α 1 2 ≤ 1, α 2 2 ≤ 1 1≤i,j≤r 1≤k,l≤n (1 + δ ij δ k,l )t R,kl i,j α k 1,i α l 2,j = sup R∈F sup α 1 , α 2 ∈ R nr α 1 2 ≤ 1, α 2 2 ≤ 1 2 n 1≤i,j≤r 1≤k≤n R[i,j]α k 1,i α k 2,j .
From this last expression, it follows that E is a supremum of L 2 operator norms

E = 2 n sup R∈F ϕ max Diag (n) (R) ,
where Diag (n) (R) is the (nr×nr) block diagonal matrix such that each diagonal block is made of the matrix R. Since the largest eigenvalue of Diag (n) (R) is exactly the largest eigenvalue of R, we get

E = 2 n sup R∈F ϕ max (R) . (5)
Applying Proposition 1.2 and gathering identities (4) and [START_REF] Moran | A Gaussian Markovian process on a square lattice[END_REF] yields

P(Z ≥ E(Z) + t) ≤ exp - t 2 L 1 E(V ) t L 2 B
, where B = E and V = D 2 .

Proof of Lemma 1.1. This result is an extension of Corollary 4 in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF]. We shall closely follow the sketch of their proof adapting a few arguments. First, we upper bound the moments of (T -E(T )) + . Then, we derive the deviation inequality from it. Here, x + = max(x, 0).

Lemma 1.3. For all real numbers q ≥ 2, (T -E(T )) + q ≤ LqE(D) + LqE , (6)
where T q q stands for the q-th moment of the random variable T . The quantities D and E are defined in Lemma 1.1.

By Lemma 1.3, for any t ≥ 0 and any q ≥ 2,

P (T ≥ E(T ) + t) ≤ E (T -E(T )) q + t q ≤ √ LqE(D) + LqE t q .
The right-hand side is at most 2 -q if √ LqE(D) ≤ t/4 and LqE ≤ t/4. Let us set

q 0 := t 2 16LE(D) 2 ∧ t 4LE
.

If q 0 ≥ 2, then P (T ≥ E(T ) + t) ≤ 2 -q0
. On the other hand if q 0 < 2, then 4 × 2 -q0 ≥ 1. It follows that

P (T ≥ E(T ) + t) ≤ 4 exp - log(2) 4L t 2 4E(D) 2 ∧ t E .
Proof of Lemma 1.3. This result is based on the entropy method developed in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF]. Let f : R N → R be a measurable function such that T = f (U 1 , . . . , U N ). In the sequel, U ′ 1 , . . . , U ′ N denote independent copies of U 1 , . . . , U N . The random variable T ′ i and V + are defined by

T ′ i := f (U 1 , . . . , U i-1 , U ′ i , U i+1 , . . . , U N ) , V + := E N i=1 (T -T ′ i ) 2 + |U N 1 ,
where U N 1 refers to the set {U 1 , . . . , U N }. Theorem 2 in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF] states that for any real q ≥ 2,

(T -E(T )) + q ≤ Lq √ V + q . (7)
To conclude, we only have bound the moments of √ V + . By definition,

T = sup t∈T {i,j} U i U j t {i,j} + N i=1 t {i} + t ∅ .
Since the set T is compact, this supremum is achieved almost surely at an element t 0 of T . For any 1

≤ i ≤ N , (T -T ′ i ) 2 + ≤ (U i -U ′ i ) j =i U j t 0 {i, j} 2 .
Gathering this bound for any i between 1 and N , we get

V + ≤ N i=1 E   (U i -U ′ i ) j =i U j t 0 {i, j} 2 U N 1   ≤ 2 N i=1 j =i U j t 0 {i, j} 2 ≤ 2 sup α∈R N , α 2 ≤1 N i=1 α i j =i t 0 {i,j} U j 2 ≤ 2 sup t∈T sup α∈R N , α 2 ≤1 N i=1 U i j =i α j t {i,j} 2 = 2D 2 .
Combining this last bound with ( 7) yields

(T -E(T )) + q ≤ Lq √ 2 D q ≤ Lq E(D) + |(D -E(D)) + q . ( 8
)
Since the random variable D defined in Lemma 1.1 is a measurable function f 2 of the variables U 1 , . . . , U N , we apply again Theorem 2 in [START_REF] Boucheron | Moment inequalities for functions of independent random variables[END_REF].

(D -E(D)) + q ≤ Lq V + 2 q
, where V + 2 is defined by

V + 2 := E N i=1 (D -D ′ i ) 2 + U N i , and
D ′ i := f 2 (U 1 , . . . , U i-1 , U ′ i , U i+1 , . . . , U N ).
As previously, the supremum in D is achieved at some random parameter (t 0 , α 0 ). We therefore upper bound V + 2 as previously.

V + 2 ≤ N i=1 E   (U i -U ′ i ) j =i α 0 j t 0 {i,j} 2 U N 1   ≤ 2 N i=1 j =i α 0 j t 0 {i,j} 2 ≤ 2 sup α (2) ∈R N , α 2≤1 N i=1 α (2) j j =i α 0 i t {i,j} 2 = 2E 2 .
Gathering this upper bound with (8) yields

(T -E(T )) + q ≤ LqE(D) + LqE .
Proof of Lemma 1.2. We shall apply the central limit theorem in order to transfer results for Rademacher chaos to Gaussian chaos. Let f be the unique function satisfying T = f (y 1 , . . . , y N ) for any (y 1 , . . . , y N ) ∈ R N . As the set T is compact, the function f is known to be continuous. Let (U (j) i ) 1≤i≤N,j≥0 an i.i.d. family of Rademacher variables. For any integer n > 0, the random variables Y (n) and T (n) are defined by

Y (n) := n j=1 U (j) 1 √ n , . . . , n j=1 U (j) N √ n , T (n) := f Y (n) .
Clearly, T (n) is a supremum of Rademacher chaos of order 2 with nN variables and a constant term. By the central limit theorem, T (n) converges in distribution towards T as n tends to infinity. Consequently, deviation inequalities for the variables T (n) transfer to T as long as the quantities E D (n) , E (n) , and E[T (n) ] converge.

We first prove that the sequence T (n) converges in expectation towards T . As T (n) converges in distribution, it is sufficient to show that the sequence T (n) is asymptotically uniformly integrable. The set T is compact, thus there exists a positive number t ∞ such that

T (n) ≤ t ∞ i,j |Y (n) i Y (n) j | + 1 ≤ t ∞ 1 + (N + 1)/2 N i=1 Y (n) i 2 .
It follows that

T (n) 2 ≤ t 2 ∞ N + 1 2 2 N + 2 2 1 + N i=1 Y (n) i 4 . (9)
The sequence Y

(n) i

does not only converge in distribution to a standard normal distribution but also in moments (see for instance [START_REF] Billingsley | Probability and measure[END_REF] p.391). It follows that limE T (n) 2 ≤ ∞ and the sequence f Y (n) is asymptotically uniformly integrable. As a consequence,

lim n→∞ E T (n) = E[T ] .
Let us turn to the limit of E D (n) . As the variable T (n) equals

T (n) = sup t∈T {i,j} t {i,j} 1≤k,l≤n U (k) i U (l) j n + i t i 1≤k≤n U (k) i √ n l =k U (l) i √ n + t ∅ + i t i ,
it follows that

D (n) = sup t∈T sup α∈R nN , α 2 ≤1 1≤i≤N 1≤k≤n U (k) i j =i t {i,j} n 1≤l≤n α (l) j + 2 l =k 2 t {i} n α (l) i ≤ sup t∈T sup α∈R nN , α 2 ≤1 i U (k) i √ n j (1 + δ i,j )t {i,j} 1≤l≤n α (l) j √ n + A (n) , (10)
where the random variable A (n) is defined by

A (n) := sup t∈T sup α∈R nN , α 2≤1 N i=1 n j=1 t {i} U (j) i n α j i .
Straightforwardly, one upper bounds

A (n) by t ∞ /n N i=1 n j=1 U (j) i 2
and its expectation satisfies

E A (n) ≤ t ∞ N n ,
which goes to 0 when n goes to infinity. Thus, we only have to upper bound the expectation of the first term in (10). Clearly, the supremum is achieved only when for all 1 ≤ j ≤ N , the sequence (α (l) j ) 1≤l≤n is constant. In such a case, the sequence (α

(1) j ) 1≤j≤N satisfies α (1) 2 ≤ 1/ √ n. it follows that E D (n) = E sup t∈T sup α∈R N α 2 ≤1 E i Y (n) i j (1 + δ i,j )α j + O 1 √ n .
Let g be the function defined by

g (y 1 , . . . , y N ) = sup t∈T sup α∈R N α 2 ≤1 i y i j (1 + δ i,j )α j ,
for any (y 1 , . . . , y N ) ∈ R N . The function g(.) is measurable and continuous as the supremum is taken over a compact set. As a consequence, g(Y (n) ) converges in distribution towards g(Y ). As previously, the sequence is asymptotically uniformly integrable since its moment of order 2 is uniformly upper bounded. It follows that lim

E D (n) = E [D].
Third, we compute the limit of E (n) . By definition,

E (n) = sup t∈T sup α1,α2∈R nN , α1 2≤1, α2 2 ≤1 N i=1 n k=1 α k 1,i j =i n l=1 α (l) 2,j t {i,j} n + 2 l =k α (l) 2,i t {i} n = sup t∈T sup α1,α2, α1 2 ≤1, α2 2≤1 N i=1 N j=1 (1 + δ i,j ) t {i,j} n n k=1 n l=1 α (k) 1,i α (l) 2,j + O 1 n .
As for the computation of D (n) , the supremum is achieved when the sequences (α k 1,i ) 1≤k≤n and (α l 2,j ) 1≤l≤n are constant for any i ∈ {1, . . . , N }. Thus, we only have to consider the supremum over the vectors α 1 and α 2 in R N .

E (n) = sup t∈T sup α1,α2∈R N αi 2≤1 N i=1 N j=1 (1 + δ ij )t i,j α 1,i α 2,j + O 1 n .
It follows that E (n) converges towards E when n tends to infinity.

The random variable

T (n) -E(T (n) ) converges in distribution towards T -E(T ).
By Lemma 1.1 ,

P(T -E(T ) ≥ x) ≤ lim exp - x 2 E[D (n) ] 2 L 1 ∧ x E (n) L 2 ,
for any x > 0. Combining this upper bound with the convergence of the sequences D (n) and E (n) allows to conclude.

Proof of Theorem 3.1

Proof of Lemma 8.3. We only consider here the anisotropic case, since the isotropic case is analogous. This result is based on the deviation inequality for suprema of Gaussian chaos of order 2 stated in Proposition 8.1. For any model m ′ belonging to M, we shall upper bound the quantities E(Z m ′ ), B m ′ , and E(W m ′ ) defined in (42) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

1. Let us first consider the expectation of Z m ′ . Let U ′ m,m ′ be the new vector space defined by

U ′ m,m ′ := U m,m ′ √ D Σ p ,
where U m,m ′ is introduced in the proof of Lemma 8.2 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. This new space allows to handle the computation with the canonical inner product in the space of matrices. Let B

(2)

m 2 ,m ′2 be the unit ball of U ′ m,m ′ with respect to the canonical inner product. If R belongs to U m,m ′ , then R H ′ = R √ D Σ /p F ,
where . F stands for the Frobenius norm.

Z m ′ = sup R∈B H ′ m 2 ,m ′2 1 p 2 tr RD Σ (YY * -I p 2 ) = sup R∈B (2) m 2 ,m ′2 tr R √ D Σ p YY * -I p 2 (11) = Π U ′ m,m ′ √ D Σ p YY * -I p 2 F ,
where Π U ′ m,m ′ refers to the orthogonal projection with respect to the canonical inner product onto the space

U ′ m,m ′ . Let F 1 , . . . , F d m 2 ,m ′2 de- note an orthonormal basis of U ′ m,m ′ . E(Z 2 m ′ ) = d m 2 ,m ′2 i=1 E tr 2 F i D Σ p 2 YY * -I p 2 = d m 2 ,m ′2 i=1 E p 2 j=1 F i [j,j] √ D Σ [j,j] p (YY * [j,j] -1) 2 = d m 2 ,m ′2 i=1 2 np 2 tr(F i D Σ F i ) ≤ d m 2 ,m ′2 i=1 2ϕ max (D Σ ) np 2 = 2d m 2 ,m ′2 ϕ max (Σ) np 2 .
Applying Cauchy-Schwarz inequality, it follows that

E(Z m ′ ) ≤ 2d m 2 ,m ′2 ϕ max (Σ) np 2 . ( 12
)
2. Using the identity (11), the quantity B m ′ equals

B m ′ = 2 n sup R∈B (2) m 2 ,m ′2 ϕ max R √ D Σ p .
As the operator norm is under-multiplicative and as it dominates the Frobenius norm, we get the following bound

B m ′ ≤ 2 ϕ max (Σ) np . (13)
3. Let us turn to bounding the quantity E(W m ′ ). Again, by introducing the ball B

(2)

m 2 ,m ′2 , we get W m ′ = 4 n sup R∈B H ′ m 2 ,m ′2 1 p 2 tr RYY * D Σ R ≤ 4ϕ max (Σ) np 2 sup R∈B (2) m 2 ,m ′2 tr RYY * R ≤ 4ϕ max (Σ) np 2 1 + sup R∈B (2) m 2 ,m ′2 tr R YY * -I p 2 R . Let F 1 , . . . F d m 2 ,m ′2 an orthonormal basis of U ′ m,m ′ and let λ be a vector in R d m 2 ,m ′2 . We write λ 2 for its L 2 norm. E sup R∈B (2) m 2 ,m ′2 tr R YY * -I p 2 R 2 = E sup λ 2 ≤1 d m 2 ,m ′2 i,j=1 λ i λ j tr F i F j (YY * /n -I p 2 ) 2 ≤ d m 2 ,m ′2 i,j=1 E tr F i F j (YY * /n -I p 2 ) 2 .
The second inequality is a consequence of Cauchy-Schwarz inequality in

R (d m 2 ,m ′2 ) 2 since the l 2 norm of the vector (λ i λ j ) 1≤i,j≤d m 2 ,m ′2 ∈ R d 2 m 2 ,m ′2
is bounded by 1. Since the matrices F i are diagonal, we get

E sup R∈B (2) m 2 ,m ′2 tr [R(YY * /n -I)R] 2 ≤ 2 n d m 2 ,m ′2 i,j=1 F i F j 2 2 .
It remains to bound the norm of the products F i F j for any i, j between 1 and d m 2 ,m ′2 .

d m 2 ,m ′2 i,j=1 F i F j 2 2 = d m 2 ,m ′2 i,j=1 p 2 k=1 F i [k,k] 2 F j [k,k] 2 = p 2 k=1   d m 2 ,m ′2 i=1 F i [k,k] 2   2 .
For any k ∈ {1, . . . , p 2 },

d m 2 ,m ′2 i=1 F i [k,k] 2 ≤ 1 since (F 1 , . . . , F d m 2 ,m ′2
) form an orthonormal family. Hence, we get

d m 2 ,m ′2 i,j=1 F i F j 2 2 ≤ p 2 k=1 d m 2 ,m ′2 i=1 F i [k,k] 2 = d m 2 ,m ′2 .
All in all, we have proved that

E(W m ′ ) ≤ 4ϕ max (Σ) np 2 1 + 2d m 2 ,m ′2 n . (14)
Gathering these three bounds and applying Proposition 8.1 allows to obtain the following deviation inequality:

P Z m ′ ≥ 2ϕmax(Σ) n 1 + α/2 d m 2 ,m ′2 + ξ ≤ exp - √ 1+α/2-1 √ d m 2 ,m ′2 +ξ 2 2L1(1+ √ 2d m 2 ,m ′2 /n) √ n √ 1+α/2-1 √ d m 2 ,m ′2 +ξ √ 2L2 ≤ exp - ω 2 m,m ′ 2L1(1+ √ 2d m 2 ,m ′2 /n) √ nω m,m ′ √ 2L2 - ξω m,m ′ L1[1+ √ 2d m 2 ,m ′2 /n] √ nξ √ 2L2
,

where ω m,m ′ = 1 + α/2 -1 d m 2 ,m ′2 .
As n and d m 2 ,m ′2 are larger than one, there exists a universal constant L ′ 2 such that

( 1 + α/2 -1) 2 d m 2 ,m ′2 2L 1 1 + 2d m 2 ,m ′2 /n √ n( 1 + α/2 -1) d m 2 ,m ′2 √ 2L 2 ≥ 4L ′ 2 d m 2 ,m ′2 1 + α/2 -1 2 ∧ 1 + α/2 -1 .
Since the vector space U m,m ′ contains all the matrices D(θ

′ ) with θ ′ belonging to m ′ , d m 2 ,m ′2 is larger than d m ′ . Besides, by concavity of the square root function, it holds that 1 + α/2 -1 ≥ α[4 1 + α/2] -1 . Setting L ′ 1 := [4L 1 (1 + √ 2)] -1 ∧ [ √ 2L 2 ] -1
and arguing as previously leads to

ξ( 1 + α/2 -1) d m 2 ,m ′2 L 1 1 + 2d m 2 ,m ′2 /n √ nξ √ 2L 2 ≥ L ′ 1 ξ α 1 + α/2 ∧ √ n .
Gathering these two inequalities allows us to conclude that

P Z m ′ ≥ 2ϕmax(Σ) n (1 + α/2) d m 2 ,m ′2 + ξ ≤ exp -L ′ 2 d m ′ α 1 + α/2 ∧ α 2 1 + α/2 -L ′ 1 ξ α 1 + α/2 ∧ √ n .
Proof of Lemma 8.4 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. The approach falls in two parts. First, we relate the dimensions d m and d m 2 to the number of nodes of the torus Λ that are closer than r m or 2r m to the origin (0, 0). We recall that the quantity r m is introduced in Definition 2.1 of [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. Second, we compute a nonasymptotic upper bound of the number of points in Z 2 that lie in the disc of radius r. This second step is quite tedious and will only give the main arguments.

Let m be a model of the collection M 1 . By definition, m is the set of points lying in the disc of radius r m centered on (0, 0). Hence,

Θ m = vect {Ψ i,j , (i, j) ∈ m} ,
where the matrices Ψ i,j are defined by Eq. ( 14) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. As Ψ i,j = Ψ -i,-j , the dimension d m of Θ m is exactly the number of orbits of m under the action of the central symmetry s.

As d m 2 is defined as the dimension of the space U m , it also corresponds to the dimension of the space which is clearly in one to one correspondence with U m . Straightforward computations lead to the following identity:

vect {C(θ), θ ∈ Θ m } + vect C(θ) 2 , θ ∈ Θ m , (15)
C(Ψ i1,j1 )C(Ψ i2,j2 ) = C(Ψ i1+i2,j1+j2 ) [1 + s i1+i2,j1+j2 ] + C(Ψ i1-i2,j1-j2 ) [1 + s i1-i2,j1-j2 ]) ,
where s x,y is the indicator function of x = -x and y = -y in the torus Λ.

Combining this property with the definition of Θ m , we embed the space (15) in the space vect {C (Ψ i1+i2,j1+j2 ) , (i 1 , j 1 ), (i 2 , j 2 ) ∈ m ∪ {(0, 0)}} , and this last space is in one to one correspondence with

vect {Ψ i1+i2,j1+j2 , (i 1 , j 1 ), (i 2 , j 2 ) ∈ m ∪ {(0, 0)}} . (16)
In the sequel, N (m) stands for the set

{(i 1 + i 2 , j 1 + j 2 ), (i 1 , j 1 ), (i 2 , j 2 ) ∈ m ∪ {(0, 0)}} .
Thus, the dimension d m 2 is smaller or equal to the number of orbits of N (m) under the action of the symmetry s.

To conclude, we have to compare the number of orbits in m and the number of orbits in N (m). We distinguish two cases depending whether 2r m + 1 ≤ p or 2r m + 1 > p. First, we assume that 2r m + 1 ≤ p. For such values the disc of radius r m centered on the points (0, 0) in not overlapping itself on the torus except on a set of null Lebesgue measure. In the sequel, ⌊x⌋ refers to the largest integer smaller than x. We represent the orbit space of m as in Figure 1. To any of these points, we associate a square of size 1. If we add 2 + 2⌊r m ⌋ squares to the d m first squares, we remark that the half disc centered on (0, 0) and with length r m is contained in the reunion of these squares. Then, we get

d m + 2 + 2⌊r m ⌋ ≥ πr 2 m 2 . ( 17
)
The points in N (m) are closer than 2r m from the origin. Consequently, all the squares associated to representants of N (m) are included in the disc of radius 2r m + √ 2.

d m 2 + 2 + 2⌊2r m ⌋ ≤ π 2 2r m + √ 2 2 .
Combining these two inequalities, we are able to upper bound

d m 2 2 + 2⌊2r m ⌋ + d m 2 ≤ 4 1 + √ 2 2r m 2 (d m + 1 + 2⌊r m ⌋) , d m 2 ≤ 4 1 + √ 2 2r m 2 d m + 4 1 + √ 2 2r m 2 (1 + 2⌊r m ⌋) .
Applying again inequality (17), we upper bound r m :

r m ≤ 2 π 1 + 1 + π 2 (1 + d m ) .
Gathering these two last bounds yields

d m 2 ≤ 4 1 + √ 2 2r m 2 1 + 1 d m 1 + 4 π 1 + 1 + π 2 (1 + d m ) d m .
This upper bound is equivalent to 4d m , when d m goes to infinity. Computing the ratio d m 2 /d m for every model m of small dimension allows to conclude.

Let us turn to the case 2r m + 1 > p. Suppose that p is larger or equal to 9. The lower bound (17) does not necessarily hold anymore. Indeed, the disc is overlapping with itself because of toroidal effects. Nevertheless, we obtain a similar lower bound by replacing r m by (p -1)/2:

d m + 2 + 2⌊ p -1 2 ⌋ ≥ π(p -1) 2 8 .
The number of orbits of Λ under the action of the symmetry s is (

p 2 + 1)/2 if p is odd and [(p + 1) 2 -1]/2 if p is even. It follows that d m 2 ≤ [(p + 1) 2 -1]/2.
Gathering these two bounds, we get

d m 2 d m ≤ (p + 1) 2 π(p -1) 2 /4 -2(p + 1)
.

This last quantity is smaller than 4 for any p ≥ 9. An exhaustive computation of the ratios when p < 9 allows to conclude. Let us turn to the isotropic case. Arguing as previously, we observe that the dimension d iso m is the number of orbits of the set m under the action of the group G introduced in in [6] Sect.1.1 whereas d m 2 is smaller or equal to the number of orbits of N iso (m) under the action of G. As for anisotropic models, we choose represent these orbits on the torus and associate squares of size 1 (see Figure 2). Assuming that r m < (p -1)/2, we bound d m and d m 2 .

d m + 1 ≥ 1 8 πr 2 m + 1 2 ⌊ √ 2r m 2 ⌋ , d m 2 ≤ 4 1 + √ 2 2r m 2 1 8 πr 2 m + 1 2 ⌊ √ 2r m ⌋ .
Gathering these two inequalities, we get

d m 2 ≤ 4 1 + √ 2 2r m 2 d m .
As a consequence, d m 2 is smaller than 4d m when d m goes to infinity. As previously, computing the ratio d m 2 /d m for models m of small dimension allows to conclude. The case r m > (p -1)/2 is handled as for the anisotropic case.

Proofs of the minimax bounds

Proof of Lemma 8.5 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. This lower bound is based on an application of Fano's approach. See [START_REF] Yu | Assouad, Fano, and Le Cam[END_REF] for a review of this method and comparisons with Le Cam's and Assouad's Lemma. The proof follows three main steps: First, we upper bound the Kullback-Leibler entropy between distributions corresponding to θ 1 and θ 2 in the hypercube. Second, we find a set of points in the hypercube well separated with respect to the Hamming distance. Finally, we conclude by applying Birgé's version of Fano's lemma.

Lemma 3.1. The Kullback-Leibler entropy between two mean zero-Gaussian vectors of size p 2 with precision matrices

I p 2 -C(θ 1 ) /σ 2 and I p 2 -C(θ 2 ) /σ 2 equals K(θ 1 , θ 2 ) = 1/2 log |I p 2 -C(θ 1 )| |I p 2 -C(θ 2 )| + tr I p 2 -C(θ 2 ) I p 2 -C(θ 1 ) -1 -p 2 ,
where for any square matrix A, |A| refers to the determinant of A.

This statement is classical and its proof is omitted. The matrices (I p 2 -C(θ 1 )) and (I p 2 -C(θ 2 )) are diagonalizable in the same basis since they are symmetric block circulant (Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]). Transforming vectors of size p 2 into p × p matrices, we respectively define λ 1 and λ 2 as the p × p matrices of eigenvalues of (I p 2 -C(θ 1 )) and (I p 2 -C(θ 2 )). It follows that

K(θ 1 , θ 2 ) = 1/2 1≤i,j≤p λ 2 [i,j] λ 1 [i,j] -log λ 2 [i,j] λ 1 [i,j] -1 .
For any x > 0, the following inequality holds

x -1 -log(x) ≤ 9 64 x - 1 x 2 .
It is easy to establish by studying the derivative of corresponding functions. As a consequence,

λ 2 [i,j] λ 1 [i,j] -log λ 2 [i,j] λ 1 [i,j] -1 ≤ 9 64 λ 2 [i,j] λ 1 [i,j] - λ 1 [i,j] λ 2 [i,j] 2 ≤ 9 64 1 λ 1 [i,j] + 1 λ 2 [i,j] 2 (λ 1 [i,j] -λ 2 [i,j]) 2 . ( 18
)
Let us first consider the anisotropic case. Let m be a model in M 1 and let θ ′ belong Θ m ∩ B 1 (0 p , 1). We also consider a positive radius r such that (1θ ′ 1 -2rd m ) is positive. For any θ 1 , θ 2 in C m (θ ′ , r) the matrices (I p 2 -C(θ 1 )) and (I p 2 -C(θ 2 )) are diagonally dominant and their eigenvalues λ 1 [i,j] and λ 2 [i,j] are larger than 1

-θ ′ 1 -2rd m . K(θ 1 , θ 2 ) ≤ 9 16(1 -θ ′ 1 -2rd m ) 2 1≤i,j≤p (λ 1 [i,j] -λ 2 [i,j]) 2 ≤ 9 16(1 -θ ′ 1 -2rd m ) 2 C(θ 1 ) -C(θ 2 ) 2 F ≤ 9d m r 2 p 2 8(1 -θ ′ 1 -2rd m ) 2 . ( 19
)
We recall that . F refers to the Frobenius norm in the space of matrices.

Let us state Birgé's version of Fano's lemma [START_REF] Birgé | A new lower bound for multiple hypothesis testing[END_REF] and a combinatorial argument known under the name of Varshamov-Gilbert's lemma. These two lemma are taken from [START_REF] Massart | Concentration inequalities and model selection[END_REF] and respectively correspond to Corollary 2.18 and Lemma 4.7.

Lemma 3.2. (Birgé's lemma) Let (S, d) be some pseudo-metric space and {P s , s ∈ S} be some statistical model. Let κ denote some absolute constant smaller than one. Then for any estimator s and any finite subset T of S, setting δ = min s,t∈T,s =t d(s, t), provided that max s,t∈T K(P s , P t ) ≤ κ log |T |, the following lower bound holds for every p ≥ 1,

sup s∈S E s [d p (s, s)] ≥ 2 -p δ p (1 -κ) .
E θ d H θ, θ ≥ d m 8 (1 -κ) , (20)
provided that

9d m r 2 p 2 n 8(1 -θ ′ 1 -2rd m ) 2 ≤ κd m 8 . (21)
Let us express (20) in terms of the Frobenius . F norm.

sup

θ∈Cm(θ ′ ,r) E θ C( θ) -C(θ) 2 F ≥ d m r 2 p 2 4 (1 -κ) .
Since for every θ in the hypercube, σ -2 (I p 2 -C(θ)) is diagonally dominant, its largest eigenvalue is smaller than 2σ -2 . The loss function l( θ, θ) equals

σ 2 /p 2 tr{[C( θ) -C(θ)](I -C(θ)) -1 [C( θ) -C(θ)]}. It follows that sup θ∈Cm(θ ′ ,r) E θ l θ, θ ≥ σ 2 d m r 2 8 (1 -κ) . (22) Condition (21) is equivalent to r 2 (1 -θ ′ 1 -2rd m ) -2 ≤ κ/(9p 2 n). If we assume that r 2 ≤ κ (1 -θ ′ 1 ) 2 18p 2 n , ( 23
) then 1-θ ′ 1 -2rd m ≥ (1 -θ ′ 1 ) 1 -2d m κ/(18np
E θ l θ, θ ≥ inf θ sup θ∈Cm θ ′ ,r∧(1-θ ′ 1) κ 18p 2 n E θ l θ, θ ≥ L r 2 ∧ (1 -θ ′ 1 ) 2 np 2 d m σ 2 .
One handles models of dimension d m between 1.5( √ 2 -1) np 2 /κ and √ np by changing the constant L in the last lower bound.

Let us turn to sets of isotropic GMRFs. The proof is similar to the nonisotropic case, except for a few arguments. Let m belongs to the collection M 1 and let θ ′ be an element of Θ iso m ∩ B 1 (0 p , 1). Let r be such that 1θ ′ 1 -8d iso m is positive. If θ 1 and θ 2 belong to the hypercube C iso m (θ ′ , r), then

K(θ 1 , θ 2 ) ≤ 9d m r 2 p 2 2(1 -θ ′ 1 -8rd iso m ) 2 .
Applying Lemma 3.2 and 3.3, it follows that

inf θ sup θ∈C iso m (θ ′ ,r) E θ d H θ, θ ≥ d iso m 8 (1 -κ) , provided that 4.5d m r 2 p 2 n(1 -θ ′ 1 -8rd iso m ) -2 ≤ κd iso m /8. As a consequence, inf θ sup θ∈C iso m (θ ′ ,r) E θ l θ, θ ≥ d iso m r 2 8 (1 -κ) , if r 2 1 -θ ′ 1 -8rd iso m -2 ≤ κ(36p 2 n) -1
. We conclude by arguing as in the isotropic case.

Proof of lemma 8.6 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. Let m be a model in M 1 , r be a positive number smaller than 1/(4d m ), and θ be an element of the convex hull of C m (0 p , r).

The covariance matrix of the vector

X v is Σ = σ 2 [I -C(θ)] -1 .
Since the field X is stationary, Var θ (X[0,0]) equals any diagonal element of Σ. In particular, Var θ (X[0,0]) corresponds to the mean of the eigenvalues of Σ. The matrix (I -C(θ)) is block circulant. As in the proof of Lemma 20, we note λ the p × p matrix of the eigenvalues of (I p 2 -C(θ)). By Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF],

λ[i,j] = 1 + (k,l)∈Λ θ[k,l] cos 2π ik p + jl p ,
for any 1 ≤ i, j ≤ p. Since θ belongs to the convex hull of

C m (0 p , r), θ[k,l] is zero if (k, l) / ∈ m and |θ[k,l]| ≤ r if (k, l) ∈ m. Thus (k,l)∈Λ |θ[k,l]
| is smaller than 1/2. Applying Taylor-Lagrange inequality, we get

1 1 + x ≤ 1 -x + x 2 (1 -|x|) 3 ,
for any x between -1 and 1. It follows that

λ[i,j] -1 ≤ 1 - k,l∈Λ θ[k,l] cos 2π ik p + jl p + 8    k,l∈Λ θ[k,l] cos 2π ik p + jl p    2 .
Summing this inequality for all (i, j) ∈ {1, . . . , p} 2 , the first order term turns out to be tr[C(θ)]/p 2 which is zero whereas the second term equals 8tr[C(θ) 2 ]/p 2 . Since there are less than 2d m non-zero terms on each line of the matrix C(θ), its Frobenius norm is smaller than 2d m p 2 r 2 . Consequently, we obtain

Var θ (X[0,0]) ≤ σ 2 1 + 16d m r 2 .
Proof of Lemma 8.7 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. This property seems straightforward but the proof is a bit tedious. Let i be a positive integer smaller than Card(M 1 ). By definition of the radius r m in Equation ( 10) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], the model m i+1 is the set of nodes in Λ\{(0, 0)} at a distance smaller or equal to r mi+1 from (0, 0), whereas the model m i only contains the points in Λ \ {(0, 0)} at a distance strictly smaller than r mi+1 from the origin.

Let us first assume that 2r mi+1 ≤ p. In such a case, the disc centered on (0, 0) with radius r mi+1 does not overlap with itself on the torus Λ. To any node in the neighborhood m i+1 and to the node (0, 0), we associate the square of size 1 centered on it. All these squares do not overlap and are included in the disc of radius r mi+1 + √ 2/2. Hence, we get the upper bound 2d mi+1 + 1 ≤ π(r mi+1 + √ 2/2) 2 . Similarly, the disc of radius r mi+1 -√ 2/2 is included in the union of the squares associated to the nodes m i ∪ {0, 0}. It follows that 2d mi + 1 is larger or equal to π r mi+1 -√ 2/2 2 . Gathering these two inequalities, we obtain

d mi+1 d mi ≤ r mi+1 + √ 2/2 2 -1 r mi+1 - √ 2/2 2 -1 , if r mi+1 is larger than 1 + √ 2/2.
If r mi+1 larger than 5, this upper bound is smaller than two. An exhaustive computation for models of small dimension allows to conclude.

If 2r mi+1 ≥ p and 2r mi < p, then the preceding lower bound of d mi and the preceding upper bound of d mi+1 still hold. Finally, let us assume that 2r mi ≥ p. Arguing as previously, we conclude that 2d mi + 1 ≥ π(p/2 -√ 2/2) 2 . The largest dimension of a model m ∈ M 1 is (p 2 -1)/2 if p is odd and ((p + 1) 2 -3)/2 if p is even. Thus, d mi+1 ≤ [(p + 1) 2 -3]/2. Gathering these two bounds yields

d mi+1 d mi ≤ 4 (p + 1) 2 -3 p - √ 2 2 ,
which is smaller than 2 if p is larger than 10. Exhaustive computations for small p allow to conclude.

Proof of Proposition 6.7 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. This result derives from the upper bound of the risk of θ ρ1 stated in Theorem 3.1 and the minimax lower bound stated in Proposition 6.6 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

Let E(a) be a pseudo-ellipsoid that satisfies Assumption (H a ) and such that a 2 1 ≥ 1/(np 2 ). For any θ in E(a) ∩ B 1 (0 p , 1) ∩ U(ρ 2 ), the penalty term satisfies pen(m) = Kσ 2 ρ 2 1 ρ 2 d m /np 2 is larger than Kd m ϕ max (Σ)/np 2 . Applying Theo-rem3.1, we upper bound the risk θ ρ1

E θ l θ ρ1 , θ ≤ L 1 (K) inf m∈M1 [l(θ m,ρ1 , θ) + pen(m)] + L 2 (K)ρ 2 σ 2 np 2 , for any θ ∈ E(a) ∩ B 1 (0 p , 1) ∩ U(ρ 2 ). It follows that sup θ∈E(a)∩B1(0p,1)∩U (ρ2) E θ l θ ρ1 , θ ≤ L(K) inf m∈M1 , dm>0 l(θ m,ρ1 , θ) + ρ 2 1 ρ 2 σ 2 d m np 2 .
Let i be a positive integer smaller or equal than Card(M 1 ). We know from Section 4.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] that the bias l(θ mi , θ) of the model m i equals Var(X[0,0]|X mi )σ 2 . Since θ belongs to the set E(a) ∩ B 1 (0 p , 1), the bias term is smaller or equal to a 2 i+1 with the convention a 2 Card(M1)+1 = 0. Hence , the previous upper bound becomes

E θ l θ ρ1 , θ ≤ L(K) inf 1≤i≤Card(M1) a 2 i+1 + ρ 2 1 ρ 2 σ 2 d mi np 2 ≤ L(K, ρ 1 , ρ 2 ) inf 1≤i≤Card(M1) a 2 i+1 + σ 2 d mi np 2 . ( 24
)
Applying Proposition 6.6 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] to the set

E(a) ∩ B 1 (0 p , 1) ∩ U(2), we get inf θ sup θ∈E(a)∩B1(0p,1)∩U (ρ2) E θ l θ, θ ≥ inf θ sup θ∈E(a)∩B1(0p,1)∩U (2) E θ l θ, θ ≥ L sup 1≤i≤Card(M1) a 2 i ∧ σ 2 d mi np 2 .
Let us define i * by

i * := sup 1 ≤ i ≤ Card(M 1 ) , a 2 i ≥ σ 2 d mi np 2 ,
with the convention sup ∅ = 0. Since a 2 1 ≥ σ 2 /np 2 , i * is larger or equal to one. It follows that inf

θ sup θ∈E(a)∩B1(0p,η) E θ l θ, θ ≥ L 2 a 2 i * +1 ∨ σ 2 d m i * np 2 .
Meanwhile, the upper bound (24) on the risk of θ ρ1 becomes

E θ l θ ρ1 , θ ≤ L(K, ρ 1 , ρ 2 ) a 2 i * +1 + σ 2 d m i * np 2 ≤ 2L(K, ρ 1 , ρ 2 ) a 2 i * +1 ∨ σ 2 d m i * np 2 ,
which allows to conclude.

Proof of the asymptotic risks bounds

Proof of Corollary 4.6 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. For the sake of simplicity, we assume that for any node (i, j) ∈ m, the nodes (i, j) and (-i, -j) are different in Λ. If this is not the case, we only have to slightly modify the proof in order to take account that Ψ i,j

2

F may equal one. The matrix V is the covariance of the vector of size

d m X i1,j1 + X -i1,-j1 , . . . , X i dm ,j dm + X -i dm ,-j dm . (25)
Since the matrix Σ of X v is positive, V is also positive. Moreover, its largest eigenvalue is larger than 2ϕ max (Σ).

Let us assume first the θ belongs to Θ + m and that Assumption (H 1 ) is fulfilled. By the first result of Proposition 4.4 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF],

lim n→+∞ np 2 E l θ m,ρ1 , θ = 2σ 4 tr IL m V -1 ≥ σ 4 ϕ max (Σ) tr[IL m ] = 2σ 4 d m ϕ max (Σ) ,
which corresponds to the first lower bound (30) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

Let us turn to the second result. We now assume that θ satisfies Assumption (H 2 ). By the identity (28) of Proposition 4.4 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], we only have to lower bound the quantity tr V W -1 .

tr V -1 W ≥ ϕ max (V ) -1 tr [W ] ≥ 1 2ϕ max (Σ) tr[W ] . Since the matrix Σ -1 = σ -2 I p 2 -C(θ) is diagonally dominant, its smallest eigenvalue is larger than σ -2 (1-θ 1 ). The matrix I p 2 -C(θ m,ρ1 ) 2 I p 2 -C(θ)
is symmetric positive. It follows that W is also symmetric positive definite. Hence, we get

tr V -1 W (26) ≥ σ -2 2 [1 -θ 1 ] dm k=1 tr C(Ψ i k ,j k ) 2 I p 2 -C(θ m,ρ1 ) 2 I p 2 -C(θ) -2 p 2 .
The largest eigenvalue of I p 2 -C(θ) is smaller than 2 and the smallest eigenvalue of

I p 2 -C(θ m,ρ1
) is larger than 1-θ m,ρ1 1 . By Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], these two matrices are jointly diagonalizable and the smallest eigenvalue of

I p 2 -C(θ m,ρ1 ) 2 I p 2 -C(θ) -2
is therefore larger than (1θ m,ρ1 1 ) 2 /4. Gathering this lower bound with (26) yields

tr V -1 W ≥ d m σ -2 2 [1 -θ 1 ] [1 -θ m,ρ1 1 ] 2 .
Lemma 4.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] states that θ m,ρ1 1 ≤ θ 1 . Combining these two lower bounds enables to conclude.

Proof of Example 4.8 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

Lemma 4.1. For any θ is the space Θ +,iso m1 , the asymptotic variance term of θ iso m1,ρ1 equals lim n→+∞ np 2 E θ l θ iso m1,ρ1 , θ = 2σ 4 tr H 2 tr (H 2 Σ) .

If θ belongs to Θ +,iso and also satisfies (H 2 ), then

lim n→+∞ np 2 E θ l θ iso m1,ρ1 , θ iso m1,ρ1 = 2 tr (I -θ iso m1,ρ1 [1, 0]H)HΣ 2 tr(H 2 Σ) , (27)
where the p 2 × p 2 matrix H is defined as H := C Ψ iso 1,0 . Proof of Lemma 4.1. Apply Proposition 4.4 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] noting that V = tr[HΣH]/p 2 and

W = tr (I -θ m iso 1 [1,0]H )HΣ 2 σ 4 p 2 .
To prove the second result, we observe that Θ +,iso m1 equals Θ +,iso m1,2 . It is stated for instance in Table 2 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

Since the matrix θ belongs to Θ +,iso m1 , we may apply the second result of Lemma 4.1. Straightforward computations lead to tr(H 2 ) = C Ψ iso 1,0 2

F = 4p 2 and tr(H 2 Σ) = 4p 2 [Var(X[0,0]) + 2cov θ (X[0,0], X[1,1]) + cov θ (X[0,0], X[2,0])] .
Since the field X is an isotropic GMRF with four nearest neighbors,

X[0,0] = θ[1,0] (X[1,0] + X[-1,0] + X[0,1] + X[0,-1]) + ǫ[0,0] ,
where ǫ[0,0] is independent from every variable X[i,j] with (i, j) = 0. Multiplying this identity by X[1,0] and taking the expectation yields

cov θ (X[0,0], X[1,0]) = θ[1,0] [Var (X[0,0]) + 2cov θ (X[0,0], X[1,1]) + cov θ (X[0,0], X[2,0])] .
Hence, we obtain tr H 2 Σ = 4cov θ (X[0,0], X[1,0])/θ [START_REF] Billingsley | Probability and measure[END_REF]0] and

tr(H 2 ) tr(H 2 Σ) = θ[1,0] cov θ (X[0,0], X[1,0]) ,
which concludes the first part of the proof.

This second part is based on the spectral representation of the field X and follows arguments which come back to Moran [START_REF] Moran | A Gaussian Markovian process on a square lattice[END_REF]. We shall compute the limit of cov θ (X[0,0], X[1,0]) when the size of Λ goes to infinity. As the field X is stationary on Λ, we may diagonalize its covariance matrix Σ applying Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. We note D Σ the corresponding diagonal matrix defined by Applying Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] to the matrix Σ -1 and noting that θ ∈ Θ iso,+ allows to get another expression of the eigenvalues of Σ

D Σ [(i-1)p+j,(i-1)p+j] = σ 2 1 -2θ[1,0] cos 2πi p + cos 2πj p .
We then combine these expression. By symmetry between i and j we get dxdy .

cov θ (X[0,0], X[1,0]) = σ 2
This last elliptic integral is asymptotically equivalent to log 16[4(1 -4θ[1,0])] -1 when θ[1,0] → 1/4 as observed for instance by Moran [START_REF] Moran | A Gaussian Markovian process on a square lattice[END_REF]. We conclude by substituting this limit in expression (33) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF].

Proof of Example 4.9 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. First, we compute [θ (p) ] iso m1 [START_REF] Billingsley | Probability and measure[END_REF]0]. By Lemma 4.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], it minimizes the function γ(.) defined in (19) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] over the whole space Θ m iso 1 . We therefore obtain

[θ (p) ] iso m1 [1,0] = tr [ΣH] tr [ΣH 2 ] .
Once again, we apply Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF] to simultaneously diagonalize the matrices H and Σ -1 . As previously, we note D Σ the corresponding diagonal matrix of Σ.

D Σ [(i-1)p+j,(i-1)p+j] = σ 2 1 -2α cos 2π pi 4p + pj 4p + cos 2π -pi 4p + pj 4p = σ 2 1 -4α cos π i 2 cos π j 2 .
Analogously, we compute the diagonal matrix D Ψ iso 1,0 .

As each term of this sum is non-negative, we may only consider the coefficients i and j which are congruent to 0 modulo 4.

tr(H

np 2 R θ (p) θ (p) iso,ρ1 m1 , [θ (p) ] iso m1 ≥ Lσ 2 1 -4α .

Miscellaneous

Proof of Lemma 1.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. Let θ be a p × p matrix that satisfies condition (3) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. For any 1 ≤ i 1 , i 2 ≤ p, we define the p × p submatrix C i1,i2 as

C i1,i2 [j1,j2] := C(θ)[(i1-1)p+j1,(i2-1)p+j2] ,
for any 1 ≤ j 1 , j 2 ≤ p. For the sake of simplicity, the subscripts (i 1 , i 2 ) are taken modulo p. By definition of C(θ), it holds that C i1,i2 = C 0,i2-i1 for any 1 ≤ i 1 , i 2 ≤ p. Besides, the matrices C 0,i are circulant for any 1 ≤ i ≤ p. In short, the matrix C(θ) is of the form

C(θ) =    C 0,1 C 0,2 • • • C 0,p . . . . . . . . . . . . C 0,p C 0,1 • • • C 0,p-1    ,
where the matrices C 0,i are circulant. Let (i 1 , i 2 , j 1 , j 2 ) be in {1, . . . , p} 4 . By definition,

C(θ)[(i1-1)p+j1,(i2-1)p+j2] = θ[i2-i1,j2-j1] .
Since the matrix θ satisfies condition (3) in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], θ[i2-i1,j2-j1] = θ[i1-i2,j1-j2]. As a consequence, C(θ)[(i1-1)p+j1,(i2-1)p+j2] = C(θ)[(i2-1)p+j2,(i1-1)p+j1] and C(θ) is symmetric.

Conversely, let B be a p 2 × p 2 symmetric block circulant matrix. Let us define the matrix θ of size p by

θ[i,j] := B[1,(i-1)p+j] ,
for any 1 ≤ i, j ≤ p. Since the matrix B is block circulant, it follows that

C(θ) = B. By definition, θ[i,j] = C(θ)[1,(i-1)p+j] and θ[-i,-j] = C(θ)[(i-1)p+j,1] for any integers 1 ≤ i, j ≤ p. Since the matrix B is symmetric, we conclude that θ[i,j] = θ[-i,-j].
Proof of Lemma 2.2 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. For any θ ′ ∈ Θ + , γ n,p (θ ′ ) is defined as

γ n,p (θ ′ ) = 1 p 2 tr (I p 2 -C(θ ′ ))X v X v * (I p 2 -C(θ ′ )) .
Applying Lemma A.1 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF], there exists an orthogonal matrix P that simultaneously diagonalizes Σ and any matrix C(θ ′ ). Let us define Y i := √ Σ -1 X i and D Σ := P ΣP * . Gathering these new notations yields

γ n,p (θ ′ ) = 1 p 2 tr I p 2 -D(θ ′ ) D Σ YY * I p 2 -D(θ ′ ) ,
where the vectors Y i are independent standard Gaussian random vectors. Except YY * , every matrix involved in this last expression is diagonal. Besides, the diagonal matrix D Σ is positive since Σ is non-singular. Thus,

tr (I p 2 -D(θ ′ ))D Σ YY * (I p 2 -D(θ ′ ))
is almost surely a positive quadratic form on the vector space generated by I p 2 and D(Θ + ). Since the function D(.) is injective and linear on Θ + , it follows that γ n,p (.) is almost surely strictly convex on Θ + .

Proof of Lemma 4.1 and Corollary 4.2 in [START_REF] Verzelen | Adaptive estimation of stationary Gaussian Markov random fields[END_REF]. The proof only uses the stationarity of the field X on Λ and the l 1 norm of θ. However, the computations are a bit cumbersome. Let θ be an element of Θ + . By standard Gaussian properties, the expectation of X[0,0] given the remaining covariates is

E θ X[0,0]|X -{0,0} = (i,j)∈Λ\(0,0) θ[i,j]X[i,j] .
By assumption (H 2 ), the l 1 norm of θ is smaller than one. We shall prove by backward induction that for any subset A of Λ\{(0, 0)} the matrix θ A uniquely defined by E θ (X[0,0]|X A ) = (i,j)∈A θ A [i,j]X [i,j] and θ A [i,j] = 0 for any (i, j) / ∈ A satisfies θ A 1 ≤ θ 1 . The property is clearly true if A = Λ\{(0, 0)}. Suppose we have proved it for any set of cardinality q larger than one. Let A be a subset of Λ\{(0, 0)} of cardinality q -1 and (i, j) be an element of Λ\(A ∪ {(0, 0)}). Let us derive the expectation of X[0,0] conditionally to X A from the expectation of X[0,0] conditionally to X A∪{(i,j)} .

E θ (X[0,0]|X A ) = E θ E(X[0,0]|X A )|X A∪{(i,j)} = (k,l)∈A θ A∪{(i,j)} [k,l]X [k,l] + θ A∪{(i,j)} [i,j]E θ [X[i,j]|X A ] . ( 28
)
Let us take the conditional expectation of X[i,j] with respect to X A∪{(0,0)} . Since the field X is stationary on Λ and by the induction hypothesis, the unique matrix θ A∪{(0,0)} (i,j)

defined by since θ A∪{i,j} [i,j]θ i,j A∪{0,0} [0,0] < 1. Then, we upper bound the l 1 norm of θ A using that θ A∪{(i,j)} 1 and θ A∪{(0,0)} (i,j)

1 are smaller or equal to θ 1 . 1θ A∪{(i,j)} [i,j]θ A∪{(0,0)} (i,j) [0,0] ≤ θ 1 + θ A∪{(i,j)} [i,j] (k,l)∈A∪{(0,0)} θ A∪{(0,0)} (i,j)

[k,l] -1θ A∪{(0,0)} (i,j) [0,0] 1θ A∪{(i,j)} [i,j]θ A∪{(0,0)} i,j [0,0] ≤ θ 1 1 + θ A∪{(i,j)} [i,j]θ A∪{(i,j)} [i,j] 1 + θ A∪{(0,0)} i,j [0,0] 1θ A∪{(i,j)} [i,j]θ A∪{(0,0)} (i,j) [0,0] ≤ θ + θ A∪{(i,j)} [i,j] ( θ 1 -1) 1 + θ A∪{(0,0)} (i,j) [0,0] 1θ A∪{(i,j)} [i,j]θ A∪{(0,0)} (i,j) [0,0]

.

Since θ 1 is smaller than one, it follows that θ A 1 ≤ θ 1 .

Let m be a model in the collection M 1 . Since m stands for a set of neighbors of (0, 0), we may define θ m as above. It follows that θ m 1 ≤ θ 1 . Since the field X is stationary on the torus, X follows the same distribution as the field X s defined by X s [i,j] = X[-i,-j]. By uniqueness of θ m , we obtain that θ m [i,j] = θ m [-i,-j]. Thus, θ m belongs to the space Θ m . Moreover, θ m minimizes the function γ(.) on Θ m . Since the l 1 norm of θ m is smaller than one, θ m belongs to Θ + m,2 . The matrices θ m and θ m,ρ1 are therefore equal, which concludes the proof in the non-isotropic case.

Let us now turn to the isotropic case. Let θ belong to Θ iso,+ and let m be a model in M 1 . As previously, the matrix θ m satisfies θ m 1 ≤ θ 1 . Since the distribution of X is invariant under the action of the group G, θ m belongs to Θ iso m . Since θ m 1 ≤ θ 1 , θ m lies in Θ +,iso m,2 . It follows that θ m = θ iso m,ρ1 .

Proof of Corollary 4.3 in [6]. Let θ be a matrix in Θ + such that (H 2 ) holds and let m be a model in M 1 . We decompose γ( θ m,ρ1 ) using the conditional

Lemma 1 . 1 .

 11 Let T be a supremum of Rademacher chaos indexed by I N of the form T := sup t∈T {i,j}

( 1 +

 1 k)} := R[i,i] n , and t R ∅ := -tr(R) , where δ i,j is the indicator function of i = j. In order to apply Lemma 1.2 with N = nr and T = t R |R ∈ F , we have to work out the quantities D and E. δ i,j δ k,l )α l j Y * R * ) .

Figure 1 .

 1 Figure 1. The black dots represent the orbit space of m and the white dots represent the remaining points of the orbit space of N (m).

Figure 2 .

 2 Figure 2. The black dots represent the orbit space of m under the action of G and the white dots represent the remaining points of the orbit space of N iso (m).

Lemma 3 . 3 .

 33 (Varshamov-Gilbert's lemma) Let {0, 1} d be equipped with Hamming distance d H . There exists some subset Φ of {0, 1} d with the following properties d H (φ, φ ′ ) > d/4 for every (φ, φ ′ ) ∈ Φ 2 with φ = φ ′ and log |Φ| ≥ d 8 . Applying Lemma 3.2 with Hamming distance d H and the set Φ introduced in Lemma 3.3 yields sup θ∈Cm(θ ′ ,r)

D

  Σ [(i-1)p+j,(i-1)p+j] = p k=1 p l=1 cov θ (X[0,0], X[k,l]) cos 2π ki p + lj p , for any 1 ≤ i, j ≤ p. Straightforwardly, we express cov θ (X[0,0], X[1,0]) as a linear combination of the eigenvalues cov θ (X[0,0], X[1,0]) (i-1)p+j,(i-1)p+j] .

. 1 -

 1 If we let p go to infinity, this sum converges to the following integral limp→+∞ cov θ (X[0,0], X[1,0]) 2θ[1,0] [cos(x) + cos(y)]

E

  θ X[i,j]|X A∪{(0,0)} = l] = 0 for any (k, l) / ∈ A ∪ {(0, 0)} satisfies θ A∪{(0,0)} (i,j) 1 ≤ θ 1 .Taking the expectation conditionally to X A of this previous expression leads toE θ (X[i,j]|X A ) = (k,l)∈A θ A∪{(0,0)} (i,j) [k,l]X [k,l] + θ A∪{(0,0)} (i,j) [0,0]E (X[0,0]|X A ) . (29)Gathering identities (28) and (29) yieldsE θ (X[0,0]|X A ) = (k,l)∈A θ A∪{i,j} [k,l] + θ A∪{(i,j)} [i,j]θ A∪{0,0} (i,j) [k,l] 1θ A∪{(i,j)} [i,j]θ A∪{0,0} (i,j) [0,0] X[k,l] ,

θ A 1 ≤

 1 (k,l)∈A θ A∪{j+1} [k,l] + (k,l)∈A θ A∪{(i,j)} [i,j]θ

.

  Let us split this sum in 16 parts depending on the congruence of i and j modulo 4. As each if of these 16 sums is shown to be zero, we conclude that tr(HΣ) =

	First, we lower bound the numerator
	tr(H 4 Σ 2 ) = σ 4	p i=1	p j=1	2 cos 2π i p + cos 2π j p 1 -4α cos π i 2 2 cos π j	2	4
	D Ψ iso 1,0 [(i-1)p+j,(i-1)p+j] = 2 cos 2π i p + cos 2π j p	.
	Combining these two last expressions, we obtain
	tr(HΣ) =	p i=1	p j=1	σ 2	2 cos 2π i p + cos 2π j p 1 -4α cos π i 2 cos π j 2
	[θ (p) ] iso m1 [1,0] = 0. By Lemma 4.1, the asymptotic risk of θ (p)	iso,ρ1 m1	therefore equals
	lim n→+∞	np 2 E θ (p) l θ (p)	iso,ρ1 m1 , [θ (p) ] iso m1	=	tr(H 4 Σ 2 ) tr(H 2 Σ)	.

  4 Σ 2 ) ≥ σ 4

					p/4-1 i=0	p/4-1 j=0	16 cos 2π i p/4 + cos 2π j p/4 (1 -4α) 2	4	.
	If we let go p to infinity, we get the lower bound
	lim p→+∞	tr(H 4 Σ 2 ) p 2	≥	σ 4 (1 -4α) 2		0	1	0	1	[cos(2πx) + cos(2πy)] 4 dxdy .
	Similarly, we upper bound tr(H 2 Σ) and let p go to infinity
	lim p→+∞	tr(H 2 Σ) p 2	≤	4σ 2 1 -4α	0	1	0	1	[cos(2πx) + cos(2πy)] 2 dxdy .
	Combining these two bounds allows to conclude
		lim p→+∞	lim n→+∞