Numerically exact correlations and sampling in the two-dimensional Ising spin glass

A powerful existing technique for evaluating statistical mechanical quantities in two-dimensional Ising models is based on constructing a matrix representing the nearest-neighbor spin couplings and then evaluating the Pfaffian of the matrix. Utilizing this technique and other more recent developments in evaluating elements of inverse matrices and exact sampling, a method and computer code for studying two-dimensional Ising models is developed. The formulation of this method is convenient and fast for computing the partition function and spin correlations. It is also useful for exact sampling, where configurations are directly generated with probability given by the Boltzmann distribution. These methods apply to Ising model samples with arbitrary nearest-neighbor couplings and can also be applied to general dimer models. Example results of computations are described, including comparisons with analytic results for the ferromagnetic Ising model, and timing information is provided.

Figure 1

Correspondence between spin state configurations and complete dimer coverings on the decorated dual graph $G_D^{*}$ for a periodic spin lattice of size $L_x\times L_y=6\times 3$. A sample configuration of Ising spins $s_i=\pm 1$ are represented by the arrows inside the large circles. These spins are coupled by horizontal and vertical bonds of strength $J_{ij}$. The (red) dashed lines indicate the bonds for one example spin. The spins and bonds are the vertices and edges, respectively, of the Ising model graph $G$. The nodes of the decorated dual graph $G_D^{*}$ are drawn as small circles and the edges are indicated by the thin and thick solid lines. The edges of $G_D^{*}$ are either internal to a Kasteleyn city (the sets of four fully connected nodes) or connect neighboring cities. Those that are internal to a city have a weight of 1 while those connecting cities have a weight $w=\exp (-2\beta J_{ij})$, where the $J_{ij}$ is the coupling strength of the bond crossing the dual edge. An example of a complete dimer covering $M$ corresponding to the displayed spin configuration is indicated by the heavy lines: Such a choice of edges includes all nodes in $G_D^{*}$ exactly once. The intercity edges belonging to the covering $M$ separate the up spins from the down spins and so compose the relative domain walls (here we are assuming that the reference configuration $S^r$ is the configuration with all spins up, $s_i=+1$). Note that for any Kasteleyn city surrounded by four spins of identical sign, there are three ways to arrange the dimers on that city. Two examples of these arrangements can be seen in the lower left and lower right cities. In other cases, the choice of Kasteleyn city edges is uniquely determined by the domain walls.

Figure 2

Depictions of the geometric dissection tree $T$ for an example decorated dual graph. The original Ising system has $L_x\times L_y=3\times 2$ spins on a periodic graph; the spins are indicated by gray circles. The decorated dual graph has $3\times 2$ Kasteleyn cities (the sets of four fully connected nodes). At each stage, a parent rectangle of Kasteleyn cities is divided into two roughly equal sibling rectangles, the children. In the algorithms described here, each rectangle $A$ has geometric information, two cluster matrices $U_A$ and $U_{\overline{A}}$, and a partial Pfaffian factor $z\left(A\right)$. (a) The nested dissection in real space. The two largest subregions with boundaries, $A^{*}$ and $B^{*}$, indicated. The region $C^{*}$ in the final stage (not shown for clarity) is the union of $A^{*}$ and $B^{*}$, $C^{*}=A^{*}\cup B^{*}$; at the end of the algorithm, it has no boundary but a matrix corresponding to this region is used as a working matrix when using periodic boundaries. (b) A diagram of the resulting nested dissection tree $T$. The leaves of the tree are Kasteleyn cities. The root of the tree has no cluster matrix associated with it, as the sample has no boundary, though it has a Pfaffian factor associated with it, which is used to find the partition function for the whole sample.

Figure 3

Depiction of merger of Kasteleyn cities $A$ and $B$ connected by an edge of weight $w=\exp (-2\beta J_{ij})$, with $w$ set to $0.5$. (a) The rows and columns of the $U_A$ and $U_B$ are indexed by nodes $\{0,1,2,3\}$ and the upper triangular part of these skew-symmetric matrices is shown with the diagonal of zeros omitted. For example, the upper left elements shown in the matrices here are in row 0 and column 1. The signs of the connections correspond to a Pfaffian orientation [3] to consistently count spin configurations. (b) The matrix $M_0$ is the result of collecting $U_A$ and $U_B$ and the edge weight $w=0.5$ and is indexed by nodes 0 through 7, with indices from $A$ for the initial part $\{0,1,2,3\}$ and indices from $B$ for the second half $\{4,5,6,7\}$. Permuting rows and columns 0 and 1 and then 1 and 7 gives the matrix $M_1$. These permutations place the nodes for the edge to be eliminated in the first two rows of $M_1$. (c) After using Pfaffian elimination to remove rows and columns 0 and 1, one is left with a superdiagonal element at $(0,1)$ of 0.5 (no pivoting is possible in this case) giving a partial Pfaffian $z\left(C\right)=0.5$ and the next generation cluster matrix $U_C$, indexed by the remaining six nodes numbered as shown.

Figure 4

Example of a higher level step of the up-sweep stage. The Ising spins on the original lattice whose couplings $J_{ij}$ are relevant to the calculations through this step are indicated by the large gray circles, while the nodes of the decorated dual lattice are indicated by the medium-sized and small circles. The solid straight and jagged lines indicate edges belonging to the decorated dual graph. Two geometric regions on the dual lattice, $A$ and $B$, each containing $3\times 3$ Kasteleyn cities, are denoted by the (blue) dashed squares. The dual edges $W$ that separate $A$ and $B$ are drawn as jagged lines and join the filled-in medium-sized nodes. The small nodes are interior to the regions $A$ and $B$ while the larger nodes (medium-sized open and solid circles) form the borders of $A$ and $B$. In the up-sweep stage (see Sec. 3c), given the cluster matrices $U_A$ and $U_B$ for $A$ and $B$, the separating edges belonging to the separator $W$ are integrated out via Pfaffian elimination, leaving a cluster matrix for $C$ that includes the entire subgraph shown. The cluster matrix $U_C$ is indexed by its border nodes, i.e., the medium-sized open circles connected by the dashed (blue) rectangle.

Figure 5

Diagram of a sample merging of cluster matrices in the down-sweep stage for the regions indicated in Fig. 4. The gray and black larger circles indicate the locations of the Ising spins in the original square grid. (The black spins are the corner spins for region $C$.) The matrix $U_{\overline{C}}$ describes the dimer correlations between the nodes that touch the outer (red) dashed line labeled $\overline{C}$. This matrix sums over correlations external to the spins shown. This cluster matrix $U_{\overline{C}}$ is merged with $U_B$ by eliminating the edges shown by the jagged solid lines. The region $B$ is indicated by the dashed (blue) square on the right of the diagram. The result of the merger is the matrix $U_{\overline{A}}$ describing correlations among the nodes on the boundary of $\overline{A}$, which is shown by the labeled square (red) dashed line.

Figure 6

Diagram of a calculation of a correlation between spins $i$ and $j$. This calculation uses information computed during the both the up-sweep and down-sweep stages. The matrices $U_{\overline{C}}$ and $U_C$ are merged using the edge weights on the jagged lines. These edges form the set $W$. Two mergings are calculated: one for all positive weights for the edges in $W$ and one where the weights are negated for the thicker (green) edges. This gives two Pfaffians for the whole sample. The first Pfaffian has positive contributions from configurations with either an odd or even number of dimer choices (i.e., domain walls) between spins $i$ and $j$. This is the partition function for the whole sample. The second Pfaffian is the difference between the partition function constrained to have an even number of domain walls between $i$ and $j$ and the partition function constrained to have an odd number of domain walls between $i$ and $j$.

Figure 7

Steps for exact sampling on a periodic Ising spin lattice of size $L_x\times L_y=8\times 8$. Empty circles indicate unknown spins. Arrows indicate assigned spins: $\uparrow$ for $s=+1$ and $\downarrow$ for $s_i=-1$. Separators $W$ are drawn as dashed light (red) lines. The regions of the graph $G_D^{*}$ used at each stage are shown by solid darker lines. (a) The first region $C^{*}$ is shown. The separator $W$ (these wrapping edges are drawn as separated half edges) connects $C^{*}$ to itself. Given a seed spin, deciding which edges in $W$ are in the matching $M$ along the bottom row fixes the spins for step (b). (b) The spins decided in (a) fill the lowest row. The remaining separator along the column is to be filled in for the start of the next step. (c) In all subsequent steps, including this step, the boundary conditions are fixed. The region $C^{*}$ is separated into the lower half $A^{*}$ and the upper region $B^{*}$. Choices are made for the dual edges connecting $A^{*}$ to $B^{*}$. (d) Three edges are chosen for each of the two separators to fix the three spins for each region pair. (e) Four separators are used to set 12 spins. (f) Eight separators are used to set eight spins. (g) In this final stage, there are 16 separators. One edge choice is made for each, fixing the remaining undecided spins. (h) The final spin assignment.

Figure 8

Plot of correlation functions for the ferromagnetic Ising model at criticality. The correlation function $\langle s(0,0)s(R,R)\rangle$ for diagonal spin separations $(R,R)$ was computed using Pfaffian methods for samples of size $L\times L$ spins, $L=8,...,128$, and compared with $L=\infty$ exact and asymptotic results.