Phil. Mag. 6 (1961), 1061-1063. MR 24:B2436

]]> 2000 American Mathematical Society Journal of the American Mathematical Society J. Amer. Math. Soc. jams 14 02 20001103 January 12, 1999 August 11, 2000 November 3, 2000 297 297-346 50 S0894-0347-00-00355-6 10.1090/S0894-0347-00-00355-6 0894-0347 1088-6834 English 82B20 82B23 82B30 2000 http://www.ams.org/jourcgi/jour-getitem?pii=S0894-0347-00-00355-6 quote The effect of boundary conditions is, however, not entirely trivial and will be discussed in more detail in a subsequent paper. flushright P. W. Kasteleyn, 1961 flushright quote Introduction introsec Description of results statementsec A domino is a 1 2 (or 2 1 ) rectangle, and a tiling of a region by dominos is a way of covering that region with dominos so that there are no gaps or overlaps. In 1961, Kasteleyn Kast1 found a formula for the number of domino tilings of an m n rectangle (with mn even), as shown in Figure fig-square for m n 68 . Temperley and Fisher TF used a different method and arrived at the same result at almost exactly the same time. Both lines of calculation showed that the logarithm of number of tilings, divided by the number of dominos in a tiling (that is, mn 2 ), converges to 2G 0.58 (here G is Catalan's constant). On the other hand, in 1992 Elkies et al. EKLP studied domino tilings of regions they called Aztec diamonds (Figure fig-aztec shows an Aztec diamond of order 48), and showed that the logarithm of the number of tilings, divided by the number of dominos, converges to the smaller number (2) 2 0.35 . Thus, even though the region in Figure fig-square has slightly smaller area than the region in Figure fig-aztec , the former has far more domino tilings. For regions with other shapes, neither of these asymptotic formulas may apply. figure scale .9 jams355el-fig-1 A random domino tiling of a square. fig-square figure figure scale .9 jams355el-fig-2 A random domino tiling of an Aztec diamond. fig-aztec figure In the present paper we consider simply-connected regions of arbitrary shape. We give an exact formula for the limiting value of the logarithm of the number of tilings per unit area, as a function of the shape of the boundary of the region, as the size of the region goes to infinity. In particular, we show that computation of this limit is intimately linked with an understanding of long-range variations in the local statistics of random domino tilings. Such variations can be seen by comparing Figures fig-square and fig-aztec . Each of the two tilings is random in the sense that the algorithm PW that was used to create it generates each of the possible tilings of the region being tiled with the same probability. Hence one can expect each tiling to be qualitatively typical of the overwhelming majority of tilings of the region in question, as is in fact the case. Figure fig-square looks more or less homogeneous, but even cursory examination of Figure fig-aztec shows that the tiling manifests different gross behavior in different parts of the region. In particular, the tiling degenerates into a non-random-looking brickwork pattern near the four corners of the region, whereas closer to the middle one sees a mixture of horizontal and vertical dominos of the sort seen everywhere in Figure fig-square (except very close to the boundary). In earlier work CEP,JPS two of us, together with other researchers, analyzed random domino tilings of Aztec diamonds in great detail, and showed how some of the properties of Figure fig-aztec could be explained and quantified. It was proved that the boundary between the four brickwork regions and the central mixed region for a randomly tiled Aztec diamond tends in probability towards the inscribed circle (the so-called arctic circle'') as the size of the diamonds becomes large JPS , and that even inside the inscribed circle, the first-order local statistics (that is, the probabilities of finding individual dominos in particular locations) fail to exhibit homogeneity on a macroscopic scale CEP . Unfortunately, the techniques of JPS and CEP do not apply to general regions. A few cases besides Aztec diamonds have been analyzed; for example, random domino tilings of square regions have been analyzed and do turn out to be homogeneous BP . (That is, if one looks at two patches at distance at least d from the boundary of an n n square, with n even, the local statistics on the two patches become more and more alike as n goes to infinity, as long as d goes to infinity with n .) However, before our research was undertaken no general analysis was known. In this paper, we will demonstrate that the behavior of random tilings of large regions is determined by a variational (or entropy maximization) principle, as was conjectured in Section 8 of CEP . We show that the logarithm of the number of tilings, divided by the number of dominos in a tiling of R , is asymptotic ( when (R) ) to equation supoverh h R ( h x , h y ) , dx : dy. equation Here the domain of integration R is a normalized version of R , the function h ranges over a certain compact set of Lipschitz functions from R to , and equation entdef (s,t) 1 (L(p a) L(p b) L(p c) L(p d)), equation where L() is the Lobachevsky function (see Mil ), defined by equation lobdef L(z) - 0 z 2 t : dt, equation and the quantities p a,p b,p c,p d are determined by the equations eqnarray 2(p a-p b) s fx , 2(p c-p d) t fy , p a p b p c p d 1 sump , (p a)(p b) (p c)(p d). sinsin eqnarray The quantities p a , p b , p c , p d can be understood in terms of properties of random domino tilings on the torus (see Subsection interpretsec ); the quantities also have attractively direct but still conjectural interpretations in terms of the local statistics for random tilings of the original planar region (see Conjecture probconj ). There exists a unique function f that achieves the maximum in supoverh . Its partial derivatives encode information about the local statistics exhibited by random tilings of the region; for example, the places (x,y) where the tilt'' ( f x , f y ) is (2,0) , (-2,0) , (0,2) , or (0,-2) correspond to the places in the tiling where the probability of seeing brickwork patterns goes to 1 in the limit, in a suitable sense. The function f need not be C , and in fact often is not. For example, in the case of Aztec diamonds f is smooth everywhere except on the arctic circle, where it is only C 1 (except at the midpoints of the sides, where it is only C 0 ). Inside this circle, f takes on a variety of values, corresponding to the fact that different local statistics are manifested at different locations. In contrast, for square regions the function f is constant, corresponding to the fact that throughout the region (except very close to the boundary), the local statistics are constant. (See Figures fig-square and fig-aztec .) In general, the locations where f is not smooth correspond to the existence of a phase transition in the two-dimensional dimer model, closely related to a phase transition first noticed by Kasteleyn Kast2 . The function f is related to combinatorial representations of states of the dimer model called height functions . Each of the different tilings of a fixed finite simply-connected region in the square grid has a height function which is a function from the vertices of the region to satisfying certain conditions (see Subsection hfsec ); there is a one-to-one correspondence between tilings and height functions (as long as we fix the height at one point, since the actual correspondence is between tilings and height functions modulo additive constants). All of the height functions agree on the boundary of the region but have different values on the interior. In an earlier paper CEP , it was shown that these height functions satisfy a law of large numbers, in the sense that for each vertex v within a very large region, the height of v in a random tiling of the region is a random variable whose standard deviation is negligible compared to the size of the region. This suggested the existence of a limit law, but did not indicate what the limit law was. We show (see Theorem maintheorem at the end of Subsection sketchsec ): for every 0 , the height function of a random tiling, when rescaled by the dimensions of the region being tiled, is, with probability tending to 1, within of f , where f is the function that maximizes the double integral in supoverh . We therefore call f an asymptotic height function . One may qualitatively summarize our main results by saying that the pattern governing local statistical behavior of uniform random tilings of a region is, in the limit, the unique pattern that maximizes the integral of the local entropy'' over the region. Moreover, the value of this maximum is an asymptotic expression for the global entropy'' of the ensemble of tilings of the region. Our main theorem (Theorem maintheorem ) thus has two intertwined components: a law of large numbers that describes the local statistics of random tilings of a large region by way of an asymptotic height function f , and an entropy estimate that tells us that the total number of tilings is determined by a certain functional of that self-same function f . Furthermore f is precisely the unique Lipschitz function that maximizes that functional. We supplement our main theorem by two supporting results: a large deviations estimate (Theorem main ) and a PDE that the entropy-maximizing function f must satisfy (Theorem pdethm ). To be technically correct, we should not assert that the function f satisfies the PDE everywhere, but only where the partial derivatives are continuous and the tilt (s,t) satisfies s t 2 . This proviso is necessary because singularities in f (corresponding to domain boundaries like the arctic circle) are an essential feature of the phenomena we are studying; in fact, we shall find (see Section partitionsection and Theorem pa ) that these singularities are related to a phase transition manifested by the dimer model as an external field is permitted to vary, and that domain boundaries can be viewed as a spatial expression of that phase transition. Our results also imply that if instead of studying the uniform distribution on the set of all the tilings of a large region R , one restricts one's attention to those tilings whose height functions (suitably normalized) approximate some asymptotic height function h (which need not be the h that maximizes the double integral supoverh ), then the entropy of this restricted ensemble (that is, the normalized logarithm of the number of tilings) approximates the value of the double integral. Note that this more general result implies that the total number of (unrestricted) tilings is at least as large as the supremum supoverh . As for the integrand itself, we remark here that (s,t) achieves its maximum at (0,0) and that it goes to zero as (s,t) goes to the boundary curve s t 2 ; that is, there are many ways to tile a patch so that its average tilt is near zero, but fewer ways to tile it as the tilt gets larger, and no ways to tile it at all if the desired tilt (s,t) fails to satisfy s t 2 . Interpretation interpretsec A dimer configuration , or perfect matching , of a graph is a set of edges in the graph such that each vertex belongs to exactly one of the selected edges. To see the equivalence between tilings by dominos and dimer configurations on a grid (the graph-theoretic dual to the graph of edges between lattice squares), one need only replace each 1 1 square by a vertex and each domino by an edge. To explain the significance of the quantities p a,p b,p c,p d , we need to digress and discuss the dimer model on a torus. Here the graph is just like the m n rectangular grid, but with m extra bonds that connect vertices on the left and right and n extra bonds that connect vertices on the top and bottom. Kasteleyn showed Kast1 that the number of dimer configurations on this graph is governed by the same asymptotic formula as for the dimer model on a rectangle. In the same article Kasteleyn considered a generalization to non-uniform distributions on the set of tilings, obtained by assigning different weights'' to horizontal and vertical bonds, and letting the probability of a particular dimer configuration be proportional to the product of the weights of its bonds. In the statistical mechanics literature, such modifications of a model are sometimes conceived of as resulting from the imposition of an external field. Here we go one step further, and impose a field that discriminates among four different kinds of bonds: a -bonds and b -bonds (both horizontal) and c -bonds and d -bonds (both vertical), staggered as in Figure evects from Section detsection (this is different from a 4-parameter external field that was considered by Kasteleyn Kast2 ). The field depends on four parameters a,b,c,d , which are the weights associated with the respective kinds of bonds; the probability of a dimer configuration on the torus graph is proportional to the product of the weights of its bonds. Taking the graph-theoretic dual, we get a non-uniform distribution on domino tilings on the torus, having less symmetry than the torus itself. One way to motivate the consideration of such asymmetrical measures is to observe that in regions like the one shown in Figure fig-aztec , the boundary conditions break the symmetry of the underlying square grid in precisely this fashion, yielding subregions in which bonds (or dominos) of one of the four types predominate. To avoid unnecessary complexity, we limit ourselves to square tori ( m n ), but the situation for more general tori is much the same. The quantities p a,p b,p c,p d have a direct (and rigorously proved) interpretation in terms of the dimer model on the n n torus. Suppose we are given positive real weights a,b,c,d as described above. Suppose that each of a,b,c,d is less than the sum of the others. Then there is a unique Euclidean quadrilateral of edge lengths a,c,b,d (in that order) which is cyclic, that is, can be inscribed in a circle. Define p a to be 1 (2) times the angle of arc of the circumscribed circle cut off by the edge a of this quadrilateral. Similarly define p b,p c , and p d . Then we shall see (see Theorem pa ) that p a is (in the large- n limit) the probability that a given a -bond belongs to a randomly-chosen dimer configuration (under the probability distribution determined by the weights a,b,c,d ), and likewise for p b,p c,p d . Technically, we only prove convergence for n in a large subset of , but we believe that convergence holds for all n . If on the other hand ab c d , then we shall see that as n , p a tends to 1 (and p b,p c,p d tend to zero). A similar phenomenon occurs when b , c , or d is greater than or equal to the sum of the others. Moreover, the quantities s 2(p a-p b) and t 2(p c-p d) have a height-function interpretation. Since the torus can be viewed as a rolled-up plane, every dimer configuration on the torus unrolls'' to give a dimer configuration on the plane, which (in the guise of a domino tiling of the plane) gives rise to a height function. The height is not in general a periodic function, but rather increases by some amounts H x and H y as one moves n vertices in the x direction or n vertices in the y direction. Here H x and H y depend on the tiling chosen, and thus are random variables; the expected values of H x n and H y n are s and t , respectively. The relationship between the uniform measure on tilings of finite simply-connected regions and the a,b,c,d -weighted measure on tilings of tori may become clearer after one has verified that the a,b,c,d -weighted measure on the tilings of a finite simply-connected region, defined in the obvious fashion, is nothing other than the uniform distribution, as long as ab cd . To see this, one may make use of the fact that any domino tiling of such a region can be obtained from any other tiling by a sequence of moves, each of which consists of applying a 90-degree rotation to a pair of dominos that form a 2 2 block . Such a move trades in an a -domino and a b -domino for a c -domino and a d -domino, or vice versa. The condition ab cd then guarantees that the two tilings have the same weight. Since such moves suffice to turn any tiling into any other, all the tilings have the same weight, and the probability distribution is uniform. This property is called conditional uniformity'', because if one conditions on the tiling outside a simply-connected region, the conditional distribution on tilings of the interior is uniform. We do not propose any concrete interpretation for the weights a,b,c,d , and it might be unreasonable to ask for one, given that the probability distribution on matchings determined by a,b,c,d is unaffected if all four are multiplied by a constant. There is a choice of scaling that makes most of our formulas relatively simple, namely, the scaling that makes the arcsines of the four weights add up to , or equivalently, the scaling that makes the product (a b c-d)(a b-c d)(a-b c d)(-a b c d) equal to the product 4(ab cd)(ac bd)(ad bc) (this is most easily seen from Theorem pa ). However, it is worth mentioning at least one instance in which the scaling -1 a -1 b -1 c -1 d does not give the simplest possible formula for a,b,c,d . Specifically, consider the normalized Aztec diamond with vertices at (1, 0) and (0, 1) . The arctangent law of CEP gives one way of expressing how p a,p b,p c,p d vary throughout the normalized diamond (corresponding to the fact that the respective densities of north-, south-, east-, and west-going dominos vary throughout large Aztec diamonds). However, an even more compact way of stating this dependence is via the formulas eqnarray a (1 y) 2 - x 2 , b (1-y) 2 - x 2 , c (1 x) 2 - y 2 , d (1-x) 2 - y 2 . eqnarray Sketch of proof and preliminary statement of results sketchsec The strategy behind our proof of the main theorem is roughly as follows. In the first few (and more qualitative) sections of the article (Sections beginpartI through endpartI ), we cover the set of all domino tilings of the region by subsets which are balls in the uniform metric on height functions (that is, each subset consists of tilings whose height functions are approximately equal); each ball is associated with an asymptotic height function h from a bounded subset of 2 into . Appealing to quantitative results proved later in more technical sections of the article (Sections partitionsection through paproof ), as well as to combinatorial arguments, we show that the logarithm of the cardinality of a ball is approximated by the double integral in supoverh times half the area of the region. However, we also show that there is a unique ball in our cover for which the corresponding h maximizes the double integral, and that the contribution that this ball makes to the total number of tilings swamps all the other contributions. In Section detsection , we compute the partition function'' for tilings (that is, the sum of the weights of all the tilings). This computation relies on Kasteleyn's original work Kast1 . The next few sections of the article, on which the first few sections depend, are close in spirit to Kasteleyn's original work on the dimer model on the torus. (Indeed, the local entropy (s,t) is equal to the asymptotic entropy for dimer covers of the nn torus as n goes to infinity, where the edges have been assigned new weights as described above, favoring those tilings that have tilt near (s,t) .) Using Kasteleyn's method of Pfaffians, we show that the dimer model in the presence of weights a,b,c,d exhibits a phase transition when any of the four weights equals the sum of the other three, and our method of analysis also requires that special attention be given to the case in which a b and c d (which was analyzed by Kasteleyn). The calculations are lengthy, but the end results, obtained in the final sections of the paper, are satisfyingly simple. For direct formulas for p a,p b,p c,p d in terms of s f x and t f y , see ( probdef ) below. Here is a loose statement of our result. See Theorems main and nbhdentropy for more precisely quantified statements, and the remainder of Section introsec for the relevant definitions. thm maintheorem Let R be a region in 2 bounded by a piecewise smooth, simple closed curve R . Let h b : R be a function which can be extended to a function on R with Lipschitz constant at most 2 in the sup norm. Let f : R be the unique such Lipschitz function maximizing the entropy functional (f) ( see equations entdef and entfuncdef ) , subject to f R h b . ( See Section beginpartI for the proof of the asserted uniqueness. ) Let R be a lattice region that approximates R when rescaled by a factor of 1 n , and whose normalized boundary height function approximates h b . Then the normalized height function of a random tiling of R approximates f , with probability tending to 1 as n . Furthermore, let g : R be any function satisfying the Lipschitz condition. Then the logarithm of the number of tilings of R whose normalized height functions are near g is n 2((g) o(1)) . thm In his thesis , Hoffe derives (in a less rigorous manner) expressions equivalent to our Zformula and pa . In DMB , Destainville, Mosseri, and Bailly set up a similar framework to our variational principle, but they use a quadratic approximation to the entropy, and they do not supply rigorous proofs. Readers of this article might also wish to read Kenyon2 and P2 . These two articles cover some of the same ground as this one, but with more of an emphasis on phenomena and less concern with proofs. Before we can continue our discussion of the variational principle and state the main result more precisely (Subsection vpsec ), we need to review some background on entropy and height functions. Height functions hfsec Height functions are a geometrical tool discovered in individual cases by Blote and Hilhorst BH and Levitov Lev and independently studied by Thurston , who situated the idea in a less ad hoc, more general context. Given a (connected and) simply-connected region R that can be tiled by dominos, domino tilings of R are in one-to-one correspondence with height functions on the set of lattice points in R , defined up to an additive constant. The height function representation is useful because the difference between the heights of two vertices encodes non-local information about a tiling. Here we will quickly summarize the basic definitions and properties of height functions. For a more extensive discussion aimed at applications to random tilings, see Subsections 6.1--6.3 and Section 8 of CEP ; for a geometrical point of view, see . Color the lattice squares in the plane alternately black and white, like a checkerboard. (It does not matter which of the two ways of doing this is used.) Define a horizontal domino to be north-going or south-going according to whether its leftmost square is white or black, and define a vertical domino to be west-going or east-going according to whether its upper square is white or black. These names come from the domino shuffling algorithm from EKLP , which is the main combinatorial tool used to study Aztec diamonds. We will not use domino shuffling, but the division of dominos into these four orientations will nevertheless be important, as will the coloring in general. (It is intended that north-going dominos correspond to a -edges, south-going to b -edges, east-going to c -edges, and west-going to d -edges. The geometrical terminology is more pleasant when we have no weights in mind.) Let R be a lattice region (i.e., a connected region composed of squares from the unit square lattice), and let V be the set of lattice points within R or on its boundary. We will always assume that R can be tiled by dominos in one or more ways. A height function h on R is a function from V to that satisfies the following two properties for adjacent lattice points u and v on which it is defined. (We consider two lattice points to be adjacent only if the edge connecting them is contained in R .) First, if the edge from u to v has a black square on its left, then h(v) equals h(u) 1 or h(u)-3 . Second, if the edge from u to v is part of the boundary of R , then h(u)-h(v) 1 . Given a height function h on R , consider the set of all pairs of adjacent lattice points u,v with h(u)-h(v) 1 ; then it is not hard to check that the set of all dominos in R that are bisected by such an edge taken together constitutes a tiling of R . Given a domino tiling of a simply-connected region R , we can reverse the process, and find a height function that leads to the tiling. Such a height function always exists, but it is not quite unique, because one can add a constant (integer) value to it everywhere to get another such height function. To avoid this ambiguity, we fix the height of one lattice point on the boundary of R . Then height functions satisfying this constraint are in one-to-one correspondence with domino tilings. From this point on, when we use the term lattice region , we will always implicitly assume that one of the lattice points on the boundary has a specified height, and when we talk about height functions we will always assume that they satisfy this condition. Notice that when the region is simply-connected, our assumption takes on an especially pleasant form, because it is then equivalent to fixing the heights along the entire boundary. (However, when the boundary is in several pieces, the situation is more complicated. In this paper, we will typically assume that our regions are simply-connected, although we will always mention that assumption when we make it.) To define a height function for a domino tiling of an n n torus, view the torus as an n n square with opposite sides identified, and view this square as sitting centered inside an (n 2) (n 2) square. A tiling of the torus determines a covering of the smaller square by dominos that lie inside the larger square. (Dominos that jump from one side of the small square to the other correspond to two dominos in the large square; the rest correspond to just one.) Then define the height function on the n n square as for a lattice region in the plane. This height function may not be well defined on the torus, since its values can be different at two identified vertices, but this will not matter for our purposes. This definition is not really natural, but it is convenient. (See spaces for a more natural definition of height functions on tori. We will not use their definitions or results.) Let R be a lattice region, and let H(R) be the set of all height functions on R . We define a partial ordering on H(R) by setting h 1 h 2 if h 1(u) h 2(u) for every lattice point u R . The set H(R) is not just a partially ordered set, but also a lattice. The join of two height functions h 1 and h 2 can be defined by (h 1 h 2)(u) (h 1(u),h 2(u)) , and their meet by (h 1 h 2)(u) (h 1(u),h 2(u)) . To show that H(R) is a lattice, we need to show that h 1 h 2 and h 1 h 2 are height functions. Notice that, at each lattice point u , the value of height functions at u is determined modulo 4 (independently of the specific tiling); this assertion is trivial when u is the point at which we have fixed the values of height functions, and if it is true at some particular lattice point, then the definition of a height function immediately implies that it is true at all adjacent lattice points. Thus, if h 1(u) is unequal to h 2(u) , then they differ by at least 4 . It now follows from the definition of a height function that if h 1(u) h 2(u) , then h 1(v) h 2(v) for all lattice points v adjacent to u . Hence, given any two adjacent lattice points, h 1 h 2 (or h 1 h 2 ) agrees with h 1 at both points, or agrees with h 2 at both of them, which is what is needed to show that h 1 h 2 and h 1 h 2 are height functions. (See P1 for details and more general results.) As Figures fig-square and fig-aztec demonstrate, the precise boundary conditions of a region can have a dramatic effect on the behavior of typical tilings. Height functions provide the proper tool for gauging the effect of boundary conditions on random tilings. (Proposition 20 of CEP is one way to make this claim precise. It says, roughly, that the height function of a typical tiling depends continuously on the heights on the boundary of the region.) For example, one can compute the rate of change of the height function in terms of the probabilities of finding dominos in given locations. If a region has statistically homogeneous random tilings, then the typical height function should be nearly planar. Because the boundary heights for Aztec diamonds are highly non-planar, typical tilings cannot be homogeneous. Define the average height function of a lattice region R to be the average of all height functions on R . (Of course, it is almost never a height function itself.) Theorem 21 of CEP implies (roughly) that if R is large, then almost all height functions on R approximate the average. Thus, the problem of describing typical tilings is reduced to the problem of describing the average height function. Entropy The entropy of a random variable that takes on n different values, with probabilities p 1,,p n , is defined as i 1 n -p i p i (with 0 0 0 ). For a uniformly distributed random variable, the entropy is simply the logarithm of the number of possible outcomes. We will nevertheless often need to deal with entropy for a non-uniform distribution. In general, we denote the entropy of a random variable X by (X) ; there should be no confusion with our use of () '' to denote local entropy depending on tilt. We will sometimes speak of the conditional entropy of a discrete random variable, conditional upon some event; this is just the entropy of the conditional distribution determined by the conditional probabilities. The only fact about entropy that we will need other than the definition is the following standard fact about conditional entropy (the proof is straightforward). lem entropy Suppose X is a random variable that takes on values x 1,,x n , and suppose that x 1,,x n is partitioned into blocks B 1,,B m . Let B be the random variable that tells which block X is in, and let X i be the random variable that takes on values in B i according to the conditional distribution of X given that X B i . Then the entropy of X is given by equation (X) (B) i 1 m (x B i) (X i). equation lem When we deal with entropy for tilings of a lattice region R , we will always normalize it by dividing by half the area of R , so that we measure the information content per domino. Thus, the normalized entropy of a set of tilings of a lattice region R (under the uniform distribution) is the logarithm of the number of tilings in the set, normalized by dividing by the number of dominos in a tiling of R . When we refer to the entropy of a region R , we mean the entropy of the set of all tilings of R , under the uniform distribution. The variational principle vpsec Let R be a region in 2 bounded by a piecewise smooth, simple closed curve R . (We will always assume our curves do not have cusps.) Suppose R is a large, simply-connected lattice region. We can normalize by scaling the coordinates by a factor of 1 n (for some appropriate choice of n ). Suppose that the normalization of R approximates R (in a sense to be clarified below). If we scale the values of height functions by dividing their values by n , then for large n the normalized height functions that one obtains approximate functions on R that satisfy a Lipschitz condition with constant 2 for the sup norm. The reason for this is that if u and v are lattice points at sup norm distance at most d from each other within R (i.e., they are connected by a path within R of length d , where the allowable steps are lattice edges or diagonal steps) and h is any height function, then one can check that h(u)-h(v) 2d 1 . (One can prove this claim directly, or deduce it from Lemmas fournier and supnorm .) Whenever we refer to a 2-Lipschitz function f from a subset of 2 to , we mean one that is locally Lipschitz with Lipschitz constant 2 for the sup norm, i.e., its domain is covered by open balls in which every pair of points (x 1,y 1) and (x 2,y 2) satisfies equation f(x 1,y 1)-f(x 2,y 2) 2 ( x 1-x 2 , y 1-y 2 ). equation Notice that where such a function f is differentiable, it must satisfy equation f x f y 2. equation Conversely, every differentiable function satisfying this condition is 2 -Lipschitz. We call (f x, f y) the tilt of the function (and use this terminology whether or not f is Lipschitz). It is important to keep in mind Rademacher's theorem (Theorem 3.16 of ), which says that every Lipschitz function is differentiable almost everywhere. Suppose that R is a simple, closed, piecewise smooth curve in 2 that bounds a region R . We say that a function h b defined on R is a boundary asymptotic height function if there exists a 2-Lipschitz function h on R such that h R h b , and we call such an h an asymptotic height function . Let (R ,h b) be the set of all asymptotic height functions on R with boundary heights h b ; notice that this set is convex and that it is compact with respect to the sup norm. Before continuing, we need to specify exactly what it means for one region bounded by a simple closed curve to approximate another. This can be defined by the Hausdorff metric on closed subsets of 2 ; specifically, two regions are defined to be within of each other if the -neighborhood of each one's boundary curve contains the other's. We also need to discuss what approximation to within means when there are boundary asymptotic height functions on the curves. Given regions R 1 and R 2 with boundary asymptotic height functions h 1 and h 2 , we say that (R 1,h 1) is within of (R 2,h 2) if for all x 1 R 1 there exists an x 2 R 2 such that d(x 1,x 2) and h 1(x 1) - h 2(x 2) , and vice versa. For technical reasons, it is most convenient to write out our arguments in the case of approximation from within, where we have a sequence R 1,R 2, of regions in the interior of R whose limit is R . In most cases, this is easily arranged. The regions R 1,R 2, will be rescaled lattice regions, and if R is a star-shaped region, then by adjusting the scaling we can assume approximation from within, without changing our asymptotic results. However, when R is not star-shaped, slightly more is required. There are two ways to deal with such domains. First, Proposition 20 of CEP tells us that average height functions behave continuously when boundary heights are perturbed, and this lets one adjust regions so that their normalizations fit within the limiting region, without changing the asymptotics. Second, one can give a direct proof by cutting a general region into star-shaped pieces. Thus, we will be able to assume approximation from within in later sections without loss of generality. In the case of approximation from within, it will be convenient to use a slightly different definition of -approximation. If R 2 is in the interior of R 1 , with boundary height functions h 1 and h 2 , then we say that (R 1,h 1) is within of (R 2,h 2) if their boundaries are within , and for every height function h on R 1 extending h 1 , the restriction of h to R 2 is within of h 2 . This clearly generates the same topology as the definitions above, but it is more convenient to work with, so we will use it in later sections. We saw in the first paragraph of this subsection that if R approximates R when rescaled, then normalized height functions on R approximate asymptotic height functions on R (because of the Lipschitz condition on height functions). Later, Proposition hfexist will tell us that every asymptotic height function is nearly equal to a normalized height function; that is, the class of asymptotic height functions has not been defined too broadly. We will show that the average height function on R is determined by finding the asymptotic height function f that maximizes the integral of local entropy, which depends on the tilt of f . Given (s,t) satisfying s t 2 (i.e., a possible tilt for an asymptotic height function), define the local entropy integrand (s,t) as follows. First define eqnarray probdef p a 4 1 2 -1 ( (s 2)-(t 2) 2), Peqns p b -4 1 2 -1 ( (s 2)-(t 2) 2), p c 4 1 2 -1 ( (t 2)-(s 2) 2), p d -4 1 2 -1 ( (t 2)-(s 2) 2), eqnarray where the values of -1 are taken from 0, . We set a (p a) , b (p b) , c (p c) , and d (p d) when s t 2 , and we define (s,t) as in entdef . We will see later that there is good reason to believe that the numbers p a,p b,p c,p d correspond to the local densities of the four orientations of dominos. However, we do not have a proof. Define the entropy functional on (R ,h b) by equation entfuncdef (h) 1 (R ) R , dx , dy. equation Suppose that R 1,R 2, is a sequence of simply-connected lattice regions with specified boundary heights, such that when R n is normalized by n (or, say, a constant times n ), as n the boundary converges to R and the boundary heights to h b (as specified above). Without loss of generality, we can assume that the normalized boundaries of R 1,R 2, all lie within R . We know from Proposition exists unique that there is a unique f (R ,h b) that maximizes (f) . We will prove (in Theorem nbhdentropy ) that the entropy of tilings of R n whose normalized height functions are close to f (that is, the logarithm of the number of such tilings, divided by half the area of R n ) is (f) o(1) as n . We can now prove a sharpened version of the claim made in the second paragraph of Theorem maintheorem : thm main For each 0 , the probability that a normalized random height function on R n differs anywhere by more than from the entropy-maximizing function f is exponentially small in n 2 ( i.e., is bounded above by r n 2 for some r 1) . thm proof Cover (R ,h b) with open sets around each asymptotic height function h , such that the entropy for normalized height functions in these sets is strictly less than (f) unless h f . By compactness, only finitely many of the sets are needed to cover (R ,h b) ; we include the one corresponding to h f . Then Theorem nbhdentropy implies that for n large, the probability that a random normalized height function on R n does not lie in the open set around f is exponentially small in n 2 . proof Theorem main provides a much stronger large deviation estimate than the best previous result (Theorem 21 of CEP ). In addition, it gives the first proof that under these conditions the normalized average height function of R n converges to f as n . (It was previously not known to converge at all, nor was a precise characterization of the entropy-maximizing function known.) It is worth pointing out that all our results generalize to tilings with unit lozenges (as studied in, for example, CLP ). The combinatorial results for lozenge tilings (or, equivalently, the dimer model on a honeycomb graph) are straightforward modifications of those we present here for domino tilings; the analytic results are special cases of ours (set one of the four weights a,b,c,d equal to 0 to move from weighted domino tilings of tori to weighted lozenge tilings). In Sections beginpartI through endpartI , we will need the following three facts proved later in the article: enumerate For every tilt (s,t) satisfying s t 2 , there exist unique weights a,b,c,d (up to scaling) that satisfy ab cd and give tilt (s,t) (see Subsection entasfnofslope ). If an n n torus has edge weights a,b,c,d such that ab cd and the tilt is (s,t) , then the normalized entropy (for the probability distribution on the tilings) is (s,t) o(1) as n (Theorem ent ). The local entropy function (,) is strictly concave and continuous as a function of tilt (Theorem concavethm ). enumerate Analytic results on the entropy functional beginpartI In this section, we will phrase all our results in terms of the local entropy integrand, but they will depend only on its continuity and strict concavity. (Thus, if desired, the results can be applied in different circumstances, for example to other sorts of tiling problems.) Everything in this section is fairly standard material from geometric measure theory and the calculus of variations; for example, Theorem 5.1.5 of is a more sophisticated version of Lemma semicontinuous . However, we will give complete proofs, partly to make this part of the paper self-contained and accessible, and partly because some of what we need does not seem to appear in the literature in quite the form we would like. In what follows, h b denotes a particular boundary asymptotic height function and h ranges over the asymptotic height functions that restrict to h b on the boundary of the region R . lem semicontinuous The functional : (R ,h b) , defined ( as above ) by equation (h) 1 (R ) R ,dx ,dy, equation is upper semicontinuous. lem For the proof of Lemma semicontinuous , as well as other applications later in the paper, we will need to know how well we can approximate an asymptotic height function h by a piecewise linear function. The simplest way to achieve such an approximation is as follows. For 0 , look at a mesh made up of equilateral triangles of side length (which we call an -mesh). Consider any piecewise linear, 2-Lipschitz function h that agrees with h at the vertices of the mesh and is linear on each triangle. (Of course, h depends on the mesh as well as on h . It is uniquely determined on each triangle that lies completely within R , but not on those that extend over the boundary.) lem approx Let h be an asymptotic height function, and let 0 . If is sufficiently small, then, on at least a 1- fraction of the triangles in the -mesh that intersect R , we have the following two properties. First, the piecewise linear approximation h agrees with h to within . Second, for at least a 1- fraction ( in measure ) of the points x of the triangle, the tilt h'(x) exists and is within of h '(x) . lem Keep in mind that h'(x) is the vector of partial derivatives of h at x ; it does not matter which norm we use to measure the distance between tilts, but for the sake of specificity we choose the 2 norm. Notice that the second property implies that (h) (h) o(1) (that is, (h) (h) as 0 ). (Keep in mind that h' 1 2 .) proof We begin with the first of the two properties. Recall that Lipschitz functions are differentiable almost everywhere (Rademacher's theorem) . For any point x at which h is differentiable, we have h(x d)-(h(x) h'(x)d) d 2 if d is sufficiently small, say d r with r 0 (where r depends on x ). If x lies within an equilateral triangle of side length with r , then on that triangle we have the approximation property we want, because there the two functions d h(x d) and d h (x d) (the unique linear function that agrees with g on the corners of the triangle) both lie within 2 of d h(x) h'(x)d . Given 0 , let S be the set of all x such that r (depending on x as above) can be taken to be at least . Take small enough that the measure of S is at least (1- 3) times (R ) . (We can do that since h is differentiable almost everywhere.) Now take . Look at any -mesh, and the piecewise linear approximation h we get from it. If is sufficiently small, then all but an o(1) fraction of the mesh triangles lie entirely within the region. At least a 1- 3-o(1) fraction of them must intersect S , which proves that the desired approximation property holds on at least a 1- 3-o(1) fraction of the triangles. For the second property, we will apply a result on metric density (see Section 7.12 of Rudin ). Let U 1,,U n be open subsets covering the set of possible tilts such that any two tilts contained within the same subset differ by at most , and for 1 i n let V i x : h'(x) U i . It follows from the theorem on metric density that (for each i ) if is sufficiently small, then for all but an 3 fraction of the points x V i , at least a 1- fraction of the ball of radius about x lies in V i (and thus the tilt at those points differs by at most from h'(x) ). Now we can take small enough that this result holds for all i from 1 to n , and then as above it follows that for , a 1- 3-o(1) fraction of the mesh triangles in any -mesh lies entirely within R and satisfies the second property. Thus, if is sufficiently small, at least a 1- fraction of the triangles satisfies both properties (since a 1- 3-o(1) fraction satisfies each). proof lem semielliptic Suppose that R is an equilateral triangle of side length , and that the asymptotic height function h satisfies h b-h on R , where h is some linear function. Then equation (h) (h) o(1) equation as 0 . lem proof Because () is concave, equation concaveineq 1 (R ) R , dx , dy (h x, ,h y, ), equation where equation h x, 1 ( R ) R h x , dx ,dy equation and equation h y, 1 ( R ) R h y ,dx ,dy. equation We have equation (h x, ,h y, ) (h x, ,h y, ) O() equation (as one can see by computing the average by integrating over cross sections and applying the fundamental theorem of calculus), and hence equation (h x, ,h y, ) (h x, , h y, ) o(1), equation since (s,t) is a continuous function of (s,t) . Now combining concaveineq with equation (h) (R )(h x, ,h y, ) equation yields the desired result. proof proof Proof of Lemma semicontinuous Let h be any asymptotic height function, and consider the neighborhood U of h consisting of all asymptotic height functions within of h . We need to show that given 0 , if is sufficiently small, then for all g U , (g) (h) . Let ' 0 . It follows from Lemma approx that if is small enough, then the piecewise linear approximation h coming from an -mesh satisfies h-h ' on all but an ' fraction of the triangles of the mesh, and (h) (h) o(1) as ' 0 . Let ' . Then if g U , g-f 2' on all but an ' fraction of the triangles, and by Lemma semielliptic the contribution to (g) from the non-exceptional triangles is at most o(1) greater than the corresponding contribution to (h) . Of course, the remaining ' fraction of the triangles contribute a total of O(') . Hence, equation (g) (1-O('))(h) o(1) O(') (h) o(1). equation By choosing ' sufficiently small, one can make the o(1) term less than . Thus, () is upper semicontinuous. proof prop exists unique There is a unique asymptotic height function f (R ,h b) that maximizes (f) . prop proof For existence, we can use a compactness argument, since (R ,h b) is compact. Because the local entropy integrand is bounded, () is bounded above, and we can choose a sequence h 1,h 2, such that (h i) approaches the least upper bound as i . By compactness, there is a subsequence that converges, and by upper semicontinuity, the limit of the subsequence must have maximal entropy. Now uniqueness is easy. Suppose that f 1 and f 2 are two asymptotic height functions that maximize entropy, with f 1 f 2 . Their derivatives cannot be equal almost everywhere, so for some 0 we must have (f 1 f 2) 2 f 1 f 2 2 on a set of positive measure, by the strict concavity of () . It follows that equation ( f 1 f 2 2) (f 1) (f 2) 2, equation which contradicts the assumption that (f 1) and (f 2) were maximal. (Notice that since f 1 and f 2 are asymptotic height functions, so is (f 1 f 2) 2 .) Therefore, only one asymptotic height function can maximize entropy. proof Combinatorial results on entropy In this section, we will compare the entropies of several regions with nearly the same shape, but slightly different boundary conditions. To deal with this sort of situation, we will use partial height functions. Let R be a simply-connected lattice region. A partial height function on R is an integer-valued function h defined on a subset of the lattice points in R , including all the lattice points on the boundary, such that h satisfies the following condition. If u and v are adjacent lattice points on which h is defined, such that the edge from u to v has a black square on its left, then h(v)-h(u) is 1 or -3 . We call h a complete height function if it is defined on all the lattice points in R , and we say that a complete height function is an extension of a partial height function if they agree where both are defined. We call a partial height function a boundary height function if the set of lattice points on which it is defined is the boundary of R and there exists a complete height function extending it. Define a free tiling (or a tiling with free boundary conditions'') of a region as a covering of the region by dominos, no two of which overlap, and each of which contains at least one cell belonging to the region. A free tiling is just like a tiling, except that the tiles are allowed to cross the boundary. Given a boundary height function h , extensions of h to R correspond to tilings where dominos cross the boundary of R at certain specified places (namely the places where the height changes by 3 ). In this way h singles out a certain subset of the free tilings of R . We could phrase all of our results in terms of free tilings, but partial height functions will be a more convenient formulation. Suppose we have a boundary height function on R , such that there is an extension to a complete height function on R . We can determine the maximal extension H (under the usual partial ordering) as follows. For any lattice point v in R , look at boundary points w and paths joining w to v such that every edge of (oriented from w towards v ) has a black square on its left. Call such a path an increasing path , since it if does not cross any dominos, then the height increases by 1 after each edge. Define a decreasing path analogously. lem Fournier fournier For each lattice point v R , define H (v) as the minimum of h(w) () , where w ranges over all boundary lattice points and ranges over all increasing paths in R from w to v , and define H (v) as the maximum of h(w)-() , where this time ranges over all decreasing paths in R from w to v . Then H is the maximal extension of h to R , and H is the minimal extension. lem For a proof, see Fournier . Notice that Lemma fournier implies that if one changes the boundary height function h by at most a constant c at each point, then H and H change by at most c as well. We can now prove a proposition we will need later that connects asymptotic height functions to actual height functions. Suppose R is a domino-tileable lattice region, such that rescaling R by a factor of 1 n gives a region whose boundary lies close to a region R , where R is a piecewise smooth, simple closed curve. Without loss of generality, we can suppose that the normalized region lies within R (as discussed in Subsection vpsec ). prop hfexist Given an asymptotic height function f which is within of the normalized boundary heights of R , there is an actual height function on R whose normalization is within O(1 n) of f . prop proof Let g be the largest height function on R whose normalization is less than or equal to f , i.e., the lattice sup of all height functions below f . (Technically, we must pick a lattice point and restrict our attention to height functions that give it height 0 modulo 4 . This ensures that all our height functions are equal modulo 4 and thus that the lattice operations are well defined.) Then all values of g are within 8 of what one would get simply by un-normalizing f , because if one takes any lattice point, assigns it a value that is correct modulo 4 and at least 4 below the un-normalized value of f , and looks at the minimal height function taking that value there, then the Lipschitz constraint on f implies that one stays below f . The height function g will typically not have the correct boundary values on R . Instead, it will have boundary heights corresponding to some different boundary height function b , although they will be within n O(1) of the actual boundary heights for R . To fix this problem, let h be the minimal height function on R , and let H be the maximal one. Then consider (g h) H , which is a height function with the correct boundary values for R . It follows from Lemma fournier that h and H differ by only n O(1) from the minimal and maximal extensions h' and H' of b , respectively. Thus, since (g h') H' g , it follows that (g h) H differs by only n O(1) from g . Because the normalization of g differs from f by only O(1 n) , we see that the normalization of (g h) H differs from f by at most O(1 n) , as desired. proof lem supnorm The minimal increasing path length ( with the path not restricted to lie in any particular region ) between two lattice points is always within 1 of twice the sup norm distance between them. lem proof Start at any lattice point and work outwards, labelling the other points with the lengths of the shortest increasing paths to them. It is easy to prove by induction that on the square at sup norm distance n from the starting point, on two opposite sides the lengths alternate between 2n and 2n 1 , and on the other two sides they alternate between 2n and 2n-1 . proof We will prove the next three results in more generality than we need, in order to make the precise hypotheses clear. We will need to apply them only to regions whose boundaries closely approximate a k k square or an equilateral triangle of side length k , such that the heights along the boundary are nearly planar (in particular are fit by a plane to within k for some fixed 0 ). prop planar Let 0 . Suppose R is a simply-connected lattice region of diameter at most n ( i.e., every lattice point in R is connected to every other by a path within R of length n or less ) , such that the heights on the boundary of R are fit to within n by a plane with tilt (s,t) satisfying s t 2 . Then the average height function is given to within O(n) by that plane ( if n is sufficiently large ) . prop Here O(n) '' simply denotes a quantity bounded in absolute value by a fixed constant times n , for sufficiently large n . proof Consider a large torus, with edge weights a,b,c,d chosen to give tilt (s,t) and satisfying ab cd (see Subsection entasfnofslope ), and view R as being contained in the torus. We will look at random free tilings of R generated according to the probability distribution on weighted tilings of the torus. If we fix the height of one point on the boundary of R , then the average height function for these tilings is given (exactly) by a linear function of the two position coordinates, so that its graph is a plane. If we select a random tiling from this distribution, then with probability differing from 1 by an exponentially small amount, the heights along the border of the patch stay within n of the plane, by Proposition 22 of CEP . (It is not hard to check that all the large deviation results from Section 6 of CEP , such as Propositions 20 and 22, apply to random tilings generated this way, not just from the uniform measures on tilings of finite regions. All that matters is conditional uniformity, in the sense that two tilings agreeing everywhere except on a subregion are equally likely to occur; conditional uniformity follows from ab cd .) Consider all possible boundary height functions on R that stay within n of a plane. The average of the corresponding average height functions, weighted by how likely they are to occur in the weighted probability distribution, is within o(1) of the plane (for n large). By Proposition 20 of CEP , all boundary height functions that stay within n of the plane have average height functions within O(n) of each other. Since the average is within o(1) of the plane, each must be within O(n) of it. This completes the proof. proof lem extremal Let 0 . Suppose R is a horizontally and vertically convex lattice region of area A with at most n rows and columns, such that n A . Assume that the boundary heights are fit to within A n by a plane with tilt (s,t) . If s t 2- , then the entropy of R is O(1 ) , if A is sufficiently large. lem proof Without loss of generality, suppose that s and t are positive. Consider any vertical line through R on which the x -coordinate is integral. If the segment that lies within R has length k , then for every tiling of R , the difference between the number of north-going dominos bisected by the line and the number of south-going dominos bisected by it is tk 4 O(A n) O(1) , as one can see by considering the total height change along the segment. Similarly, on a horizontal segment of length k , the number of west-going dominos bisected minus the number of east-going ones bisected is sk 4 O(A n) O(1) . If we add up all of these quantities, the error term becomes O(A) O(n), which is O(A) since n A . The total number of dominos in any tiling of R is A 2 . By the results of the previous paragraph, this quantity is also twice the number of east-going or south-going dominos plus (s 4 t 4)A O(A) . We have (s 4 t 4)A ((2-) 4)A . Hence, the total number of east-going or south-going dominos is O(A) . Every tiling is determined by the locations of its east-going and south-going dominos. (To see why, recall that superimposing the matchings corresponding to two tilings gives a collection of cycles. Any disagreement between the two tilings yields a cycle of length at least 4 , which must contain an east-going or south-going edge that is in one tiling but not the other.) Hence, there are at most equation O(A) O(A) 2 O(A) equation tilings, since there are O(A) possibilities for the number of east-going or south-going dominos, at most equation O(A) equation places to put them, and at most equation 2 O(A) equation ways to choose which are east-going. It is not hard to check that Stirling's formula implies that equation O(A) O((1 ) A). equation Therefore, the entropy of R is O(1 ) . proof prop continuity Fix k 0 . Let 0 , and let R be a horizontally and vertically convex region of area A with at most n rows and columns, such that kn 2 A and n A . Suppose h is a boundary height function on R that is fit to within A n by a fixed plane. Then for A sufficiently large and sufficiently small, the entropy of extensions of h to R is independent of the precise boundary conditions, up to an error of O( 1 2 1 ) . ( The entropy does depend on the tilt of the plane, and this proposition says nothing about whether it depends on the shape of R . Notice also that the dependence on k is hidden within the implicit constant in the big- O term. ) prop proof Let (s,t) be the tilt of the plane. If s t 2- 1 2 , then the conclusion follows from Lemma extremal . Thus, we can assume that the tilt satisfies s t 2- 1 2 . Suppose g is another boundary height function on R , which agrees with the same plane as h , to within O(A n) . We need to show that extensions of g and h have nearly the same entropy. Without loss of generality we can assume that g h , since otherwise we can go from g to g h and from h to g h . Given any extension f of g , let H(f) be the infimum of f and the maximal extension of h , so that H(f) is an extension of h . Similarly, given any extension f of h , let G(f) be the supremum of f and the minimal extension of g , so that G(f) is an extension of g . The maps H and G are not inverses of each other, but they come fairly close to being inverses. Given an extension f of h , H(G(f)) agrees with f at every lattice point except those with heights less than or equal to their heights in the minimal extension of g . By Lemma fournier , the minimal extension of g is within O(A n) of the minimal extension of h , so these points have heights within O(A n) of their minimal heights. Similarly, given an extension f of g , G(H(f)) agrees with f at all lattice points that are not within O(A n) of their maximal heights. By assumption, the tilt (s,t) of the plane satisfies s t 2- 1 2 . Thus, the height difference between any two points in the plane is bounded by 2- 1 2 times the sup norm distance between them. Therefore, Lemma fournier and Lemma supnorm imply that the extreme heights at a point differ from the heights on the plane by at least 1 2 times the sup norm distance to the boundary. Look at all lattice points not within sup norm distance ( 1 2 1 ) A n of the boundary. By Proposition planar , at any such point the average height for extensions of g or h is within O(A n) of the plane (notice that because kn 2 A n 2 , being O(n) is the same as being O(A n) ), and the extreme heights differ from the plane by at least (1 )A n , so the extreme heights differ from the average heights by an amount on the order of (1 )A n , for small. By Proposition 22 of CEP , the probability that any such point will have height within O(A n) of its extreme heights is exponentially small in n . Thus, given any point not within sup norm distance ( 1 2 1 )A n of the boundary, the probability that H G or G H will not be the identity at that point is exponentially small. It follows that with probability nearly 1 , H G and G H are the identity except on at most O(( 1 2 1 )A) lattice points. Thus, the numbers of extensions of g and h differ by at most a factor of equation 4 O(( 1 2 1 )A) , equation so the entropy of extensions of g to R differs from that of extensions of h by O( 1 2 1 ) . proof Proof of the variational principle endpartI thm precise Let 0 . Suppose R is an n n square, with a boundary height function h fit to within n by a plane with tilt (s,t) satisfying s t 2 . Then for n sufficiently large, the entropy of extensions of h to R is equation (s,t) O( 1 2 1 ), equation as is the entropy for free boundary conditions staying within n of the plane ( i.e., the entropy for the set of all free tilings of R whose boundary heights stay within n of the plane ) . thm proof We know from Proposition continuity that the entropy is independent of the precise boundary conditions, but we still need to prove that it equals (s,t) . To do so, we will compare with an n n torus that has edge weights a,b,c,d satisfying ab cd and yielding tilt (s,t) . (We can suppose that s t 2 , since otherwise the result follows from Lemma extremal and Proposition continuity .) The torus is obtained by identifying opposite sides of R , so that tilings of R give tilings of the torus, but not vice versa. Keep in mind that because of the weighted edges in the torus, the probability distribution on its tilings will not be uniform. However, the equation ab cd implies conditional uniformity (as mentioned earlier), so if we fix the behavior on the boundary of R , then the conditional distribution on extensions to the interior will be uniform. In Section edgeprobsection , we will define a set W depending on (s,t) . By Lemma halfin , for sufficiently large even n , either n or n 2 is in W . In Proposition ent1 , we show that the entropy of n n tori with tilt (s,t) converges to (s,t) as n in W . First we will suppose that n W , so that the entropy of the n n torus is (s,t) o(1) . It follows from Proposition 22 of CEP that with probability exponentially close to 1 , in a random tiling of the torus, the heights on the boundary of the square will be fit to within n by the average height function, which is a linear function with tilt (s,t) . The number of toroidal boundary conditions is exponential in n , and by Proposition continuity each has about the same entropy (except ones that are not nearly planar, but they are very unlikely to appear). Lemma entropy tells us that the entropy of the torus equals the average of the entropies for the different boundary conditions, plus a negligible quantity for large n (since we have normalized by dividing by half of the area n 2 ). Because all the nearly planar boundary conditions have the same entropy, the torus must as well, so since we know it has entropy (s,t) o(1) as n , each of the nearly planar boundary conditions must have entropy (s,t) O( 1 2 1 ) . (Since is fixed as n , we can absorb the o(1) into the big- O term.) Now it is easy to deal with the case of n W . If n W , then n 2 W and n-2 W . The entropies for tilings of (n-2) (n-2) , n n , and (n 2) (n 2) squares with nearly planar boundary conditions are nearly the same. (To prove that, embed an (n-2) (n-2) square into an n n one, and an n n one into an (n 2) (n 2) one, extending the boundary conditions arbitrarily. Then the number of tilings increases with each embedding, so the entropy of the n n square is caught between the other two, to within a 1 o(1) factor coming from the differing areas. By Proposition continuity , this result holds for all nearly planar boundary conditions.) Note that the same argument as above goes back from squares to the torus, thus proving entropy convergence for all n (not just n W ). The claim about free boundary conditions follows easily (since the number of boundary conditions that stay within n of the plane is only exponential in n , and all of them have about the same entropy). proof cor triangle Theorem precise holds if R is an equilateral triangle of side length n , instead of a square. cor proof It is not hard to check that equilateral triangles satisfy the hypotheses of Proposition planar , Lemma extremal , and Proposition continuity . Thus, we just need to deal with the case of an equilateral triangle whose boundary heights are within O(1) of being planar. To prove that the entropy of the triangle is at least what we expect, tile the triangle with smaller squares, such that their boundary heights are within O(1) of being planar. (Of course, those near the edges will stick out over the boundary, but if n is large enough, we can make the squares small enough compared to the triangle that only an fraction of the squares will cross the boundary.) Except for an error of O() from the squares that cross the boundary, the entropy of the triangle is at least the entropy of the squares, which is what we wanted. An analogous argument (involving tiling a square with triangles) shows that the entropy of the triangle is at most what we expect. proof For the next theorem, we will use the same setup as in Proposition hfexist . thm nbhdentropy Let R be the region bounded by a simple closed curve, and let h b be a boundary height function on R . Suppose R is a simply-connected lattice region such that when R is normalized by a factor of 1 n , it approximates R to within , and its normalized boundary heights approximate the region R with boundary heights h b to within ; we assume that the normalization of R lies within R . Given an asymptotic height function h (R ,h b) , the logarithm of the number of tilings of R whose normalized height functions are within O() of h is the area of R times equation (h) o(1) equation as 0 ( for n sufficiently large ) . thm proof Notice that the set of tilings whose normalized height functions are within O() of h is non-empty, by Proposition hfexist . Call the set of such tilings U . Fix 0 . Choose small enough that we can apply Lemma approx to the piecewise linear approximation h to h derived from an -mesh (with approximation tolerance ). Then, as is pointed out after the statement of Lemma approx , (h) (h) o(1) . We will take , and show that the entropy we want to compute is (h) o(1) . We know that h-h on all but at most an fraction of the triangles in the mesh. Those triangles can change the entropy by only O() , so we can ignore them. We can also ignore the O() fraction of the triangles that do not lie within the normalization of R (which change the entropy by O() O() ). We will call the triangles within the normalization of R on which h-h the included triangles, and the others the excluded triangles. Let g be any element of U . The entropy of U is bounded below by the sum over all included triangles of the entropy of g restricted to that triangle, plus the O() contribution from the excluded triangles. It is bounded above by the same sum (including the O() ), but with free boundary conditions on the included triangles (subject to the condition of staying within of f ). We can now apply Corollary triangle . It tells us that each included triangle's contribution to the entropy of U is approximately equal to its contribution towards the entropy of h . It follows that our upper and lower bounds for the entropy of U both equal (h) O() O( 1 2 1 ) . This gives us the desired conclusion. proof Note that Theorem nbhdentropy implies Theorem maintheorem . 3 0in 3 psfile 1 2 Overview of the remaining sections In Section partitionsection we compute the partition function Z n(a,b,c,d) for matchings of the toroidal graph G n 2 2n 2 with 4n 2 vertices. In Subsection limitsection we compute the limit, as n , of Z n 1 (2n 2) . In Section edgeprobsection we compute the limit of the edge-inclusion probabilities for edges of each type, with respect to the measures n , and also a bound on the variance of the number of edges of a fixed type in G n . This computation is only done for n in an infinite subset W . Since the variance is o(n 4) , the measure is concentrating near tilings with the mean number of edges of each type. This fact allows us, in Section entropysection , to compute the limit for nW of the entropies. As explained in the proof of Theorem precise , it follows that the limit for arbitrary n must be the same as the limit for nW . In Section concavesection we show that the entropy is strictly concave as a function of the tilt (s,t) . In Section PDEsection we present the PDE which a C 1 entropy-maximizing Lipschitz function must satisfy (at least in the distributional sense). The partition function partitionsection Let the graph G be the infinite square grid. Define an a -edge to be a horizontal edge whose left vertex has even coordinate sum, a b -edge to be a horizontal edge whose left vertex has odd coordinate sum, a c -edge to be a vertical edge whose lower vertex has even coordinate sum, and a d -edge to be a vertical edge whose lower vertex has odd coordinate sum. Let a,b,c,d be four non-negative real numbers. Weight the a -edges with weight a , the b -edges with weight b , and so on. (For comparison with our earlier, more geometrical terminology, a -edges are north-going, b -edges south-going, c -edges east-going, and d -edges west-going; in other words, points with even coordinate sum correspond to white squares.) We assume without loss of generality that ab , cd , and ac . For n an even positive integer let G n denote the quotient of G by the action of translation by (2n,0) and (0,2n) . Then G n is a graph on the torus and has 4n 2 vertices and 8n 2 edges ( 2n 2 edges of each type). A vertex (a,b) of G (where a,b 0,2n-1 ) will be denoted in what follows not by an ordered pair (a,b) but rather by an ordered triple (x,y,t) , where x a 2 , y b 2 , and t 1,2,3, or 4 corresponding to (a,b) being congruent to (0,0),(1,1),(1,0), or (0,1) modulo 2 . The partition function Z n is by definition the sum, over all perfect matchings of G n , of the product of the edge weights in the matching: equation Z n(a,b,c,d) matchings a N a b N b c N c d N d , equation where N a is the number of matched a -edges, etc. There is a natural probability measure n n(a,b,c,d) on the set of all matchings, where the probability of a matching that has N a a -edges, etc., is equation a N a b N b c N c d N d Z n. equation The physical interpretation of is the following. Let E a , E b , E c , and E d be energies associated with dimers'' on an a -edge, b -edge, c -edge, and d -edge, respectively. Define weights a,b,c,d (also called activities'' in the statistical mechanics literature) by equation a e -E a , b e -E b , c e -E c , d e -E d , equation where is a constant depending on the temperature. Then is the Boltzmann measure associated to these energies; specifically, the probability of a configuration of energy E is proportional to e -E , where E is the sum of the energies of the individual dimers. In what follows we will take a,b,c,d as the fundamental quantities and will not deal with or temperature as such. Our concern will be with the situation in which n goes to infinity with the field-parameters a,b,c,d fixed: the so-called thermodynamic limit''. Even though for each finite n , Z n(a,b,c,d) is a smooth function (indeed a polynomial function) of the 4-tuple (a,b,c,d) , the limit Z(a,b,c,d) n Z n 1 (2n 2) and other thermodynamic quantities need not be. We will see that Z(a,b,c,d) is C 1 everywhere but not C 2 in the vicinity of the locus of (a b c-d)(a b-c d)(a-b c d)(-a b c d) 0 . It's worth pointing out that the 4-parameter field determined by (a,b,c,d) actually only has two meaningful degrees of freedom: one degree of freedom drops out because of the imposition of the constraint ab cd , and the other drops out by virtue of the fact that multiplying a,b,c,d by a constant has no effect on any of the quantities of interest. These two degrees of freedom correspond to the two degrees of freedom associated with (s,t) (the tilt). Determinants detsection The goal of this section is to compute Z n : we use Proposition Kas and equations ( detA )--( detA4 ) below. Given an enumeration of the 4n 2 vertices of G n , we define the weighted adjacency matrix of G n as the 4n 24n 2 matrix whose i,j th entry is the weight of the edge connecting vertex i to vertex j (interpreted as 0 if there is no such edge). Define a matrix A 1 by multiplying the weights on vertical edges of the weighted adjacency matrix of G n by i -1 . The matrix A 2 is obtained from A 1 by multiplying by -1 the weights on the vertical edges from vertices (j,n-1,4) to (j,0,1) and edges from vertices (j,n-1,2) to (j,0,3) for all j 0,n-1 . The matrix A 3 is obtained from A 1 by multiplying by -1 the weights on horizontal edges from vertices (n-1,k,2) to (0,k,4) and horizontal edges from vertices (n-1,k,3) to (0,k,1) , for k 0,n-1 . The matrix A 4 is obtained from A 1 by multiplying by -1 the weights on both these sets of edges. By the method of Kasteleyn Kast1 , we have the following proposition. prop Kas For i 1,2,3,4 the quantities A i are non-negative, satisfy Z n A i , and satisfy equation Z n(a,b,c,d ) 12(- A 1 A 2 A 3 A 4 ). equation prop Let V denote the set of vertices of G n . The matrix A 1 operates on V in the following way: for f V and wV , equation (A 1f) w vV a vw f v. equation Let T (j,k) be the linear operator on V corresponding to the translation by (j,k) on G n . The operators T (2,0) and T (0,2) commute with each other and with A 1 . The eigenvalues of T (2,0) are e 2ij n for integers j 0,n-1 . The eigenspace of T (2,0) for eigenvalue e 2ij n is 4n -dimensional: a vector v in this eigenspace is determined by its coordinates in two consecutive columns of G n . Similarly T (0,2) has eigenvalues e 2ik n and a vector in the e 2ik n -eigenspace is determined by its coordinates in two consecutive rows of G n . The intersection of a maximal eigenspace of T (2,0) and one of T (0,2) is 4 -dimensional: a vector in the intersection is determined by its coordinates on a 22 square of vertices (as in Figure evects , vertices v 1,v 2,v 3,v 4 ). Let z e 2i n . For (j,k) 0,n-1 2 and s 1,2,3,4 define a vector e j,k (s) by equation e j,k (s) (x,y,t) cases z jx ky if t s , and 0 otherwise . cases equation The e (s) j,k are in the intersection of the z j -eigenspace of T (2,0) and the z k -eigenspace of T (0,2) . Let S be the 4n 24n 2 matrix whose columns are these eigenvectors: equation S (e 0,0 (1) ,,e 0,0 (4) , e 1,0 (1) ,,e 1,0 (4) ,,e n-1,n-1 (4) ). equation Note that S satisfies S -1 1 n 2 S t , where t denotes the transpose. Because both T (2,0) and T (0,2) commute with A 1 , S -1 A 1S has the block-diagonal form equation SAS S -1 A 1S ( array ccccc B 0,0 0 0 B 1,0 0 0 0 0 0 0 B n-1,n-1 array ), equation with 44 blocks B j,k for j,k 0,n-1 . The block B j,k is the action of A 1 on the intersection of the z j -eigenspace of T (2,0) and the z k -eigenspace of T (0,2) . We have equation Bjk B j,k ( array cccc 0 0 a bz -j i(c dz -k ) 0 0 i(d cz k) b az j a bz j i(d cz -k ) 0 0 i(c dz k) b az -j 0 0 array ) equation for the ordering e j,k (1) ,,e j,k (4) (see Figure evects ). figure htbp equation z kv 1 - d d z kv 3 - d c z -j v 2 - r a v 4 - r b - d c v 2 - d d - l a z jv 4 z -j v 3 - r b v 1 - r a v 3 - l b z jv 1 z -k v 4 - u d z -k v 2 - u c equation evects A vector v in the intersection of the z j -eigenspace of T (2,0) and the z k -eigenspace of T (0,2) . figure Since the upper right and lower left 2 -by- 2 subdeterminants of B j,k are complex conjugates, the determinant of A 1 is equation A 1 j 0 n-1 k 0 n-1 2ab a 2z -j b 2z j 2cd c 2z -k d 2z k 2, equation so equation detA A 1 j 0 n-1 k 0 n-1 (a bz j) 2 z j (c dz k) 2 z k 2. equation The matrices A 2,A 3,A 4 have similar determinants. For example for A 2 , let R (0,2) act on V by translation by (0,2) followed by negation of the coordinates in the first two rows, that is, vertices (j,0,s) for 0jn-1 and s 1,2,3,4 . Then R (0,2) n -I , and R (0,2) and T (2,0) commute with A 2 . One thus finds that equation detA2 (A 2) j 0 n-1 k 0 n-1 (a be 2ij n ) 2 e 2ij n (c de i(2k 1) n ) 2 e i(2k 1) n 2, equation and similarly eqnarray detA3 (A 3) j 0 n-1 k 0 n-1 (a be i(2j 1) n ) 2 e i(2j 1) n (c de 2ik n ) 2 e 2ik n 2, detA4 (A 4) j 0 n-1 k 0 n-1 (a be i(2j 1) n ) 2 e i(2j 1) n (c de i(2k 1) n ) 2 e i(2k 1) n 2. eqnarray We show that, with the exception of A 1 , the square roots of A i are polynomials in a,b,c,d . Define eqnarray P2 P 2 j 0 n-1 k 0 n-1 ( (a be 2ij n ) 2 e 2ij n (c de i(2k 1) n ) 2 e i(2k 1) n ), P3 P 3 j 0 n-1 k 0 n-1 ( (a be i(2j 1) n ) 2 e i(2j 1) n (c de 2ik n ) 2 e 2ik n ), P4 P 4 j 0 n-1 k 0 n-1 ( (a be i(2j 1) n ) 2 e i(2j 1) n (c de i(2k 1) n ) 2 e i(2k 1) n ). eqnarray The functions P 2,P 3 , and P 4 are polynomials in a,b,c,d . Note that in ( P2 ), the involution (j,k)(-j,-k 1) maps each term to its complex conjugate (and this involution has no fixed point). Therefore P 2 A 2 . Similarly P 3 A 3 and P 4 A 4 . Define equation P1def P 1 j 0 n-1 k 0 n-1 ( (a be 2ij n ) 2 e 2ij n (c de 2ik n ) 2 e 2ik n ), equation where the sign holds when a b c d and the - sign holds when a b c d . The involution (j,k)(-j,-k) maps each term in ( P1def ) to its complex conjugate; the product of the terms which are fixed under this involution is ((a b) 2 (c d) 2)(-(a-b) 2 (c d) 2)((a b) 2-(c-d) 2)((a-b) 2 (c-d) 2). This product is positive or negative depending on whether a b c d or a b c d . Thus on the domain a b c d , P 1 is a polynomial (taking non-negative values), and on the domain a b c d , P 1 is the negative of this polynomial (and this polynomial also takes non-negative values). Note that P 1 A 1 . (The quantity P 1 is defined similarly for all positive values of a,b,c,d ; when one of a,b,c,d is greater than the sum of the other three, the - sign is used, otherwise the sign is used.) From Proposition Kas we have equation Ps Z n 12(-P 1 P 2 P 3 P 4). equation Eigenvalues and roots roots Here we study in greater detail the function equation q(z,w) (a bz) 2 z (c dw) 2 w, equation where a,b,c,d are non-negative reals. Let equation r(z) cd (a bz) 2 2z . equation Then equation q(z,w) c 2 2r(z)w d 2 w 2 w (dw-c)(dw-c) w, equation where equation alphabeta (z),(z) -r(z) cd ( r(z) cd ) 2-1 . equation We will choose (z) to be the larger root (in modulus). Let (z) (a bz) 2 z and (w) -(c dw) 2 w , so that (z)-(w) q(z,w) . The critical points of (z),(w) are respectively z a b and w c d . Recall our assumption that ab and cd . lem eta The map maps the punctured unit disk z: , z 1, z0 injectively onto the exterior of the ellipse E 1 whose major axis has endpoints -(a-b) 2,(a b) 2 and whose minor axis has endpoints 2ab i(a 2-b 2),2ab-i(a 2-b 2) . The ellipse E 1 has foci 0 and 4ab and center 2ab . In the case a b , this ellipse degenerates to the line segment 0,4ab . The map maps the punctured unit disk w 1, w0 injectively onto the exterior of the ellipse E 2 with major axis -(c d) 2,(c-d) 2 , minor axis -2cdi(c 2-d 2) , and foci at 0 and -4cd . lem See Figures ellipses and missellips . proof Recall that ab by hypothesis. Suppose a b . We have (z) a 2 z b 2z 2ab . For z 1 write z t i . Then eqnarray (z) a 2(-i) b 2( i) 2ab 2ab (a 2 b 2)- i(a 2-b 2). eqnarray Thus as z runs counterclockwise around S 1 , (z) runs clockwise around the ellipse E 1 . Since has no critical points in the unit disk, it is injective on the unit disk, mapping it to the exterior of E 1 (with 0 mapping to ). When a b , E 1 degenerates to a segment; still, maps the punctured open unit disk injectively onto the exterior of the segment. The case of is similar. proof figure scale .9 jams355el-fig-4 ellipses The ellipses E 1 and E 2 in the case a b c d . figure figure scale .9 jams355el-fig-5 missellips The ellipses E 1 and E 2 in the case a b c d . figure The function q can also be written equation mod1 q ( b i( d )) ( b-i( d )). equation The coefficients of a,b,c,d in either term of ( mod1 ) have modulus 1 when z 1 w . Suppose a b c d . Then q(e i ,e i ) 0 for all and , so E 1E 2 , and hence E 2 is contained inside the bounded region delimited by E 1 . See Figure missellips . For each fixed z satisfying z 1 , (z) is on E 1 and so by Lemma eta the roots w of 0 q(z,w) (z)-(w) are situated one outside and one inside the unit disk. Since is the larger of and , (z) c d is the root inside the disk. On the other hand if a b c d but not both a b and c d , then (a-b) 2 (c d) 2 (and (c-d) 2 (a b) 2 by the hypothesis that a is the largest of a,b,c,d ) so the two ellipses intersect as in Figure ellipses (because the places where they cross the x -axis are interlaced). Let (z 0,w 0) (e i 0 ,e i 0 ) and ( z 0 , w 0 ) (e -i 0 ,e -i 0 ) be the roots (satisfying z w 1 ) of q(z,w) 0 , where 0(0,) (the angle cannot be 0 or since, as we noted, the ellipses are interlaced on the x -axis). Again by Lemma eta , for each z with z 1 , exactly one of the roots w c d,c d is inside the (closed) unit disk when - 0 0 , and for - 0, 0 , both roots are outside. We will again take (z)c d to be the smaller root. In the case a b c d , the two ellipses are tangent, and their single intersection point is at z -1,w 1 , that is, ( 0, 0) (,0) . In the case when both a b and c d , the two degenerate ellipses intersect only when z w -1 , so that the single intersection point is at ( 0, 0) (,) . An important fact about the above four cases is that (z)c d 1 always, and (z)c d 1 unless - 0, 0 . Going back to the case ab c d , at the point ( 0, 0) , the four quantities equation z 0 , b z 0 , ic w 0 ,id w 0 equation sum to zero by ( mod1 ), assuming we choose the correct signs for z 0 and w 0 . They therefore form the edge vectors of a quadrilateral. When taken in the order as in Figure cycquad , the quadrilateral is in fact cyclic since opposite angles sum to : The (interior) angle between sides a z 0 and ic w 0 is equation -( z 0 w 0 ic ) 3 2- w 0 z 0 , equation and the (interior) angle between sides b z 0 and id w 0 is equation -( b z 0 id w 0 ) 3 2- z 0 w 0 . equation Summing these two gives angle (modulo 2 , of course). figure scale .9 jams355el-fig-6 cycquad The cyclic quadrilateral formed at ( 0, 0) . figure The limit of the partition functions limitsection thm The limit equation Z Z(a,b,c,d) n Z n 1 (2n 2) equation exists, and equation Z 1 8 2 0 2 0 2 (a be i ) 2 e i (c de i ) 2 e i ,d , d. equation thm Kasteleyn Kast1 computed Z in the special case a b and c d . proof From Proposition Kas , for positive reals a,b,c,d we have P Z for all . Combined with Proposition Kas , this gives equation P (a,b,c,d) Z n(a,b,c,d)32 P (a,b,c,d). equation This implies that equation n Z n 1 (2n 2) n ( P ) 1 (2n 2) equation (assuming these limits exist). Here note that the j for which P j P may also depend on n . Let I denote the integral in the statement of the theorem. We show that ( 2n 2 ) -1 P converges to I . Let equation F(,) q(e i ,e i ) a 2e -i 2ab b 2e i c 2e -i 2cd d 2e i . equation We then have equation P 1 j,k F( 2j n, 2k n), equation and similar expressions for P 2,P 3,P 4 (see ( P2 )--( P1def )). The quantity n -2 (P ) is a Riemann sum for the integral of F . The function F is continuous, hence Riemann integrable, on the complement of a small neighborhood of its two (possible) singularities ( 0, 0) and (- 0,- 0) , defined in Subsection roots . We break the proof into four cases: the case a b c d , the case a b c d but not both a b and c d , the case a b and c d , and finally the case a b c d . If a b c d , then F is continuous everywhere on 0,2 0,2 and so the Riemann sums converge to I . In the case a b c d but not both a b and c d , there are two singularities. It suffices to show that the Riemann sums n -2 P for each are small on a small neighborhood of the singularities. This is not quite true since for some the product P may have a factor in which (,) lands close to a singularity; as a consequence this P may be very small. However, we will show that this can happen for at most one of the four products P . For each term in the four products ( P2 )--( P1def ), (,) is of the form ( j' n, k' n) for integers j',k' . Furthermore for each pair of integers (j',k'), exactly one of the four products has a term with (,) ( j' n, k' n) . Thus at most one of the four products has a term closer than ( 2n ) to the singularity ( 0, 0) . The same P will have the term closest to the other singularity (- 0,- 0) . Fix a small constant 0 and let U 0,2 0,2 be the -neighborhood of a singularity. In U we use the Taylor expansion equation taylor F(,) F ( 0, 0)(- 0) F ( 0, 0) (- 0) O((- 0) 2, (- 0) 2), equation and by Lemma notpar below, the ratio of F ( 0, 0) and F ( 0, 0) is not real. The sum of those terms in n -2 P for which (,)U is 1 n 2 (,)U F(,) 1 n 2 (,)U C 1(- 0) C 2(- 0) O( 3), for constants C 1 F ( 0, 0),C 2 F ( 0, 0) with C 1 C 2 . (We used here the fact that (x O(x 2)) (x) O(x) for small x ; we then could take the big-O term out of the summation because of the 1 n 2 factor.) To bound this sum, note that C 1x C 2y C 3 x C 4 y C 5 x iy for some positive constants C 3,C 4,C 5 since C 1 C 2 . We use polar coordinates around the singularity. In the annulus around ( 0, 0) of inner radius K ( 2n ) and outer radius (K 1) (2n) , there are at most constant K points (,) which contribute to the sum, and each such point contributes n -2 ( constant K ( 2n )) to the sum. Therefore the sum on U (for those P j without terms within (2n) of the singularity) is bounded by equation constant n 2 1Kn K Kn O( 2). equation For the P j which does have a term closer than (2n) to the singularity, the above calculations give an upper bound on ( 2n 2 ) -1 P j . (Including in the factor close to the singularity only decreases the product.) Thus we have shown that ( 2n 2 ) -1 P converges to I . This proves the convergence of ( 2n 2 ) -1 Z to I . In the case a b and c d we have only one singularity (,) , and P 1 0 . For the other P 's, when (,) is close to (,) we have eqnarray F(,) a 2(2 2 ) c 2(2 2 ) a 2(-) 2 c 2(-) 2 O((-) 4,(-) 4). eqnarray An argument similar to the previous case holds: on U we have equation 1 n 2 (,)U F 1 n 2 (,)U a 2(-) 2 c 2(-) 2 O( 3). equation Now a 2x 2 c 2y 2 C(x 2 y 2) . Summing over annuli as before gives the bound. Finally, suppose a b c d . For any 0 we have equation Z n(a ,b,c,d)Z n(a,b,c,d)Z n(a-,b,c,d) equation (recall that coefficients of the polynomial Z n are non-negative). For fixed , the limits equation n Z n(a,b,c,d) 1 (2n 2) equation both exist, and (as we shall see in ( Zformula ) below) converge to the same value as 0 . Thus equation n Z n(a,b,c,d) 1 (2n 2) equation exists and converges to this same value. This completes the proof. proof The integral in Theorem can be written more usefully as follows. As before set equation r r(z) cd (a bz) 2 2z . equation We can evaluate the first integral (the integral with respect to ) in I as follows. We have (with w e i ) equation 1 2 0 2 c 2 2rw d 2w 2 w d 1 2 0 2 (dw-c)(dw-c) :d, equation where , are chosen as in Subsection roots . Using the identity equation 1 2 0 2 t s e i :d cases t if t s , and s if s t cases equation (note that the logarithmic singularity in the case s t makes no contribution to the integral), we find (with z e i ) align 4Z c(z) d d :d c(z) d c(z) :d c(z) d d :d c(z) d c(z) :d. align From Lemma eta we know that c (z) d for all z on the unit circle. If ab c d , then c(z) d for all z 1 and so equation 4Z 2d 2c - (z) :d equation or equation a bcd Z 12d 12c 1 4 - (z) :d. equation If a b c d , then recall c(z) d if and only if (- 0, 0) . Thus we have equation 4Z - 0 0 d :d 0 2- 0 c(z) :d - c(z) :d, equation which gives (using 1 ) equation Zformula Z 0 2 d (1- 0 2 )c 1 4 - 0 0 (z) :d. equation Comparing ( a bcd ) and ( Zformula ), we see that ( Zformula ) holds for all a,b,c,d , as long as we define 0 when ab c d . lem notpar For the function F of , the ratio F ( 0, 0) F ( 0, 0) unless a b c d or both a b and c d . lem proof We must show equation quotient -c 2e -i 0 d 2e i 0 -a 2e -i 0 b 2e i 0 . equation For this proof only, let e i and e i be square roots of e i 0 ,e i 0 , respectively, with signs chosen so that equation linear ae -i be i i(ce -i de i ) 0 equation (cf. ( mod1 )). We can then factor ( quotient ) as equation (de i -ce -i ) be i -ae -i (de i ce -i ) be i ae -i , equation and the second quotient is i . That is, we are left to show that equation Im ( (d-c)i -(d c) (b-a) i(b a) )0 . equation Separating the real and imaginary parts of ( linear ) we have eqnarray (a b) (c-d) 0, and (-a b) (c d) 0. eqnarray Solving these for , and plugging in to the above gives equation Im ( (d-c)i a-b c d -(d c) a b d-c (b-a) i(b a) ) equation which is zero only if the real and imaginary parts of the numerator and denominator are in proportion: either 0 or (clearing denominators) equation (d-c) 2(a-b) 2 (d c) 2(a b) 2. equation Neither of these is possible (recall that 0 only when a b and c d ), so the proof is complete. proof The edge-inclusion probabilities edgeprobsection Let n denote the measure on matchings of G n , where each matching has weight which is the product of the edge weights of its matched edges. The expected number of a -edges occurring in a n -random matching is simply equation (N a) a Z n Z n a equation (this follows from the definition of n ). From ( Ps ), the probability p a(n) of a particular a -edge occurring in a n -randomly chosen matching is therefore equation bigfrac p a(n) 2n 2 a (-P 1 P 2 P 3 P 4) -P 1 P 2 P 3 P 4 . equation From ( P1def ) we obtain (assuming P 1(a,b,c,d)0 ) equation ddaA a P 1 P 1 j,k 2(b ae -2ij n ) (a be 2ij n ) 2e -2ij n (c de 2ik n ) 2e -2ik n . equation Note that this holds for all values of a,b,c,d for which P 1(a,b,c,d)0 , independently of whether or not a b c d . Similar expressions hold for P 2,P 3,P 4 . In what follows we may no longer assume that a is the largest of b,c,d since we are computing a non-symmetric function p a . Recall that the quadruple (a,b,c,d) determines two possible singularities ( 0, 0) of the function F of (). We will define a set W , depending on a,b,c,d , on which our remaining convergence arguments work. If one of a,b,c,d is greater than the sum of the others, or if both a b and c d , take W . If one of a,b,c,d equals the sum of the others, we define W below in the proof of Proposition pamn . In the remaining case, there are two distinct singularities ( 0, 0) , where 0(0,) . If 0 2, let W be the set of n for which 0 is not well approximated by rationals of denominator n , in the following sense: for all integers j we have equation 0- j n 1 n 3 2 . equation If 0 2 , define W as above using 0 instead: note that 0 and 0 cannot both equal 2 , for F( 2, 2) -i(a bi) 2-i(c di) 2 cannot be zero (its real part is 2ab 2cd ). lem halfin When 0(0,) , for any sufficiently large even n , one of n,n 2 is in W . lem proof Without loss of generality 02 . Suppose equation 0- j n 1 n 3 2 equation and equation 0- j' n 2 1 (n 2) 3 2 . equation Then j' must be equal to one of j,j 1 , or j 2 ; but then equation j n- j' n 2 2j, 2j-n , 2j-2n n(n 2) constant n, equation a contradiction. proof prop pamn If none of a,b,c,d equals the sum of the others, then for n tending to in W , the edge-inclusion probability p a(n) converges to equation sumintegral p a 4 2 0 2 0 2 (b ae -i ) :d :d (a be i ) 2 e -i (c de i ) 2e -i . equation If one of a,b,c,d equals the sum of the other three, then p a(n) converges to the above integral along a subsequence W containing at least one of each pair n,n 2 with n even. prop proof We will show that the sum equation psum 2n 2 a P 1 2n 2 , 2(b ae -i ) F(,) equation converges to the desired integral ( sumintegral ), where the sum is over (,) ( 2j n, 2k n) . Similar arguments hold for P 2,P 3 , and P 4 . The value ( bigfrac ) is then a weighted average of these sums, with weights P ( 2Z n ) . Since Z nP 0 (see Proposition Kas ), the weights are bounded in absolute value (less than 1 2 ) and sum to 1 . Therefore the weighted average also converges to ( sumintegral ). We separate the proof into four cases. In the first case, where one of a,b,c,d is greater than the sum of the other three, there are no singularities ( F(,) is never zero), so the summand is a continuous function on 0,2 2 . Therefore ( psum ) converges to the integral in ( sumintegral ). For the second case, suppose each of a,b,c,d is strictly less than the sum of the others, but we are not in the case where both a b and c d . Since nW , none of the four P can have a term with (,) within n -3 2 of a singularity. We claim that the sum ( psum ) converges to the integral ( sumintegral ). This is proved in the same manner as in Theorem : one only needs to check that the contribution on a small neighborhood of the singularities is small. As before, let U be a -neighborhood of a singularity. Ignore for a moment the single term closest to the singularity. Lemma notpar and ( taylor ) give us the estimate equation est 2n 2 (,) U 2(b ae -i ) F(,) equation equation 2a 2n 2 U b ae -i C 1(- 0) C 2(- 0) O((- 0) 2, (- 0) 2) equation for constants C 1 F ( 0, 0),C 2 F ( 0, 0) where C 1 C 2 . Now C 1(- 0) C 2(- 0) C 3 (- 0) i(- 0) C 4 - 0 , - 0 for some positive constants C 3,C 4 . Using 1 ( x O(x 2) ) 1 x O(1) for small x , where x C 3 (- 0) i(- 0) , we get the bound equation 2a 2n 2 U b ae -i C 1(- 0) C 2(- 0) O(1) . equation Taking the O(1) term out of the summation turns it into a O( 2) . Summing over annuli concentric about the singularity, we may bound the left-hand side of ( est ) by equation O( 2) constant n 2 1Kn K O(). equation The single term closest to the singularity contributes a negligible amount equation constant n 2 n 3 2 . equation In the third case, when both a b and c d , we have ( 0, 0) (,) . Then P 1 0 and the other P are non-zero. Furthermore the pairs (,) appearing in the products for P 2,P 3,P 4 do not come within distance n of the singularity. Since near (a,b,c,d) , P 1 is a polynomial taking non-negative values which is zero at (a,b,c,d) , it must have a double root there. Thus its derivative with respect to a is zero at (a,b,c,d) . We can therefore remove P 1 from the expression ( bigfrac ) and just deal with the remaining three P . When a b and c d , ( est ) becomes eqnarray 2a 2n 2 U a(1 e -i ) a 2(2 2) c 2(2 2) 2a 2 2n 2 U i(-) a 2(-) 2 c 2(-) 2 O( 3), eqnarray where we used 1 ( x 2 O(x 4) ) 1 x 2 O(1) . This is O() (sum over annuli as before). Finally we consider the case when one of a,b,c,d is equal to the sum of the other three. Suppose first that a b c d . Note that p a(n)(a,b,c,d) is a monotonic increasing function of a : the expected number of a -edges increases with their relative weight. For each 0 sufficiently small, choose n so that p a(n)(a-,b,c,d)-p a(a-,b,c,d) . Such an n exists because (a-,b,c,d) is in the domain of case two, above. By monotonicity equation p a(n)(a,b,c,d)p a(n)(a-,b,c,d)p a(a-,b,c,d)- equation for this n . Take a sequence of 's tending to 0 . On the corresponding sequence of n 's, p a(n)(a,b,c,d) tends to 1 , which is equal to the value of the integral p a(a,b,c,d) (see Theorem pa below). The set W in this case is obtained from the concatenation of the appropriate subintervals of the sets W that are defined for each quadruple (a-,b,c,d) . Since p a(n)(a,b,c,d)1 we must have that eqnarray p b(n)(a,b,c,d) 0, p c(n)(a,b,c,d) 0, p d(n)(a,b,c,d) 0. eqnarray This (and symmetry) takes care of the remaining cases. proof A rather lengthy calculation yields the following result: thm pa If ab c d , then p a 1 . If one of b,c,d is larger than the sum of the other three of a,b,c,d , then p a 0 . Otherwise, let Q be a cyclic quadrilateral with edge lengths a,c,b,d in cyclic order. Then p a is 1 (2) times the angle of the arc cut off by the edge a of Q . That is, equation paeqn p a 1 sin -1 ( a (a b c-d)(a b-c d)(a-b c d)(-a b c d) 2 (ab cd)(ac bd)(ad bc) ). equation ( Here the branch of the arcsine that is needed is the one given by the preceding geometrical condition. ) thm The proof is in Section paproof . Note that by this theorem, in the case of non-extremal tilt, the 4 -tuple equation ((p a),(p b), (p c),(p d)) equation is proportional to (a,b,c,d) . Since a constant of proportionality has no effect on the measures n , we can assume that a (p a) , etc. We also note a simple relation between the singularity ( 0, 0) and the edge-inclusion probabilities (hereafter simply called edge probabilities), which will be useful later. From Figure cycquad and Theorem pa , when ab and cd we find eqnarray thetapcpd 0 -(p c-p d), and 0 -(p a-p b). eqnarray We now bound the variance of N a , the number of a -edges in a matching. prop sigmasquared For all a,b,c,d and for nW , 2(N a) o(n 4) . varbd prop proof The graph G n has 2n 2 a -edges. For k 1,2n 2 let q k be the 0,1 -valued random variable indicating the presence of the k th a -edge in a random matching. Then N a q 1 q 2n 2 , and so equation var 2(N a) k 2(q k) k ((q kq )- (q k)(q )). equation We have (q k) p a(n) and so 2(q k) p a(n)-p a(n) 2 . In the case when one of a,b,c,d is greater than the sum of the others, we know from Proposition pamn that p a converges to 1 or 0 ; as a consequence 2(q k)0 , and the covariances converge to 0 also, so 2(N a) o(n 4) as well. Similarly when one of a,b,c,d equals the sum of the others; then p a converges to 1 or 0 along W , and so 2(N a) is o(n 4) on this same subsequence. The remaining cases require more work. By (a straightforward extension of) Theorem 6 of Kenyon1 , we have equation detavg (q kq ) - (A 1 -1 ) q k,q P 1 (A 2 -1 ) q k,q P 2 (A 3 -1 ) q k,q P 3 (A 4 -1 ) q k,q P 4 -P 1 P 2 P 3 P 4 , equation where equation (A i -1 ) q k,q ( array cc A i -1 (v k,w k) A i -1 (v k,w ) i -1 (v ,w k) A -1 (v ,w ) array ), equation and where v k,w k are the vertices of the edge associated with q k ( v k being the left vertex) and v ,w the vertices of the edge associated with q ( v being the right vertex). The inverses A i -1 are only defined when the corresponding P i are non-zero. When nW tends to the diagonal entries A i -1 (v k,w k), A i -1 (v ,w ) tend to p a (see Proposition pamn ). Writing each 22 determinant in ( detavg ) as the product of the two diagonal entries minus the product of the off-diagonal entries, the diagonal entries can be taken out of the quotient in ( detavg ) and contribute p a 2 o(1) . It remains to estimate the contribution of the off-diagonal entries. We can compute the inverses of the A i as follows. Note that from ( Bjk ) we have equation B j,k -1 ( array cc 0 D 1 2 0 array ), equation where equation D 1 1 (a bz j) 2z -j (c dw k)w -k ( array cc b az -j -i(d cw -k ) -i(c dw k) a bz j array ) equation (we won't need the expression for D 2 ). Recall the definition of the matrix S of (). Let x,y,s be the vector equation x,y,s (j,k,t) cases 1 if (j,k,t) (x,y,s) , and 0 otherwise. cases equation We have equation S -1 ( x,y,s )(j,k,t) cases 1 n 2 e -2ijx n e -2iky n if t s , and 0 otherwise. cases equation From ( SAS ) we therefore find, for example, equation Ainverse split A 1 -1 ((0,0,1),(x,y,3)) 1 n 2 j 0 n-1 k 0 n-1 e -2i(jx ky) n (b ae -2ij n ) (a be 2ij n ) 2e -2ij n (c de 2ik n ) 2e -2ik n split equation and equation split A 1 -1 ((0,0,1),(x,y,4)) 1 n 2 j 0 n-1 k 0 n-1 e -2i(jx ky) n (-i)(d ce -2ik n ) (a be 2ij n ) 2z -2ij n (c de 2ik n ) 2e -2ik n . split equation We also have A 1 -1 ((0,0,1),(x,y,t)) 0 when t 1 or t 2 . Similar expressions hold for inverses of A 2,A 3,A 4 . An argument identical to the proof of Proposition pamn (the only difference is the factor z -(jx ky) , which has modulus 1 ) shows that the parts of the sums ( Ainverse ) over a -neighborhood U of the singularities are O() . We will show that for all y with (1-)n y n (later we will set n -1 4 ) the value A 1 -1 ((0,0,1),(x,y,3)) tends to zero as n in W . Similar results hold for A 2 , A 3 , and A 4 . For simplicity of notation let 0,v denote the vertices (0,0,1) and (x,y,3) . The equation ( Ainverse ) has the form equation A 1 -1 (0,v) 1 n 2 j,k e -2i(jx ky) n G 1(j n,k n), equation where G 1 is a smooth function on the complement of the region U . We already know that the sum over U is O() , so let us replace G 1 by a new function G 2 which agrees with G 1 outside U and is zero on U . We sum by parts over the variable k to get eqnarray A 1 -1 (0,v) 1 n 2 j 0 n-1 k 0 n-1 ( 0 k e -2iy n ) (G 2(,)- G 2(, k 1 n)) 1 n 2 j 0 n-1 ( 0 n-1 e -2iy n )G 2( jn,0) O(). eqnarray Since (1-)n y n , the sum over of the exponentials is equation 1-e -2i(k 1)y n 1-e -2iy n O(1) equation for each k . The difference G 2(j n,k n)-G 2(j n,( k 1 ) n) is bounded by 1 n times the supremum of G 1 y on the complement of U , except that at the points adjacent to the boundary of U , the difference is bounded by the supremum of G 1 near the boundary. One can check (see the proof of Proposition pamn ) that the sup of G 1 y on the complement of U is O( -2 ) , and the supremum of G 1 on the boundary of U is O( -1 ) . Only O(n) pairs (j,k) correspond to points adjacent to the boundary of U , so we have eqnarray A 1 -1 (0,v) 1 n 2 ( j kO( -1 ) O( -2 )1n) 1 n 2 O(n)O( -1 )O( -1 ) 1 n 2 j O( -1 )O( -1 ) O() O(1 2n ) O(1 n ) O(1 n ) O(). eqnarray Choosing n -1 4 , we have A -1 (0,v) O(n -1 4 ) . A similar argument holds in the case where both a b and c d . From ( var ) we have equation 2(N a)n 2O(1) n 4 o(1) edges k i 1 4 A i -1 (v k,w )A i -1 (v ,w k) , equation but A i -1 (v k,w ) and A i -1 (v ,w k) are O(n -1 4 ) as soon as the edges k and are separated by at least n n 3 4 in their y -coordinates. Therefore (using n -1 4 ) equation 2(N a) O(n 2) o(n 4) O(n 2n 2) O(n 4 n -1 2 ) o(n 4). equation Here the term O(n 2n 2) comes from edges whose y coordinate is less than n or greater than (1-)n , and the term O(n 4 n -1 2 ) consists of the remaining pairs of edges. As a consequence 2(N a) o(n 4) . proof One can show from this argument (using the result of Kenyon1 ) that for nW the measures n converge weakly to a measure on the set of matchings on 2 . The result is as follows. Define equation P(2x 1,2y) 1 4 2 0 2 0 2 e -i(x y) (b a e -i ) d d (a be i ) 2e -i (c d e i ) 2e -i equation and equation P(2x,2y 1) 1 4 2 0 2 0 2 e -i(x y) (-i)(d c e -i ) d d (a be i ) 2e -i (c d e i ) 2e -i . equation Also, define a colored configuration of dominos as a configuration of dominos with a checkerboard coloring of the underlying square grid. prop couplingfn As n within W , the probability of finding a certain colored configuration of dominos in an (a,b,c,d) -weighted n n torus converges to w M , where w is the product of the weights of those dominos and M i,j P(v i,j ) , with v i,j 2 the displacement from the i -th white square to the j -th black one. prop Ben Wieland has pointed out that when ab cd , one can write this more simply. Define P'(2x 1,2y) c(a b) x(c d) yP(2x 1,2y) and P'(2x,2y 1) c(a b) x(c d) yP(2x,2y 1) . One can check that if ab cd , then the probability we seek is simply the determinant of the matrix with entries P'(v i,j ) , with no need to multiply by the products of the weights of the included dominos. This formulation makes the conditional uniformity immediately apparent. The entropy entropysection Entropy as a function of edge-inclusion probabilities Since the set of matchings on G n is finite, the entropy of a measure on the set of matchings is simply equation H() M -(M)(M), equation where the sum is over all matchings M and (M) is the probability of M occurring for the measure . The entropy per dimer is by definition equation () 1 2n 2 H() equation (recall that a matching of G n has 2n 2 matched edges). Recall that for real z , L(z) is the Lobachevsky function, defined by lobdef . prop ent1 As n in W , ( n) converges to equation volume (a,b,c,d) 1 (L(p a) L(p b) L(p c) L(p d)), equation where p a,p b,p c,p d are given by paeqn . prop proof Let C(N a,N b,N c,N d) denote the coefficient of equation a N a b N b c N c d N d equation in Z n . As we computed earlier, on the toroidal graph G n the n -probability of an a -edge (resp. b -, c -, d -edge) is given by p a(n) (resp. p b(n) , p c(n) , p d(n) ). The expected number of a -edges is a def (N a) 2n 2p a(n) . Let eqnarray U (N a,N b,N c,N d): N a- a a, N b- b b, N c- c c, N d- d d . eqnarray Let V be the corresponding set of matchings, i.e., those where the corresponding quadruples (N a,N b,N c,N d) are in U . Because the variance is o(n 4) , for all ,' 0 there exists n 0 such that for nW greater than n 0 we have equation U C(N a,N b,N c,N d)a N a b N b c N c d N d (1-')Z n(a,b,c,d). equation Note that if p j are probabilities and p 1 p k ' 1 , then equation -p jp j-'( ' k) equation (since the left-hand side is maximized when the p j 's are equal). Thus for the entropy we may write eqnarray H() - MV (M)(M) - MV (M)(M) - MV (M)(M) - MV (M)( a N a b N b c N c d N d Z n ) O('( ' constant n 2 )) - MV (M)( a a b b c c d d Z n a N a- a d N d- d ), eqnarray but for MV , we have (a N a- a ) (a a ) , and similarly for b,c,d , so H() MV (M) (-(a a b b c c d d ) (1 O()) Z n) O(n 2''). Note also that MV (M)1-' . Letting ,'0 as n we have finally that the limiting entropy per dimer is eqnarray (a,b,c,d) n 1 2n 2 (Z n(a,b,c,d)- (a a b b c c d d )) Z(a,b,c,d)- p a (a)-p b(b)-p c(c)-p d(d). eqnarray Without loss of generality we may assume that ab and cd ; then from ( thetapcpd ) we have 0 -(p c-p d) . Plugging in from ( Zformula ) now gives equation entent (a,b,c,d) 1 4 - 0 0 (z) :d 1-p c-p d 2(cd)-p a(a)-p b(b). equation To prove the equivalence of this formula and ( volume ), we show that they agree when a b c d , and show that their partial derivatives are equal for all a,b,c,d . Formula ( entent ) gives eqnarray (1,1,1,1) 1 4 - (2 (2 ) 2-1 )d 1 4 - 2(( 2) 2( 2) 1 )d. eqnarray This is two times the value of the entropy per site given in formula (17) of Kast1 , as it should be since the entropy (1,1,1,1) as we defined it is the entropy per dimer. Kasteleyn also shows that this value equals 2G , where G is Catalan's constant equation G 1-1 3 2 1 5 2 -1 7 2 . equation From the expansion equation L(x) 12 k 1 (2kx) k 2 equation (see Mil ) we have 2G (4 )L( 4) , so the two formulas agree when a b c d . It remains to compute the derivatives. Equation volume is symmetric under the full symmetry group S 4 , and entent is by definition symmetric under the operations of exchanging under the operations of exchanging a and b , exchanging c and d , and exchanging a,b with c,d . These operations are transitive on a, b, c, d , so it suffices to show equality of the derivatives with respect to a . We have eqnarray a (a,b,c,d) a (Z-p a(a)-p b(b)-p c(c)-p d(d)) -(a) p a a -(b) p b a - (c) p c a -(d) p d a eqnarray (recall that p a (a Z) Z a ). On the other hand when x is one of a,b,c,d we have equation a 1 L(p x) -(2(p x)) p x a . equation Taking a of ( volume ) and recalling that a (p a) , etc., gives equation -(2a) p a a -(2b) p b a - (2c) p c a - (2d) p d a . equation Since equation (2) a (p a p b p c p d) 0, equation the proof is complete. proof As explained in the proof of Theorem precise , the entropy converges for all n , not just for nW : thm ent As n the entropy per edge of matchings on G n of the measure n(a,b,c,d) converges to equation (a,b,c,d) 1 (L(p a) L(p b) L(p c) L(p d)), equation where p a,p b,p c,p d are given by paeqn . thm Intriguingly, this formula can be used to show that when each of a,b,c,d is less than the sum of the others (the only case in which the entropy (s,t) is non-zero), (s,t) is equal to 1 times the volume of a three-dimensional ideal hyperbolic pyramid, whose vertices in the upper-half-space model are the vertex at infinity and the four vertices of the cyclic quadrilateral of Euclidean edge lengths a,c,b,d in cyclic order (otherwise the entropy is 0 ). We have no conceptual explanation for this coincidence. We do not even fully understand why the limiting behavior of the measures n , viewed as a function of a,b,c,d , turns out to be symmetrical in its four arguments. This symmetry is not merely combinatorial, since it emerges only in the limit as the size of the torus goes to infinity. R. Baxter suggests (in personal communication) that it is almost certainly the same symmetry that occurs in the checkerboard Ising model. In that setting it can be proved using the Yang-Baxter relation (see JM,MR ). Entropy as a function of tilt entasfnofslope Given a tilt (s,t) satisfying s t 2 , we claim that there is a unique (up to scale) 4 -tuple of weights a,b,c,d satisfying the conditional uniformity property ab cd and such that the average tilt of a,b,c,d is (s,t) . To determine a,b,c,d , we note that p a,p b,p c,p d are determined by the equations ( fx )--( sinsin ). To solve these equations, note that ( sinsin ) can be written equation ((p a-p b))-((p a p b)) ((p c-p d))-((p c p d)), equation and so equation (s 2)-((p a p b)) (t 2)-(-(p a p b)), equation giving equation 2((p a p b)) (s 2)-(t 2). equation This combined with () gives ( probdef ), where the values of -1 are taken from 0, (to see why, notice that -1 (((s 2) - (t 2)) 2) (p a p b) , which is between 0 and ). Finally we can determine a,b,c,d as functions of the tilt by equation a (p a), b (p b), c (p c), d (p d). equation Concavity of the entropy concavesection thm concavethm The entropy per edge (s,t) is a strictly concave function of s,t over the range s t 2 . thm proof We show that the Hessian (matrix of second derivatives) is negative definite, that is, ss (s,t) 0, tt (s,t) 0 , and ss (s,t) tt (s,t)- st (s,t) 2 0 , for all s,t , except at the four points (s,t) (2,0) or (0,2) . A computation using Theorem ent and equations ( Peqns ) gives eqnarray s(s,t) p a p a s p d p d s -14( (p a) (p b) ). eqnarray A second differentiation yields eqnarray 2 (s,t) s 2 - 32((p a p b))(p a)(p b) ( 2( t 2) (( s 2) ( t 2)) 2 2), eqnarray and this quantity is strictly negative except at the points (s,t) (2,0),(0,2) . A similar calculation holds for tt (s,t) . We have equation st (s,t) - 32 (s 2)(t 2) (p a) (p b)((p a p b)) . equation Finally, multline ss tt - st 2 ( 32((p a p b))(p a)(p b) ) 2 ( 2( t 2) (( s 2) ( t 2)) 2 2) ( 2( s 2) (( s 2) ( t 2)) 2 2). . - 2( s 2) 2( t 2) , multline which is clearly positive. proof Proof of Theorem pa paproof In what follows, we must be careful to distinguish the differential dw from the product dw (we will write the product as wd to avoid confusion). Let z e i , w e i , and r(z) cd (a bz) 2 ( 2z ) as before. Then (see ()) eqnarray p a 4 2 S 1 S 1 w(b a z) c 2 2rw w 2d 2 dw iw dz iz -a 4 2 S 1 (b a z)dz z S 1 dw (wd-c)( wd-c) . eqnarray Recall (see the second-to-last paragraph of Subsection roots ) that when z 1 , (z)c d always, and (z) c d if and only if (- 0, 0) (remember that we defined 0 in the case ab c d ). If (z) c d , then the residue of ((wd-c)(wd-c)) -1 is equation 1 (-)cd . equation Thus we have eqnarray p a -a 4 2 2i cd - 0 0 b a z z(-) ,dz -ai 2cd - 0 0 (bz a) ,dz 2z 2 ( cd -1)( cd 1) , eqnarray and recalling the definition of r and simplifying yields equation p a -ai 2 - 0 0 dz z (a bz) 2 4zcd . equation We don't have to worry about keeping track of the sign of the square root since we know that we want p a0 . In fact we only need to be careful about the sign when we get to ( sign ), below. This integral can be explicitly evaluated, giving equation p a 2 ( a 2 (ab 2cd)z a (a bz) 2 4zcd z) e -i 0 e i 0 . equation This expression can be simplified using the variable (or rather, one of the two variables) w w(z) such that (a bz) 2 z (c wd) 2 w 0 : the expression under the square root is then equation (a bz) 2 4zcd z(- (c wd) 2 w 4cd) -z w(c-wd) 2. equation Plugging this in yields equation split p a 2 ( a 2 z ab 2cd a(c-wd) wz ) e -i 0 e i 0 2 (2cd b z i( d) -2ida ) e -i 0 e i 0 . split sign equation Up until ( sign ), changing the sign of the square root will only change the sign of the integral. In ( sign ), we choose the sign of w z (or, what is the same, the sign of wz ) so that the expression in curly brackets is zero (cf. ( mod1 )). We then have equation p a 2 (-2d(-c ia )) e -i 0 e i 0 . equation Had we chosen the other sign we would have gotten a similar expression with c and d interchanged. But now Figure 2quads (whose lower quadrilateral is a rotation of Figure cycquad , and whose upper quadrilateral is the reflection of the lower across the edge c ) shows that, as runs from - 0 to 0 , the quantity -c ia w z sweeps out an angle of a , the angle of arc cut off by edge a in a cyclic quadrilateral of edge lengths c,a,d,b . Thus p a a ( 2 ) . In the case a b c d , one can similarly show that -c ia w z sweeps out an angle of 2 , and thus p a 1 . figure scale .9 jams355el-fig-7 2quads The integral defining p a . Note that the angle between the two dotted rays at the origin is a . (Recall that a 2p a is the angle of arc cut off by the edge a .) figure Finally, to prove the formula ( paeqn ), recall the well-known formula for the radius r of the circumcircle of a triangle of sides a,b,c : equation circum r 2 a 2b 2c 2 (a b c)(a b-c)(a-b c)(-a b c) . equation From this one can compute, for a cyclic quadrilateral of sides a,c,b,d , the length s of the diagonal having the b and c edges on the same side: equation s 2 (ab cd)(ac bd) ad bc . equation Plugging this value back into the formula ( circum ) with a,b,c replaced with a,d,s gives equation r 2 (ab cd)(ac bd)(ad bc) (a b c-d)(a b-c d)(a-b c d)(-a b c d) , equation from which ( paeqn ) follows by p a -1 sin -1 (a ( 2r )) . The PDE PDEsection Under the assumption that the entropy-maximizing function f is C 2 , the Euler-Lagrange equation for f is equation dx ( s(f x,f y)) dy ( t(f x,f y)) 0, equation where s, t are the partial derivatives with respect to the first and second variable, respectively. This equation holds only at points where the tilt (f x,f y) is non-extremal, i.e., satisfies f x f y 2 ; otherwise, perturbing f will mess up the Lipschitz condition. In case f is only C 1 , this equation still holds in a distributional sense: it is true when integrated against any smooth test function g vanishing on the boundary (and such that f g is 2-Lipschitz for sufficiently small 0 ). In such a case f is called a weak solution GT . We computed in the proof of Theorem concavethm that eqnarray Epq s(s,t) -14( (p a) (p b) ), and t(s,t) -14( (p c) (p d) ). eqnarray Plugging in from ( Peqns ) and simplifying yields the following PDE: thm pdethm At the points where the entropy-maximizing function f is C 2 and has non-extremal tilt, it satisfies the PDE eqnarray (2(1-D 2)- 2( f x 2))f xx 2( f x 2) ( f y 2)f xy (2(1-D 2)- 2( f y 2))f yy 0, eqnarray where D 12((f x 2)-(f y 2)) . thm Conjectures and open problems conjs In this article, p a , p b , p c , and p d were defined in relation to the dimer model on an n n torus, in the thermodynamic limit as n . We have proved no corresponding interpretation of these quantities for the thermodynamic limit of planar regions. However, our results on the asymptotic height function associated with large regions imply that, in a patch of a large region where the associated asymptotic height function f satisfies ( f x , f y ) (s,t) , the local density of a -edges minus the local density of b -edges equals p a-p b , and likewise for the c - and d -edges. To be precise here, one would define local densities as averages over mesoscopic patches within the tiling (which we recall are defined as patches whose absolute size goes to infinity but whose relative size goes to zero). conj denseconj The local densities of a -edges, b -edges, c -edges, and d -edges are given by p a , p b , p c , and p d , respectively, in the thermodynamic limit. conj Here our use of the phrase thermodynamic limit'' carries along with it the supposition that we are dealing with an infinite sequence of ever-larger planar regions whose normalized boundary height functions converge to some particular boundary asymptotic height function, and that the subregions we are studying stay away from the boundary by at least some mesoscopic distance. The local density of a -edges is just the average of the inclusion probabilities of all the a -edges within a mesoscopic patch. The authors have empirically observed that such averages do not arise from the smoothing out of genuine fluctuations. Rather, it seems that apart from degenerate cases, all the a -edges within a patch have roughly the same inclusion probability. These degenerate cases occur in patches where the asymptotic height function is not smooth (so that (s,t) ( f x , f y ) is undefined), or where the tilt is nearly extremal (i.e., s t is close to 2). Hence we believe: conj In the thermodynamic limit, the probability of seeing a domino in a particular location is given by the suitable member of the 4-tuple (p a,p b,p c,p d) , wherever the tilt (s,t) ( given by the partial derivatives of the entropy-maximizing height function ) is defined and satisfies s t 2 . probconj conj This is the conjectural interpretation of p a,p b,p c,p d that was alluded to in Subsection statementsec . We are two removes from being able to prove it in the sense that we do not even know how to prove Conjecture denseconj . The restriction s t 2 deserves some comment. It is not hard to devise a large region composed of long diagonal herringbones'' that has only one tiling (see the final section of CEP ). For such regions, the edge-inclusion probabilities can be made to fluctuate erratically between 0 and 1. Such regions have asymptotic height functions in which the tilt is extremal everywhere it is defined, so we can rule out such behavior on the basis of tilt. Incidentally, a conjecture analogous to Conjecture probconj can be made for the case of lozenge tilings (which can be studied using the methods of this paper, by setting one edge weight equal to 0 ). Here the analogue of Conjecture denseconj is actually a theorem; that is, for lozenges there are only three kinds of orientations of tiles, so that their relative frequencies, which jointly have two degrees of freedom, both determine and are determined by the local tilt (s,t) of the asymptotic height function. Straying further into the unknown, we might inquire about the probabilities of finite (colored) configurations of tiles. Here again we are guided by the Ansatz that what is true for large tori should be true for large finite regions as well. Given any tilt (s,t) satisfying s t 2 , choose weights a,b,c,d satisfying ab cd and giving tilt (s,t) in accordance with our earlier formulas, and use the formula in Proposition couplingfn to define a measure s,t on the space of domino tilings of the plane. This measure is invariant under color-preserving translations and satisfies the property of conditional uniformity''---given any finite region, the conditional distribution upon fixing a tiling of the rest of the plane is uniform. Proposition couplingfn invites us to surmise: conj torusconj Let s t 2 , and let a,b,c,d be weights satisfying ab cd such that the average height function for weighted torus tilings has tilt (s,t) . Then for any colored configuration of dominos, the probability of finding it in a specified location in a random n n torus tiling converges as n to the value given by the measure s,t . conj Proposition couplingfn approaches this claim, but it restricts n to lie in a large subset W of the integers. The measures s,t have positive entropy whenever s t 2 . Furthermore, the coupling function calculations in the proof of Proposition sigmasquared show that these measures are mixing (correlations between distant cylinder sets tend to 0 ) and hence ergodic. We believe that these measures can be characterized uniquely by the properties mentioned so far, although we cannot prove it. conj Every ergodic, conditionally uniform measure on the set of tilings of the plane that is invariant under color-preserving translations and has positive entropy is of the form s,t for some (s,t) satisfying s t 2 . conj Assuming the truth of Conjecture torusconj , it would be natural to advance a further claim that would come close to being the final, definitive answer to the question Kasteleyn raised nearly four decades ago: conj lastconj In the thermodynamic limit for a sequence of finite regions converging to a fixed shape, the probability of seeing any colored configuration of dominos is given by the measure s,t , wherever the tilt (s,t) given by the variational principle is defined and satisfies s t 2 . conj Finally, we discuss the seemingly miraculous geometric interpretation of our formula for the entropy. As was pointed out in Section introsec , when each of a,b,c,d is less than the sums of the others, the (asymptotic) entropy of torus tilings with weights a,b,c,d equals 1 times the volume of a three-dimensional ideal hyperbolic pyramid, whose vertices in the upper-half-space model are the vertex at infinity and the four vertices of the cyclic quadrilateral of Euclidean edge lengths a,c,b,d (in cyclic order). Notice that we needn't assume any particular scaling for the weights, because homotheties (x,y,z) (rx,ry,rz) are isometries of hyperbolic space. prob What do tilings have to do with hyperbolic geometry prob Explaining that connection is one of the most intriguing open problems in this area. Acknowledgements We thank Matthew Blum, Ben Raphael, Ben Wieland, and David Wilson for helpful discussions.