The Arctic Circle Revisited

TL;DR

This paper analyzes limit shapes in the six-vertex model with domain wall boundary conditions using a correlation function called the emptiness formation probability, deriving new results including the Arctic Circle and Ellipses theorems.

Contribution

It provides a closed-form expression for the emptiness formation probability and links it to a matrix model, extending the Arctic shape results to a broader class of vertex weights.

Findings

Derived a determinant and integral expression for the correlation function.

Established the connection to a Triple Penner matrix model at the free-fermion point.

Proved the Arctic Circle and Ellipses theorems via saddle-point analysis.

Abstract

The problem of limit shapes in the six-vertex model with domain wall boundary conditions is addressed by considering a specially tailored bulk correlation function, the emptiness formation probability. A closed expression of this correlation function is given, both in terms of certain determinant and multiple integral, which allows for a systematic treatment of the limit shapes of the model for full range of values of vertex weights. Specifically, we show that for vertex weights corresponding to the free-fermion line on the phase diagram, the emptiness formation probability is related to a one-matrix model with a triple logarithmic singularity, or Triple Penner model. The saddle-point analysis of this model leads to the Arctic Circle Theorem, and its generalization to the Arctic Ellipses, known previously from domino tilings.

Equations

$$Z_N=\sum _{\begin{array}{c}\text{arrow configurations}\\ \text{with DWBC}\end{array}}{w}_1^{n_1}{w}_2^{n_2}\dots w_6^{n_6}.$$

$$w_1=w_2=:a,w_3=w_4=:b,w_5=w_6=:c.$$

$$\Delta =\frac{a^2+b^2-c^2}{2ab}.$$

$$a=\sin (\lambda +\eta ),b=\sin (\lambda -\eta ),c=\sin 2\eta .$$

$$F_N(r,s)=\frac{1}{Z_N}\sum _{\begin{array}{c}\text{‘constrained’}\\ \text{ arrow configurations}\\ \text{with DWBC}\end{array}}{w}_1^{n_1}{w}_2^{n_2}\dots w_6^{n_6}.$$

$$\phi :=\frac{c}{ab}=\frac{\sin 2\eta }{\sin (\lambda +\eta )\sin (\lambda -\eta )},$$

$$\mu (x):={\mathrm{e}}^{x(\lambda -\pi /2)}\frac{\sinh (\eta x)}{\sinh (\pi x/2)},$$

$$\phi ={\int }_{-\infty }^{\infty }\mu (x)dx.$$

$${\int }_{-\infty }^{\infty }{P}_n(x)P_m(x)\mu (x)dx=h_n{\delta }_{nm}.$$

$$K_n(x)=n!{\phi }^{n+1}\frac{1}{h_n}{P}_n(x).$$

$$\omega (ϵ):=\frac{a}{b}\frac{\sin (ϵ)}{\sin (ϵ-2\eta )},\tilde{\omega }(ϵ):=\frac{b}{a}\frac{\sin (ϵ)}{\sin (ϵ+2\eta )}.$$

$$a^2\tilde{\omega }-2\Delta ab\tilde{\omega }\omega +b^2\omega =0,$$

$$F_{N}(r,s)=(-1)^{s}\det_{1\leq j,k\leq s}\left[K_{N-k}(\partial_{\epsilon_{j}})\right]\,\prod_{j=1}^{s}\frac{\left[\omega(\epsilon_{j})\right]^{N-r}}{\left[\omega(\epsilon_{j})-1\right]^{N}}\\ \times\prod_{1\leq j<k\leq s}\frac{\left[{\tilde{\omega}}(\epsilon_{j})-1\right]\left[\omega(\epsilon_{k})-1\right]}{{\tilde{\omega}}(\epsilon_{j})\omega(\epsilon_{k})-1}\bigg{|}_{\epsilon_{1}=0,\dots,\epsilon_{s}=0}\,.$$

$$H_N(r):=F_N(r,1)-F_N(r-1,1).$$

$$H_{N}(r)=K_{N-1}(\partial_{\epsilon})\frac{[\omega(\epsilon)]^{N-r}}{[\omega(\epsilon)-1]^{N-1}}\bigg{|}_{\epsilon=0}\,.$$

$$h_N(z):=\sum _{r=1}^N{H}_N(r)z^{r-1}.$$

$$K_{N-1}(\partial_{\epsilon})f(\omega(\epsilon))\bigg{|}_{\epsilon=0}=\frac{1}{2\pi i}\oint_{C_{0}}\frac{(z-1)^{N-1}}{z^{N}}\;h_{N}(z)f(z)\,dz\,.$$

$$\begin{gathered} F_N(r,s)={\left(-\frac{1}{2\pi i}\right)}^s{\oint }_{C_0}\cdots {\oint }_{C_0}{d}^s\omega \underset{1\le j,k\le s}{\det }\left[h_{N-k+1}({\omega }_j){\left(\frac{{\omega }_j-1}{{\omega }_j}\right)}^{N-k}\right] \\ \times \prod _{j=1}^s\frac{{\omega }_j^{N-r-1}}{({\omega }_j-1{)}^N}\prod _{1\le j<k\le s}\frac{({\tilde{\omega }}_j-1)({\omega }_k-1)}{{\tilde{\omega }}_j{\omega }_k-1}. \end{gathered}$$

$$\tau ={\tan }^2(\lambda -\pi /4),0<\tau <\infty .$$

$$\tilde{\omega }=-\tau \omega .$$

$$h_N(z)={\left(\frac{1+\tau z}{1+\tau }\right)}^{N-1}.$$

$$\begin{gathered} F_N(r,s)={\left(-\frac{1}{2\pi i}\right)}^s{\oint }_{C_0}\cdots {\oint }_{C_0}{d}^s\omega \underset{1\le j,k\le s}{\det }\left\{{\left[\frac{(1+\tau {\omega }_j)({\omega }_j-1)}{(1+\tau ){\omega }_j}\right]}^{N-k}\right\} \\ \times \prod _{j=1}^s\frac{{\omega }_j^{N-r-1}}{({\omega }_j-1{)}^N}\prod _{1\le j<k\le s}\frac{(1+\tau {\omega }_j)({\omega }_k-1)}{1+\tau {\omega }_j{\omega }_k}. \end{gathered}$$

$$\prod _{j=1}^s{\left[\frac{(1+\tau {\omega }_j)({\omega }_j-1)}{(1+\tau ){\omega }_j}\right]}^{N-s}$$

$$\prod _{1\le j<k\le s}\left[\frac{(1+\tau {\omega }_j)({\omega }_j-1)}{(1+\tau ){\omega }_j}-\frac{(1+\tau {\omega }_k)({\omega }_k-1)}{(1+\tau ){\omega }_k}\right]\frac{(1+\tau {\omega }_j)({\omega }_k-1)}{1+\tau {\omega }_j{\omega }_k}.$$

$$\frac{1}{s!}(-1{)}^{s(s-1)/2}\prod _{j=1}^s\frac{1}{{\omega }_j^{s-1}}\prod _{1\le j<k\le s}({\omega }_j-{\omega }_k{)}^2.$$

$$\begin{gathered} F_N(r,s)=\frac{(-1{)}^{s(s+1)/2}}{s!(1+\tau {)}^{s(N-s)}(2\pi i{)}^s} \\ \times {\oint }_{C_0}\cdots {\oint }_{C_0}{d}^s\omega \prod _{1\le j<k\le s}({\omega }_j-{\omega }_k{)}^2\prod _{j=1}^s\frac{(1+\tau {\omega }_j{)}^{N-s}}{({\omega }_j-1{)}^s{\omega }_j^r}. \end{gathered}$$

$$r/N=x,s/N=y,$$

$$F(x,y):=\underset{N\to \infty }{\lim }{F}_N(xN,yN),x,y\in [0,1].$$

$$\begin{gathered} I_N(r,s):=\frac{(-1{)}^{s(s+1)/2}}{s!(1+\tau {)}^{s(N-s)}(2\pi i{)}^s} \\ \times {\oint }_{C_1}\cdots {\oint }_{C_1}{d}^s\omega \prod _{1\le j<k\le s}({\omega }_j-{\omega }_k{)}^2\prod _{j=1}^s\frac{(1+\tau {\omega }_j{)}^{N-s}}{({\omega }_j-1{)}^s{\omega }_j^r}, \end{gathered}$$

$$I_N(r,s)=1$$

Full text

1. Introduction

The Arctic Circle has first appeared in the study of domino tilings of large Aztec diamonds [EKLP, JPS]. The name originates from the fact that in most configurations the dominoes are ‘frozen’ outside the circle inscribed into the diamond, while the interior of the circle is a disordered, or ‘temperate’, zone. Further investigations of the domino tilings of Aztec diamonds, such as details of statistics near the circle, can be found in [CEP, J1, J2]. Here we mention that the Arctic Circle is a particular example of a limit shape in dimer models, in the sense that it describes the shape of a spatial phase separation of order and disorder. Apart from domino tilings, many more examples have been discussed recently, see, among others, papers [CKP, CLP, KO, KOS, OR].

As long as only dimer models are considered, this amounts to restrict to discrete free-fermionic models, although with nontrivial boundary conditions. Indeed, many of them can be viewed as a six-vertex model at its Free Fermion point (the correspondence being however usually not bijective), with suitably chosen fixed boundary conditions. In particular, this is the case of domino tilings of Aztec diamonds [EKLP], and the corresponding boundary conditions of the six-vertex model are the so-called Domain Wall Boundary Conditions (DWBC). Hence the problem of limit shapes extends to the six-vertex model with generic weights, and with fixed boundary conditions, among which the case of DWBC is the most interesting.

Historically, the six-vertex model with DWBC was first considered in paper [K] within the framework of Quantum Inverse Scattering Method [KBI] to prove the Gaudin hypothesis for norms of Bethe states. The model was subsequently solved in paper [I] where a determinant formula for the partition function was given; see also [ICK] for a detailed exposition. Quite independently, the model was later found, under certain restrictions on the vertex weights, to be deeply related with enumerations of alternating sign matrices (see, e.g., [Br] for a review) and, as already mentioned, to domino tilings of Aztec diamonds [EKLP].

Concerning the problem of limit shapes for the six-vertex model with DWBC, as far as the Free Fermion point is considered, the relation with domino tilings provided apparently an indirect proof of the corresponding Arctic Circle. The non-bijective nature of the correspondence between the two models asked for more direct results, purposely for the free-fermion six-vertex model, see [Zi1, FS, KP]. Out of the Free Fermion point, however, only very few analytical results are available, such as exact expressions for boundary one-point [BPZ] and two-point [FP, CP1] correlation functions. The present knowledge on the subject is based mainly on numerics [E, SZ, AR]; some steps towards finding the limit shapes of the model have been done recently in [PR].

In the present note we propose a rather direct strategy to address the problem: after briefly reviewing the six-vertex model with DWBC, we define a bulk correlation function, the Emptiness Formation Probability (EFP), which discriminates the ordered and disordered phase regions. We give for this correlation function two equivalent representations, in terms of a determinant and of a multiple integral. The core derivation of EFP is heavily based on the Quantum Inverse Scattering Method [KBI], along the lines of papers [BPZ, CP1]; it is out of the scope of the present paper, corresponding details being given in a separate publication [CP4]. Here our aim is to demonstrate how the limit shapes for the considered model can be extracted from EFP in a suitable scaling limit, by making use of ideas and techniques of Random Matrix Models.

To be more specific, and to establish a contact with previous results, we specialize here our further discussion to the case of free-fermion six-vertex model. We show that the asymptotic analysis of multiple integral formula for EFP in the scaling limit reduces to a saddle-point problem for a one-matrix model with a triple logarithmic singularity, or triple Penner model. We argue that the limit shape corresponds to condensation of all saddle-point solutions to a single point. This allows us to recover the known Arctic Circle and Ellipses.

As a comment to our approach, it is to be stressed that it is directly tailored on the six-vertex model, rather than domino tilings. For this reason it is not restricted to the free-fermion models, even if, of course, further significant efforts might be necessary, essentially from the point of view of Random Matrix Model reformulation, for application to more general situations. On the basis of our previous results in [CP2], however, the application of the method to the particular case of the so-called Ice Point of the model should be straightforward. This would provide the limit shape of alternating sign matrices.

2. The model

2.1. The six-vertex model

The six-vertex model (for reviews, see [LW, Ba]) is formulated on a square lattice with arrows lying on edges, and obeying the so-called ‘ice-rule’, namely, the only admitted configurations are such that there are always two arrows pointing away from, and two arrows pointing into, each lattice vertex. An equivalent and graphically simpler description of the configurations of the model can be given in terms of lines flowing through the vertices: for each arrow pointing downward or to the left, draw a thick line on the corresponding edge. This line picture implements the ‘ice-rule’ in an automated way. The six possible vertex states and the Boltzmann weights $w_{1},w_{2},\dots,w_{6}$ assigned to each vertex according to its state are shown in Figure 1.

2.2. Domain Wall Boundary Conditions

The Domain Wall Boundary Conditions (DWBC) are imposed on the $N\times N$ square lattice by fixing the direction of all arrows on the boundaries in a specific way. Namely, the vertical arrows on the top and bottom of the lattice point inward, while the horizontal arrows on the left and right sides point outward. Equivalently, a generic configuration of the model with DWBC can be depicted by $N$ lines flowing from the upper boundary to the left one. A possible state of the model both in terms of arrows and of lines is shown in Figure 2.

2.3. Partition function

The partition function is defined, as usual, as a sum over all possible arrow configurations, compatible with the imposed DWBC, each configuration being assigned its Boltzmann weight, given as the product of all the corresponding vertex weights,

[TABLE]

Here $n_{1},n_{2},\dots,n_{6}$ denote the numbers of vertices with weights $w_{1},w_{2},\dots,w_{6}$, respectively, in each arrow configuration ($n_{1}+n_{2}+\dots+n_{6}=N^{2}$).

2.4. Anisotropy parameter and phases of the model

The six-vertex model with DWBC can be considered, with no loss of generality, with its weights invariant under the simultaneous reversal of all arrows,

[TABLE]

Under different choices of Boltzmann weights the six-vertex model exhibits different behaviours, according to the value of the parameter $\Delta$, defined as

[TABLE]

It is well known that there are three physical regions or phases for the six-vertex model: the ferroelectric phase, $\Delta>1$; the anti-ferroelectric phase, $\Delta<-1$; the disordered phase, $-1<\Delta<1$. Here we restrict ourselves to the disordered phase, where the Boltzmann weights are conveniently parameterized as

[TABLE]

With this choice one has $\Delta=\cos 2\eta$. The parameter $\lambda$ is the so-called spectral parameter and $\eta$ is the crossing parameter. The physical requirement of positive Boltzmann weights, in the disordered regime, restricts the values of the crossing and spectral parameters to $0<\eta<\pi/2$ and $\eta<\lambda<\pi-\eta$.

The special case $\eta=\pi/4$ (or $\Delta=0$) is related to free fermions on a lattice, and there is a well-known correspondence with dimers and domino tilings. In particular, at $\lambda=\pi/2$, the $\Delta=0$ six-vertex model with DWBC is related to the domino tilings of Aztec diamond. For arbitrary $\lambda\in[\pi/4,3\pi/4]$, we shall refer to the $\Delta=0$ case as the Free Fermion line.

The case $\eta=\pi/6$ (i.e. $\Delta=1/2$) and $\lambda=\pi/2$, where all weights are equal, $a=b=c$, is known as the Ice Point; all configurations are given the same weight. In this case there is a one to one correspondence between configurations of the model with DWBC and $N\times N$ alternating sign matrices.

2.5. Phase separation and limit shapes

The six-vertex model exhibits spatial separation of phases for a wide choice of fixed boundary conditions, and, in particular, in the case of DWBC. Roughly speaking, the effect is related to the fact that ordered configurations on the boundary can induce, through the ice-rule, a macroscopic order inside the lattice.

The notion of phase separation acquires a precise meaning in the scaling limit, that is the thermodynamic/continuum limit, performed by sending the number of sites $N$ to infinity and the lattice spacing to zero, while keeping the total size of the lattice fixed, e.g., to $1$. On a finite lattice, several macroscopic regions may appear, which in the scaling limit are expected to be sharply separated by some curves, the so-called Arctic curves.

For the six-vertex model with DWBC the shape of the Artic curve, or limit shape, has been found rigorously only on the Free Fermion line, and for the closely related domino tilings of Aztec diamond [JPS, CEP, Zi1, FS, KP]. For generic values of weights the limit shapes are not known, but the whole picture is strongly supported both numerically [E, SZ, AR] and analytically [KZ, Zi2, BF, PR].

3. Emptiness Formation Probability

3.1. Definition

We shall use the following coordinates on the lattice: $r=1,\dots,N$ labels the vertical lines from right to left; $s=1,\dots,N$ labels the horizontal lines from top to bottom. We may now introduce the correlation function $F_{N}(r,s)$, measuring the probability for the first $s$ horizontal edges between the $r$-th and $r+1$-th line to be all ‘full’ (i.e. thick in the line picture, or with a left arrow in the standard picture of the six-vertex model):

[TABLE]

Here the sum is performed over all arrow configurations on the $N\times N$ lattice, subjected to the restriction of DWBC, and to the condition that all arrows on the first $s$ edges between the $r$-th and $r+1$-th line should point left, see Figure 3.

Although this correlation function may appear rather sophisticated, it is computable in some closed form by means of the Quantum Inverse Scattering Method, on which DWBC are indeed tailored. It is the natural adaptation of the Emptiness Formation Probability of quantum spin chains to the present model. For this reason, and to link to the common practice in the quantum integrable models community, even if $F_{N}(r,s)$ actually describes ‘fullness’ formation probability, we shall call it Emptiness Formation Probability (EFP).

3.2. Qualitative discussion of $F_{N}(r,s)$

Let us restrict ourselves to the disordered regime, $-1<\Delta<1$, for definiteness. From previous analytical and numerical work, in the large $N$ limit the emergence of a limit shape, in the form of a continuous closed curve touching once each of the four sides of the lattice, is expected. It follows that five regions emerge in the lattice: a central region, enclosed by the curve, and four corner regions, lying outside the closed curve and delimited by the sides of the lattice. The central region is disordered, while the four corners are frozen, with mainly vertices of type 1, 3, 2, 4 (see Figure 1) appearing in the top-left, top-right, bottom-right and bottom-left corner, respectively.

By construction, EFP is expected to be almost one in frozen regions of type $1$, or $3$, bordering the top side of the lattice, and to be rather small otherwise. DWBC exclude a region of type $3$ to emerge in the upper part of the lattice. Hence $F_{N}(r,s)$ describes, at a given value of $r$, as $s$ increases, a transition from a frozen region of vertices of type $1$, where $F_{N}(r,s)\sim 1$, to a generic region where $F_{N}(r,s)\sim 0$.

It follows that $F_{N}(r,s)$ can describe only the upper left portion of the closed curve, between its top and left contact points. Nevertheless, it should be mentioned that the full curve can be built from the knowledge of its top left portion, just exploiting the crossing symmetry of the six-vertex model. Hence EFP, $F_{N}(r,s)$, is well suited to describe limit shapes.

3.3. Some notations

For a given choice of parameters $\lambda$, $\eta$ we define

[TABLE]

and the integration measure on the real line

[TABLE]

related to $\varphi$ as follows:

[TABLE]

Let us introduce the complete set of monic orthogonal polynomial $\left\{P_{n}(x)\right\}_{n=0,1,\dots}$ associated to the integration measure $\mu(x)$, with the orthogonality relation

[TABLE]

The square norms $h_{n}$ are completely determined by the measure $\mu(x)$, and may be expressed, in principle, in terms of its moments. In the following we shall be interested in the complete set of orthogonal polynomials $\left\{K_{n}(x)\right\}_{n=0,1,\dots}$ defined as

[TABLE]

We moreover define

[TABLE]

Note that the following relation holds

[TABLE]

allowing to express $\tilde{\omega}$ in terms of $\omega$.

3.4. Determinant representation

For EFP in the six-vertex model with DWBC, the following representation holds:

[TABLE]

This representation has been obtained in the framework of the Quantum Inverse Scattering Method [KBI], along the lines of analogous derivations worked out for one-point and two-point boundary correlation functions of the model [BPZ, CP1]. The details of the derivation can be found in [CP4].

3.5. The boundary correlation function

If we consider expression (3.3) when $s=1$, we recover the boundary polarization, introduced and computed in [BPZ]. It is convenient to consider the closely related boundary correlation function

[TABLE]

As shown in [BPZ, CP1], the following representation holds:

[TABLE]

We define the corresponding generating function

[TABLE]

Noticing that $\omega(\epsilon)\to 0$ as $\epsilon\to 0$, it can be shown that, given any arbitrary function $f(z)$ regular in a neighbourhood of the origin, the following inverse representation holds

[TABLE]

Here $C_{0}$ is a closed counterclockwise contour in the complex plane, enclosing the origin, and no other singularity of the integrand.

3.6. Multiple integral representation

Plugging (3.5) into representation (3.3), we readily obtain the following multiple integral representation for EFP:

[TABLE]

Here $\tilde{\omega}_{j}$’s should be expressed in terms of $\omega_{j}$’s through (3.2). Indeed, due to (3.5), relation (3.2) for functions $\omega(\epsilon)$, $\tilde{\omega}(\epsilon)$, translates directly into the same relation between $\omega_{j}$ and $\tilde{\omega}_{j}$, $j=1,\dots,s$.

Representation (3.6), and all results in this Section hold for any choice of parameters $\lambda$ and $\eta$ within the disordered regime. Moreover, by analytical continuation in parameters $\lambda$ and $\eta$, these results can be easily extended to all other regimes.

The determinant in expression (3.6) is a particular representation of the partition function of the six-vertex model with DWBC, when the homogeneous limit is performed only on a subset of the spectral parameters [CP3]. The structure of the previous multiple integral representation therefore closely recalls analogous ones for the Heisenberg XXZ quantum spin chain correlation functions [JM, KMT].

For generic values of $\lambda$ and $\eta$, the orthogonal polynomials $K_{n}(x)$, or the generating function $h_{N}(z)$, are known only in terms of rather implicit representations. Fortunately, there are three notable exceptions [CP2]: the Free Fermion line ($\eta=\pi/4$, $-\pi/4<\lambda<\pi/4$, $\Delta=0$), the Ice Point ($\eta=\pi/6$, $\lambda=\pi/2$, $\Delta=1/2$), and the Dual Ice Point ($\eta=\pi/3$, $\lambda=\pi/2$, $\Delta=-1/2$). In these three cases, the $K_{n}(x)$ turn out to be classical orthogonal polynomials, namely Meixner-Pollaczek, Continuous Hahn and Continuous Dual Hahn polynomial, respectively. Correspondingly, the generating function can be represented explicitly in terms of terminating hypergeometric functions that may simplify considerably further evaluation of EFP. In the next Section we shall focus on the case of Free Fermion line.

4. Multiple integral representation at $\Delta=0$

4.1. Specialization to $\eta=\pi/4$

We shall now restrict ourselves to the case $\eta=\pi/4$. We have $\Delta=0$, and the six-vertex model reduces to a model of free fermions on the lattice. The parameter $\lambda$ can still assume any value in the interval $(-\pi/4,\pi/4)$. It is convenient to trade $\lambda$ for the new parameter

[TABLE]

The symmetric point (related to the domino tiling of Aztec Diamond) corresponds now to $\tau=1$. For generic values of $\tau$ we have:

[TABLE]

The generating function (3.4) is known explicitely (see [CP2] for details):

[TABLE]

Plugging this expression into (3.6), we get

[TABLE]

4.2. Symmetrization

After extracting a common factor

[TABLE]

from the determinant in (4.1), we recognize it to be of Vandermonde type. We can therefore collect from the integrand of (4.1) the double product

[TABLE]

Noticing that the integration and the remaining of integrand are fully symmetric under permutation of variables $\omega_{1},\dots,\omega_{j}$, we can perform total symmetrization of the previous double product over all its variables, with the result

[TABLE]

Hence, we finally obtain the following representation for EFP on the Free Fermion line:

[TABLE]

The appearance of a squared Vandermonde determinant in this expression naturally recalls the partition functions of $s\times s$ Random Matrix Models.

5. Triple Penner model and Arctic Ellipses

5.1. Scaling limit

We shall now address the asymptotic behaviour of expression (4.2) for EFP in the $\Delta=0$ case. We are interested in the limit $N,r,s\to\infty$, while keeping the ratios

[TABLE]

fixed. In this limit, $x,y\in[0,1]$ will parameterize the unit square to which the lattice is rescaled. Correspondingly EFP is expected to approach a limit function

[TABLE]

We shall exploit the standard approach developed for instance in the investigation of asymptotic behaviour for Random Matrix Models. Before this let us however point out some facts which holds already for any finite value of $s$.

5.2. A useful identity

Consider the quantity

[TABLE]

which differs from (4.2) only in the integration contours. Here $C_{1}$ is a closed, clockwise oriented contour (note the change in orientation) in the complex plane enclosing point $\omega=1$, and no other singularity of the integrand. We have the identity

[TABLE]

for any integer $r,s=1,\dots,N$. The simplest way to prove the previous identity is by shifting $\omega_{j}\to\omega_{j}+1$, and rewriting $I_{N}(r,s)$ as an Hankel determinant; indeed we have

[TABLE]

The entries of the Hankel matrix vanish for $j+k>s+1$, and hence the determinant is simply given by the product of the antidiagonal entries, $j+k=s+1$ (modulo a sign $(-1)^{s(s-1)/2}$ emerging from the permutation of all columns). Identity (5.1) follows immediately.

5.3. Saddle-point evaluation for large $N$ and finite $s$

When using the saddle-point method in variables $\omega_{1},\dots,\omega_{s}$ to evaluate the behaviour of $F_{N}(r,s)$ for large $N$ and $r$, and finite $s$, it is rather easy to see that the saddle-point equations decouple at leading order, and that each saddle-point will be on the real axis, contributing with a factor $\mathrm{e}^{-NS_{j}}$ with $S_{j}$ positive.

If a given saddle-point is smaller than $1$, the contour $C_{0}$ can be deformed through the saddle-point without encountering any pole, and its contribution will vanish as $\mathrm{e}^{-NS_{j}}$ in the large $N$ limit. If however the saddle-point, still on the real axis, happens to be larger than $1$, the deformation of the contour $C_{0}$ through the saddle-point will pick up the contribution of the pole at $\omega=1$ (with a reversed orientation of the contour), and the $j$-th integral will behave as $1+\mathrm{e}^{-NS_{j}}$. Hence, in the large $N$ limit (at fixed $s$) the quantity $F_{N}(r,s$) will vanish unless all the saddle-points are greater than $1$, in which case $F_{N}(r,s)\sim I_{N}(r,s)=1$. Note that in the present situation the $s$ saddle-points coincide. A detailed analysis shows that in this case the position of the $s$ saddle-points depends on the value $x=r/N$ as $\omega_{0}=\frac{x}{\tau(1-x)}$. In correspondence to the value $x_{0}=\frac{\tau}{1+\tau}$, for which these saddle-points are exactly $1$, the function $F(x,0)$ has a step discontinuity. More precisely, it is easy to show that for $x\in[0,1]$, $F(x,0)=\theta(x-x_{0})$, where $\theta(x)$ is Heaviside step function. From a physical point of view $x_{0}$ is the contact point between the limit shape and the boundary. What have been discussed here can easily be verified in the case $s=1$. The extension to finite $s>1$ is rather direct as well.

5.4. Saddle-point equation

Having in mind the analogy with $s\times s$ Random Matrix Models, and the scaling limit specified in Section 5.1, we rewrite our expression for $F_{N}(r,s)$ at $\Delta=0$ as follows:

[TABLE]

Both sums in the exponent are $O(s^{2})$. The corresponding (coupled) saddle-point equations read

[TABLE]

A standard physical picture reinterprets the saddle-point equations as the equilibrium condition for the positions of $s$ charged particle confined to the real axis, with logarithmic electrostatic repulsion, in an external potential. In the present case the latter can be seen as generated by three external charges, $1$, $x/y$, and $-(1/y-1)$ at positions $1$, [math], and $-1/\tau$, respectively. It is natural to refer to this model as the triple Penner model. Although the simple Penner [P] matrix model has been widely investigated, not so much is known about the much more complicate double Penner model [M, PW]. We have not been able to trace any previous study concerning the triple Penner model.

5.5. The exact Green function at finite $s$

To investigate the structure of solutions of the saddle-point equations (5.3) for large $s$ we first introduce the Green function

[TABLE]

which, if the $\omega_{j}$’s solves (5.3), has to satisfy the differential equation:

[TABLE]

The coefficients $\alpha$, $\beta$ and $\gamma$ are readily obtained as the coefficients of the second order polynomial appearing in the numerator, when setting to common denominator the left hand side of (5.3). We give them explicitly for later convenience:

[TABLE]

The derivation of the differential equation is very standard (see, e.g., [SD]). The left hand side is built by suitably combining the explicit definition of the Green function and its derivative. The result has to be a polynomial of the first degree in $z$, whose coefficients are constructed by matching the leading and first subleading behaviour of the left hand side as $|z|\to\infty$.

5.6. The first moment $\Omega$

The quantity ${\Omega}$ appearing in (5.4) is defined as the first moment of the solutions of the saddle-point equations:

[TABLE]

It is related in a obvious way to the first subleading coefficient of $G_{s}(z)$; indeed, from the definition of the Green function, it is evident that

[TABLE]

It is worth to emphasize that ${\Omega}$ is still unknown, and that in principle its value should be determined self consistently by first working out the explicit solution of $G_{s}(z)$ (which will depend implicitly on ${\Omega}$), from (5.4) and then demanding that $\frac{1}{s}\sum_{j=1}^{s}\omega_{j}$ evaluated from this solution coincides with $\Omega$. The appearance of the undetermined parameter ${\Omega}$ is a manifestation of the ‘two-cuts’ nature of the Random Matrix Model related to (5.2), see, e.g., par. 6.7 of [D1].

5.7. The asymptotic Green function

We are now in condition to perform the large $s$ (and large $N$, $r$) limit at fixed $x,y$. In the limit, we can neglect terms of order $O(s)$ in the differential equation (5.4), which therefore reduces to an algebraic equation for the limiting Green function $G(z)$:

[TABLE]

The previous algebraic equation has to be supplemented by the normalization condition

[TABLE]

Hence the Green function describing the large $s$ asymptotic distribution of solutions for the saddle equation (5.3) reads:

[TABLE]

We have selected the positive branch of the square root, to satisfy the normalization condition (note that the coefficient of $z^{4}$ under the square root is $(\alpha-2\tau)^{2}$, and $\alpha-2\tau$ is negative for any $x,y\in[0,1]$). However, the expression for $G(z)$ is not completely specified yet, because $\Omega$ is still undetermined.

5.8. Limit shape and condensation of roots

The polynomial under the square root is of fourth order, hence $G(z)$ will have in general two cuts in the complex plane. The emergence of a two-cut problem was already expected from the appearance of the undetermined first moment $\Omega$ in (5.4). The discontinuity of $G(z)$ across these cuts defines, when positive, the density of solutions of the saddle-point equations (5.3) when $s\to\infty$. The problem of explicitly finding this density, for arbitrary $\alpha$, $\beta$, $\gamma$ (or $x$, $y$), is a formidable one, not to mention the evaluation of the corresponding ‘free energy’, and of the saddle-point contribution to the integral in (5.2). But our aim is much more modest, since we are presently interested only in the expression of the limit shape, i.e. in the curve in the square $x,y\in[0,1]$, delimiting regions where $F(x,y)=0$ from regions where $F(x,y)=1$. Of course we are here somehow assuming that the transition of $F(x,y)$ from [math] to $1$ is stepwise in the scaling limit, but this is supported both by the physical interpretation of EFP (in the disordered region, by definition, the number of ‘thin’ lines is macroscopic, and the probability of finding no ‘thin’ horizontal edges immediately vanishes in the scaling limit) and by the discussion of Section 5.3.

As explained in the discussion of the double Penner model in paper [PW], the logarithmic wells in the potential can behave as condensation germs for the saddle-point solutions. In our case, this can role can be played only by the ‘charge’ at $\omega=1$ in the Penner potential since the charge at $\omega=-1/\tau$ is always repulsive, while the one at $\omega=0$ is larger than $1$, at least in the region of interest. [PW] have shown that condensation can occur only for charges less than or equal to $1$, since this will be the fraction of condensed solutions. This consideration, together with the expected stepwise behaviour and the discussion in Section 5.3, suggest the following picture for the evolution of saddle-point solution density from the disordered region, $F(x,y)\sim 0$, to the upper left frozen region, $F(x,y)\sim 1$: in the disordered region there is a macroscopic fraction of solutions which are real and smaller than $1$, while in the upper left frozen region this fraction vanish. On the basis of the discussion here and in Sections 3.2 and 5.3, we shall assume that at the transition between the two regions all saddle-point solutions have condensed at $\omega=1$.

5.9. Main assumption

We claim that the Arctic curve in the square $x,y\in[0,1]$ separating the disordered phase from the upper left frozen phase is defined by the condition that all solutions of the saddle-point equation lies at $\omega=1$.

In the derivation of the limit shape, this is indeed the only assumption to which we are unable to provide a proof. There is in fact no guarantee, at this level, for this possibility to occur, and limit shapes could in principle emerge from a different condition. But if for some values of $x,y\in[0,1]$ we have all solutions of the saddle-point equation condensing at $\omega=1$, then this provides a transition mechanism from [math] to $1$ for $F(x,y)$, and this might correspondingly define some limit shape.

If all saddle-point solutions condensate at $\omega=1$, then we obviously have:

[TABLE]

Moreover, the complicate expression (5.7) for $G(z)$ should simply reduce to

[TABLE]

since we expect to have no cuts, and only one pole at $z=1$ with unit residue.

5.10. Arctic Ellipses

Consider the quartic polynomial under the square root in (5.7). It is convenient to rewrite it in terms of

[TABLE]

Note that ${\tilde{\alpha}}$ and ${\tilde{\gamma}}$ are always positive for $x,y\in[0,1]$. When $\Omega=1$, our quartic polynomial reads

[TABLE]

which may be equivalently rewritten as

[TABLE]

We see that the quartic polynomial reduces to a perfect square, and hence, when $\Omega=1$, the two cuts of $G(z)$ disappear, as expected.

Now, when $\Omega=1$, in our new notations, the Green function reads:

[TABLE]

We now require the coefficients ${\tilde{\alpha}}$, ${\tilde{\beta}}$, ${\tilde{\gamma}}$ to be such that the polynomial under the square root combines with the first part of the numerator in (5.10) to give $2z(\tau z+1)$ and simplify the Green function according to (5.8). Once we have chosen a given branch of the square root (the positive one, in order to satisfy normalization condition (5.6)), it is obvious that the required simplification can occur for any $z$ in the complex plane only if the second order polynomial ${\tilde{\alpha}}z^{2}+{\tilde{\beta}}z+{\tilde{\gamma}}$ does not change its sign, i.e. only if its two roots coincide, implying:

[TABLE]

Rewriting the last relation in terms of $x$, $y$, through (5.9), we readily get

[TABLE]

We have therefore recovered the limit shape, which in this Free Fermion case is the well-known Arctic Ellipse (Arctic Circle for $\tau=1$) [JPS, CEP]. We recall that, as discussed in Section 3.2, $F(x,y)$ is non-vanishing only in the upper left region of the unit square. Therefore, concerning EFP, only the upper left portion of the Arctic curve, between the two contact points at $(\frac{\tau}{1+\tau},0)$ and $(1,\frac{1}{1+\tau})$, is relevant.

6. Concluding remarks

Our starting point has been the definition of a relatively simple but relevant correlation function for the six-vertex model with DWBC, the Emptiness Formation Probability. We have provided both a determinant representation and a multiple integral representation for the proposed correlation function. This is the first example in literature of a bulk (as opposed to boundary) correlation function for the considered model, for generic weights.

The multiple integral representation, specialized to the Free Fermion case, has been studied in the scaling limit. In the standard picture of Random Matrix Models, we recognize the emergence of a triple Penner model. Assuming condensation of the roots of saddle point equations in correspondence to a limit shape, we recover the well-known Arctic Circle and Ellipse. It would be interesting to investigate whether universality considerations of Random Matrix Models (see, e.g., [D2]) can be extended to the Penner model in the neighbourhood of its logarithmic singularities. This would imply directly the results of [CEP, J1, J2] on the Tracy-Widom distribution and the Airy process, emerging in a suitably rescaled neighbourhood of the Arctic Ellipse.

It is worth to stress that the multiple integral representation for EFP presented in Section 3 can be studied beyond the usual Free Fermion situation. We expect that condensation of roots of the saddle point equation in correspondence of the limit shape is a general phenomenon. We believe that this assumption could be of importance in addressing the problem of limit shapes in the six-vertex model with DWBC.

Our derivation of the limit shape in the Free Fermion case uses the explicit knowledge of function $h_{N}(z)$, standing in the multiple integral representation (3.6). It is worth mentioning that function $h_{N}(z)$ is also known explicitly at Ice Point, ($\Delta=1/2$), and Dual Ice Point, ($\Delta=-1/2$), being expressible in terms of (polynomial) Gauss hypergeometric function [Ze, CP2]. For instance, at Ice Point the triple Penner model discussed above generalizes to a two-matrix Penner model. This model can be studied along the lines presented here, thus providing a solution to the longstanding problem of limit shape for Alternating Sign Matrices.

Acknowledgements

We thank Nicolai Reshetikhin for useful discussion, and for giving us a draft of [PR] before completion. FC is grateful to Percy Deift, and Courant Institute of Mathematical Science, for warm hospitality. AGP thanks INFN, Sezione di Firenze, where part of this work was done. We acknowledge financial support from MIUR PRIN program (SINTESI 2004). One of us (AGP) is also supported in part by Civilian Research and Development Foundation (grant RUM1-2622-ST-04), by Russian Foundation for Basic Research (grant 04-01-00825), and by the program Mathematical Methods in Nonlinear Dynamics of Russian Academy of Sciences. This work is partially done within the European Community network EUCLID (HPRN-CT-2002-00325), and the European Science Foundation program INSTANS.