arXiv is now an independent nonprofit! Learn more
License: arXiv.org perpetual non-exclusive license
arXiv:1809.02703v1 [cs.DS] 07 Sep 2018

Torpid Mixing of Markov Chains for the Six-vertex Model on 2\mathbb{Z}^{2}

Tianyu Liu Thanks: Department of Computer Sciences, University of Wisconsin-Madison. Supported by NSF CCF-1714275. Email: tl@cs.wisc.edu
Abstract

In this paper, we study the mixing time of two widely used Markov chain algorithms for the six-vertex model, Glauber dynamics and the directed-loop algorithm, on the square lattice 2\mathbb{Z}^{2}. We prove, for the first time that, on finite regions of the square lattice these Markov chains are torpidly mixing under parameter settings in the ferroelectric phase and the anti-ferroelectric phase.

1 Introduction

Introduced by Linus Pauling [Pau35] in 1935 to describe the properties of ice, the six-vertex model or the ice-type model was originally studied in statistical mechanics as an abstraction of crystal lattices with hydrogen bonds. During the following decades, it has attracted enormous interest in many disciplines of science, and become one of the most fundamental models defined on the square lattice. In particular, the discovery of integrability of the six-vertex models with periodic boundary conditions was considered a milestone in statistical physics [Lie67c, Lie67a, Lie67b, Sut67, FW70].

For computational expediency and modeling purposes, physicists almost entirely focused on planar lattice models. On the square lattice 2\mathbb{Z}^{2}, every vertex is connected by an edge to four “nearest neighbors”. States of the six-vertex model on 2\mathbb{Z}^{2} are orientations of the edges on the lattice satisfying the ice-rule — every vertex has two incoming edges and two outgoing edges, i.e., they are Eulerian orientations. The name of six-vertex model comes from the fact that there are six ways of arranging directions of the edges around a vertex (see Figure 1).

Refer to caption
(a) 11
Refer to caption
(b) 22
Refer to caption
(c) 33
Refer to caption
(d) 44
Refer to caption
(e) 55
Refer to caption
(f) 66
Figure 1: Valid configurations of the six-vertex model.

In general, each of the six local arrangements will have a weight, denoted by w1,,w6w_{1},\dots,w_{6}, using the ordering of Figure 1. The total weight of a state is the product of all vertex weights in the state. If there is no ambient electric field, by physical considerations, then the total weight of a state should remain unchanged when flipping all arrows [Bax82]. Thus one may assume without loss of generality that w1=w2=a,w3=w4=b,w5=w6=cw_{1}=w_{2}=a,w_{3}=w_{4}=b,w_{5}=w_{6}=c. This complementary invariance is known as arrow reversal symmetry or zero field assumption. In this paper, we assume a,b,c>0a,b,c>0, as is the case in classical physics. We study the six-vertex model restricted to a finite region of the square lattice with various boundary conditions customarily studied in statistical physics literature. On a finite subset Λ2\Lambda\subset\mathbb{Z}^{2}, denote the set of valid configurations (i.e. Eulerian orientations) by Ω\Omega. The probability that the system is in a state τΩ\tau\in\Omega is given by the Gibbs distribution

μ(τ)=1Z(an1+n2bn3+n4cn5+n6),\mu\left(\tau\right)=\frac{1}{Z}\left(a^{n_{1}+n_{2}}b^{n_{3}+n_{4}}c^{n_{5}+n_{6}}\right),

where nin_{i} is the number of vertices in type ii (1i6)(1\leq i\leq 6) on Λ\Lambda in the state τ\tau, and the partition function ZZ is a normalizing constant which is the sum of the weights of all states.

In 1967, Elliot Lieb [Lie67c] famously showed that, for parameters (a,b,c)=(1,1,1)(a,b,c)=(1,1,1) on the square lattice graph, as the side NN of the square approaches \infty, the value of the “partition function per vertex” W=Z1/N2W=Z^{1/N^{2}} approaches (43)3/21.5396007\left(\frac{4}{3}\right)^{3/2}\approx 1.5396007\ldots (this is called Lieb’s square ice constant). This result is called an exact solution of the model, and is considered a triumph. After that, exact solutions for other parameter settings have been obtained in the limiting sense [Lie67a, Lie67b, Sut67, FW70]. Readers are referred to [CLL17] for known results in the computational complexity of (both exactly and approximately) computing the partition function ZZ of the six-vertex model on general 4-regular graphs.

In statistical physics, Markov chain Monte Carlo (MCMC) is the most popular tool to numerically study the properties of the six-vertex model. A partial list includes [RS72, YN79, BN98, Elo99, SZ04, AR05, LKV17]. In the literature, two Markov chain algorithms are mainly used. The first one is Glauber dynamics. It can be shown that there is a correspondence between Eulerian orientations of the edges and proper three-colorings of the faces on a rectangle region of the square lattice. (See Chapter 8 of [Bax82] for a proof). Therefore, the Glauber dynamics for the three-coloring problem on square lattice regions (which changes a local color at each step) can be employed to sample Eulerian orientations. In fact, this simple Markov chain is used in numerical studies (e.g. in [Elo99, AR05, LKV17] for the density profile) of the six-vertex model under various boundary conditions. The second one is the directed-loop algorithm. Invented by Rahman and Stillinger [RS72] and widely adopted in the literature (e.g., [YN79, BN98, SZ04]), the transitions of this algorithm are composed of creating, shifting, and merging of two “defects” on the edges. An interesting aspect is that this process depicts the Bjerrum defects happening in real ice [BN98]. More detailed descriptions of the two Markov chain algorithms can be found in Section 2.

With the heavy usage of MCMC in statistical mechanics for the six-vertex model, the efficiency of Markov chain algorithms was inevitably brought into focus by physicists. Many of them (e.g. [BN98, SZ04, LKV17]) reported that Glauber dynamics and the directed-loop algorithms of the six-vertex model experienced significant slowdown and are even “impractical” for simulation purposes when the parameter settings are in the ordered phases (see Figure 2(a), in the regions FE & AFE). Despite the concern and numerical experience for the convergence rate of these algorithms, there is no previous provable result except for one point (that corresponds to the unweighted case) in the parameter space. This is in stark contrast to the popular studies on the mixing rate of Markov chains for the ferromagnetic Ising model [MO94a, MO94b, CGMS96, LS12] and hardcore gas model on lattice regions [BCK+99, Ran06, BGRT13].

Prior to [CLL17], to our best knowledge, the only provable result in the complexity of approximate sampling and counting for the six-vertex model is at the single, unweighted, parameter setting (a,b,c)=(1,1,1)(a,b,c)=(1,1,1) where the partition function counts Eulerian orientations. In the unweighted case, all known results are positive. Mihail and Winkler’s pioneering work [MW96] gave the first fully polynomial randomized approximation scheme (FPRAS) for the number of Eulerian orientations on a general graph (not necessarily 4-regular). Luby, Randall, and Sinclair showed that Glauber dynamics with extra moves is rapidly mixing on rectangular regions of the square lattice with fixed boundary conditions [LRS01]. Randall and Tetali proved the rapid mixing of the Glauber dynamics (without extra moves) with fixed boundary conditions by a comparison technique applied to this Markov chain and the Luby-Randall-Sinclair chain [RT00]. Goldberg, Martin, and Paterson extended further the rapid mixing of Glauber dynamics to the free-boundary case [GMP04]. The unweighted setting is the single green point depicted in the blue region of Figure 2(b).

In [CLL17], Cai, Liu, and Lu showed that under parameter settings (a,b,c)(a,b,c) with a2b2+c2a^{2}\leq b^{2}+c^{2}, b2a2+c2b^{2}\leq a^{2}+c^{2}, and c2a2+b2c^{2}\leq a^{2}+b^{2} (the blue region in Figure 2(b)), the directed-loop algorithm mixes in polynomial time with regard to the size of input for any general 4-regular graph, resulting in an FPRAS for the partition function of the six-vertex model. Moreover, it is shown that in the ordered phases (FE & AFE in Figure 2(a)), the partition function on a general graph is not efficiently approximable unless NP=RP. Although the rapid mixing property for the directed-loop algorithm on general 4-regular graphs implies the same on the lattice region, the hardness result for general 4-regular graphs has no implications on the mixing rate of Markov chains for the six-vertex model on the square lattice in the ordered phases (FE & AFE).

In this paper, we give the first provable negative results on mixing rates of the two Markov chains for the six-vertex model under parameter settings in the ferroelectric phases and the anti-ferroelectric phase. Our results conform to the phase transition phenomena in physics. Here we briefly describe the phenomenon of phase transition of the zero-field six-vertex model (see Baxter’s book [Bax82] for more details). On the square lattice in the thermodynamic limit: (1) When a>b+ca>b+c (FE: ferroelectric phase) any finite region tends to be frozen into one of the two configurations where either all arrows point up or to the right (Figure 1-1), or all point down or to the left (Figure 1-2). (2) Symmetrically when b>a+cb>a+c (also FE) all arrows point down or to the right (Figure 1-3), or all point up or to the left (Figure 1-4). (3) When c>a+bc>a+b (AFE: anti-ferroelectric phase) configurations in Figure 1-5 and Figure 1-6 alternate. (4) When c<a+bc<a+b, b<a+cb<a+c, and a<b+ca<b+c, the system is disordered (DO: disordered phase) in the sense that all correlations decay to zero with increasing distance; in particular on the dashed curve c2=a2+b2c^{2}=a^{2}+b^{2} the model can be solved by Pfaffians exactly [FW70], and the correlations decay inverse polynomially, rather than exponentially, in distance. See Figure 2(a).

Refer to caption
(a) Phase diagram of the six-vertex model.
Refer to caption
(b) Mixing time of Markov chains for the six-vertex model on 2\mathbb{Z}^{2}.
Figure 2:

Let Λ\Lambda be a square region on the square lattice. We show the following two theorems.

Theorem 1.1 (Ferroelectric phase).

The directed-loop algorithm for the six-vertex model under parameter settings (a,b,c)(a,b,c) with a>b+ca>b+c or b>a+cb>a+c (i.e. the whole FE) mixes torpidly on Λ\Lambda with periodic boundary conditions.

Remark 1.1.

We note that for periodic boundary conditions Glauber dynamics is not irreducible, so we do not consider that.

Theorem 1.2 (Anti-ferroelectric phase).

Both Glauber dynamics and the directed-loop algorithm for the six-vertex model under parameter settings (a,b,c)(a,b,c) with c2.639max(a,b)c\geq 2.639\max(a,b) (in AFE) mix torpidly on Λ\Lambda with free boundary conditions and periodic boundary conditions.

Parameter settings covered by the above two theorems are depicted as the grey region in Figure 2(b). Given that the FF model in statistical mechanics is a special case of the six-vertex model when a=b=1a=b=1 [Lie67a], Theorem 1.2 holds for the FF model with c2.639c\geq 2.639.

Our proofs build on the equivalence between small conductance and torpid mixing by Jerrum and Sinclair [SJ89]. When arguing Markov chains for the six-vertex model in the anti-ferroelectric phase have small conductance, we switch our view between finite regions of the square lattice and their medial graphs. This transposition allows us to adopt a Peierls argument which has been used in statistical physics to prove the existence of phase transitions (e.g., [Pei36, BKW73]), and in theoretical computer science to prove the torpid mixing of Markov chains (e.g., [Ran06, BGRT13]).

In the proof of Theorem 1.2, we introduce a version of the fault line argument for the six-vertex model. Fault line arguments are introduced by Dana Randall [Ran06] for the lattice hardcore gas and latter adapted in [LPW06] for the lattice ferromagnetic Ising, which proves torpid mixing of Markov chains via topological obstructions. The constant 2.639 comes from an upper bound for the connective constant for the square lattice self-avoiding walks [GC01].

2 Preliminaries

2.1 Markov chains

2.1.1 Glauber dynamics

Denote by Λn\Lambda_{n} a square lattice region where there are nn vertices of degree 4 on each row and each column. Λn\Lambda_{n} is in periodic boundary condition if it forms a two-dimensional torus; the free boundary condition can be formulated in the following way: there are n+2n+2 vertices on each row and each column, where the “boundary vertices” are of degree 1 and don’t need to satisfy the ice-rule (and don’t take weights) in a valid six-vertex configuration. For convenience, we assume there are “virtual edges” connecting every two boundary vertices with unit distance on 2\mathbb{Z}^{2}. A virtual edge does not have orientations, serving only the purpose that every unit square inside the (n+1)×(n+1)(n+1)\times(n+1) region is closed.

Let Ω\Omega be the set of all valid configurations of the six-vertex model (Eulerian orientations) on Λn\Lambda_{n}. The Glauber-dynamics Markov chain, which we will denote by G\mathcal{M}_{G}, has state space Ω\Omega. To move from one configuration to another, this chain selects a unit square (a face) ss on Λn\Lambda_{n} (together with the virtual edges) uniformly at random. If all the non-virtual edges along the unit square ss are oriented consistently (clockwise or counter-clockwise), the chain picks a direction dd (clockwise or counter-clockwise) and reorients the non-virtual edges along ss according to the Gibbs measure.

One can easily check that such transitions take valid configurations to valid configurations. Actually, this Markov chain is equivalent to that in [GMP04] for sampling three-colorings on the faces of Λn\Lambda_{n}. The ergodicity of that chain translates straightforwardly to the ergodicity of G\mathcal{M}_{G} (with free boundary conditions) thanks to the equivalence between Eulerian orientations and three-colorings on 2\mathbb{Z}^{2}. Besides, the heat-bath move indicates that the stationary distribution of G\mathcal{M}_{G} is the Gibbs distribution for the six-vertex model.

2.1.2 Directed-loop algorithm

The directed-loop algorithm Markov chain, denoted by D\mathcal{M}_{D}, is formally defined in [CLL17] for general 4-regular graphs, so here we only describe D\mathcal{M}_{D} at a high level.

The state space of D\mathcal{M}_{D} is not only Ω\Omega, the “perfect” Eulerian orientations, but also the set of all “near-perfect” Eulerian orientations, denoted by Ω\Omega^{\prime}. For example, in Figure 3 the state τur\tau_{ur} is in Ω\Omega and all other five states are in Ω\Omega^{\prime}. We think of each edge in Λn\Lambda_{n} as the two half-edges cut in the middle, and each of the half edge can be oriented independently. We say an orientation of all the half-edges is perfect (in Ω\Omega) if every pair of half-edges is oriented consistently and the ice-rule is satisfied at every vertex (except for boundary vertices under free boundary conditions); an orientation is near-perfect (in Ω\Omega^{\prime}) if there are exactly two pairs of half-edges p1p_{1} and p2p_{2} not oriented consistently and the ice-rule is satisfied at every vertex (except for boundary vertices under free boundary conditions), with the restriction that if two half-edges in p1p_{1} are oriented toward each other then in p2p_{2} the two half-edges must be oriented against each other and vice versa.

The transitions in D\mathcal{M}_{D} are Metropolis moves among “neighboring” states. An Ω\Omega state τ\tau and an Ω\Omega^{\prime} state τ\tau^{\prime} are neighboring if τ\tau^{\prime} can be transformed from τ\tau by picking two half-edges e1,e2e_{1},e_{2} incident to a vertex vv with one pointing inwards vv and the other pointing outwards vv (or two half-edges e1,e2e_{1},e_{2} on the boundary with one pointing towards the boundary and the other pointing against the boundary), and reverse the direction of e1e_{1} and e2e_{2} together. For instance, in Figure 3 {τur,τ1b}\{\tau_{ur},\tau_{1b}\} and {τur,τ1c}\{\tau_{ur},\tau_{1c}\} are two pairs of neighboring states. An Ω\Omega^{\prime} state τ1\tau^{\prime}_{1} and another Ω\Omega^{\prime} state τ2\tau^{\prime}_{2} are neighboring if τ2\tau^{\prime}_{2} can be transformed from τ1\tau^{\prime}_{1} by “shifting” one pair of conflicting half-edges one step away, while fixing the other pair of conflicting half-edges. For example, in Figure 3 τ1c\tau_{1c} and τ2\tau_{2} are neighboring to each other. D\mathcal{M}_{D} can be proved to be ergodic and converges to the Gibbs measure on ΩΩ\Omega\cup\Omega^{\prime} with both free boundary conditions and periodic boundary conditions [CLL17].

2.2 Mixing time

The mixing time tmixt_{\operatorname{mix}} measures the time required by a Markov chain to evolve to be close to its stationary distribution, in terms of total variation distance. (The definition of mixing time can be found in [LPW06].) We say a Markov chain is torpid mixing if the mixing time is exponentially large in the input size. A common technique to bound the mixing time is via bounding conductance, defined by Jerrum and Sinclair [SJ89].

Let π\pi denote the stationary distribution of an ergodic and time reversible (π(x)P(x,y)=π(y)P(y,x)\pi(x)P(x,y)=\pi(y)P(y,x) for any x,yΩx,y\in\Omega) Markov chain \mathcal{M} on a finite state space Ω\Omega, with transition probabilities P(x,y)P(x,y), x,yΩx,y\in\Omega. The conductance of \mathcal{M} is defined by

Φ=Φ()=minSΩ0<π(S)12Q(S,S¯)π(S),\Phi=\Phi(\mathcal{M})=\min_{\begin{subarray}{c}S\subset\Omega\\ 0<\pi(S)\leq\frac{1}{2}\end{subarray}}\frac{Q(S,\overline{S})}{\pi(S)},

where Q(S,S¯)Q(S,\overline{S}) denotes the sum of Q(x,y)=π(x)P(x,y)Q(x,y)=\pi(x)P(x,y) over edges in the transition graph of \mathcal{M} with xSx\in S, and yS¯=ΩSy\in\overline{S}=\Omega\setminus S.

In order to show a Markov chain mixes torpidly, we only need to prove that the conductance is (inverse) exponentially small due to the following bound [LPW06]:

tmix=tmix(14)14Φ.t_{\operatorname{mix}}=t_{\operatorname{mix}}\left(\frac{1}{4}\right)\geq\frac{1}{4\Phi}.

As is usually assumed, Markov chains studied in this paper are all lazy (P(x,x)=12P(x,x)=\frac{1}{2} for any xΩx\in\Omega) and transition probabilities (P(x,y)P(x,y) for x,yΩx,y\in\Omega) between neighboring states (where P(x,y)>0P(x,y)>0) are at least inverse polynomially large. Therefore, armed with the above bound, we can prove the torpid mixing of a Markov chain if we can establish the following:

  1. 1.

    Partition the state space Ω\Omega into three subsets ΩLEFTΩMIDDLEΩRIGHT\Omega_{\text{LEFT}}\cup\Omega_{\text{MIDDLE}}\cup\Omega_{\text{RIGHT}} as a disjoint union.

  2. 2.

    Show that for any state τlΩLEFT\tau_{l}\in\Omega_{\text{LEFT}} and τrΩRIGHT\tau_{r}\in\Omega_{\text{RIGHT}}, P(τl,τr)=0P(\tau_{l},\tau_{r})=0. Under the assumption that the Markov chain is irreducible (i.e., the transition graph is strongly connected), this indicates that in order to go from states in ΩLEFT\Omega_{\text{LEFT}} to states in ΩRIGHT\Omega_{\text{RIGHT}}, the Markov process has to go through the “middle states” ΩMIDDLE\Omega_{\text{MIDDLE}}.

  3. 3.

    Demonstrate that π(ΩMIDDLE)\pi(\Omega_{\text{MIDDLE}}) is exponentially small (compared with min(π(ΩLEFT),π(ΩRIGHT))\min(\pi(\Omega_{\text{LEFT}}),\pi(\Omega_{\text{RIGHT}}))) in the input size. This means that starting from any state in ΩLEFT\Omega_{\text{LEFT}}, the probability of going through ΩMIDDLE\Omega_{\text{MIDDLE}} (and consequently to any state in ΩRIGHT\Omega_{\text{RIGHT}} and reach stationarity) is exponentially small. Hence the conclusion of torpid mixing.

3 Ferroelectric phase

In this section we prove Theorem 1.1 that D\mathcal{M}_{D} in the directed-loop algorithm for the six-vertex model in the ferroelectric phase is torpid mixing on Λn\Lambda_{n} with periodic boundary conditions.

For any parameter setting (a,b,c)(a,b,c) in the ferroelectric phase, either a>b+ca>b+c or b>a+cb>a+c. By symmetry, without loss of generality, suppose a>b+ca>b+c. This implies that vertex configurations as shown in Figure 1-1 and Figure 1-2 have higher weights than others. Under the periodic boundary condition, there is a state τur\tau_{ur} in which every vertical edge points upwards and every horizontal edge points to the right (Figure 3(a)), i.e., every vertex on Λn\Lambda_{n} is in local configuration shown in Figure 1-1. The total weight of τur\tau_{ur} is an2a^{n^{2}} as there are n2n^{2} vertices on Λn\Lambda_{n}.

Refer to caption
(a) τur\tau_{ur}
Refer to caption
(b) τ1b\tau_{1b}
Refer to caption
(c) τ1c\tau_{1c}
Refer to caption
(d) τ2\tau_{2}
Refer to caption
(e) τsaw\tau_{\text{saw}}
Refer to caption
(f) τcircle\tau_{\text{circle}}
Figure 3: Some states in the state space of D\mathcal{M}_{D}.

For D\mathcal{M}_{D}, the three-way partition of the state space ΩΩ\Omega\cup\Omega^{\prime} is as follows. Denote by TiT_{i} the states that can be reached from τur\tau_{ur} in at most ii steps of transitions where ii is a nonnegative integer. Write Ti=TiTi1\partial T_{i}=T_{i}\setminus T_{i-1} for i1i\geq 1. Let ΩLEFT=Tn1\Omega_{\text{LEFT}}=T_{n-1}, ΩMIDDLE=Tn\Omega_{\text{MIDDLE}}=\partial T_{n}, and ΩRIGHT=(ΩΩ)(ΩLEFTΩMIDDLE)\Omega_{\text{RIGHT}}=(\Omega\cup\Omega^{\prime})\setminus(\Omega_{\text{LEFT}}\cup\Omega_{\text{MIDDLE}}). It is obvious that ΩΩ=ΩLEFTΩMIDDLEΩRIGHT\Omega\cup\Omega^{\prime}=\Omega_{\text{LEFT}}\cup\Omega_{\text{MIDDLE}}\cup\Omega_{\text{RIGHT}} is a partition of the state space. Clearly τurΩLEFT\tau_{ur}\in\Omega_{\text{LEFT}}, thus the total weight of ΩLEFT\Omega_{\text{LEFT}} is no less than an2a^{n^{2}}, the weight of τur\tau_{ur}.

Before proving the total weight of ΩMIDDLE\Omega_{\text{MIDDLE}} is exponentially small compared with that of ΩLEFT\Omega_{\text{LEFT}} or ΩRIGHT\Omega_{\text{RIGHT}}, let us look at what is in TiT_{i} with 0in0\leq i\leq n. T0T_{0} is just {τur}\{\tau_{ur}\}. T1\partial T_{1} consists of all the states evolved from τur\tau_{ur} by picking a vertex vv on Λn\Lambda_{n} and two incident half-edges (one pointing towards vv and the other away from vv), and then reversing the orientations on these two edges. After such a transition, two pairs of conflicting half-edges are created, so T1Ω\partial T_{1}\subseteq\Omega^{\prime}.

For example, the states shown in Figure 3(b) (state τ1b\tau_{1b}) and Figure 3(c) (state τ1c\tau_{1c}) are in T1\partial T_{1}. The weight of τ1b\tau_{1b} is an21ba^{n^{2}-1}b and that of τ1c\tau_{1c} is an21ca^{n^{2}-1}c. For every state in T1\partial T_{1} obtained by transitions from τur\tau_{ur}, there is exactly one vertex vv^{*} on Λn\Lambda_{n} no longer in the local configuration Figure 1-1. (Of course no vertex can be in state Figure 1-2.) Actually, depending on whether the two pairs of conflicting half-edges are: (1) both vertical, (2) both horizontal, or (3) one horizontal and the other vertical, the vertex vv^{*} is in configuration shown in (1) Figure 1-3, (2) Figure 1-4, or (3) Figure 1-5/6, respectively. Therefore, every state in T1\partial T_{1} has weight an21ba^{n^{2}-1}b in case (1) and case (2) or an21ca^{n^{2}-1}c in case (3).

Transitions from states in T1\partial T_{1} to states in T2\partial T_{2} are composed of “shifting” one of the two conflicting pairs of half-edges to a neighboring edge on Λn\Lambda_{n}. For example, the state in Figure 3(d) is in T2\partial T_{2}. This process will result in exactly two vertices on Λn\Lambda_{n} not in local configuration Figure 1-1 (nor in Figure 1-2). As a consequence, the weight of any state in T2\partial T_{2} is among an22b2a^{n^{2}-2}b^{2}, an22bca^{n^{2}-2}bc, and an22c2a^{n^{2}-2}c^{2}. The state shown in Figure 3(d) has weight an22c2a^{n^{2}-2}c^{2}.

This line of argument can be extended for Ti\partial T_{i} for 1in1\leq i\leq n, in any state of which there are exactly ii vertices on Λn\Lambda_{n} not in local configuration Figure 1-1 (nor in Figure 1-2). When two conflicting pairs of half-edges are created in T1\partial T_{1}, one of them is above or to the right of another (or both). Denote the former by purp_{ur} (the green pair in Figure 3) and the latter by pdlp_{dl} (the red pair in Figure 3). Observe that as the Markov chain evolves, by a single step, from a state in Ti\partial T_{i} to another in Ti+1\partial T_{i+1} (where 1in11\leq i\leq n-1), either purp_{ur} is “pushed” up or to the right, or pdlp_{dl} down or to the left. For example, from Figure 3(c) to Figure 3(d), purp_{ur} is pushed to the right. By induction, purp_{ur} is always above or to the right of pdlp_{dl} (when i<ni<n). A direct consequence is that there can be no state containing a closed circuit formed by the reversed edges (with regard to τur\tau_{ur}) in Ti\partial T_{i} until i=ni=n. Therefore, the edges reversed in any state in Ti(1in)\partial T_{i}\ (1\leq i\leq n) can be seen as either a self-avoiding walk between the middle points of the two pairs of conflicting half-edges (e.g. Figure 3(e)) or a self-avoiding circuit (e.g. Figure 3(f)). In fact, when the reversed edges form a circuit, the circuit must “go straightforward” at each step. This circuit is a circle parallel or perpendicular to the torus equatorial plane. The weight of any state in Ti\partial T_{i} is an2ibjcka^{n^{2}-i}b^{j}c^{k} with i=j+ki=j+k, where the values of jj and kk depend on how many “turnarounds” are there in the self-avoiding walk.

Therefore, the total weight of states in Tn\partial T_{n} is at most n2an2n(b+c)nn^{2}\cdot a^{n^{2}-n}(b+c)^{n}, where n2n^{2} is an upper bound on all the possible starting points for self-avoiding walks, and each monomial in (b+c)n(b+c)^{n} is from a unique self-avoiding walk. Combining with the fact that total weight of ΩLEFT\Omega_{\text{LEFT}} is at least an2a^{n^{2}} (the weight of τur\tau_{ur}) and is at most that of ΩRIGHT\Omega_{\text{RIGHT}} (because there is a weight-preserving injective map from ΩLEFT\Omega_{\text{LEFT}} to ΩRIGHT\Omega_{\text{RIGHT}} by reversing orientations of all the edges), we know that the conductance of D\mathcal{M}_{D} is at most n2an2n(b+c)nan2=n2(b+ca)n\frac{n^{2}\cdot a^{n^{2}-n}(b+c)^{n}}{a^{n^{2}}}=n^{2}\left(\frac{b+c}{a}\right)^{n}. This is exponentially small in nn since a>b+ca>b+c are fixed constants in the ferroelectric phase.

4 Anti-ferroelectric phase

In this section, we prove the following theorem which is part of Theorem 1.2. After proving Theorem 4.1, we state the ideas needed to extend it to Theorem 1.2, the full proof of which is omitted due to space limit.

As we did in the ferroelectric phase, the intuition behind our proof for the anti-ferroelectric phase is to find a partition Ω=ΩLEFTΩMIDDLEΩRIGHT\Omega=\Omega_{\text{LEFT}}\cup\Omega_{\text{MIDDLE}}\cup\Omega_{\text{RIGHT}} of the state space of G\mathcal{M}_{G}, i.e., all the Eulerian orientations on Λn\Lambda_{n}. However, the strategy is different from that used in Section 3 — here the subset ΩMIDDLE\Omega_{\text{MIDDLE}} is determined in terms of a topological obstruction.

Theorem 4.1 (Anti-ferroelectric phase).

Glauber dynamics for the six-vertex model under parameter settings (a, b, c) with c2.639max(a,b)c\geq 2.639\max(a,b) mix torpidly on Λn\Lambda_{n} with free boundary conditions.

Observe that there are two states in Ω\Omega with maximum weights: τG\tau_{\text{G}} (Figure 4(a)) and τR\tau_{\text{R}} (Figure 4(b)) where every vertex is in local configuration Figure 1-5 or Figure 1-6, and thus has vertex weight cc. Since τG\tau_{\text{G}} and τR\tau_{\text{R}} are total reversals of each other in edge orientations, for any edge in any state τΩ\tau\in\Omega, it is oriented either as in τG\tau_{\text{G}} or as in τR\tau_{\text{R}}. Let us call an edge to be green if it is oriented as is in τG\tau_{\text{G}} and red otherwise. Observe that in order to satisfy the ice-rule (2-in-2-out), the number of green (and thus also two red) edges incident to any vertex (except for the boundary vertices) is always even (0,2, or 4)(0,2,\text{ or }4), and if there are two green (or red) edges they must be rotationally adjacent to each other. See Figure 4(c) for an example. Also note that the four edges along a unit square on 2\mathbb{Z}^{2} are all red edges or all green edges if and only if they are oriented consistently, hence flippable by a single move of G\mathcal{M}_{G}.

Refer to caption
(a) τG\tau_{\text{G}}
Refer to caption
(b) τR\tau_{\text{R}}
Refer to caption
(c) τ\tau
Figure 4: Some states in the state space of G\mathcal{M}_{G}.

We say a simple path from a horizontal edge on the left boundary of Λn\Lambda_{n} to a horizontal edge on the right boundary of Λn\Lambda_{n} is a horizontal green (or red) bridge if the path consists of only green (or red, respectively) edges; a vertical green (or red) bridge is defined similarly. A state τΩ\tau\in\Omega has a green cross if it has both a green horizontal bridge and a green vertical bridge; a red cross is defined similarly. Let CGΩC_{\text{G}}\subset\Omega denote the states having a green cross and CRC_{\text{R}} the states having a red cross. In the following lemma, we prove that CGCR=C_{\text{G}}\cap C_{\text{R}}=\emptyset.

Lemma 4.2.

A green cross and a red cross cannot coexist.

Proof.

It suffices to show that a green horizontal bridge precludes a red vertical bridge. Consider a virtual point vLv_{\text{L}} sitting to the left of Λn\Lambda_{n} connected by an edge to every (external) vertices of Λn\Lambda_{n} on the left boundary, and another virtual point vRv_{\text{R}} connected by an edge to every vertex of Λn\Lambda_{n} on the right boundary. Connect vLv_{\text{L}} and vRv_{\text{R}} by an edge below Λn\Lambda_{n}.

If there is a green horizontal bridge, then by definition there is a continuous closed curve 𝒞\mathcal{C} formed by the bridge and some edges we added (Figure 5(a)). According to the Jordan Curve Theorem, 𝒞\mathcal{C} separates the plane into two disjoint regions, the inside and the outside. Vertices of Λn\Lambda_{n} that are on the bottom boundary are inside; vertices on the top boundary are outside. Therefore, in order to have a red vertical bridge, there must be a simple red path going across 𝒞\mathcal{C}. That is to say, a red vertical bridge must cross the green horizontal bridge.

Refer to caption
(a) A closed curve 𝒞\mathcal{C}.
Refer to caption
(b) An impossible configuration.
Figure 5:

However, this is impossible. Clearly, being of different colors, a red bridge and a green bridge cannot share any edge. Since the local configuration shown in Figure 5(b) means that the four edges incident to a vertex (i,j)(i,j) on Λn\Lambda_{n} are all pointing inwards (when i+ji+j is even) or all pointing outwards (when i+ji+j is odd), it is not allowed in any valid states of six-vertex configurations. Similarly, the local configuration of a vertex surrounded by two red horizontal edges and two green vertical edges (a 90 degree rotation of Figure 5(b)) is also not allowed. ∎

Next we characterize the states in Ω(CGCR)\Omega\setminus(C_{\text{G}}\cup C_{\text{R}}). Define a shifted lattice11 1 Strictly speaking, a lattice is a discrete subgroup of n\mathbb{R}^{n}. A shifted copy of a lattice does not contain 0. LL to be 2+(12,12)\mathbb{Z}^{2}+\left(\frac{1}{2},\frac{1}{2}\right) where two points (a,b)(a,b) and (c,d)(c,d) in LL are neighbors if |ac|=|bd|=1|a-c|=|b-d|=1 (aa, bb, cc, and dd are all half integers), i.e., they are at the center of a square in 2\mathbb{Z}^{2} and are connected by “diagonal” edges of length 2\sqrt{2}. An example of LL and its relationship with 2\mathbb{Z}^{2} is shown in Figure 6(a). LL is not connected — it is composed of two sub-lattices L0L_{0} and L1L_{1} (depicted with different colors in Figure 6(b)). Denote by LnL_{n} the restriction of LL on the finite region inside Λn\Lambda_{n}. Note that in graph theoretical terms, the square lattice Λn\Lambda_{n} is planar and 4-regular, and thus can be seen as the medial graph of two planar graphs. In fact, they are L0L_{0} and L1L_{1} (restricted onto LnL_{n}).

From now on, we use Λn\Lambda_{n}-vertices/edges as an abbreviation for vertices/edges in Λn\Lambda_{n}; and we use LnL_{n}-vertices/edges and other similar notations whenever it has a clear meaning in the context. For any state τΩ\tau\in\Omega, there is a subset LτL_{\tau} of LnL_{n}-edges associated with τ\tau. Each LnL_{n}-edge ee “goes diagonally through” exactly one Λn\Lambda_{n}-vertex, denoted by vev_{e}. We say ee is in LτL_{\tau} if and only if the four Λn\Lambda_{n}-edges incident to vev_{e} are 2-green-2-red and ee separates the two green edges from the two red edges (Remember that in this case edges in the same color must be rotationally adjacent to each other). See Figure 6(c) for an instance of a state τ\tau and its associated LτL_{\tau}. In the following we abuse the notation and use LτL_{\tau} as its induced subgraph of LL. This view was adopted by [BKW73] for establishing the existence of the spontaneous staggered polarization in the anti-ferroelectric phase of the six-vertex model.

For any τΩ\tau\in\Omega and LτL_{\tau}, we make the following observations:

  • There is always an even number of LτL_{\tau}-edges meeting at any LnL_{n}-vertex, except for the LnL_{n}-vertices on the boundary. Because this number is equal to the number of times for the color change on the four Λn\Lambda_{n}-edges surrounding the LnL_{n}-vertex, if we start from any one of the four Λn\Lambda_{n}-edges and go rotationally over the four Λn\Lambda_{n}-edges, which is even.

  • For any Λn\Lambda_{n} vertex, there can be at most one LτL_{\tau}-edge going through which is either in L0L_{0} or in L1L_{1}.

  • If τ¯Ω\bar{\tau}\in\Omega is the state by a total edge reversal of τ\tau, then Lτ¯=LτL_{\bar{\tau}}=L_{\tau}.

Refer to caption
(a) LL
Refer to caption
(b) L0L_{0} and L1L_{1}
Refer to caption
(c) τ\tau and LτL_{\tau}
Figure 6:

For a state τΩ\tau\in\Omega, we say τ\tau has a horizontal (or vertical) fault line if there is a self-avoiding path in LτL_{\tau} connecting a LnL_{n}-vertex on the left (top, respectively) boundary of LnL_{n} to a LnL_{n}-vertex on the right (bottom, respectively) boundary of LnL_{n}. See Figure 7(c) for an example where a state has both a horizontal fault line and a vertical fault line. Denote by CFLC_{\text{FL}} the set of states containing a horizontal fault line or a vertical fault line. Since a fault line separates green edges from red edges, a vertical (horizontal) fault line precludes any horizontal (vertical, respectively) monochromatic bridge (the proof is basically the same as Lemma 4.2). This is to say, CGC_{\text{G}}, CFLC_{\text{FL}}, and CRC_{\text{R}} are pairwise disjoint. Next we show the following lemma and its direct implication (Corollary 4.4).

Lemma 4.3.

If in a state τ\tau there is no monochromatic cross, then there is a fault line.

Proof.

If there is no monochromatic cross (i.e., a green cross or a red cross) in τ\tau, we can assume that

there is no green horizontal bridge and there is no red horizontal bridge. (*)
  • Suppose τ\tau has a green horizontal bridge. There is no red vertical bridge since it cannot “cross” the green horizontal bridge; there is no green vertical bridge since there is no green cross. Therefore, this case is symmetric to (4), switching horizontal for vertical.

  • Suppose τ\tau has a red horizontal bridge. This case is similar to the above case, and thus is also symmetric to (4).

Next we show there is a vertical fault line if there is no monochromatic horizontal bridge. We introduce another graph MnM_{n} that is the medial graph of Λn\Lambda_{n} where every vertex of MnM_{n} corresponds to an edge of Λn\Lambda_{n}, i.e., two MnM_{n}-vertices are neighboring if the two corresponding Λn\Lambda_{n}-edges are rotationally adjacent. Note that MnM_{n} is part of another shifted square lattice. An example of MnM_{n} and its relationship with 2\mathbb{Z}^{2} is shown in Figure 7(a). For any state τΩ\tau\in\Omega, there is a subset MτM_{\tau} of MnM_{n}-edges associated with τ\tau. An MnM_{n}-edge ee is in MτM_{\tau} if the two vertices that ee is incident to (as Λn\Lambda_{n}-edges) have the same color (both green or both red). See Figure 7(d) for an example.

Observe that the correspondence between Λn\Lambda_{n}-edges and MnM_{n}-vertices translates into the correspondence between simple monochromatic paths in Λn\Lambda_{n} to simple connected paths in MτM_{\tau}. In fact, the connected components in MτM_{\tau}, capturing the notion of monochromatic regions of Λn\Lambda_{n}-edges, and connected components in LτL_{\tau}, capturing the notion of separation between regions of Λn\Lambda_{n}-edges of different colors, are in a dual relationship.

This duality is depicted in Figure 7(e) and helps us find a fault line. Let VMV_{M} be the collection of MnM_{n}-vertices that can be reached from the left boundary of MnM_{n} by a simple path in MτM_{\tau}. Since there is no monochromatic horizontal bridge, VMV_{M} does not contain any MnM_{n}-vertex on the right boundary. As a consequence, there is a cutset in MnM_{n} separating VMV_{M} from the right boundary. This cutset, composed of MnM_{n}-edges, corresponds to a vertical fault line. For instance, in Figure 7(e) the blue solid LτL_{\tau}-path γ\gamma is a fault line defined by the above argument. ∎

Corollary 4.4.

Ω=CGCFLCR\Omega=C_{\text{{G}}}\cup C_{\text{{FL}}}\cup C_{\text{{R}}} is a partition of the state space.

Refer to caption
(a) MnM_{n}
Refer to caption
(b) τ\tau
Refer to caption
(c) LτL_{\tau}
Refer to caption
(d) MτM_{\tau}
Refer to caption
(e) Duality
Refer to caption
(f) The injective map
Figure 7:

Before moving on to prove Theorem 4.1, we introduce the notion of almost fault lines. A horizontal (or vertical) almost fault line is a self-avoiding LnL_{n}-path connecting a LnL_{n}-vertex on the left boundary of LnL_{n} to a LnL_{n}-vertex on the right boundary of LnL_{n} where all edges except for one are in LτL_{\tau}. Denote by CAFLC_{\text{AFL}} the set of states containing an almost fault line. Let CG\partial C_{\text{G}} be the set of states outside CGC_{\text{G}} which are one-flip away from CGC_{\text{G}} in the state space of G\mathcal{M}_{G}.

Lemma 4.5.

CGCFLCAFL\partial C_{\text{{G}}}\subset C_{\text{{FL}}}\cup C_{\text{{AFL}}}.

Proof.

If a state τCG\tau_{\partial}\in\partial C_{\text{G}} is not in CFLC_{\text{FL}} (does not have a fault line), then by Corollary 4.4 τCR\tau_{\partial}\in C_{\text{R}} since by definition it is outside of CGC_{\text{G}}. Because τ\tau_{\partial} is one move away from CGC_{\text{G}}, there exists a state τGCG\tau_{\text{G}}\in C_{\text{G}} such that flipping four monochromatic edges along a unit square ss on 2\mathbb{Z}^{2} yields τ\tau_{\partial}. We know that in τG\tau_{\text{G}} there exists a green cross and no red cross; in τ\tau_{\partial} there exists a red cross and no green cross. Then it must be true that the four edges along ss are all green in τG\tau_{\text{G}} and all red in τ\tau_{\partial}. See Figure 8 for a pictorial illustration. Moreover, any green cross in τG\tau_{\text{G}} must contain edges along ss and so is any red cross in τ\tau_{\partial}; otherwise, a green horizontal (vertical) bridge in τG\tau_{\text{G}} must go across a red vertical (or horizontal, respectively) bridge in τ\tau_{\partial} at some other place on Λn\Lambda_{n}, which is impossible.

Therefore, there exists a simple green path ΓG\Gamma_{\text{G}} from some vertices on ss to the top boundary of Λn\Lambda_{n} and a simple red path ΓR\Gamma_{\text{R}} from some vertices on ss to the top boundary of Λn\Lambda_{n}. In the following, we prove that the above conditions suffice to show that there exists an LτL_{\tau}-path from ss to the top boundary of LnL_{n}. Similar conclusions can be made for the existence of LτL_{\tau}-paths from ss to the bottom boundary of LnL_{n}. By adding at most one LnL_{n}-edge, we can concatenate these paths to obtain a vertical almost fault line.

ΓG\Gamma_{\text{G}} and ΓR\Gamma_{\text{R}} cannot cross each other (Figure 5). Without loss of generality, suppose ΓG\Gamma_{\text{G}} is to the left of ΓR\Gamma_{\text{R}}. Then we use the medial lattice view MτM_{\tau_{\partial}} as in the proof of Lemma 4.3. ΓG\Gamma_{\text{G}} corresponds to a connected component mGm_{\text{G}} in MτM_{\tau_{\partial}}. Denote by VGV_{\text{G}} the set of MnM_{n}-vertices which can be reached by mGm_{\text{G}} in MτM_{\tau_{\partial}}, and let VGV^{\prime}_{\text{G}} be VGV_{\text{G}} together with the set of MnM_{n}-vertices which are separated from the right boundary of MnM_{n} by VGV_{\text{G}}. Then there is a cutset in MnM_{n} separating VGV^{\prime}_{\text{G}} from the right boundary. This cutset, composed of MnM_{n}-edges, corresponds to a dual LτL_{\tau}-path from ss to the top boundary of LnL_{n}. ∎

Refer to caption
(a) τG\tau_{\text{G}}
Refer to caption
(b) τ\tau_{\partial}
Figure 8: A step in G\mathcal{M}_{G}.

Now we are ready to show Theorem 4.1. Let ΩLEFT=CG\Omega_{\text{LEFT}}=C_{\text{G}}, ΩMIDDLE=CFL(CAFLCR)\Omega_{\text{MIDDLE}}=C_{\text{FL}}\cup(C_{\text{AFL}}\cap C_{\text{R}}), and ΩRIGHT=CRCAFL\Omega_{\text{RIGHT}}=C_{\text{R}}\setminus C_{\text{AFL}}. Theorem 4.1 is a consequence of the following lemma which uses a Peierls argument to show that π(CFLCAFL)\pi(C_{\text{FL}}\cup C_{\text{AFL}}) is exponentially small.

Lemma 4.6.

π(CFLCAFL)O(n)(2.639max(a,b)c)n\pi(C_{\text{{FL}}}\cup C_{\text{{AFL}}})\leq O(n)\left(\frac{2.639\max(a,b)}{c}\right)^{n}.

Proof.

For a self-avoiding path γ\gamma in LnL_{n} connecting a vertex on the top boundary to a vertex on the bottom boundary, denote by FγF_{\gamma} the set of states in Ω\Omega that contain γ\gamma as vertical fault line or almost fault line. Reversing directions of all the edges to the left side of γ\gamma defines an injective mapping from FγF_{\gamma} to ΩFγ\Omega\setminus F_{\gamma} that magnifies probability by a factor of at least min(a,b)c(cmax(a,b))|γ|1\frac{\min(a,b)}{c}\cdot\left(\frac{c}{\max(a,b)}\right)^{|\gamma|-1}. This is because: if γ\gamma is a fault line in a state, every Λn\Lambda_{n}-vertex sitting on γ\gamma would have four incident edges in the same color after the map, which increase its weight to cc (from aa or bb); if γ\gamma is an almost fault line in a state, the above is true except that for one Λn\Lambda_{n}-vertex, its weight decrease from cc to aa or bb after the map, as orientations on half of the four monochromatic edges are reversed. This indicates that π(Fγ)cmin(a,b)(max(a,b)c)|γ|1\pi(F_{\gamma})\leq\frac{c}{\min(a,b)}\cdot\left(\frac{\max(a,b)}{c}\right)^{|\gamma|-1}. See Figure 7(f) for an example. The same goes for horizontal (almost) fault lines.

Since every (almost) fault line is a self-avoiding walk on LL, the number of fault lines of length ll is upper bounded by 2n2n times the number of self-avoiding walks of that length starting at a vertex on the left or bottom boundary. The latter can be bounded by an well-studied estimate μl\mu^{l} on the number of self-avoiding walks of length ll on 2\mathbb{Z}^{2}, where μ\mu is called the connective constant. The best proved bound is μ2.638158\mu\approx 2.638158\cdots [GC01]. Summing this over fault lines of length from nn to n2n^{2} completes the proof. ∎

We have proved the torpid mixing of Glauber dynamics for the six-vertex model on the lattice region Λn\Lambda_{n} with free boundary conditions. Next we state the idea to extend the proof for the case when the Markov chain is D\mathcal{M}_{D} and the case when the boundary of Λ~n\widetilde{\Lambda}_{n} is periodic. Theorem 1.2 is a combination of Theorem 4.1 and the extensions.

Refer to caption
Figure 9: A state in D\mathcal{M}_{D}.

To extend Theorem 4.1 to hold for the directed-loop algorithm 𝒟\mathcal{M_{D}} whose state space is ΩΩ\Omega\cup\Omega^{\prime}, we need to pay extra attention for the states in Ω\Omega^{\prime}, the near-perfect Eulerian orientations. For a state τΩ\tau^{\prime}\in\Omega^{\prime}, there are two “defects” on the edges (Figure 9). Apart from the diagonal LnL_{n}-edges, there are two (2+(12,12))\left(\mathbb{Z}^{2}+\left(\frac{1}{2},\frac{1}{2}\right)\right)-edges separating green (half-)edges from red (half-)edges. The adaption we make is to put such (2+(12,12))\left(\mathbb{Z}^{2}+\left(\frac{1}{2},\frac{1}{2}\right)\right)-edges also into the set LτL_{\tau^{\prime}}. Notice that a connected component in LτL_{\tau^{\prime}} could possibly lie on L0L_{0}-edges as well as L1L_{1}-edges. Everything we prove is still correct if we allow fault lines to have (2+(12,12))\left(\mathbb{Z}^{2}+\left(\frac{1}{2},\frac{1}{2}\right)\right)-edges. Due to the possible positions of such two defects, the weight of CFLCAFLC_{\text{FL}}\cup C_{\text{AFL}} only increases by a polynomial factor in nn, hence not affecting the fact of torpid mixing.

To extend Theorem 4.1 to hold for Λ~n\widetilde{\Lambda}_{n} with periodic boundary condition (i.e., a 2-dimensional torus), we make the following modification. When nn is even, there still are two states with maximum weights cn2c^{n^{2}} (similar to τG\tau_{\text{G}} and τR\tau_{\text{R}} in Figure 4). Again, for any state τ\tau, Λ~n\widetilde{\Lambda}_{n}-edges can be classified as green or red, and its associated LτL_{\tau} separates Λ~n\widetilde{\Lambda}_{n}-edges of different colors. For the 2-dimensional torus T2T^{2}, the homology group H1(T2)×H_{1}(T^{2})\cong\mathbb{Z}\times\mathbb{Z}. We say a state τ\tau has a green (or red) cross if there are two non-contractable cycles of green (or red, respectively) edges of homology classes (a1,b1)(a_{1},b_{1}) and (a2,b2)(a_{2},b_{2}) with det[a1a2b1b2]0\operatorname{det}\left[\begin{smallmatrix}a_{1}a_{2}\\ b_{1}b_{2}\end{smallmatrix}\right]\neq 0; τ\tau has a pair of fault lines if there is a pair of non-contractable cycles of LτL_{\tau}-edges. (By parity, if there is one LτL_{\tau}-cycle there must be two.) Then the proofs in this section can be naturally adapted for the torus case, and the torpid mixing result follows.

Acknowledgement

The author thanks Professor Jin-Yi Cai for the valuable suggestions on the preliminary version of this paper.

References

  • [AR05] David Allison and Nicolai Reshetikhin. Numerical study of the 6-vertex model with domain wall boundary conditions. Annales de l’Institut Fourier, 55(6):1847–1869, 2005.
  • [Bax82] R. J. Baxter. Exactly Solved Models in Statistical Mechanics. Academic Press, 1982.
  • [BCK+99] Christian Borgs, Jennifer T. Chayes, Jeong Han Kim, Alan Frieze, Prasad Tetali, Eric Vigoda, and Van Ha Vu. Torpid mixing of some Monte Carlo Markov chain algorithms in statistical physics. In Proceedings of the 40th Annual Symposium on Foundations of Computer Science (FOCS), pages 218–229, 1999.
  • [BGRT13] Antonio Blanca, David Galvin, Dana Randall, and Prasad Tetali. Phase coexistence and slow mixing for the hard-core model on 2\mathbb{Z}^{2}. In Approximation, Randomization, and Combinatorial Optimization. Algorithms and Techniques (APPROX/RANDOM), pages 379–394, 2013.
  • [BKW73] H. J. Brascamp, H. Kunz, and F. Y. Wu. Some rigorous results for the vertex model in statistical mechanics. Journal of Mathematical Physics, 14(12):1927–1932, 1973.
  • [BN98] G. T. Barkema and M. E. J. Newman. Monte Carlo simulation of ice models. Phys. Rev. E, 57:1155–1166, Jan 1998.
  • [CGMS96] F. Cesi, G. Guadagni, F. Martinelli, and R. H. Schonmann. On the two-dimensional stochastic Ising model in the phase coexistence region near the critical point. Journal of Statistical Physics, 85(1):55–102, Oct 1996.
  • [CLL17] Jin-Yi Cai, Tianyu Liu, and Pinyan Lu. Approximability of the six-vertex model. CoRR, abs/1712.05880, 2017.
  • [Elo99] K. Eloranta. Diamond Ice. Journal of Statistical Physics, 96:1091–1109, September 1999.
  • [FW70] Chungpeng Fan and F. Y. Wu. General lattice model of phase transitions. Phys. Rev. B, 2:723–733, Aug 1970.
  • [GC01] A.J. Guttmann and A.R. Conway. Square lattice self-avoiding walks and polygons. Annals of Combinatorics, 5(3):319–345, Dec 2001.
  • [GMP04] Leslie Ann Goldberg, Russell Martin, and Mike Paterson. Random sampling of 3-colorings in 2\mathbb{Z}^{2}. Random Structures & Algorithms, 24(3):279–302, 2004.
  • [Lie67a] Elliott H. Lieb. Exact solution of the FF model of an antiferroelectric. Phys. Rev. Lett., 18:1046–1048, Jun 1967.
  • [Lie67b] Elliott H. Lieb. Exact solution of the two-dimensional Slater KDP model of a ferroelectric. Phys. Rev. Lett., 19:108–110, Jul 1967.
  • [Lie67c] Elliott H. Lieb. Residual entropy of square ice. Phys. Rev., 162:162–172, Oct 1967.
  • [LKV17] I Lyberg, V Korepin, and J Viti. The density profile of the six vertex model with domain wall boundary conditions. Journal of Statistical Mechanics: Theory and Experiment, 2017(5):053103, 2017.
  • [LPW06] David A. Levin, Yuval Peres, and Elizabeth L. Wilmer. Markov chains and mixing times. American Mathematical Society, 2006.
  • [LRS01] Michael Luby, Dana Randall, and Alistair Sinclair. Markov chain algorithms for planar lattice structures. SIAM Journal on Computing, 31(1):167–192, 2001.
  • [LS12] Eyal Lubetzky and Allan Sly. Critical Ising on the square lattice mixes in polynomial time. Communications in Mathematical Physics, 313(3):815–836, Aug 2012.
  • [MO94a] F. Martinelli and E. Olivieri. Approach to equilibrium of Glauber dynamics in the one phase region. I. the attractive case. Communications in Mathematical Physics, 161(3):447–486, Apr 1994.
  • [MO94b] F. Martinelli and E. Olivieri. Approach to equilibrium of Glauber dynamics in the one phase region. II. the general case. Communications in Mathematical Physics, 161(3):487–514, Apr 1994.
  • [MW96] M. Mihail and P. Winkler. On the number of Eulerian orientations of a graph. Algorithmica, 16(4):402–414, Oct 1996.
  • [Pau35] Linus Pauling. The structure and entropy of ice and of other crystals with some randomness of atomic arrangement. Journal of the American Chemical Society, 57(12):2680–2684, 1935.
  • [Pei36] R. Peierls. Statistical theory of adsorption with interaction between the adsorbed atoms. Mathematical Proceedings of the Cambridge Philosophical Society, 32(3):471–476, 1936.
  • [Ran06] Dana Randall. Slow mixing of Glauber dynamics via topological obstructions. In Proceedings of the Seventeenth Annual ACM-SIAM Symposium on Discrete Algorithm (SODA), pages 870–879, 2006.
  • [RS72] Aneesur Rahman and Frank H. Stillinger. Proton distribution in ice and the Kirkwood correlation factor. The Journal of Chemical Physics, 57(9):4009–4017, 1972.
  • [RT00] Dana Randall and Prasad Tetali. Analyzing Glauber dynamics by comparison of Markov chains. Journal of Mathematical Physics, 41(3):1598–1615, 2000.
  • [SJ89] Alistair Sinclair and Mark Jerrum. Approximate counting, uniform generation and rapidly mixing Markov chains. Information and Computation, 82(1):93–133, 1989.
  • [Sut67] Bill Sutherland. Exact solution of a two-dimensional model for hydrogen-bonded crystals. Phys. Rev. Lett., 19:103–104, Jul 1967.
  • [SZ04] Olav F. Syljuåsen and M. B. Zvonarev. Directed-loop Monte Carlo simulations of vertex models. Phys. Rev. E, 70:016118, Jul 2004.
  • [YN79] A. Yanagawa and J.F. Nagle. Calculations of correlation functions for two-dimensional square ice. Chemical Physics, 43(3):329–339, 1979.