Cédric Boutillier

Béatrice

Béatrice De Tilière

Statistical mechanics on isoradial graphs

teaching and research institutions in France or abroad, or from public or private research centers.

Introduction

Statistical mechanics aims at describing large scale properties of physics systems based on models which specify interactions on a microscopic level. In this setting, physics systems are modelled by random configurations of graphs embedded in a d-dimensional space. Since vertices typically represent atoms, the goal is to let the mesh size tend to zero and rigorously understand the limiting behavior of the system.

Although real world suggests that we should focus on the case of 3-dimensional systems, we restrict ourselves to the case of dimension 2. There are two main reasons guiding this choice: first, for models we consider, only very few rigorous results Cédric Boutillier Laboratoire de Probabilités et Modèles Aléatoires, Université Pierre et Marie Curie, 4 place Jussieu, F-75005 Paris. e-mail: cedric.boutillier@upmc.fr Béatrice de Tilière Laboratoire de Probabilités et Modèles Aléatoires, Université Pierre et Marie Curie, 4 place Jussieu, F-75005 Paris. e-mail: beatrice.de tiliere@upmc.fr exist in dimension 3, and more importantly, 2-dimensional systems exhibit beautiful and rich behaviors which are strongly related to this choice of dimension.

Historically, the most studied graph is certainly Z 2 , followed by the honeycomb and triangular lattices. Since solving statistical mechanics models involves dealing with combinatorial and geometric features of the underlying graph, it has been most convenient to handle the simplest and most regular ones. However in this paper, we deal with models defined on a more general class of graphs, called isoradial graphs.

The motivation behind this generalization is to find the most appropriate setting, which exhibits some essential features required to solve the questions addressed, thus allowing for proofs revealing the true nature of the problems, and hopefully a full understanding of the issues at stake. Before going any further, let us define isoradial graphs and the corresponding rhombus graph. Definition 1. A graph G = (V, E) is said to be isoradial [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF], if it has an embedding in the plane such that every face is inscribed in a circle of radius 1, and such that all circumcenters of the faces are in the closure of the faces.

From now on, when we speak of the graph G, we mean the graph together with a particular isoradial embedding in the plane. Examples of isoradial graphs are given in Figure 1 below, in particular Z 2 , the honeycomb and triangular lattices are isoradial.

Fig. 1 Examples of isoradial graphs: the square lattice (left), the honeycomb lattice (center), and a more generic one (right). Every face is inscribed in a circle of radius 1, represented in dashed lines.

To such a graph is naturally associated the diamond graph, denoted by G ⋄ , defined as follows. Vertices of G ⋄ consist in the vertices of G, and the circumcenters of the faces of G. The circumcenter of each face is then joined to all vertices which are on the boundary of this face, see Figure 2. Since G is isoradial, all faces of G ⋄ are side-length-1 rhombi. Moreover, each edge e of G is the diagonal of exactly one rhombus R e of G ⋄ ; we let θ e be the half-angle of the rhombus at the vertex it has in common with e. For later purposes, we label vertices of the rhombus R e by v 1 , v 2 , v 3 , v 4 , as in Figure 2 (right).

Isoradial graphs naturally arise in two different approaches to statistical mechanics. We describe them in the next two sections.

θe e G R e G v 1 v 2 3 v v 4
Fig. 2 An isoradial graph (left). The white dots are the circumcenters of the faces, which are also the vertices of the dual G * . Its diamond graph is represented in the center. On the right is the rhombus R e and the half-rhombus angle θ e assigned to an edge e.

Transfer matrices, star transformations and Z-invariance

The first question addressed when solving a model of statistical mechanics is to compute the free energy, which is the growth rate of the partition function, counting the weighted number of configurations. An important technique introduced by physicists to solve two-dimensional models is the use of transfer matrices, which appears for the first time in the work of Kramers and Wannier [KW41a,[START_REF] Kramers | Statistics of the two-dimensional ferromagnet. part ii[END_REF]. If a model of statistical mechanics with local interactions is defined on a torus of size n × m, then the partition function Z m,n can be expressed as the trace of the n-fold product of a transfer matrix T , whose rows and columns are indexed by the configurations of the model in a strip m × 1, and if C and C ′ are two such configurations, then the matrix element T C ,C ′ is the contribution to the Boltzmann weight of the interactions in C and C ′ are in consecutive strips. The free energy of the model f can be expressed in terms of the largest eigenvalue of T denoted by λ 1 :

f = -lim m,n→∞ 1 mn log(Z m,n ) = -lim m→∞ 1 m log(λ 1 ).
If the interaction constants on the graph are not homogeneous, then the transfer matrices used to pass from one strip to another are different. If we want to be able to compute the free energy, as in the homogeneous case, in terms of spectral characteristics of the transfer matrices, then these should commute. Since the interactions are local, the transfer matrices are themselves sums of local operators R acting on configurations at neighboring sites, and diagrammatically represented by one rhombus of a strip. A sufficient condition to ensure the commutation of transfer matrices is to demand that these R-matrices satisfy some algebraic relation, called the Yang-Baxter equations, which can be loosely formulated as follows. Consider the two ways to tile a hexagon with the same three rhombi on Figure 4. The transformation from one to the other is called a star-star transformation. The model satisfies the Yang-Baxter equations if the sum of Boltzmann weights over local configurations of the left hand-side equals those of the right hand side. This very strong constraint gives a set of equations which needs to be satisfied by Boltzmann weights, see [START_REF] Perk | Yang-baxter equations[END_REF][START_REF] Baxter | Exactly solved models in statistical mechanics[END_REF] for more details. If a non trivial solution is found, we say that the model is Z-invariant.

The condition of Z-invariance yields a natural generalization. Suppose that the graph embedded on the torus of size m × n is a periodic tiling by rhombi (where the period is independent of m and n). Then by performing a sequence of star-star moves, this graph can be transformed into big pieces of tilted copies of Z 2 [START_REF] Kenyon | An introduction to the dimer model[END_REF]. If boundary effects can be neglected, then it suffices to solve the model on each of the copy of pieces of Z 2 using transfer matrices. Summarizing, if the model is Zinvariant, then it can be solved on any periodic graph consisting of rhombi.

Consider an isoradial graph G and recall that a unique rhombus of the underlying diamond graph G ⋄ can be assigned to every edge of G. Then, looking at edges rather than rhombi, the star-star transformation becomes a star-triangle transformation, see Figure 5 below. The definition of Z-invariance naturally extends to this setting, and if the model is Z-invariant the transfer matrix approach can be performed on G ⋄ , thus explaining the occurrence of isoradial graphs and Z-invariance in the context of transfer matrices.

Conformal field theory and discrete complex analysis

Conformal field theory (CFT), introduced by Belavin, Polyakov and Zamolodchikov [START_REF] Belavin | Infinite conformal symmetry in two-dimensional quantum field theory[END_REF] (see also [START_REF] Di | Conformal field theory[END_REF]) is a theory which aims at describing models at criticality, supposed not to depend on specificities of the graph, and to be conformally invariant. This very strong statement remained largely inaccessible to the mathematics community (except for a few models, such as dimers on the square lattice [START_REF] Kenyon | Conformal invariance of domino tiling[END_REF]) until the introduction of the SLE process by Schramm in 1999. The SLE is conformally invariant and conjectured to be the limiting process for many models of statistical mechanics, thus filling a huge gap with CFT. For relations between CFT and SLE, see for example [START_REF] Friedrich | Conformal restriction, highest-weight representations and SLE[END_REF][START_REF] Bauer | CFTs of SLEs: the radial case[END_REF]. Several of these conjectures are now solved: loop erased random walk and uniform spanning tree (Lawler, Schramm, Werner [START_REF] Lawler | Conformal invariance of planar loop-erased random walks and uniform spanning trees[END_REF]), site percolation on the triangular lattice (Smirnov [Smi01]), Ising (Smirnov [Smi10], Chelkak and Smirnov [CSb]). . .

It is thus of key importance to have a setting suitable for doing discrete complex analysis, and proving convergence to its continuous counterpart. To this purpose, isoradial graphs are a perfectly suited object. Indeed, when considering an edge and its dual, they consist in the two diagonals of a rhombus, and are thus orthogonal. This allows for a natural discretization of the Cauchy-Riemman equation. Since a lot of developments have happened in this direction, we dedicate it the next section.

Discrete complex analysis on isoradial graphs

As mentioned above, isoradial graphs provide a natural framework for a generalization of the construction in the case of the square lattice. The ideas presented here go back to Duffin [START_REF] Duffin | Potential theory on a rhombic lattice[END_REF], and have been developed later by Mercat [Mer01], Kenyon [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF], Chelkak and Smirnov [CSa], Cimasoni [Cim] and others. There are other possible discretizations of complex analysis on isoradial graphs with remarkable properties. An example using the notion of cross-ratios is given in [START_REF] Bobenko | Linear and nonlinear theories of discrete analytic functions. Integrable structure and isomonodromic Green's function[END_REF].

Discrete holomorphic and discrete harmonic functions

Let us start with the notion of discrete holomorphy. In the continuous setting, the Cauchy-Riemann equations satisfied by a holomorphic function imply that its partial derivatives along two orthogonal unit vectors differ by a factor i. In the setting of isoradial graphs, one takes advantage of the fact that faces of G ⋄ are rhombi, and thus have perpendicular diagonals, to write down finite difference equations which are discretizations of Cauchy-Riemann equations.

Let f be a C-valued function on vertices of G ⋄ , let e be an edge of G, and R e be the corresponding rhombus with the labeling of its vertices v 1 , v 2 , v 3 , v 4 , given in Figure 2. We say that the function f is discrete holomorphic at the rhombus R e , if:

f (v 3 ) -f (v 1 ) v 3 -v 1 = f (v 4 ) -f (v 2 ) v 4 -v 2 .
More generally we define a discrete operator, ∂ : G) , by:

C V (G ⋄ ) → C E(
∂ f (e) = i(v 4 -v 2 )(v 3 -v 1 ) f (v 4 ) -f (v 2 ) v 4 -v 2 - f (v 3 ) -f (v 1 ) v 3 -v 1 = 4 ∑ j=1 ∂e,v j f (v j ), where ∂e,v j = i(v j-1 -v j+1 ), (1)
with indices written mod 4. A function f is then said to be discrete holomorphic on G ⋄ , if ∂ f ≡ 0. Similarly, we define an operator ∂ by replacing the coefficients ∂e,v j by their complex conjugate.

The operators ∂ and ∂ can be extended to 1-forms (from

C E(G) to C V (G ⋄ ) ): let g be a C-valued function on edges of G, then for all v ∈ V (G ⋄ ), ∂ g(v) = ∑ e=(v,v ′ )∼v 1 (v ′ -v) g(e).
The extension of the operator ∂ is again obtained by replacing coefficients by their complex conjugate. Note that the normalization of the operators ∂ and ∂ we adopt here differs from that chosen by Mercat, or by Chelkak and Smirnov. One can recover their normalization by multiplying our operators by diagonal operators acting on edges. So these variations do correspond to the same notion of discrete holomorphy.

An notion intimately connected to holomorphy is that of harmonicity. The Laplace operator ∆ on G can be written as the restriction of ∂ ∂ to functions supported on G1 . Its action on a function f defined on vertices of G is given by

∀v ∈ G, ∆ f (v) = ∂ ( ∂ f )(v) = ∑ u∼v tan θ uv ( f (u) -f (v)).
A function f is said to be discrete harmonic if ∆ f ≡ 0. This corresponds to choosing for every edge e, a conductance equal to tan θ e .

From the factorisation ∆ = ∂ ∂ , one sees that, like in the continuous case, the restriction of a discrete holomorphic function to G (resp. G * ) is harmonic. Conversely, if H is a harmonic function on G, there exists a harmonic function

H * on G * such that H + iH * is discrete holomorphic.
From this choice of Laplacian, one can define the discrete analogues of classical quantities in continuous potential theory (Green function, Poisson kernel, harmonic measure,. . . ). A surprising amount of estimates have a discrete version, and can be found in [CSa].

Discrete complex analysis can serve as an approximation of the continuous theory. Mercat [START_REF] Mercat | Discrete Riemann surfaces and the Ising model[END_REF] proved that the pointwise limit of a sequence of discrete holo-morphic functions on isoradial graphs in a domain with a mesh going to zero is holomorphic. See [Cim] for a generalisation of this result to isoradial discretization of compact manifolds. These local convergence results have been supplemented by global convergence theorems [CSb], implying in particular the uniform convergence for the discrete potential theory objects to their continuous counterparts.

Discrete exponential functions

There is a class of functions playing a special rôle in this theory. These functions, called discrete exponential functions, are defined recursively on G ⋄ . For any given vertex v 0 of G ⋄ and any λ ∈ C, the function Exp v 0 (•; λ ) is defined as follows: its value at v 0 is 1, and if v 1 and v 2 are neighbors in G ⋄ , let e iα = v 2v 1 , then:

Exp v 0 (v 2 ; λ ) = Exp v 0 (v 1 ; λ ) 1 + λ e iα 1 -λ e iα .
The name comes from the fact that these functions satisfy the following identity: for any pair of neighbouring vertices v 1 and

v 2 of G ⋄ , Exp v 0 (v 2 ; λ ) -Exp v 0 (v 1 ; λ ) = λ Exp v 0 (v 2 ; λ ) + Exp v 0 (v 1 ; λ ) 2 (v 2 -v 1 ), which is a discrete version of the differential equation d exp(λ x) = λ exp(λ x)dx.
satisfied by the usual exponential function. It is straightforward to check from the definition that Exp v 0 (•; λ ) is discrete holomorphic. In fact, Bobenko, Mercat and Suris [START_REF] Bobenko | Linear and nonlinear theories of discrete analytic functions. Integrable structure and isomonodromic Green's function[END_REF] show that any discrete holomorphic function can be written as a (generalized) linear combination of discrete exponential functions, at least under the quasicrystallic assumption (a finite number of slopes ±e iθ 1 , . . . , ±e iθ d for the rhombi chains): if f is a discrete holomorphic function on G ⋄ , then there exists a function g(λ ) defined on a neighborhood in C of {±e iθ 1 , . . . ,

±e iθ d } such that ∀v ∈ G ⋄ , f (v) = Γ Exp v 0 (v; λ )g(λ )dλ ,
where Γ is a collection of disjoint small positively oriented contours around the possible pôles ±e iθ j of the integrand.

Geometric integrability of discrete Cauchy-Riemann equations

An important feature of the discrete Cauchy-Riemann is the property of 3D consistency [START_REF] Bobenko | Linear and nonlinear theories of discrete analytic functions. Integrable structure and isomonodromic Green's function[END_REF]. Let f be a discrete holomorphic function on G ⋄ . Given its value on three vertices of a rhombus, the other can be computed using the discrete Cauchy-Riemann equation. If a star-star transformation is realised on the graph G ⋄ (see Figure 4), can we determine the value of the function f at the new vertex in the center, so that f is again discrete holomorphic? A priori, there are three different ways to compute a possible value: one for each rhombus. It turns out that the three values obtained are equal.

Under the quasicrystallic assumption, G ⋄ can be seen as a monotone surface in Z d , projected back to a properly chose plane. The star-star move is a local displacement of this surface along the faces of a cube. If given two monotonic surfaces Σ 1 and Σ 2 in Z d that can be deduced one from the other by a sequence of star-star moves, then there is a canonical way to extend a discrete holomorphic function f defined on Σ 1 to Σ 2 by "pushing" its values along the deformation making use of the Cauchy-Riemann equations. In this sense, this 3D consistency can be considered as a notion of integrability [START_REF] Bobenko | Discrete differential geometry[END_REF]. This property is very much related to the Z-invariance in statistical mechanics.

Generalization of the operator ∂

The double graph G D of G is a planar bipartite graph constructed as follows, see also Figure 6 The operator ∂ : G) of Section 2.1 can be interpreted now as an operator ∂ :

C V (G ⋄ ) → C E(
C B(G D ) → C W (G D )
, where B(G D ) (resp. W (G D )) denotes the set of black (resp. white) vertices of G D . Let f be a C-valued function defined on black vertices of G D , then for every white vertex of

G D : ∂ f (w) = ∑ b∼w ∂w,b f (b).
Let w, x, b, y be the vertices of a rhombus enumerated in cclw order, so that w ∈ W (G D ) and b ∈ B(G D ). Then by Formula (1), (2)

G D
This definition can be extended to any bipartite isoradial graph, see Section 3.1.

For a detailed study of this operator on isoradial graphs on compact surfaces and its connections with discretization of geometric structures (especially spin structures), see [Cim].

Dimer model

The dimer model represents the adsorption of diatomic molecules on the surface of a crystal. The surface of the crystal is modeled by a planar graph G = (V, E), which we assume to be finite for the moment. A dimer configuration of G is a perfect matching of G, that is a subset of edges M such that every vertex of G is incident to a unique edge of M. Denote by M (G) the set of dimer configurations of G. An example of dimer configuration is given in Figure 8. Assume that a positive weight function ν is assigned to edges of G, that is every edge e of G has weight ν e . Then, when the graph G is finite, the probability of occurrence of a dimer configuration M is given by the dimer Boltzmann measure:

P dimer (M) = ∏ e∈M ν e Z dimer ,
where Z dimer = ∑ M∈M (G) ∏ e∈M ν e is the normalizing constant known as the parti- tion function.

The dimer model has the attractive feature of being exactly solvable, meaning that there is an explicit expression for the partition function. This fundamental result is due to Kasteleyn [Kas61] and independently to Temperley and Fisher [START_REF] Hnv Temperley | Dimer problem in statistical mechanics-an exact result[END_REF]. It relies on the Kasteleyn matrix, which we now define. A Kasteleyn orientation of the graph G is an orientation of the edges such that every face is clockwise odd meaning that, when traveling clockwise around edges of a face, an odd number of them are co-oriented. The corresponding Kasteleyn matrix K is the associated signed, weighted, adjacency matrix of the graph: lines and columns of K are indexed by vertices of G, and for any two vertices x, y of G, the coefficient K x,y is:

K x,y =      ν xy
if x ∼ y and the edge xy is oriented from x to y -ν xy if x ∼ y and the edge xy is oriented from y to x 0 otherwise.

Theorem 1 ([Kas61, TF61]). The dimer partition function is,

Z dimer = | Pf(K)| = det(K).
Remark 1. The Pfaffian, denoted by Pf, of a skew-symmetric matrix is a polynomial in the entries of the matrix, which is a square root of the determinant.

Proof (sketch). Refer to [START_REF] Kasteleyn | Graph theory and crystal physics[END_REF] for details. When writing out the Pfaffian of an adjacency matrix as a sum over permutations, there is exactly one term per dimer configuration. The issue is that each term comes with a signature. The purpose of the Kasteleyn orientation is to compensate the signature with signs of coefficients.

As a consequence of Theorem 1, Kenyon derives an explicit expression for the dimer Boltzmann measure.

Theorem 2 ([Ken97]). Let E

= {e 1 = x 1 y 1 , • • • , e k = x k y k } be
a subset of edges of the graph G. Then, the probability that these edges are in a dimer configuration chosen with respect to the dimer Boltzmann measure is;

P dimer ({e 1 , • • • , e k }) = k ∏ i=1 K x i ,y i Pf((K -1 ) E ) ,
where (K -1 ) E is the submatrix of K -1 whose lines and columns are indexed by vertices of E.

Proof. Since the proof is very short, we repeat it here. The probability of dimer configurations containing edges of E is given by the weighted number of these configurations divided by the weighted number of all configurations. By expanding the Pfaffian along lines and columns, this yields:

P dimer ({e 1 , • • • , e k }) = k ∏ i=1 K x i ,y i Pf(K E c ) Pf(K) . Using Jacobi's formula for Pfaffians [IW00], | Pf(K E c )| = | Pf(K) Pf((K -1 ) E )| yields the result.
In Section 3.1, we state Kenyon's results for the ∂ operator on infinite, bipartite, isoradial graphs. Then, in Section 3.2, we relate them to the corresponding dimer model.

Dirac operator and its inverse

As mentioned in Section 2.4, the ∂ operator, also referred to as the Dirac operator in [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF], can be extended to infinite bipartite isoradial graphs. Let G be such a graph, the set of vertices can be divided into two subsets B ∪ W , where vertices in B (the black ones) are only incident to vertices in W (the white ones). Then ∂ maps C B to C W : let f be a C-valued function on black vertices, then

( ∂ f )(w) = ∑ b∼w ∂w,b f (b),
where, if w, x, b, y, are the vertices of the rhombus R wb in cclw order, ∂w,b = i(xy) is given by Formula (2).

One of the main results of the paper [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF] is an explicit expression for the coefficients of the inverse ∂ , as a contour integral of an integrand which only depends on a path joining the two vertices. In order to state the result, let us introduce the following functions defined on vertices of G. Fix a reference white vertex w 0 and let λ be a complex parameter. Let v be a vertex of G, and

w 0 = v 0 , v 1 , • • • , v k = v be a
path in G ⋄ from w 0 to v. Each edge v j v j+1 has exactly one edge of G (the other is an edge of G * ). Direct it away from this vertex if it is white, and towards it if it is black, and let e iα j be the corresponding vector. The function f v is defined inductively along the path, starting from f v 0 ≡ 1, and

f v j+1 (λ ) =        f v j (λ )(λ -e iα j
) if the edge v j v j+1 leads away from a black or towards a white vertex,

f v j (λ ) (λ -e iα j )
otherwise.

It is easy to see that the function is well defined, i.e. independent of the choice of path from w 0 to v. An important point is that, by using a branched cover of the plane over w 0 , the angles α j are defined in R and not [0, 2π[. Then, Kenyon has the following theorem.

Theorem 3. Coefficients of the inverse Dirac operator can explicitly be expressed as:

( ∂ -1 ) b,w 0 = 1 4π 2 i C f b (λ ) log(λ )dλ , ( 3
)
where C is a closed contour oriented counterclockwise, containing all poles of the integrand and with the origin in its exterior.

Remark 2. The remarkable features of this theorem are the following. Explicit computations become tractable, since they only involve residues of rational functions. Moreover, the Formula (3) has the surprising feature of being "local", meaning that if the graph is changed away from the vertices w 0 , b and a path joining them, the value of the inverse Dirac operator stays the same.

Dirac operator and dimer model

The Dirac operator ∂ can be represented by an infinite matrix K, whose lines (resp. columns) are indexed by white (resp. black) vertices of G, and the coefficient K w,b is given by:

K w,b = ∂w,b = i(x -y) = 2 sin(θ wb )e iα wb ,
where θ wb is the half-angle of the rhombus R wb , and e iα wb is the unit-vector in the direction from w to b. The matrix K resembles a Kasteleyn matrix of the graph G, but differs in two instances.

1. Rows and columns are indexed by only "half" the vertices. To overcome this, one can define a weighted adjacency matrix K whose lines and columns are indexed by all vertices of the graph:

K = 0 K -K T 0 . Then, we would have Pf( K) = ± det( K) = ±| det(K)|.
2. The weights of the edges have an extra complex factor of modulus one, instead of a sign given by a Kasteleyn orientation. In [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF], Kenyon shows that the matrix is Kasteleyn flat, meaning that for each face of G whose vertices are labeled by u 1 , v 1 , • • • , u m , v m , the coefficients of the matrix K satisfy:

arg[K u 1 ,v 1 • • • K u m ,v m ] = arg[(-1) m-1 K v 1 ,u 2 • • • K v m ,u 1 ].
Using a result of Kuperberg [START_REF] Kuperberg | An exploration of the permanent-determinant method[END_REF], this implies that when the graph is finite and planar, K behaves as a usual Kasteleyn matrix, thus yielding an explicit formula for the partition function.

One of the main conjectures of the paper [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF] is to show that the inverse ∂ operator, which we write K -1 using the matrix notation, is related to the dimer model on the infinite, bipartite, isoradial, graph G. Its relation to the Gibbs measure is given in the following.

Theorem 4 ([dT07c]). There is a unique Gibbs measure P dimer on dimer configurations of G, such that the probability of occurrence of a subset of edges E{e

1 = w 1 b 1 , • • • , e k = w k b k } in a dimer
configuration of G chosen with respect to the Gibbs measure P dimer is:

P dimer ({e 1 , • • • , e k }) = k ∏ i=1 K w i ,b i det 1≤i, j≤k [(K -1 ) b i ,w j ], (4)
where (K -1 ) b i ,w j = ( ∂ -1 ) b i ,w j is given by (3).

Proof (sketch). Fix the edge set E, and cut out a finite piece of the rhombus graph G ⋄ containing E and paths joining vertices of E. Use a result of [dT07c] by which any finite piece of the graph G ⋄ can be filled with rhombi in order to become part of a periodic rhombus tiling of the plane. Define the probability of the edge set E as the weak limit of the Boltzmann measures on the natural toroidal exhaustion of the periodic graph. Use the uniqueness of the inverse Dirac operator and its locality property to deduce that this expression coincides with (4). Use Kolmogorov's extension theorem to show existence of a unique measure on G which specifies as (4) on cylinder events.

Other results

Free energy. Assume that the infinite bipartite isoradial graph is also periodic, and let G 1 = G/Z 2 . In [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF], Kenyon proves an explicit expression for the "logarithm of the normalized log of the Dirac operator", which has the surprising feature of only depending on rhombus angles of G 1 :

log(det 1 ∂ ) = 1 |V (G 1 )| ∑ e∈E(G 1 ) 1 π L(θ ) + θ π log(2 sin θ e ) .
Kenyon conjectures it to be related to the free energy of the corresponding dimer model, result that we prove this conjecture in [dT07b].

Gaussian free field.

There is a natural way of assigning a height function to dimer configurations of G. In [dT07a], we show that when dimer configurations are chosen with respect to the Gibbs measure of Theorem 4, the fluctuations of the height function are described by a Gaussian Free Field. Note that this proof holds for all dimer models defined on bipartite isoradial graphs, in particular for Z 2 and the honeycomb lattice.

Ising Model

The Ising model, first considered by Lenz [Len20], has been introduced in the physics literature as a model for ferromagnetism. The vertices of a graph G = (V, E) represent atoms in a crystal with magnetic moment reduced to one component that can take two values ±1. A spin configuration is thus a function σ on V , with values in {+1, -1}. Magnets with opposite moments tend to repel each other, which has a cost in terms of energy. From these considerations, the energy of a configuration σ on a finite graph is defined by:

E (σ ) = -∑ e=uv∈E J e σ u σ v ,
where the positive numbers (J e ) are called interaction constants. The probability of a occurrence of a spin-configuration is then defined using the Ising Boltzmann measure:

P Ising (σ ) = exp(-E (σ )) Z Ising .
The normalising factor, Z Ising = ∑ σ exp(-E (σ )) is the Ising partition function on G.

If the graph is drawn on an orientable surface with no boundary2 , then the partition function can be written as a combinatorial sum over contours. There are in fact two well-known such expansions: the so-called high-temperature expansion (contours on G) and the low-temperature expansion (contours on G * ). The contours in the low-temperature expansion have a nice interpretation in terms of the dual closed curves separating zones of different spins.

The most studied case is the Ising model on (pieces of) regular lattices where the interaction constants are taken to be constant equal to β , the inverse temperature. Kramers and Wannier [START_REF] Kramers | Statistics of the two-dimensional ferromagnet. part i[END_REF] discovered a duality for the Ising model on the square lattice: the measure on contours on G * coming from the low-temperature expansion is equal to that obtained by considering the high temperature expansion of the Ising model on G * at another temperature β * satisfying sinh(2β ) sinh(2β * ) = 1. This showed that the self-dual temperature for which β = β * = log 1 + √ 2 is the critical point if this one is unique. Later, Lebowitz and Martin-Löf [START_REF] Lebowitz | On the uniqueness of the equilibrium state for Ising spin systems[END_REF] proved that this is indeed the case.

Note that the case of Z 2 is singular since the graph is isomorphic to its dual. This duality can be generalized to any graph and its dual. For example the lowtemperature contour measure for the Ising model on the honeycomb lattice with inverse temperature β is the same as the high-temperature measure for the Ising model on its dual which is the triangular lattice. But using star-triangle transformations, one can go back to the Ising model on the honeycomb lattice. This generalizes the duality we had for the square lattice.

On isoradial graphs, instead of taking the same interaction constants for all edges, it is natural to take J e to be a function of the angle θ e . In order to ensure some integrability at the discrete level, we can impose that the Boltzmann weights satisfy the star-triangle relations, yielding a one-parameter family of Z-invariant interaction constants. There is a generalized duality relation: the high temperature expansion on G corresponds to the low temperature expansion on G * for a dual value of parameter. We qualify as critical the interaction constants corresponding to the self-dual value of the parameter, given by the following formula,

J(θ e ) = 1 2 log 1 + sinθ e cos θ e ,
and refer to the corresponding Ising model as critical Z-invariant. See [Bax89], chapter 7.13 for the parametrization of the star-triangle relation. They coincide with the critical inverse temperatures of the homogeneous triangular, square, and honeycomb lattices for θ e identically equal to π 6 , π 4 , π 3 respectively. The Ising model on two-dimensional graphs is in correspondence with other well-known models of statistical mechanics: the dimer model and the q-random cluster model with q = 2. We differ the discussion of the correspondence with dimers to Subsection 4.2.

There are essentially two ways of taking the large graph limit of the Ising model. Either you let the mesh tend to zero at the same time as the number of vertices of the graph goes to infinity, in order to get continuous objects in a bounded region of the plane (macroscopic level); or you let the graph tend to infinity, keeping the mesh size fixed, to get a model of statistical mechanics on an infinite graph.

The first approach has been adopted by Chelkak and Smirnov [CSb] to prove conformal invariance of the Ising model on bounded domains of the plane with Dobrushin boundary conditions. The second approach was used by the authors in [START_REF] Boutillier | The critical Z-invariant ising model via dimers: the periodic case[END_REF][START_REF] Boutillier | The critical Z-invariant ising model via dimers: Locality property[END_REF] to construct a Gibbs measure for the critical Ising model on infinite isoradial graphs.

Before explaining these two approaches, let us mention some works by physicists on the Z-invariant model at and off criticality by Au-Yang and Perk [AYP87, AP04, AP07] and by Martinez [START_REF] Reyes Martìnez | Correlation functions for the Z-invariant ising model[END_REF][START_REF] Reyes Martìnez | Multi-spin correlation functions for the Z-invariant ising model[END_REF].

Conformal invariance

Let Ω be a bounded domain of R 2 , and fix two points a and b on the boundary. For every ε, let G ⋄ ε be an approximation of Ω by a rhombus-graph3 with rhombi of side length ε. Let also a ε and b ε be approximations of a and b, located at the center of half-rhombi on the boundary.

Consider the Ising model on the isoradial graph G ε corresponding to black vertices of G ⋄ ε , with Dobrushin boundary conditions, i.e fixed spins +1 on one arc from a ε to b ε , and -1 on the other. Introduce, following Chelkak and Smirnov [CSb]

F G ε ,a ε ,b ε (z) = F ε (z) ∝ Zε (z) Zε (b ε ) ,
up to a multiplicative factor, introduced for technical reasons.

Chelkak and Smirnov prove that the function F ε and is discrete holomorphic by comparing configurations differing by local arrangements near z, and that it solves some discrete Riemann-Hilbert boundary value problem. They then deduce using approximation results [CSa] that, as ε → 0, the function F ε converges to the function solving the analogue continuous Riemann-Hilbert boundary value problem, which is conformally covariant.

Moreover, they prove that this observable satisfies a martingale property with respect to the growth of the curve from a ε to b ε : given the first n steps of the curve (γ 0 = a ε , γ 1 , . . . , γ n ), the expected value of the observable

F G ε \[γ 0 ,γ n+1 ],γ n+1 ,b ε (z) is equal to F G ε \[γ 0 ,γ n ],γ n ,b ε (z).
This martingale property together with the convergence to a conformally covariant object is sufficient to imply the convergence of the interface between a and b to a chordal SLE, with parameter κ = 3.

In the same paper, Chelkak and Smirnov construct another discrete holomorphic observable using not the loops separating clusters of spins, but those separating the clusters in the corresponding random cluster model4 with q = 2. This observable is a direct generalization of the one introduced by Smirnov in [START_REF] Smirnov | Conformal invariance in random cluster models. i. holomorphic fermions in the ising model[END_REF] for the square lattice. Again, this observable is discrete holomorphic and satisfies a martingale property. In this case, the scaling limit of the interface is a chordal SLE(16/3).

The two-dimensional Ising model as a dimer model

It turns out that a general Ising model on a graph drawn on a surface without boundary can be reformulated, through its contour expansion, as a dimer model on a decorated graph. This correspondence is due to Fisher [START_REF] Fisher | On the dimer solution of planar ising models[END_REF]. Since a lot of exact computations can be carried out on the dimer model, this correspondence has been proven to be a useful tool to study the Ising model [START_REF] Mccoy | The two-dimensional Ising model[END_REF] ([Kas67], . . . )

The idea is to construct a decorated version G of the graph G with a measure preserving mapping between polygonal contours of G and dimer configurations of G . This is an example of holographic reduction [START_REF] Valiant | Holographic algorithms[END_REF]. We now present the version used by the authors in [START_REF] Boutillier | The critical Z-invariant ising model via dimers: the periodic case[END_REF][START_REF] Boutillier | The critical Z-invariant ising model via dimers: Locality property[END_REF], which has the advantage of using decoration with cyclic symmetry. For other versions of Fisher's correspondence, see [START_REF] Fisher | On the dimer solution of planar ising models[END_REF], [START_REF] Cimasoni | A generalized kacward formula[END_REF].

The decorated graph, on which the dimer configurations live, is constructed from G as follows. Every vertex of degree k of G is replaced by a decoration consisting of 3k vertices: a triangle is attached to every edge incident to this vertex, and these triangles are linked by edges in a circular way, see Figure 7 below. This new graph, denoted by G , is also embedded on the surface without boundary and has vertices of degree 3. It is referred to as the Fisher graph of G.

Here comes the correspondence: to any contour configuration C coming from the high-temperature expansion of the Ising model on G, we associate 2 |V (G)| dimer G G Fig. 7 Left: a vertex of G with its incoming edges. Right: corresponding decoration in G . configurations on G : edges present (resp. absent) in C are absent (resp. present) in the corresponding dimer configuration of G . Once the state of these edges is fixed, there is, for every decorated vertex, exactly two ways to complete the configuration into a dimer configuration. Figure 8 below gives an example in the case where G is the square lattice Z 2 . In order to have a measure preserving mapping with dimers, the dimer weight function ν is given for the special case of the Z-invariant critical Ising model by, ν e = cot θ 2 if e comes from a rhombus of G with half-angle θ , ν e = 1 if e belongs to a decoration.

A Kasteleyn matrix K is constructed on G . One of the main result of [BdT10a] is that its inverse has a 'local' expression in the spirit of that obtained by Kenyon in [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF] for the inverse Dirac operator, see also Theorem 3.

Theorem 5 ([BdT10a]

). If x (resp. y) is a vertex of G and belongs to the decoration corresponding to the vertex x (resp. y) of G, then (K -1 ) x,y has the following expression:

(K -1 ) x,y = 1 4π 2 C x,y g x (λ )g y (-λ ) Exp x (y; λ ) log(λ )dλ ,
where g x (resp. g y ) is a simple rational fraction of λ that depends only on the geometry of the decorated graph in an immediate neighborhood of x (resp. of y), and the contour of integration C x,y is a simple closed curve containing all poles of the integrand, and avoiding the half-line starting from x in the direction of y.

Then, as for the dimer model, we obtain an explicit expression for a Gibbs measure, and recover the explicit formula for the partition function, obtained by Baxter [START_REF] Baxter | Free-fermion, checkerboard and Z-invariant lattice models in statistical mechanics[END_REF], at the critical point, see [START_REF] Boutillier | The critical Z-invariant ising model via dimers: Locality property[END_REF] for precise statements.

Other models

We briefly discuss now some aspects of other models of statistical mechanics: first some models related to the Laplacian, then the q-Potts model and its random cluster representation, and finally the 6-vertex and 8-vertex models.

Random walk and the Green function

We already mentioned the Laplace operator on functions of vertices, corresponding to conductances on edges given by c(e) = c(θ e ) = tan θ e . For this particular choice of conductances, the associated random walk is a martingale and the covariance matrix associated to a step is scalar. Up to a time reparametrization, its scaling limit is a standard two-dimensional Brownian motion [START_REF] Chelkak | Discrete complex analysis on isoradial graphs[END_REF][START_REF] Beffara | Is critical 2D percolation universal? In In and out of equilibrium[END_REF]. Kenyon [START_REF] Kenyon | The laplacian and dirac operators on critical planar graphs[END_REF] proves a local formula for the inverse of the Laplacian, the Green function ∆ -1 , on an infinite isoradial graph: ∀x, y ∈ V (G), ∆ -1

x,y = -

1 8iπ 2 C Exp x (y; λ ) log λ λ dλ .
where C is a positively oriented contour containing in its interior all the poles of the discrete exponential function, and the cut of the determination of the logarithm. He also proves precise asymptotics of G, improved in [B 08]. This Green function gives information on the random walk on G with conductances tan θ , but can also be used to gather properties of other models from statistical mechanics with an interpretation in terms of electric network, such as random spanning trees.

A spanning tree of G is an acyclic connected subgraph of G whose vertex set contains all vertices of G. On a finite graph, if the weight of a spanning tree equals the product of the conductances of its edges, the partition function is given, via Kirchhoff's matrix-tree theorem [START_REF] Kirchhoff | über die Auflösung der Gleichungen, auf welche man bei der Untersuchung der linearen Vertheilung galvanischer Ströme geführt wird[END_REF], by the determinant of any principal minor of the Laplacian. One can construct a measure on spanning trees of G as limit of measures on finite graphs where a spanning tree would have a probability proportional to its weight. The statistics of edges present in the random spanning tree are given by a determinantal process on edges: the probability that edges {e 1 = v 1 w 1 , . . . , e k = v k w k } are present is given by P tree ({e 1 , ..,

e k }) = k ∏ i=1 tan θ i det 1≤i, j≤k
H(e i , e j ) ,

where H is the impedance transfer matrix [START_REF] Burton | Local characteristics, entropy and limit theorems for spanning trees and domino tilings via transfer-impedances[END_REF] defined by H(e i , e j ) = ∆ -1 (v i , v j ) -∆ -1 (v i , w j ) -∆ -1 (w i , v j ) + ∆ -1 (w i , w j ).

q-Potts models and the random cluster model

A natural generalization of the Ising model is the Potts model with q states (or colors). Spins at neighboring vertices interact only if they have different colors.

The energy of a spin configuration σ is E (σ ) = -∑ e=uv∈E J e I σ u =σ v . For q = 2, we recover the Ising model. This model can be mapped to the (inhomogeneous) random cluster model [START_REF] Grimmett | The random-cluster model, volume 333 of Grundlehren der Mathematischen Wissenschaften[END_REF] as follows:

• spins of different colors are in different percolation clusters,

• neighboring spins sharing an edge e of the same color are connected with probability, p e = 1e -2J e .

Once again, we want the interaction constants to be functions of the half-angle of the corresponding rhombus. It is possible to take them to satisfy the Yang-Baxter equation. Moreover, if we impose a generalized self-duality, then there is a unique choice for the function p [Ken02a]:

p(θ ) = √ q sin 2r π θ sin 2r π ( π 2 -θ ) + √ q sin 2r π θ

, with r = cos -1 ( √ q/2) Although a formula for the free energy of that model is known [START_REF] Baxter | Exactly solved models in statistical mechanics[END_REF], little is known except for the case q = 2 corresponding to the Ising model.

6-vertex and 8-vertex models

To conclude, let us briefly mention that the 6-and 8-vertex models, introduced originally on the square lattice, have natural generalization on isoradial graphs, or more precisely on the dual of its rhombic graph G ⋄ * , which is a 4-regular graph. A configuration of the 6-vertex model (resp. of the 8-vertex model) is an orientation of the edges of G ⋄ * such that the number of ingoing edge (and thus of outgoing edges) at each vertex is equal to 2 (resp. is even). These models can be solved through Bethe Ansatz on the square lattice (in the sense that their free energy can be computed [START_REF] Lieb | Exact solution of the problem of the entropy of two-dimensional ice[END_REF][START_REF] Baxter | Eight-vertex model in lattice statistics[END_REF]), and possess a very rich algebraic structure. There are some studies on their Z-invariant generalizations [START_REF] Baxter | Perimeter bethe ansatz[END_REF][START_REF] Baxter | Free-fermion, checkerboard and Z-invariant lattice models in statistical mechanics[END_REF], but many questions still need to solved.

Fig. 3

 3 Fig. 3 Transfer matrices for an inhomogeneous model on the m × n torus.

Fig. 4

 4 Fig. 4 Star-star transformation (the name star referring to the 3-branches stars).

Fig. 5

 5 Fig. 5 Star-triangle transformation for isoradial graphs.

  below. Vertices of G D are decomposed into two classes: black vertices, which are vertices of G and of G * , and white vertices, which are the centers of the rhombi of G ⋄ , seen either as edges of G or as those of G * . A black vertex b and a white vertex w are connected by an edge in G D if the vertex corresponding to b is incident to the edge corresponding to w in either G or G * . That is, edges of G D correspond to half diagonals of rhombi of G ⋄ . It turns out that G D has an isoradial embedding (with rhombi of side-length 12 ) obtained by splitting each rhombus of G ⋄ into four identical smaller rhombi.

Fig. 6

 6 Fig.6The double graph G D of the graph on the right of Figure1: in thick full lines is the original graph G, in full lines is the dual graph G * , in dotted lines is the underlying G ⋄ graph.

  , a twisted version of the partition function, considered as a function on edges of G ε , or equivalently on inner rhombi of G ⋄ ε . Let z denote a generic rhombus of G ⋄ ε , ZG ε ,a ε ,z = Zε (z) = sin θ z γ) is the winding of the curve γ from a to z. If we remove the prefactor and the contribution of the winding, this would be for z = b ε the low temperature expansion of the partition function of the Ising model on Ω ε with Dobrushin boundary conditions between a ε and b ε .Now define for all inner rhombi z of G ⋄ ε ,

Fig. 8

 8 Fig. 8 Polygonal contour of Z 2 , and corresponding dimer configurations of the Fisher graph.

One can define in a similar way a Laplacian on G * by restricting the same operator ∂ ∂ to G *

The correspondence can be extended to surfaces with boundary by including in addition to the closed contours a certain number of paths connected to boundary.

On the boundary, we put only half-rhombi such that only "black" vertices are exposed on the boundary.

also known as the Fortuin-Kasteleyn percolation.


Dan Beltoft

Cédric Boutillier

Nathanaël Enriquez

Random Young Diagrams in a Rectangular Box

Keywords: Young diagrams, Gauss polynomials, Ornstein-Uhlenbeck process. Math. Subject Classification: 60C05, 60F17, 05A10, 05A15, 05A30

We exhibit the limit shape of random Young diagrams having a distribution proportional to the exponential of their area (grand-canonical ensemble), and confined in a rectangular box. The Ornstein-Uhlenbeck bridge arises from the fluctuations around the limit shape. The fluctuations for the unconfined case lead to a two-sided stationary Ornstein-Uhlenbeck process.

Introduction

A partition of an integer n is a finite non-increasing sequence of integers π = (π 1 , π 2 , . . . , π k ), with π 1 + • • • + π k = n. A convenient graphical representation in (Z + ) 2 of the partition π is its Young diagram, also denoted by π. The Young diagram consists of a stack of π i unit squares on the i-th column. See Figure 1. The area |π| of the Young diagram π is the total number of squares, which is equal to n.

If the number of summands k is less than some integer a and all the summands π j are less than b, then the Young diagram stays in the box [0, a] × [0, b].

In this paper we study the limit shape and fluctuations of large random Young diagrams assigned to stay inside a (large) rectangular box, when the probability measure is proportional to q |π| , and q is a suitable parameter. It is easy to see that if one doesn't want the system to degenerate in the limit, one has to make q depend on the size of the box. Namely 1q has to be taken of the order of the inverse of the side of the box.

The study of the combinatorics of partitions of integers goes back to Hardy and Ramanujan and then Erdős in the 40's. The statistical physics point of view was introduced by Vershik [START_REF] Vershik | Statistical mechanics of combinatorial partitions, and their limit configurations[END_REF] in his study of the typical shape of the partition of a large integer and obtained what we will call Vershik's curve. Recently, Funaki and Sasada obtained Vershik's result, in [START_REF] Funaki | Hydrodynamic limit for an evolutional model on two-dimensional young diagrams[END_REF], as a by-product of an hydrodynamic result for the corresponding particle system model.

Our study lies at the intersection of classical topics of probability theory, combinatorics, and statistical physics.

Along the paper, we make an extensive use of the classical Gauss polynomials, wellknown in combinatorics as the generating functions for the number of Young diagrams with given area (see e.g. [START_REF] Andrews | Integer partitions[END_REF]). In order to get to the asymptotic regime, we are led to state a q-analogue for Stirling's formula.

For all values of the parameters of the problem (limiting aspect ratio of the box, and c = lim n(1q)), the limit shape obtained turns out to be a restriction of Vershik's curve as it was recently noticed by Petrov in [START_REF] Petrov | Limit shapes of young diagrams. two elementary approaches[END_REF].

The main part of the paper deals with the fluctuations around this limit shape and requires a fine understanding of the boundary of the random diagram. From a probabilistic point of view, the boundary is described by a Markov which can be considered as the q-analogue of the bridge of the simple random walk. In the classical case (c = 0), the properly rescaled interface converges to the Brownian bridge. In the general case, the limiting process is an Ornstein-Uhlenbeck bridge. Let us mention that, in the framework of unlimited partitions, and under some specifications on the summands, Vershik and Yakubovich [START_REF] Vershik | The limit shape and fluctuations of random partitions of naturals with fixed number of summands[END_REF][START_REF] Yu | The central limit theorem for normalized Young diagrams of partitions into different summands[END_REF] already observed Gaussian fluctuations. For strict partitions, see [START_REF] Vershik | A local limit theorem for random partitions of natural numbers[END_REF]. The fluctuations for the unconfined case lead to a two-sided stationary Ornstein-Uhlenbeck process. However, in contrast with the limit shape problem, there is no direct argument to deduce fluctuations for partitions in a box from the unlimited case.

The paper is organized as follows. After a presentation of the combinatorics of the problem (Section 2), we perform its asymptotic analysis through a q-version of Stirling's formula (Section 3). Section 4 is devoted to the derivation of the limit shape phenomenon. In Section 5, we study the fluctuations around the limit shape: we first compute the limiting 2-correlation function and deduce the convergence of finite-dimensional marginals. This section ends with the delicate proof of the tightness of the fluctuations.

Presentation of the model

We study the asymptotic distribution of random Young diagrams fitting in a large rectangular box with dimensions a×b. Given a real number q > 0, we assign to each diagram π the probability

P q a,b (π) = q |π| Z a,b (q) ,
where |π| is the number of boxes of π and Z a,b (q) is the partition function, the sum of all q |π| . For the sake of clarity, the probability measure P q a,b will be simply denoted, by P.

We fix a parameter ρ ∈ (0, 1), and choose the dimensions of the box for each n to be a n × b n , where (a n ) and (b n ) are sequences of positive integers satisfying:

a n + b n = 2n, lim n→∞ a n b n = ρ 1 -ρ .
We are interested in the limiting behavior when n goes to ∞. To obtain a non-degenerate limit, q must go to 1 as n goes to infinity. We fix a real parameter c and pose q = e -c n . Note that c = 0 corresponds to a uniform probability measure. We will assume that c = 0. The results for c = 0 can be obtained by taking limits. The physical meaning of the parameter c is that of a pressure, since it is the variable conjugated to the (two-dimensional) volume.

Combinatorics of the partitions

Let us start with a fixed box with dimensions a × b. The partition function Z a,b (q) is expressed in terms of Gaussian polynomials or q-binomial coefficients, where integers j are replaced by their q-analogues (j) q = 1-q j 1-q . Definition 1. Let q > 0,

n! q = n j=1 (j) q = n j=1 1 -q j 1 -q n m q = n! q (n -m)! q m! q = m-1 j=0 1 -q n-j 1 -q m-j . Lemma 1. For all (a, b) ∈ (N * ) 2
, and all q > 0, the partition function Z a,b (q) is equal to

Z a,b (q) = a + b a q = a + b b q .
Proof. a+b a q and Z a,b (q) both follow the recursion relation

Z a,b (q) = Z a,b-1 (q) + q b Z a-1,b (q),
where the first term corresponds to those diagrams with all parts strictly smaller than b, and the second term to the diagrams with at least one part of size b.

We use the following coordinates. The bounding rectangle is the rectangle in the plane with corners at the points (0, 0) and (a + b, ba) and sides with slopes ±1. The boundary of a diagram is encoded as a lattice path X k 0≤k≤n from the origin to the point (a + b, ba) and such that Computing the probability P (X k = ℓ) is straightforward. A lattice path passing through (k, ℓ) is composed of a path from (0, 0) to (k, ℓ) and a path from (k, ℓ) to (a + b, ba). Since X k has the same parity as k, the exact formulas also depend on this parity.

X k+1 = X k ± 1 for all k. 2n -a n b n 2k X 2k
Proposition 1. The 1-dimensional marginal of X under P is given by

P q a,b (X 2k = 2i) = q (k+i)(a-k+i) Z a,b (q) Z k-i,k+i (q)Z a-k+i,b-k-i (q) (1)
and

P q a,b (X 2k+1 = 2i + 1) = q (k+i+1)(a-k+i) Z a,b (q) Z k-i,k+i+1 (q)Z a-k+i,b-k-i-1 (q) (2)
We focus our attention to the behaviour of the process at even times, which is sufficient to study its scaling limit, since the process has bounded steps (+1/-1).

We now mention a unimodality result for the distribution above which turns out to be useful later.

k -i a n k + i b n a n -k + i b n -k -i 2k 2i Figure 3: Illustration of Equation (1). Lemma 2. The function ℓ → P q a,b (X 2k = 2ℓ) is unimodal: there exists an integer L ♯ (k) = L ♯ a,b,q (k) such that P q a,b (X 2k = 2(ℓ + 1)) P q a,b (X 2k = 2ℓ) ≤ 1 ⇐⇒ ℓ ≥ L ♯ (k). ( 3
)
Proof. Writing the ratio

P q a,b (X 2k = 2(ℓ + 1)) P q a,b (X 2k = 2ℓ) = Z k-ℓ-1,k+ℓ+1 (q)Z a-k+ℓ+1,b-k-ℓ-1 (q) Z k-ℓ,k+ℓ (q)Z a-k+ℓ,b-k-ℓ (q) q (a-k+ℓ+1)(k+ℓ+1)-(a-k+ℓ)(k+ℓ) = (1 -q k-ℓ )(1 -q b-k-ℓ ) (1 -q k+ℓ+1) (1 -q a-k+ℓ+1 ) q a+2ℓ+1
This ratio is smaller than 1 if and only if

(1 -q k-ℓ )(1 -q b-k-ℓ )q a+2ℓ+1 ≤ (1 -q k+ℓ+1 )(1 -q a-k+ℓ+1 )
⇐⇒ (q a+b+1 -1) ≤ q a+1 (q -1)q 2ℓ + q k+1 (q a -1) + q a-k+1 (q b -1) q ℓ

Both terms on the right hand side are increasing continuous functions of ℓ, proving the existence of an integer L ♯ (k) with the asserted property.

Proposition 2. The 2-dimensional marginal of (X 2k ) 0≤k≤n is given by

P q a,b (X 2k = 2i, X 2l = 2j) = Z k+i,k-i (q)Z j+l-i-k,l+i-k-j (q)Z a-l+j,b-l-j (q) Z a,b (q) q (a-(k-i))(k+i)+(l+j-k-i)(a-l+j) .
Now we have exact expressions for probabilistic quantities, we let the dimensions of the box depend on n and investigate the limiting behaviour of the random Young diagram under P q an,bn in the regime

a n + b n = 2n → +∞, a n 2n → ρ ∈ (0, 1), -n log q → c ∈ R.
3 Asymptotics of q-factorials In order to study the limit of the process as the size of the box goes to infinity, we first need the asymptotic behavior of the q-factorial.

Proposition 3. Let ε > 0. For ℓ > (1 -q) -1 ε, (ℓ!) q = 2πℓ q exp ℓ 0 log t q dt (1 + O(1 -q))
beginequation* Consider log ℓ! q == ℓ k=1 log k q . By the Euler-Maclaurin formula, we have

log ℓ! q = ℓ 0 log(t q ) dt + 1 2 (log ℓ q ) + ℓ 1 B 1 (t) d log t q dt - 1 0 log(t q ) dt.
where B 1 (t) = t -⌊t⌋ -1 2 is the first periodic Bernoulli polynomial. The integral 1 0 log t q dt = 1 + O(1q). For the second integral in (3), we add and subtract 1 t to d log tq dt in order to get

ℓ 1 B 1 (t) d log t q dt dt = ℓ 1 B 1 (t) 1 t dt + ℓ 1 B 1 (t) d log t q dt - 1 t dt (4)
The function

f (t) = d log tq dt -1 t is continuous on R + and its derivative is bounded by an absolute constant M times (1 -q) 2 , i.e. |f ′ (t)| ≤ M (1 -q) 2 for all t ∈ R + . On the interval [k -1, k], the function B 1 (t) is given by B 1 (t) = t -k + 1 2 and has the antiderivative Bn,k (t) = 1 2 t 2 -(k -1 2 )t + k(k-1)
2 with Bn,k (k -1) = Bn,k (k) = 0. By partial integration we get

ℓ 1 B 1 (t)f (t) dt ≤ M ℓ-1 k=1 k+1 k | Bn,k (t)| dt = M ℓ(1 -q) 2 12
proving that this term is O(1q) as n → ∞. The first term of ( 4) is found to be

ℓ 1 B 1 (t) 1 t dt = ℓ -1 -(ℓ + 1 2 ) log ℓ + log ℓ! = -1 + 1 2 log 2π + O( 1 ℓ ) (5)
by the classical Stirling approximation ℓ! = √ 2πℓ ℓ e ℓ 1 + O( 1 ℓ ) . Combining (3), ( 4) and ( 5), we get the result.

Using the notations of Section 1, we restate the q-Stirling's formula in a convenient form for this particular context.

Corollary 1 (q-Stirling's Formula). Let c ∈ R * + and fix ε > 0 In the limit when n goes to infinity, with q = e -c n , the following asymptotics hold for all ℓ > εn

ℓ! q = √ 2πn e cℓ n -1 c n ℓ exp(nS c ( ℓ n )) 1 + O( 1 ℓ ) ,
where

S c (α) = α 0 log 1 -e -cx c dx = αS αc (1).
In particular, there exist positive constants m and M such that for all n, and all ℓ between 1 and n,

m √ 2πn e cℓ n -1 c n ℓ exp(nS c ( ℓ n )) < ℓ! q < M √ 2πn e cℓ n -1 c n ℓ exp(nS c ( ℓ n )).
We turn now to the the asymptotics of the 1-and 2-dimensional marginals of X under P given by Propositions 1 and 2. For this purpose, define f c (x, y) = S c (x + y) -S c (x) -S c (y) and h c (x, y) = c(e c(x+y) -1) (e cx -1)(e cy -1)

.

If x = ℓ n and y = k n , the asymptotics of the q-binomial coefficient is given by

Z (q) l,k = ℓ + k ℓ q = 1 2πn h c (x, y) exp(nf c (x, y)) 1 + O( 1 n ) .
As a consequence, Corollary 2. For all ε > 0, for all k between εn and (1-ε)n, with s = k/n and x = i/n, we have

P(X 2k = 2i) = 1 2πn H (1) ρ,c (s, x) exp(nF (1)
an 2n ,c (s, x)) 1 + O( 1 n ) ,
and there exists a constant M > 0 such that for all n, i and k,

P(X 2k = 2i) ≤ M 1 2πn
H (1) ρ,c (s, x) exp(nF (1) ρ,c (s, x)),

where

F (1) ρ,c (s, x) = -c(2ρ -s + x)(s + x) + f c (s -x, s + x) + f c (2ρ -s + x, 2 -2ρ -s -x) -f c (2ρ, 2 -2ρ), H (1) ρ,c (s, x) = h c (s + x, s -x)h c (2ρ -s + x, 2 -2ρ -s -x) h c (2ρ, 2 -2ρ) .
In the first equality, when replacing an 2n by ρ in the indices of H (1) , we make an error of order O( 1n ), which is absorbed in the factor (1 + O( 1 n )). Making the same substitution in the indices of F (1) would change the multiplicative constant in the asymptotics.

Corollary 3. For all ε > 0, for all k < l between εn and (1ε)n then, with s = k/n, t = l/n and x = i/n, y = j/n

P(X 2k = 2i, X 2l = 2j) = 1 2πn H (2) ρ,c (s, t, x, y) exp(nF (2)
an 2n ,c (s, t, x, y))(1 + O( 1 n )),
where

F (2) ρ,c (s, t, x, y) = -c (2ρ -s + x)(s + x) + (t -s + y -x)(2ρ -t + y) + f c (s -x, s + x) + f c (t -s + y -x, t -s -y + x) + f c (2ρ -t + y, 2 -2ρ -t -y) -f c (2ρ, 2 -2ρ), H (2) ρ,c (s, t, x, y) = h c (s + x, s -x)h c (t -s + y -x, t -s -y + x)h c (2ρ -t + y, 2 -2ρ -t -y) h c (2ρ, 2 -2ρ) .
We give now a slight refinement of Corollary 3 that will be useful in the fluctuations.

Corollary 4. For all ε > 0, for all k < l between εn and (1ε)n then, with s = k/n, t = l/n and x = i/n, y = j/n, for all

s ′ = s + o(1), t ′ = t + o(1), x ′ = x + o(1), y ′ = y + o(1): P n (X 2k = 2i, X 2l = 2j) = 1 2πn H (2) ρ,c (s ′ , t ′ , x ′ , y ′ ) exp(nF (2)
an 2n ,c (s, t, x, y))(1 + o(1)).

Limit Shape

We associate to the lattice path (X k ) the continuous piecewise linear function, X : s → X s defined on [0, 2n], which coincides with X k when s = k. The graph of the function X is the boundary of the random Young diagram we consider. Let L ρ,c be the function on [0, 1] defined by

∀t ∈ [0, 1], L ρ,c (t) = 1 -2ρ + 1 c log sinh(ct) + e c(2ρ-1) sinh(c(1 -t)) sinh c = 1 c log e -ct -e ct + e -c(2-2ρ-t) -e -c(t-2ρ) e -c(2-2ρ) -e 2cρ . ( 6
)
This is the limit shape of the rescaled random process, in the sense made precise in Theorem 1 below. First we need to estimate the proximity between the value of the function L ρ,c and the most probable value for X 2k denoted by L ♯ n (k) = L an,bn,e -c/n (k) in Lemma 2.

Lemma 3. For all n and all

0 ≤ k ≤ n, 1 n L ♯ n (k) -L ρ,c ( k n ) ≤ 1 n + an 2n -ρ , where q = e -c
n . Consequently, in the limit when n goes to infinity, and a n /2n goes to ρ,

1 n L ♯ n (k) -L ρ,c ( k n ) → 0. ( 7
)
Proof. The ratio of probabilities (3) used to define

L ♯ (k) = L ♯ n (k) can be rewritten in terms of the function (ρ, c, t, x) → R ρ,c (t, x) = e -c sinh( c 2 (t -x)) sinh( c 2 (2 -2ρ -t -x)) sinh( c 2 (t + x)) sinh( c 2 (2ρ -t + x))
,

namely as R an 2n ,c ( k n , ℓ n ). The function x → R ρ,c (t, x)
is decreasing, and by definition

L ♯ n (k) is the smallest integer ℓ such that R an 2n ,c ( k n , ℓ n ) ≤ 1.
On the other hand, a computation shows that for each t ∈ [0, 1], the equation

R ρ,c (t, x) = 1
has the (unique) solution x = L ρ,c (t). Substituting an 2n for ρ, we conclude that

1 n L ♯ n (k) -L an 2n ,c ( k n ) ≤ 1 n (8)
for all n and all k ≤ n.

Differentiating L ρ,c with respect to ρ, we find that

∂L ρ,c ∂ρ (t) ≤ 1
for all t, and so by the mean value theorem,

L an 2n ,c ( k n ) -L ρ,c ( k n ) ≤ a n 2n -ρ (9)
for all k ≤ n. From ( 8), ( 9) and the triangle inequality, we get the expected result.

Remark 1. It turns out that log R ρ,c (t, x) coincides with the derivative of F

ρ,c (s, x) with respect to x. As a consequence, L ρ,c (s) can also be viewed as the argmax of F (1) ρ,c (s, •), which is nonpositive and vanishes at L ρ,c (s). See Corollary 6.

Theorem 1. The boundary of the rescaled random Young diagram converges in probability, for the uniform topology, to the curve of

t → L ρ,c (t). ∀ε > 0, lim n→0 P e -c/n an,bn sup t∈[0,1] 1 2n X 2tn -L ρ,c (t) > ε = 0, where L ρ,c (t) = 1 c log e -ct -e ct + e -c(2-2ρ-t) -e -c(t-2ρ) e -c(2-2ρ) -e 2cρ .
Remark that in the case of a square box (ρ = 1/2), the expression for L ρ,c boils down to:

L 1 2 ,c (t) = 1 c log cosh(c(t -1 2 )) cosh c 2 Proof. Fix ε > 0. For t < ε 2 or t > 1 -ε 2 , the difference 1 2n X 2tn -Lρ, c(t) is always smaller than ε. We have to control what happens for t ∈ ( ε 2 , 1 -ε 2 )
. Using the fact that L ρ,c is differentiable, and its derivative with respect to t is bounded by 1, we have that

L ρ,c 1 n ⌊tn⌋ -L ρ,c (t) ≤ 1 n ,
which is smaller than ε 3 for n sufficiently large. The same is true of 1 2n X 2tn -1 2n X 2⌊tn⌋ . Thus, by an ε 3 -argument, to control the sup over ( ε 2 , 1 -ε 2 ), it is sufficient to control what happens at points of the form t = k n :

P sup t∈( ε 2 ,1-ε 2 ) 1 2n X 2tn -L ρ,c (t) > ε ≤ k∈N∩n( ε 2 ,1-ε 2 ) P X 2k > 2n(L ρ,c ( k n ) + ε) or X 2k < 2n(L ρ,c ( k n ) -ε) (10)
But it follows from [START_REF] Pittel | On a likely shape of the random Ferrers diagram[END_REF], that for n sufficiently large, 2n

L( k n ) + ε ≥ 2L ♯ n (k) for all k ∈ N ∩ n( ε 2 , 1 -ε 2 )
, and thus using the unimodality of the law of X 2k

P X 2k > 2n(L ρ,c ( k n ) + ε) = l>n(Lρ,c( k n )+ε) P(X 2k = 2l) ≤ n × P X 2k = 2 n(L ρ,c ( k n ) + ε) .
which by Corollary 2 is exponentially small, uniformly in k, as n goes to infinity. Thus for n sufficiently large, the sum on the RHS of (10) can be made smaller than any positive number.

Remark 2. From the proof and Lemma 3, we see that the convergence statement in the theorem holds uniformly for a family of sequences (a n ) and (b n ), so long as the convergence an 2n → ρ is uniform for the family.

Fluctuations

We now study the fluctuations of the interface around the limit shape. We define for n ∈ N * a new rescaled process:

∀t ∈ [0, 1] , Xt = X(n) t = √ n 1 2n X 2nt -L ρ,c (t) . (11)
We place ourselves in the space D of càd-làg paths on [0, 1] endowed with its usual topology. We state now our second main result for the convergence of the fluctuations of the interface to the Ornstein Uhlenbeck bridge. See the appendix for the definition and some properties of this process. Define

f (s) = 1 √ 2 (e 2cρ -1)(1 -e -2c(1-ρ) ) sinh cs + e c(2ρ-1) sinh(c(1 -s)) = 2 sinh(cρ) sinh(c(1 -ρ)) e c( 1 2 -ρ) sinh(cs) + e c(ρ-1 2 ) sinh(c(1 -s)) . ( 12
) Theorem 2. The sequence ( X(n) s /f (s)) n , converges weakly in D to the Ornstein-Uhlenbeck bridge (Y s ) s∈[0,1] , which is the Gaussian process on [0, 1] with covariance E[Y s Y t ] = sinh(cs) sinh(c(1 -t)) c sinh(c) , for 0 ≤ s < t ≤ 1.

Two-point correlations

To prove convergence of the two-dimensional marginal to a Gaussian process, we apply a saddle-point method, ie. we need to show that the function F (2) , which governs the exponential decay of 2-dimensional marginals (see Corollary 3), has a critical point 'on the limit shape', and that it takes its maximal value of 0 at this point. We will prove that this is indeed the case, but first we need a lemma describing the limit shape in a subrectangle.

Lemma 4. The limit shape L satisfies the relations

1 1 -s L ρ,c (t) -L ρ,c (s) = L ρ ′ ,c ′ t-s 1-s (13
)

and 1 t L ρ,c (s) = L ρ ′′ ,c ′′ s t ( 14
)
for 0 < s < t < 1, with ρ ′ = 2ρ -s + L ρ,c (s) 2(1 -s) , c ′ = c(1 -s)
and

ρ ′′ = t -L ρ,c (t) 2t , c ′′ = tc.
Proof. These relations can be checked analytically. We sketch now a less computational argument for Equation (14). The same strategy applies for Equation (13). Fix t ∈ (0, 1), and take a sequence of boxes with sidelengths [0, a n ] × [0, b n ] (with a n + b n = 2n, a n /2n → ρ). Theorem 1 states that, when n goes to infinity, the probability under P e -c/n an,bn that X nt /n converges to L ρ,c (t) goes to 1. As a consequence, the restriction of the limit shape L ρ,c to the, rescaled by n, lower left subbox [0,

t-Lρ,c(t) 2 ] × [0, t+Lρ,c(t) 2
] is the rescaled limit shape for boxes with ratio

t-Lρ,c(t) t-Lρ,c(t)+t+Lρ,c(t) = t-Lρ,c(t) 2t
. The rescaling factor is the perimeter of the subbox t. The value of the parameter q remains the same.

Since q = e -c/n = e -c ′′ nt , we get c ′′ = tc. This lemma basically describes the limit shape of the process restricted to a subrectangle defined by a point (s, L ρ,c (s)) and either the right or left corner of the original rectangle.

Lemma 5. For all 0 < s < t < 1, the functions F (1) and F (2) satisfy the relations

F (2) ρ,c (s, t, x, y) = F (1) ρ,c (s, x) + (1 -s)F (1)
ρ ′ ,c ′ (t ′ , y ′ ) ( 15
)
and

F (2) ρ,c (s, t, x, y) = F (1) ρ,c (t, y) + tF (1)
ρ ′′ ,c ′′ (s ′′ , x ′′ ) ( 16
)
with

ρ ′ = 2ρ -s + x 2(1 -s) , t ′ = t -s 1 -s , y ′ = y -x 1 -s , c ′ = c(1 -s)
and

ρ ′′ = t -y 2t , s ′′ = s t , x ′′ = x t , c ′′ = tc.
Proof. The relation (15) follows from the identity

P(X 2k = 2i, X 2ℓ = 2j) = P(X 2k = 2i)P(X 2ℓ = 2j|X 2k = 2i).
allied to scaling arguments similar to those used in the proof of the previous lemma. Equation ( 16) follows from the conditioning on X 2l instead of X 2k .

Lemma 6. The partial derivatives of F (1) with respect to ρ, s and x all vanish at x = L ρ,c (s).

Corollary 5. The partial derivatives of F (2) with respect to ρ, s, t, x and y all vanish at x = L ρ,c (s) and y = L ρ,c (t).

Proof. For the partial derivatives with respect to ρ, use (16). For ∂F (2) ∂s and ∂F (2) ∂x , use (16) together with (14). Similarly, for ∂F (2) ∂t and ∂F (2) ∂y , use (15) together with (13).

Corollary 6. F

ρ,c (s, L ρ,c (s)) = 0 and F

ρ,c (s, t, L ρ,c (s), L ρ,c (t)) = 0.

Proof. By the chain rule and Lemma 6 we get d ds (F (1) ρ,c (s, L ρ,c (s))) = 0.

Taking the limit s → 0 yields zero, proving the claim. For the second identity, use the first one together with ( 16) and (14).

Proposition 4. Assume an

2n = ρ + o 1 √ n . Let s < t in [0, 1]. The joint law of ( X(n) s , X(n) t )
converges to the law of the centered 2-dimensional Gaussian vector

( Xs , Xt ) = (f (s)Y s , f (t)Y t ),
where f (s) is defined in as in Equation (12), and (Y t ) t∈[0,1] is an Ornstein-Uhlenbeck bridge on the interval [0, 1] with parameter c (see Appendix).

Remark 3. In the case of a square (ρ = 1 2 ), the expression of f simplifies drastically to

f (x) = 1 √ 2 cosh c s -1 2 .
Proof. Let a 1 , b 1 , a 2 , b 2 be four real numbers such that a 1 < b 1 and a 2 < b 2 . In the following, we write L for L ρ,c . Clearly, X(n

) s ∈ [a 1 , b 1 ] if and only if X 2ns is in the interval 2nL(s) + 2a 1 √ n, 2b 1 √ n .
Since X 2ns -X 2⌊ns⌋ ≤ 2 and X 2nt -X 2⌊nt⌋ ≤ 2, we just need to compute the limit of

P X 2⌊ns⌋ , X 2⌊nt⌋ ∈ A 1 × A 2 .
as n goes to infinity. Set k n = ⌊ns⌋ and ℓ n = ⌊nt⌋. One has:

P X 2⌊ns⌋ , X 2⌊nt⌋ ∈ A 1 × A 2 = (i,j)∈ 1 2 A 1 × 1 2 A 2 P (X 2kn = 2i, X 2ℓn = 2j) = (i,j)∈ 1 2 A 1 × 1 2 A 2 1 2πn
H (2) ρ,c (s, t, L(s), L(t)) exp nF Here, F (2) ≤ 0 everywhere and the number of terms is O(n), so the error terms o(1) can be replaced with an o(1) term outside the sum. Thus, we are left with a Riemann sum for the double integral

1 2πn H (2) ρ,c (s, t, L(s), L(t)) nL(s)+b √ n nL(s)+a √ n nL(t)+d √ n nL(t)+c √ n exp nF (2) an 2n ,c kn n , ℓn n , u n , v n dv du (17)
in which we make the substitution

u = x √ n + nL(s), v = y √ n + nL(t) to get 1 2π H (2) ρ,c (s, t, L(s), L(t)) b a d c exp nF (2) an 2n ,c kn n , ℓn n , x √ n + L(s), y √ n + L(t) dy dx
We make a second degree Taylor expansion of the function F

(2) an 2n ,c kn n , ℓn n , •, • at the point (L(s), L(t)). Let z n be the point ( kn n , ℓn n , L(s), L(t)) in R 4 .

nF

(2)

an 2n ,c kn n , ℓn n , x √ n + L(s), y √ n + L(t) = nF (2) an 2n ,c (z n ) + ∂F (2) an 2n ,c ∂x (z n )x √ n + ∂F (2) an 2n ,c ∂y (z n )y √ n + 1 2 ∂ 2 F (2) an 2n ,c ∂x 2 (z n )x 2 + 1 2 ∂ 2 F (2) an 2n ,c ∂y 2 (z n )y 2 + ∂ 2 F (2) an 2n ,c ∂x∂y (z n )xy + O n -1 2 (18)
We must find the limit of this expression as n → ∞. The first term on the right has limit zero, which can be seen by Taylor expansion around (s, t, L(s), L(t)), using Corollaries 5 and 6.

Since F (2) is an analytic function, we have in general

F (2) ρ+o n -1 2 ,c s + o n -1 2 , t + o n -1 2 , x, y = F (2) ρ,c (s, t, x, y) + o n -1 2 (19)
and similarly for all partial derivatives of F (2) . Since

an 2n = ρ + o n -1 2 , kn n = s + o n -1 2 and ℓn n = t + o n -1 2
, this shows that the second and third terms on the right in (18) tend to zero as n → ∞. Thus, if z = (s, t, L(s), L(t)), nF

an 2n ,c kn n , ℓn n , x √ n + L(s), y √ n + L(t) → 1 2 ∂ 2 F (2) ρ,c ∂x 2 (z)x 2 + 1 2 ∂ 2 F (2) ρ,c ∂y 2 (z)y 2 + ∂ 2 F (2)
as n → ∞. From (19) it also follows that the sequence F

an 2n ,c ( kn n , ℓn n , •, •) is equicontinuous, hence by Ascoli's theorem that the sequence converges uniformly, hence that the double integral (17) converges to

1 2π H (2) ρ,c (s, t, L(s), L(t)) b a d c exp 1 2 ∂ 2 F (2) ρ,c ∂x 2 (z)x 2 + 1 2 ∂ 2 F (2) ρ,c ∂y 2 (z)y 2 + ∂ 2 F (2) ρ,c

∂x∂y (z)xy dy dx

To find these double derivatives, we use (15) and ( 16). Hence we need

∂ 2 F (1) ρ,c ∂x 2 (s, L(s)) = c(1 -e -2c )e 2cs (1 -e -c(2-2ρ) -e -2cs + e -c(2s-2ρ) ) 2 (1 -e 2cρ )(1 -e -2cs )(1 -e -c(2-2ρ) )(1 -e -c(2-2s) )
and exploit the fact that

∂ 2 F (2) ρ,c ∂x 2 (s, t, L(s), L(t)) = 1 t ∂ 2 F (1) ρ ′′ ,c ′′ ∂x 2 (s ′′ , L ρ ′′ ,c ′′ (s ′′ )) ∂ 2 F (2) ρ,c ∂y 2 (s, t, L(s), L(t)) = 1 1 -s ∂ 2 F (1) ρ ′ ,c ′ ∂x 2 (t ′ , L ρ ′ ,c ′ (t ′ ))
This way, we find the double derivatives of F

(2) ρ,c (s, t, •, •), evaluated at the critical point (x, y) = (L(s), L(t)), to be

∂ 2 F (2) ρ,c ∂x 2 (s, t, L(s), L(t)) = c(1 -e 2ct )(1 -e 2cs + e -c(2-2ρ-2s) -e 2cρ ) 2 (1 -e 2cs )(e 2ct -e 2cs )(1 -e 2cρ )(1 -e -c(2-2ρ) ) = -c sinh ct sinh cs sinh c(t -s)f (s) 2 , ∂ 2 F (2) ρ,c ∂y 2 (s, t, L(s), L(t)) = c(e 2cs -e c )(1 -e 2ct + e -c(2-2ρ-2t) -e 2cρ ) 2 (e 2cs -e 2ct )(e c -e 2ct )(1 -e 2cρ )(1 -e -c(2-2ρ) ) = -c sinh c(1 -s) sinh c(1 -t) sinh c(t -s)f (t) 2 ∂ 2 F (2) ρ,c ∂x∂y (s, t, L(s), L(t)) = c(1 -e 2cs + e -c(2-2ρ-2s) -e 2cρ )(1 -e 2ct + e -c(2-2ρ-2t) -e 2cρ ) (e 2cs -e 2ct )(1 -e 2cρ )(1 -e -c(2-2ρ) ) = c sinh c(t -s)f (s)f (t)
The last one, the mixed derivative, was calculated from scratch, i.e. by calculating

∂ 2 F (2) ρ,c
∂x∂y (s, t, x, y) = c(e 2cte 2cs ) (e c(t+y)e c(s+x) )(e c(t-y)e c(s-x) ) and evaluating at the critical point (x, y) = (L(s), L(t)). Thus, the Hessian of

F (2) ρ,c (s, t, •, •) at the critical point is H(F (2) ρ,c ) = -c sinh cs sinh c(1 -t) sinh c(t -s) × 1 f (s) 0 0 1 f (t) sinh ct sinh c(1 -t) -sinh cs sinh c(1 -t) -sinh cs sinh c(1 -t) sinh cs sinh c(1 -s) 1 f (s) 0 0 1 f (t)
.

The covariance matrix of the limiting Gaussian distribution is the negative of the inverse of H(F (2) ), which we compute to be

Σ = 1 c sinh c f (s) 0 0 f (t) • sinh cs sinh c(1 -s) sinh cs sinh c(1 -t) sinh cs sinh c(1 -t) sinh ct sinh c(1 -t) • f (s) 0 0 f (t) (21)
The matrix in the middle together with the factor 1 c sinh c is the covariance matrix of the Ornstein-Uhlenbeck bridge on the interval [0, 1] with parameter c (see Appendix). Further computations reveal that

1 det(Σ) = H (2) ρ,c (s, t, L(s), L(t)) since both sides equal c f (s)f (t) sinh c sinh cs sinh c(t -s) sinh c(1 -t) . A B C D s 0 g - ε (s 0 ) g + ε (s 0 )
Figure 4: The ε-parallelogram of a rectangle.

Lemma 7. Let δ ∈ (0, 1) and A > 0. Take ε as in Fact 1. Then there exists a positive constant C such that for all ρ ∈ [δ, 1δ], and all c ∈ [-A, A],

∀(s, x) ∈ B ε ρ , F (1) ρ,c (s, x) ≤ -Cs(1 -s),
where B ρ is the (macroscopic) box with perimeter 1 and aspect ratio ρ.

Proof. For a fixed s, the function x → F

ρ,c (s, x) is concave, and its maximum is reached at x = L ρ,c (s), and the point (s, L ρ,c (s)) is inside the ε-interior of the box. Therefore

F (1) ρ,c (s, x) ≤ max{F (1) ρ,c (s, g + ε (s)), F (1)
ρ,c (s, g - ε (s))}. In a neighborhood of 0, g + ε (s) = (1ε)s, and s → F

ρ,c (s, (1ε)s) ≍ s because it is differentiable and vanishes at s = 0 with a non zero derivative. The same argument applied to a neighborhood of s = 1, and for g - ε gives the result. The uniform bound on C results from a compactness argument.

Lemma 8 (L.U.F.). Let δ ∈ (0, 1) and A > 0. Take ε as in Fact 1. For this ε, there exist two constants M and κ such that if a + b = 2n and a a+b ∈ (δ, 1δ) and for all c ∈ [-A, A], then

∀ 0 ≤ k ≤ n, P e -c/n a,b (X 2k = 2 nL ρ,c ( k n ) + y √ n ) ≤ M √ n exp -y 2 2κs(1-s)
2πκs(1s) .

as soon as (2k, nL ρ,c

( k n ) + y √ n ) is in the ε-interior of the box.
Proof. Use Corollary 2. Notice that h c (x, y) ≍ x+y xy uniformly in c. Therefore

H ρ,c (s, x) ≍ ρ(1 -ρ) s (s -x)(s + x) 1 -s (2ρ -s + x)(2 -2ρ -s -x)
.

The first term is bounded. The second term is of order 1 √ s in a neighborhood of (0, 0) in the ε-interior of the box. A simple change of variable (s, x) → (1s, 2ρ -1x) exchanges the second and the third term, which shows that the third term is of order

1 √ 1-s in a neigborhood of (1, 1 -ρ) in the ε-interior of the box. Therefore, H ρ,c (s, x) ≍ 1 s(1 -s) ,
as long as (s, x) is in the ε-interior of the box.

In order to bound the exponential term in Corollary 2, we bound from below the absolute value of

∂ 2 F (1) ρ,c (s,x) ∂x 2
. Using that S ′′ c (u) = c e cu -1 , we get

∂ 2 F (1) ρ,c (s, x) ∂x 2 = -2c - c e c(s+x) -1 - c e c(s-x) -1 - c e c(2ρ-s+x) -1 - c e c(2-2ρ-s-x) -1 .
For s close to 0, the main contribution comes from the second and the third term, that are both negative, and of order s -1 . For s close to 1, the main contribution comes from the fourth and the fifth term, that are both negative, and of order (1s) -1 .

Lemma 9. Let δ ∈ (0, 1) and A > 0, and ε ∈ (0, 1) given by Fact 1. Then there exists a constant C > 0 such that: for all n ≥ 1, for all c ∈ [-A, A], for all sequences of boxes

(B n ) with sides a n , b n such that a n + b n = 2n, an 2n ∈ (δ, 1 -δ), ∀λ > 0, ∀s ∈ [0, 1], P( Xs ∈ [-λ, λ]; (X 2ns , 2ns) ∈ B ε n ) ≤ C (s(1 -s)) 2 λ 4
Proof. Lemma 8 gives the following bound:

P( Xs ∈ [-λ, λ]; (X 2ns , 2ns) ∈ B ε n ) ≤ |y|≥λ M √ n exp (-y 2 2κs(1-s) ) 2πκs(1 -s) ,
where the index of the sum y runs on the set 1

√ n Z - √ nL ρ,c (s).
Comparing the sum on the right-hand side and the corresponding integral, we get

P( Xs ∈ [-λ, λ]; (X 2ns , 2ns) ∈ B ε n ) ≤ M ′ R\(-λ,λ) exp (-y 2 2κs(1-s) ) 2πκs(1 -s) ds.
This last integral is equal to

P(|N | ≥ λ √ κs(1-s)
), where N is a standard Gaussian variable.

Conclude by using Markov inequality for the fourth moment.

We can now verify Billingsley's condition for tightness [START_REF] Billingsley | Convergence of probability measures[END_REF].

Lemma 10. Let ((a n , b n )) n be a sequence such that a n + b n = 2n, and an 2n = ρ + o( 1 √ n ).
Then there exists a constant C such that for all n > 0, for all 0 ≤ r ≤ s ≤ t ≤ 1, for all λ > 0,

P e -c/n an,bn (| Xs -Xr | ≥ λ; | Xt -Xs | ≥ λ) ≤ C(t -r) 2 λ 4 . ( 22
)
Proof. Inequality ( 22) is automatically satisfied, as soon as |r -s| or |t -s| is less than 1/n, or λ is greater than √ n(|r -s| ∧ |t -s|). We suppose now that none of these conditions are satisfied. Introduce now the three following ε-interiors: B ε , the ε-interior of the box

B := [0, a n ] × [0, b n ], B ε L the ε-interior of the box B L := [0, ns -1 2 X 2sn ] × [0, sn + 1 2 X 2sn ] and B ε R the ε-interior of the box B R := [ns -1 2 X 2sn , a n ] × [sn + 1 2 X 2sn , b n ]. See Fig. 5. ns X ns B η(ε) L B η(ε) R Figure 5: The ε-interiors B ε (grey), B η(ε) L (blue) and B η(ε) R (red).
Let us first control the probability that the random interface exits the ε-interior of the box. We use unimodality of the distribution of X 2ns , Corollary 2, and Lemma 7 to get:

P((2ns, X 2ns ) ∈ B ε ) ≤ 4εns(1 -s)P(X 2ns = ng ε ± (s)) ≤ M ns(1 -s) √ n H ρ,c (s, g ε ± (s)) exp(-nCs(1 -s)).
Using the uniform bound of

H ρ,c ≍ 1 √ s(1-s)
obtained in the beginning of the proof of Lemma 8, we get

P((2ns, X 2ns ) ∈ B ε ) ≤ M ns(1 -s) exp(-nC(s(1 -s))) ≤ M (nCs(1 -s)) 2 .
Using our assumption that n ≥ λ 2 |s-r| 2 ∧|t-s| 2 , we obtain:

P((2ns, X 2ns ) ∈ B ε ) ≤ M |r -s| 4 ∧ |t -s| 4 λ 4 s 2 (1 -s) 2 ≤ M |r -s| 2 ∧ |t -s| 2 λ 4 ≤ M |t -r| 2 λ 4 . (23) Remark 4. The weaker condition n ≥ λ 2 s 2 ∧(1-s) 2 yields P((2ns, X 2ns ) ∈ B ε ) ≤ M s 2 ∧ (1 -s) 2 λ 4 . (24)
Now we consider three cases depending on the values of r and t.

Case 1: r < λ 4 √ n . Then | Xr | ≤ 2r √ n < λ 2 . Thus, P(| Xs -Xr | ≥ λ; | Xt -Xs | ≥ λ) ≤ P e -c/n an,bn (| Xs -Xr | ≥ λ) ≤ P(| Xs | ≥ λ 2 ) ≤ P(| Xs | ≥ λ 2 , (ns, X ns ) ∈ B ε ) + P((ns, X ns ) ∈ B ε ).
But we just saw in (23) that

P((2ns, X 2ns ) ∈ B ε ) ≤ M |t -r| 2 λ 4 .
Moreover, to Lemma 9,

P(| Xs | ≥ λ 2 , (2ns, X 2ns ) ∈ B ε ) ≤ C s 2 (1 -s) 2 λ 4 ≤ Cs 2 λ 4 .
Since r < λ 4 √ n < s 4 , we get that s < 4 3 (sr), and thus

P(| Xs | ≥ λ 2 , (2ns, X 2ns ) ∈ B ε ) ≤ C(s -r) 2 λ 4 ≤ C(t -r) 2 λ 4 .
And the tightness condition ( 22) is satisfied in that case.

Case 2: 1t < λ

4

√ n . This case is treated similarly as the previous one, and ( 22) is again satisfied.

Case 3: r ≥ λ

4

√ n and 1t ≥ λ

4

√ n . This is the generic situation. Conditional on {(2ns, X 2ns ) ∈ B ε }, the boxes B L and B R have an aspect ratio in the interval (ε, 1ε). Applying Lemma 9 with δ = ε in the subboxes B L and B R , we get the existence of an η(ε) such that with sufficient probability (2nr, X 2nr ) (resp. (2nt,

X 2nt )) is in B η(ε) L (resp. B η(ε) R ).
More precisely, the assumptions on r and t allow to repeat the argument to derive (24) with the proper scaling in the subboxes B L and B R and get for some constant M :

P((2nr, X 2nr ) ∈ B η(ε) L |(2ns, X 2ns ) ∈ B ε ) ≤ M ( r s ) 2 ∧ ( s-r s ) 2 (λ/ √ s) 4 ≤ M (s -r) 2 λ 4 , (25)
P((2nt, X 2nt ) ∈ B η(ε) R |(2ns, X 2ns ) ∈ B ε ) ≤ M ( t-s 1-s ) 2 ∧ ( 1-t 1-s ) 2 (λ/ √ 1 -s) 4 ≤ M (t -s) 2 λ 4 . ( 26
)
By the Markov property of X, the variables X r and X t are independent conditional on the value of X s . Hence we can write √ n).

P(| Xs -Xr | ≥ λ; | Xt -Xs | ≥ λ) = = j,
Let E ε r,s,t = {(2ns, X 2ns ) ∈ B ε }∩{(2nr, X 2nr ) ∈ B η(ε) L }∩{(2nt, X 2nt ) ∈ B η(ε) R }.
From Inequalities (23), ( 25), ( 26), there exists a constant M such that:

1 -P(E ε r,s,t ) ≤ M (t -r) 2 λ 4 .
As a consequence, it suffices to bound the sum in ( 27) for values of (j, k, l) such that

{ Xr = j, Xs = k, Xt = l} ⊂ {| Xs -Xr | ≥ λ} ∩ {| Xt -Xs | ≥ λ} ∩ E ε r,s,t .
By a scaling argument, when looking at what happens on the left of s,

P e -c/n an,bn ( Xr = j|X ns ) = P e -c/n 2ns-X 2ns 2 , 2ns+X 2ns 2 X r s = j √ s + n s L ρ,c (r) - √ nsL 2ns-X 2ns 4ns ,sc ( r s
) .

Note that on the right hand side, Xr/s is defined with respect to the limit shape inside the box [0, 2ns-X 2ns 4ns

] × [0, 2ns+X 2ns

4ns

]. See Fig. 6. When summing over the values of j such that |j -Xs | ≥ λ, one can get an upper bound for P an,bn (| Xr -Xs | ≥ λ| Xs ) using Lemma 9 as long as the limit shape in the original box B and the one in B L are sufficiently close to each other. We express this proximity by the following sticking condition for B L :

1 √ s [ Xs -λ, Xs + λ] + n s L ρ,c (r) - √ nsL 2ns-X 2ns 2ns ,sc ( r s ) ⊃ - λ √ 2s , λ √ 2s . (28)
Indeed, when the sticking condition (28) is satisfied, P e -c/n an,bn (| Xr -Xs | ≥ λ| Xs ) takes the form:

P e -c/n 2ns-X 2ns 2 , 2ns+X 2ns 2 X r s ∈ 1 √ 2s [ Xs -λ, Xs + λ] + n s L ρ,c (r) - √ nsL 2ns-X 2ns 4ns
,sc ( r s )

≤ P e -c/n 2ns-X 2ns 2 , 2ns+X 2ns 2 X r s ∈ - λ √ 2s , λ √ 2s ≤ M ( r s s-r s ) 2 (λ/ √ 2s) 4 ≤ 4M (t -r) 2 λ 4
by Lemma 9 in the small box. By the mean value theorem applied to ρ → L ρ,sc (r/s), there exists a ρ in the interval

( s-Lρ,c(s) 2s , s-Lρ,c(s)-Xsn -1/2 2s ) such that s √ n L s-Lρ,c(s) 2s ,sc r s -L s-Lρ,c(s)-Xsn -1/2 2s ,sc r s = 1 2 Xs ∂L ρ,cs (r/s) ∂ρ ρ=ρ .
But differentiating (6) with respect to ρ, for generic values of ρ, c, t yields

0 ≤ 1 + 1 2 ∂L ρ,c (t) ∂ρ = e c(2ρ-1) sinh(c(1 -t)) sinh(ct) + e c(2ρ-1) sinh(c(1 -t)) ≤ 1 ∧ e c(2ρ-1) sinh(c(1 -t)) sinh(ct) ≤ 1 ∧ K 1 -t t
with a constant K which works for all c ∈ [-A, A], all ρ ∈ (0, 1), and all t ∈ [0, 1]. As a consequence, the sticking condition (28) is satisfied as soon as

| Xs |(1 ∧ K s -r r ) ≤ λ 2 .
Therefore the probability that (28) is not verified is less than

P (| Xs |(1 ∧ K s -r r ) > λ 2 , (2sn, X 2sn ) ∈ B ε ), which, by Lemma 9 is bounded by λ -4 s 2 (1 -s) 2 (1 ∧ K 4 ( s-r r ) 4
). This bound is less than some constant times λ -4 (sr) 2 ≤ λ -4 (tr) 2 , as one can check in both regimes s-r r ≶ K. This finishes the proof of Lemma 10.

The unbounded case

In this section, we relax the constraint of remaining in a box. For any q ∈ (0, 1), we define a probability measure P q on all partitions, by P q (λ) = 1 Z(q) q |λ| , for all partitions λ, where Z(q) = ∞ i=1 (1q i ) -1 is the generating function of all partitions. When q goes to 1, a law of large number for the shape of the random partition occurs. Namely, when rescaled by (1q), the boundary of λ converges to a deterministic curve.

Theorem 3 (Vershik). Let L ∞ (s) = log(2 cosh s), and (X s ) s∈R be the piecewise linear function describing the boundary of λ in the (s, x)-coordinates. Then, for all ε > 0,

lim q→1 P q (sup s∈R |(1 -q)X s(1-q) -1 -L ∞ (s)| > ε) = 0.

Link between L ∞ and L ρ,c

Petrov recently discussed in [START_REF] Petrov | Limit shapes of young diagrams. two elementary approaches[END_REF] the link between the limit shapes obtained above and the infinite shape of the unbounded problem. We explain it in this short section.

The family of limit shapes obtained above has the following property: if we fix s 0 ∈ (0, 1) and take the point (s 0 , L ρ,c (s 0 )) as the right corner of a new family of bounding boxes, then the limit shape of this new problem is simply the restriction of L ρ,c to the interval [0, s 0 ], rescaled by s 0 (See Lemma 4). We are therefore led to the idea that there exists some 'inverse limit' of these curves, from which they can all be recovered by restriction. This universal curve is precisely L ∞ from Theorem 3. Indeed, for every ρ ∈ (0, 1) and c ∈ R * + , one can find two real numbers s 0 < s 1 such that the restriction of L ∞ to [s 0 , s 1 ] is (up to an affine transformation) L ρ,c (see Figure 8): ∀t

∈ [0, 1], L ρ,c (t) = L ∞ (s 0 + t(s 1 -s 0 )) -L ∞ (s 0 ) 1 -s 0 . (29)
For given values of s 0 < s 1 , the parameter ρ is given by

ρ = 1 2 - L ∞ (s 1 ) -L ∞ (s 0 ) 2(s 1 -s 0 ) ,
and c is found for example by comparing the slopes of the limit shape in the corner and identifying ∂L∞ ∂s s=s 0 (s 1s 0 ) and ∂Lρ,c ∂s s=0

= e c(2ρ-1) -cosh c sinh c

, which is an increasing function of c.

The limit case c → 0 correspond to s 1 → s 0 . For negative values of c one can take advantage of the symmetry of the model.

Fluctuations around L ∞

As in the boxed case, we can compute explicitely the two point distribution. The analogue of Proposition 2 in the unbounded case is Proposition 5. The 2-dimensional marginal of (X 2k ) k∈Z is given by P q (X 2k = 2i, X 2l = 2j) = Z k+i,∞ (q)Z j+l-i-k,l+i-k-j (q)Z j-l,∞ (q) Z(q) q (i-k)(k+i)+(l+j-k-i)(j-l) .

where Z a,∞ (q) = lim b→∞ Z a,b (q) = a j=1 (1q j ) -1 . The asymptotic analysis of this formula as q → 1 relies on q-Stirling's Formula (Corollary 1) and goes as in Section 5.1.

One can then deduce the following theorem:

Theorem 4. As q goes to 1 -, the random function When conditioning the process Y s / √ 2 cosh s to have zero integral over R, one obtains a new centered Gaussian process whose covariance is equal to e -|t-s| 2 cosh s cosh t -π 2 6 h(s)h(t),

s → √ 2 cosh s √ 1 -q (1 -q)X s(1-q) -1 -L ∞ (s)
where h(t) = 6 π 2 t tanh tlog(2 cosh t). Expressed in the original coordinates (by rotating back the picture by 45 degrees), this covariance is the one obtained by Pittel [START_REF] Pittel | On a likely shape of the random Ferrers diagram[END_REF] who deals with the microcanonical ensemble of partitions (with a fixed area |λ| = n → ∞). The fluctuations in that case are thus given by this Gaussian, but non-Markov process.

Appendix: The Ornstein-Uhlenbeck Bridge

In this short appendix, we give a description of the Ornstein-Uhlenbeck bridge. Let (B t ) be a Brownian motion on [0, ∞). The Ornstein-Uhlenbeck process (Z t ) is the Gaussian random process defined by (30)

Using the right-hand side of (30), we can also represent Z t as:

Z t =
Observe that, for 0 ≤ t ≤ 1,

E [Y t Z 1 ] = E Z t Z 1 -h(t)Z 2 1 = e -c sinh(ct) c - sinh(ct) sinh c e -c sinh c c = 0,
showing that Y t is independent of Z 1 . Now let f be any test function and consider

E [f (X t ) | X 1 = 0] = E [f (Y t + h(t)X 1 ) | X 1 = 0] = E [f (Y t ) | X 1 = 0] = E[f (Y t )],
where the last equality is a consequence of the independence just proved. Thus, the process (Y t ) has the probability law of the Ornstein-Uhlenbeck process conditional on X 1 = 0. Its covariance is Notice that this covariance is also the Green function of a Brownian motion on [0, 1] killed an exponential rate c 2 [START_REF] Chan | On polymer conformations in elongational flows[END_REF]. Let us mention also, that this process was used by C. Donati in her solution of Buffon-Synge problem concerning the typical distance between the extremities of a string of given length, thrown at random [START_REF] Donati-Martin | Le problème de Buffon-Synge pour une corde[END_REF].

E[Y s Y t ] = E [Y s (Z t -h(t)Z 1 )] = E[Y s Z t ] = E[Z s Z t ] -h(s)E[Z 1 Z t ] = sinh (

Figure 1 :

 1 Figure 1: The Young diagram for the partition π = (8, 7, 7, 5, 4, 3, 3, 2, 2, 1) of n = 42.

Figure 2 :

 2 Figure 2: A partition as a lattice path.

  k,l |j-k|≥λ |l-k|≥λ P( Xs = k)P( Xr = j| Xs = k)P( Xt = l| Xs = k), (27) where j, k, l runs through all possible values for Xr , Xs , Xt (two successive values of j, k, l differ by 1/

Figure 6 :

 6 Figure 6: Limit shapes in the box B and the subboxes B L and B R . The double arrow represents the difference of the limit shapes in B and B L at time r, whose expression, 2nL ρ,c (r) -2nsL 2ns-X 2ns 4ns ,sc

Figure 7 :

 7 Figure 7: The limit shape in a subbox.

s 0 s 1 Figure 8 :

 18 Figure 8: Limit shapes are restrictions of the universal curve.

  converges weakly in D to the two-sided stationary Ornstein-Uhlenbeck process (Y s ) s∈R , which is the Gaussian process on R with covariance∀(s, t) ∈ R 2 , E[Y s Y t ] = e -|t-s| .

Z t = B t -c t 0 Z 0 e

 00 s ds = e -ct t cs dB s .

  cs) c sinh c) e -ct sinh ce -c sinh(ct) = sinh(cs) sinh(c(1t)) c sinh c , cf. the middle matrix in (21).

  e -ct B t 0 e cs ds = e -ct B e 2ct -1 Tying down' this process at t = 1 gives the Ornstein-Uhlenbeck (Y t ) t∈[0,1] bridge of length 1. Set

			.		(31)
			2c
	From (31) we derive the covariance. For 0 ≤ s ≤ t,
	E [Z s Z t ] = e -c(s+t) e 2cs -1 2c	= e -ct sinh(cs) c	.	(32)
	'h(t) =	sinh(ct) sinh c		(33)
	and define
	Y

t = Z th(t)Z 1 .

Acknowledgements. The authors are pleased to thank Nikolai Reshetikhin for having proposed the question solved in this paper, Marc Yor for having recognized the covariance we obtained as that of the Ornstein-Uhlenbeck bridge, and Richard Kenyon for showing us the way to get the limit shapes we found from the infinite limit shape. The work of the third author was partially supported by the ANR MEMEMO grant.

This completes the proof.

We explain now that for m points 0 < t 1 < t 2 < • • • < t m < 1, the limit of the corresponding m-dimensional marginal X(n)

tm is Gaussian, with covariance matrix Σ = (σ i,j ) defined by

for all i ≤ j.

A first approach would be to repeat similar computations as for the 2-dimensional case, and observe that the matrix σ is the inverse of the m-dimensional Hessian of

Another way takes advantage of the Markov property satisfied by the process ( Xn ). Denote by P (n) s,t (x, y) the transition kernel of this Markov chain. The law of the m-tuple

tm is given by

Proposition 4 implies the convergence of each √ nP

(n) t i ,t i+1 (x n , y n ) to the density of the kernel of the Ornstein-Uhlenbeck bridge as n → ∞, and (x n ), (y n ) converge. The convergence in law follows from Lebesgue dominated convergence theorem. Domination is ensured by Lemma 8 below and unimodality stated in Lemma 2.

Tightness and proof of Theorem 2

We already proved the convergence of the finite-dimensional distributions. We need to show that the sequence of the distribution is tight. For this we use Theorem 13.5 p.142 of [START_REF] Billingsley | Convergence of probability measures[END_REF]. The criterion is checked in Lemma 10.

Before entering the proof, we need a geometric definition.

Definition 2. Let ABCD be a rectangle with sides having slopes ±1 in the (s, x) coordinates. For ε ∈ (0, 1), the ε-parallelogram of the rectangle ABCD is the unique parallelogram with diagonal AC and sides with slopes 1ε and -1 + ε in the (s, x) coordinates. The ε-interior of the rectangle is the interior of the ε-parallelogram. The complement of the ε-interior in the rectangle is called the ε-boundary of the rectangle. See Figure 5.2. For s 0 ∈ (0, 1), the sides of the parallelogram intersect the straight line s = s 0 in two points. We denote the ordinates of these two points by g + ε (s 0 ) and g - ε (s 0 ). We first state the following useful fact. Fact 1. For all δ ∈ (0, 1), for all A > 0, there exists ε, such that for all rectangular boxes of side lengths a, b satisfying ρ = a a+b ∈ (δ, 1δ), and all c ∈ [-A, A], the limit shape L ρ,c is entirely included in the ε-interior of the box.

The following lemma controls the function F