History and significance

Perfect matchings illustrate the distinction between finding a combinatorial structure and counting its realizations. Edmonds gave a polynomial-time algorithm for finding a maximum matching in a general graph [5]. Valiant proved that exact counting is #P-complete even for bipartite graphs, where it is the permanent of a zero–one matrix [19]. Pfaffian methods give polynomial-time exact counting for planar graphs [13].

For approximation, three problems must be distinguished: counting all matchings, counting bipartite perfect matchings, and counting perfect matchings in general graphs. Jerrum and Sinclair gave an FPRAS for the total number of matchings in an arbitrary graph, summed over all cardinalities [10]. Their perfect-matching algorithm applied when the ratio of near-perfect matchings, which leave exactly two vertices unmatched, to perfect matchings was bounded by a fixed polynomial [10]. Without such a bound, perfect matchings can carry too little stationary mass for this sampling approach. Approximating the total number of matchings therefore does not isolate its perfect-matching contribution with relative accuracy.

Jerrum, Sinclair, and Vigoda removed the ratio restriction in the bipartite setting and obtained an FPRAS for the permanent of every nonnegative matrix [11]. Their method assigns weights to near-perfect states according to their two unmatched vertices and learns suitable weights while gradually changing edge activities. The general-graph question was explicitly raised by Jerrum and Sinclair [10] and remained a benchmark for approximate counting. Cai and Liu relate it to approximation for the eight-vertex model [2]; Fei, Goldberg, and Lu distinguish the spin-system parameter line equivalent to general perfect-matching approximation from the region covered by their FPRAS [7]. 1 answers the general-graph perfect-matching approximation problem affirmatively.

A direct extension of the bipartite chain faces a precise obstruction. Štefankovič, Vigoda, and Wilmes constructed graphs on which every JSV-type chain has either exponentially small stationary probability of perfect matchings or exponentially slow mixing, even with arbitrary weights depending on the hole pattern [16]. Here JSV type means the perfect-and-near-perfect state space with the specified Broder proposals and a Metropolis acceptance rule. The sampler below uses a different state space: a product of perfect-matching spaces on an enlarged colored graph, with exchanges of alternating cycles between coordinates. Its energy comparison is proved for these states and moves.

Subsequent permanent algorithms sharpen the bipartite result. Bezáková, Štefankovič, Vazirani, and Vigoda accelerated the cooling schedule and improved the analysis of the matching chain [1]. Chen, Vigoda, and Yang obtain a further improvement through restricted Poincaré inequalities and coupled flows [4]; Chen, Guo, Vigoda, and Yang develop an energy-weighted multicommodity-flow comparison and improve the running time again [3]. These results concern bipartite perfect matchings.

For nonbipartite graphs, Yi gives deterministic relative approximation for hafnians, the weighted sums over perfect matchings, when the \(n\)-vertex support has minimum degree at least \((1/2+\gamma)n\) and nonzero weights lie in \([\theta,1]\) [21]. Here \(0<\gamma<1/2\) and \(0<\theta\le1\) are fixed constants, and the polynomial running-time exponent may depend on them. 1 has no density promise on its simple unweighted input.

The proof below also builds on the methodological framework of weighted matching chains. Alternating overlays and weight-preserving switches underlie the encodings of Jerrum and Sinclair and of Jerrum, Sinclair, and Vigoda [10, 11]. Their sampling-to-counting strategy learns local partition ratios while gradually changing activities. Here these ratios are exposed by marked edges in an enlarged graph, and vertex scaling maintains a row-wise balance condition. The capacity bound, subdivision and marked-edge identities needed for this construction are proved explicitly.

Energy comparison through weighted demands has classical antecedents in multicommodity-flow bounds for mixing [15]. Our comparison differs from routing every demand along a path of chain transitions: long-arc conversions enter a signed identity and need not be legal single updates. The resulting estimate initially controls only additive functions of two matching coordinates. The separate product-space argument is essential to obtain a spectral gap for arbitrary functions. It uses the classical orthogonal decomposition into coordinate interactions [6]; the ordered assignment of residual terms and the replication estimate are proved below.

Proof architecture

The counting algorithm starts with the complete graph, whose perfect matchings are easy to count, and gradually suppresses the nonedges of \(G\) by decreasing their positive activities. It estimates the ratios of successive partition functions and multiplies them. To keep the observations statistically useful, it updates vertex scales from estimates of the two-hole partition ratios: partition functions after deleting two vertices, divided by the current partition function. Accurate estimates keep each row maximum between two absolute constants. A vertex scale multiplies every complete matching weight by the same factor, preserving the probability law used to estimate the current counting ratio.

The purpose of balance is to control the extra mass introduced by an auxiliary graph. We replace each logical edge by a weighted path, with a bijection between logical perfect matchings and matchings using only these real path edges. We then add weak clique edges in overlapping bags indexed by a tree. Under balance, the added edges increase the partition function by at most a factor of two, so a sample uses only real edges with probability at least one half. A specially marked auxiliary edge exposes the corresponding two-hole partition ratio. These are exact event identities; they provide both the counting estimate and the data for the next scale update.

To prove the inflation bound, 3 assigns a pairing capacity to each even hole set using the strongest bottleneck paths between its vertices. The capacity pays for the weights of possible added edges on those holes. Subdivision transfers this bound to the real graph while making incident effective strengths comparable. The resulting geometry supplies the tree of clique bags. Thus the hole estimates control the probability of the observable events, while the tree-bag structure supplies a sampler for their law.

That sampler uses local switches in one matching and exchanges alternating cycles between two matching coordinates. Its energy comparison initially controls only additive functions of the two coordinates. A quadrangulation adapted to the bag labels gives an exact signed identity for a whole-cycle difference ([lem:quadrangulation,thm:pair-energy]). The identity is not a route through chain transitions: long arc conversions may change many edges. Two-matching encodings instead compare its terms with local switches and cycle exchanges, with polynomially bounded multiplicity.

A ladder of nearby weight distributions then places many independent coordinates at each tier. An ordered orthogonal decomposition ensures that each coordinate-interaction term is charged to at most one pair. Replication reduces the coefficient of their combined contribution enough to absorb it. This yields a spectral gap for arbitrary functions of the product state (20). The bottom tier has unit activities and admits exact counting and ideal sampling by a tree dynamic program. The sampler requires no balance assumption. It therefore remains defined, with the same time bound, after an inaccurate empirical estimate; balance is used only to prove that accurate estimates propagate through the counting stages.

1 summarizes the changes of law. Finally, [sec:counting,sec:complexity] account for all arithmetic and finite-bit random choices. Every polynomial has an absolute degree; the constants are chosen for transparent estimates rather than practical efficiency.

The specialized combinatorial and sampling arguments are proved here. We use Edmonds’s algorithm for the initial matching-existence test and Hoeffding’s bounded-variable inequality for concentration and confidence amplification [5, 9].

Weighted matchings and holes

We allow parallel edges distinguished by colors. A matching specifies the actual colored edges it uses. For positive rational activities \(a=(a_e)\), put \[\operatorname{wt}_a(M)=\prod_{e\in M}a_e,\qquad Z_a(-U)=\sum_{M\in\mathcal M(H-U)}\operatorname{wt}_a(M), \qquad Z_a=Z_a(-\varnothing),\] where \(H\) is the underlying graph and \(\mathcal M(H-U)\) is its set of perfect matchings after deleting \(U\). The empty graph has one perfect matching, of weight one. Whenever \(Z_a>0\), define \[g_a(U)=\frac{Z_a(-U)}{Z_a},\qquad \pi_a(M)=\frac{\operatorname{wt}_a(M)}{Z_a}.\] These weighted colored graphs are intermediate spaces constructed from the simple input graph. Their representation is explicit: vertices, colored edges, and rational activities are stored individually. The main theorem retains the finite simple undirected input model; it does not assert a counting algorithm for compressed edge multiplicities or for singleton-covering loops.

We call the vertices in \(U\) holes. For distinct vertices we abbreviate \(g_a(\{i_1,\ldots,i_k\})\) to \(g_a(i_1\cdots i_k)\). Collapsing parallel activities into their sum preserves every partition function. Conditional on a collapsed edge, its original type can be recovered with probabilities proportional to the summands.

The union of two matchings, with common edges retained twice, has maximum degree two. Its nontrivial components are alternating paths and even alternating cycles. Swapping the layers along a component preserves the product of their weights. In the union of two perfect matchings, deleting common identical edges leaves disjoint alternating cycles; different colors on the same pair of vertices form a cycle of length two. We fix a vertex order and a color order once and for all. They resolve every combinatorial tie in the proof.

Bottleneck similarities

For a positive weighted connected graph define \[ B_{ij}=\max\left\{1,\ \max_{P:i\leadsto j}\ \min_{e\in P}a_e\right\} \quad (i\ne j). \tag{1}\] Simple paths suffice. This symmetric similarity satisfies \(B_{ij}\ge\min\{B_{ik},B_{kj}\}\): concatenating two paths and erasing loops cannot lower their minimum activity. For \(t>1\), its equivalence classes are the connected components of the graph retaining edges of activity at least \(t\). Singletons are included. At \(t=1\) we take one class containing every vertex, consistently with the floor at one. For an even set \(U\) define its pairing capacity \[ \Phi_B(U)=\max_{\mathcal P}\prod_{ij\in\mathcal P}B_{ij}, \qquad \Phi_B(\varnothing)=1, \tag{2}\] where \(\mathcal P\) runs over pairings of \(U\). We use this quantity only in proofs; the counting algorithm never needs to enumerate pairings or compute \(\Phi_B(U)\) for a general hole set.

Tree bags

Every vertex pair has at most two edge colors. At any bag, the interfaces with its parent and its different children are disjoint: a vertex in two such interfaces would have at least three memberships. Empty bags and repeated underlying vertex sets are permitted. Clamping is monotone and does not increase the ratio of a larger activity to a smaller one. In particular, reference activities of nearby labels remain comparable after clamping.

Dirichlet energy and approximation

For a finite reversible Markov kernel \(P\) with stationary law \(\mu\), its Dirichlet energy is \[ \mathcal E_{\mu,P}(f) =\frac12\sum_{x,y}\mu(x)P(x,y)(f(x)-f(y))^2. \tag{4}\] A Poincaré inequality \(\mathop{\mathrm{Var}}_\mu f\le\gamma^{-1}\mathcal E_{\mu,P}(f)\) means that the spectral gap is at least \(\gamma\). A symmetric proposal \(Q\) has Metropolis acceptance probability [14, 8] \(\min\{1,\mu(y)/\mu(x)\}\) for a proposed move from \(x\) to \(y\). Its off-diagonal stationary capacity is \[ \mu(x)P(x,y)=Q(x,y)\min\{\mu(x),\mu(y)\}. \tag{5}\] We regard inapplicable proposals as holding moves. All state spaces used below contain at least one perfect matching; thus their positive activities define probability laws with full support on that space.

Total variation is \(\|\mu-\rho\|_{\mathop{\mathrm{TV}}}=\frac12\sum_x|\mu(x)-\rho(x)|\). A deterministic map, or a common randomized postprocessing kernel, does not increase total variation, by the triangle inequality. Two finite laws at total variation distance \(\eta\) admit a coupling that disagrees with probability \(\eta\): first match the common masses \(\min(\mu(x),\rho(x))\), then couple the remainders. Applying this construction successively bounds the discrepancy probability of an adaptive simulation by the sum of its conditional one-step errors.

The path profile

Orient each logical edge \(e=uv\) by the fixed vertex order, and replace it by a path \[x_0=u,x_1,\ldots,x_{4p},x_{4p+1}=v.\] Internal vertices of different paths are distinct. Denote the activity of \(x_{j-1}x_j\) by \(t_j=t_j^e\). For \(1\le r\le p\), set \[\begin{align*} t_{2r-1}=t_{2r} &=\max\{B_{uv},H_u2^{-(r-1)}\},\\ t_{4p+3-2r}=t_{4p+2-2r} &=\max\{B_{uv},H_v2^{-(r-1)}\},\tag{16}\\ t_{2p+1}&=\lambda_{uv}. \end{align*}\] We call these the real edges and write \(Z'\) and \(g'\) for their partition function and hole ratios. Let \(B'\) be their bottleneck similarity. The two innermost equal pairs have activity \(B_{uv}\), because \(H_u2^{-(p-1)},H_v2^{-(p-1)}\le1/2\).

The enlarged matching law

For each integer \(0\le k\le\lfloor\log_2\max_uH_u\rfloor\), take one label for every component of the enlarged graph whose retained edges have effective strength at least \(2^k\), including singleton components. A label at level \(k>0\) has as parent its containing component at level \(k-1\). At level zero there is one component; these labels therefore form a rooted tree. Components recurring at different levels are still different labels.

Label each real edge \(xy\) by the component containing it at level \(\lfloor\log_2B'_{xy}\rfloor\). A vertex belongs precisely to the labels of its incident real edges. Equal-level incident labels are identical, and the factor-two bound in 6 implies that there are at most two incident levels. If there are two, they are consecutive and their labels are parent and child. Thus these memberships satisfy 2. In particular, a bag is the set of vertices with that membership; it need not be the entire threshold component defining its label.

Within every bag of label \(l\) and height \(h_l\), add a virtual edge of activity \(h_l/D\) for each vertex pair. Collapse any real edge of color \(l\) into the virtual edge with the same endpoints and color, adding its activity. This gives one colored edge per common label, of activity \[W_{xy,l}=h_l/D+ \begin{cases} t_{xy},&xy\text{ is a real edge with label }l,\\ 0,&\text{otherwise}. \end{cases}\] Since a real edge assigned to \(l\) has \(t_{xy}\le B'_{xy}<2h_l\), these weights satisfy [eq:bag-weight-range]. Also \(1/D\le W_{xy,l}\le3\cdot4^K\).

For each logical edge \(e\), split off activity \(1/D\) from the virtual activity on its central colored edge, and mark that portion as the probe for \(e\). This is possible because \(h_l\ge1\). After sampling a colored matching, independently reveal the portion of each used edge in proportion to its real, probe, and remaining virtual activities; zero portions are omitted.

1 shows the two real tilings and the probe on one illustrative path.

Let \(Z_{\mathrm{bag}}\) be the colored partition function of the enlarged graph, and put \(I=Z_{\mathrm{bag}}/Z'\). This ratio measures the mass introduced by the auxiliary edges.

The identities show what the enlarged graph provides to the counting algorithm. To use its all-real event with polynomially many samples, we also need a uniform lower bound on that event’s probability. The next lower-bound estimate requires balance; the construction and the exact laws above do not.

We prove this bound in the remainder of the section. A collection of virtual edges leaves a set \(U\) of endpoints to be removed from a real matching. Its added weight must therefore be compared with \(g'(U)\). The next two subsections transfer the logical hole estimate to these enlarged hole sets.

Exact accounting for deleted vertices

We next describe exactly what deleting internal vertices does to the logical graph. A path is broken if at least one of its internal vertices is deleted. Let \(E_{\mathrm{br}}(U)\) be the set of logical edges whose paths are broken by a hole set \(U\).

Consider the ordered internal holes on a broken path. If a real completion exists, successive hole indices alternate in parity, because the intervening internal run must have even order. A first hole of even index forces the prefix to match the terminal \(u\); we say that this hole consumes \(u\) and map it to \(u\). A last hole of odd index similarly consumes and maps to \(v\). Remove all consumed holes (zero, one, or two) from the list. The remaining holes split, in their order, into pairs \((j,k)\) with \(j\) odd and \(k\) even; call these neutral pairs. There is no ambiguity when the list has one hole, since it cannot both have even and odd index.

Let \(R\) be the set of explicitly deleted terminals together with all consumed terminals. Feasibility requires these terminals to be distinct. To each consumed hole \(x_j\) assign the factor \[f_j= \begin{cases} 1,&\text{if it maps to its home terminal},\\ \lambda_e/T_j,&\text{if it maps to the other terminal}. \end{cases}\] For a neutral pair \((j,k)\) assign \[ f_{jk}= \begin{cases} T_j^{-1},&j,k\le2p,\\ T_k^{-1},&j,k\ge2p+1,\\ \lambda_e/(T_jT_k),&j\le2p<k. \end{cases} \tag{20}\] Let \(F=F(U)\) be the product of all these factors on all broken paths; an empty product is one. Both \(R\) and \(F\) depend only on \(U\) and the fixed path profiles, not on a choice of completion.

For a set \(E_0\) of logical edges, write \(Z_\lambda(-R;E_0\text{ forbidden})\) for the partition function after deleting \(R\) and forbidding those edges.

Charging the pairing capacity

The deletion identity gives a factor \(F(U)\) and a logical hole set \(R\). To transfer the logical bound, we must compare the pairing capacity of \(U\) with that of \(R\). The next lemma expresses both capacities through threshold classes. Throughout the comparison, these classes belong to the original undeleted graphs, so they are independent of the completion being counted.

The pair kernel

On \((A,B)\in\mathcal M^2\), choose each of the following proposal types with probability \(1/3\).

1. Two-edge switch. Choose an unordered pair of distinct edges of \(A\) uniformly, and choose either alternative pairing of their four endpoints with probability \(1/2\). For each proposed edge, choose index \(1\) or \(2\) with probability \(1/2\) from the ordered list of common labels of its endpoints. If an index is absent, stay put. If \(A\) has fewer than two edges, stay put.

2. Color change. Choose an edge of \(A\) uniformly and choose index \(1\) or \(2\) from the common-label list with probability \(1/2\). An absent index gives a holding move.

3. Cycle swap. Choose uniformly one discrepancy cycle of \(A,B\) and interchange its two layers. If there is no discrepancy cycle, stay put.

Accept each proposal by the Metropolis rule for \(\pi\otimes\pi\). Denote the resulting kernel by \(P\) and its energy by \(\mathcal E\). The proposal is symmetric. For the first two types, the reverse choice has the same edge, pairing, and color-index probabilities; the possible holding choices do not change this fact. For the third type, the discrepancy cycles and their number are unchanged by the swap. In particular, a cycle swap always has acceptance probability one, since its product weight is unchanged.

In the proof write \(F(A,B)=f_1(A)+f_2(B)\) and \(E=\mathcal E(F)\). A two-edge switch is called safe if its four removed or inserted colored edges have activities within a factor \(A_0\) of one another. Let \(\mathcal S(A)\) be its possible safe destinations, and let \(\mathcal K(A)\) be its possible destinations by changing one edge color. Consecutive labels have comparable heights, so every color change has activity ratio in \([A_0^{-1},A_0]\).

An exact identity on a cycle

Fix a pair of demand matchings \((I_1,I_2)\) and a simple discrepancy cycle \(C\) of length at least four. All orientations and choices below are deterministic, using the orders fixed in 2. Let \(R\) be the edges of \(I_1\) outside \(C\), and let \(O_0=I_1\) and \(O_1\) be the matching obtained by toggling \(C\) in \(I_1\). Thus both orientations use \(R\) outside \(C\). Put \[D_C=f_1(O_0)-f_1(O_1).\] Apply 12 and fix the colors of its diagonals consistently. A cell’s four sides always mean the four odd arcs of the original whole cycle, not paths formed from earlier diagonals.

Use the through, interior, and closed patterns \(T_J,N_J,C_J\) from 5, with the chosen chord colors. A side \(J\) is through in \(O_\varepsilon\) if \(T_J\subseteq O_\varepsilon\). Exactly one orientation \(O_\varepsilon\) contains \(T_J\); define \[ H_J=(O_\varepsilon\setminus T_J)\cup C_J, \qquad g_J=f_1(O_\varepsilon)-f_1(H_J). \tag{30}\] For a short side, the chord is its actual cycle edge, so \(g_J=0\).

The reason for using the original cycle on every cell side is already visible on a six-cycle \(C=(v_0,\ldots,v_5,v_0)\). Take \[O_0=R\cup\{v_0v_1,v_2v_3,v_4v_5\},\qquad O_1=R\cup\{v_1v_2,v_3v_4,v_5v_0\},\] and suppose the diagonal \(e=v_0v_3\) is admissible. In the two cells it represents the complementary arcs \(J=(v_3,v_4,v_5,v_0)\) and \(\bar J=(v_0,v_1,v_2,v_3)\), respectively. Both conversions yield the same matching \[H_J=H_{\bar J}=H=R\cup\{v_1v_2,e,v_4v_5\}.\] Thus \(g_{\bar J}-g_J=f_1(O_0)-f_1(O_1)=D_C\); the shared diagonal contributes one whole-cycle difference, not zero. 3 displays this common converted matching. This example isolates the complementary-gain identity. Its cells have only one long side each, so the context errors introduced below vanish; cells with two opposite long sides require a further estimate.

In a cell, let \(P\) be its opposite pair of sides through in \(O_0\), and \(Q\) its opposite pair through in \(O_1\). Let \(A_P\) use \(N_J\) on every side and the chords of the two sides of \(P\), together with \(R\). Define \(A_Q\) analogously and put \[\Delta=f_1(A_P)-f_1(A_Q).\] These are legal perfect matchings: the interior vertices are matched by the \(N_J\)’s, and the two chosen chords match the four corners. Their difference is a safe two-edge switch by 13.

Order each opposite pair as \((Y,X)\). For the orientation in which both are through, first convert \(Y\) and then \(X\). Write \(A^{(0)}=O_\varepsilon\) and \(A^{(1)}=H_Y\) for the two pre-conversion contexts of \(X\). For a matching \(R'\) on the vertices outside \(X\), define \[ \partial_X f(R')=f(T_X\cup R')-f(C_X\cup R'). \tag{31}\] The context error of the ordered opposite pair is \[ E_{Y,X} =\partial_X f_1(A^{(1)}\setminus T_X) -\partial_X f_1(A^{(0)}\setminus T_X). \tag{32}\] It is zero if either side is short. Denote these errors by \(E_P,E_Q\) for the two ordered pairs of the cell.

Encoding the demands

A switch demand specifies \((I_1,I_2)\), a simple discrepancy cycle \(C\), and one cell of its fixed quadrangulation. Its value is \(\Delta\) above. An error demand specifies these data and one opposite pair with two long sides, in its fixed order \((Y,X)\). Its value is \(E_{Y,X}\). The weight of either demand is \(\rho(d)=\pi(I_1)\pi(I_2)\). Repeated demands from different cells or different cycles are included separately in all sums.

The encodings provide two kinds of control. Their bounded preimage count converts a sum over demands to a sum over pairs of matchings. For context errors, the guide also contains a prescribed closed arc; this additional condition must survive that summation. We use the bounded preimage count for switches below, then retain the closed-arc condition when treating context errors in the next subsection.

We apply this lemma first to switches. For a fixed encoded pair \((A,B)\) there are at most \(\mathcal T\) demands, and each has its associated destination \(A_Q\) in \(\mathcal S(A)\). Therefore \[\begin{align*} \sum_{d\text{ switch}}\rho(d)\Delta_d^2 &\le A_0^4\mathcal T\sum_{A,B}\pi(A)\pi(B) \sum_{A'\in\mathcal S(A)}(f_1(A)-f_1(A'))^2 \\ &\le60N^2A_0^6\mathcal TE. \tag{35}\end{align*}\]

Controlling context errors by cycle swaps

Here an identity \(X\) means a directed, colored, simple odd path of length at least three together with a specified colored edge joining its endpoints. It determines \(T_X,C_X\) and (31). We may sum over all such identities: this is a finite proof device, and the algorithm does not enumerate them. Put \[p_X=\pi\{B:C_X\subseteq B\}.\] We omit identities with \(p_X=0\). Every identity arising from an error demand has \(p_X>0\), by the guide constructed in 17. For any first matching \(A\supseteq T_X\) and either \(j=0,1\), that lemma gives the following conditional demand-mass bound: \[ \sum_{\substack{d\text{ error}:X(d)=X\\ A^{(j)}(d)=A}} \rho(d) \le A_0^4\mathcal T\,\pi(A)p_X. \tag{36}\] Indeed, for each guide \(B\supseteq C_X\) there are at most \(\mathcal T\) preimages of \((A,B)\); now sum its product weight over these guides. The factor \(p_X\) is essential and will cancel the conditioning below.

The next comparison averages over precisely those guides counted by \(p_X\). That conditional law can be concentrated on a very rare event. We therefore retain \(p_X\) in the demand bound instead of estimating it from below.

For a fresh guide \(B\sim\pi\), define \[\mu_X=\mathbb E\bigl[\partial_X f_2(B\setminus C_X) \mid C_X\subseteq B\bigr].\] For an error demand, (32) and \((a-b)^2\le2(a-c)^2+2(b-c)^2\) imply \[E_d^2\le 2\sum_{j=0}^1 \bigl(\partial_X f_1(A^{(j)}\setminus T_X)-\mu_X\bigr)^2.\] Applying (36), and then Jensen’s inequality to the conditional expectation defining \(\mu_X\), gives \[\begin{align*} \sum_{d\text{ error}}\rho(d)E_d^2 &\le4A_0^4\mathcal T \sum_X\sum_{A\supseteq T_X}\pi(A)p_X \bigl(\partial_X f_1(A\setminus T_X)-\mu_X\bigr)^2 \\ &\le4A_0^4\mathcal T\,\Sigma_{\rm swap}, \tag{37}\end{align*}\] where \[ \Sigma_{\rm swap}= \sum_X\sum_{\substack{A\supseteq T_X\\ B\supseteq C_X}} \pi(A)\pi(B) \bigl(\partial_X f_1(A\setminus T_X) -\partial_X f_2(B\setminus C_X)\bigr)^2. \tag{38}\] In detail, for a fixed \(X,A\) Jensen bounds the squared difference from \(\mu_X\) by \[\frac1{p_X}\sum_{B\supseteq C_X}\pi(B) \bigl(\partial_X f_1(A\setminus T_X) -\partial_X f_2(B\setminus C_X)\bigr)^2.\] Multiplication by the \(p_X\) in (36) removes this denominator exactly. No lower bound on \(p_X\) is used.

Completion of the energy proof

An exact refresh at unit activities

The product chain

Set \[ c=1+N^{-2},\qquad T=\min\{t\in\mathbb Z_{\ge0}:c^t\ge\max(D,3\cdot4^K)\}, \qquad q=32L, \tag{42}\] where \(L=(10^4ND)^{100}\) is the integer in 14. At tier \(t\in\{0,\ldots,T\}\), clamp every reference activity to \([c^{-t},c^t]\) and denote the matching law by \(\pi_t\). Thus \(\pi_0\) is the uniform colored matching law and \(\pi_T\) is the desired law. Moreover, \[ \frac12\le\frac{\pi_{t-1}(A)}{\pi_t(A)}\le2 \qquad(A\in\mathcal M,\ 1\le t\le T). \tag{43}\] Indeed, an individual activity changes by a factor in \([c^{-1},c]\), so unnormalized matching weights and partition functions each change by factors in \([c^{-N/2},c^{N/2}]\). The probability ratio is therefore in \([c^{-N},c^N]\), and \(c^N\le\exp(1/N)<2\). Also \(\log(1+N^{-2})\ge(2N^2)^{-1}\), whence \[ T\le1+2N^2\log\max(D,3\cdot4^K) =O\bigl(N^2(K+\log D+1)\bigr). \tag{44}\] Here and below, unbased logarithms are natural.

Use \(q\) slots at each tier, with product stationary law \[\nu=\bigotimes_{t=0}^T\pi_t^{\otimes q}.\] An allowed pair \((x,y)\) consists of a slot \(x\) at tier \(t\ge1\) and any slot \(y\) at tier \(t-1\). On this pair use exactly the symmetric proposal of 14, with local proposals acting on \(x\), and accept by Metropolis for \(\pi_t\otimes\pi_{t-1}\). Across unequal tiers a cycle swap uses this Metropolis acceptance ratio as well; the acceptance-one observation for two identical laws in 6 no longer applies. For each base slot there is instead the update that refreshes it from \(\pi_0\), using 19. There are \[ M=Tq^2+q \tag{45}\] updates. The chain \(P\) holds with probability \(1/2\); otherwise it chooses one of these \(M\) updates uniformly. Every update preserves \(\nu\) and is reversible, so the same is true of \(P\).

Write \(\mathcal E_{xy}\) for the energy of one allowed pair update, with all other coordinates integrated against \(\nu\), and \(\mathcal E_i\) for the energy of a refresh of base slot \(i\). In this notation \[ \mathcal E_{\nu,P}(f) =\frac1{2M}\left(\sum_{xy}\mathcal E_{xy}(f) +\sum_{i\text{ base}}\mathcal E_i(f)\right). \tag{46}\]

For later use, the two-coordinate inequality also applies across adjacent tiers. The pair state space and symmetric proposal are identical for the targets \(\pi_t\otimes\pi_t\) and \(\pi_t\otimes\pi_{t-1}\). By [eq:adjacent-density], the latter target density is at least half the former at every pair state. The Metropolis capacity formula [eq:metropolis-capacity] thus shows that the unequal-tier energy dominates half the equal-tier energy. Consequently, for any functions \(a\) and \(b\) of the two slots, \[ \mathop{\mathrm{Var}}_{\pi_t}a \le 2L\,\mathcal E_{xy}\bigl(a(x)+b(y)\bigr). \tag{47}\] In this display the energy is on the two slots alone; it can also be applied conditionally on any fixed values of other slots.

A quantitative sampling guarantee

The schedule and vertex scales

If \(n=|V|=0\), return one. If a binary vertex count is used and \(n>2|E|\), return zero: the listed edges cannot cover all vertices. On the remaining nonempty input, perform a deterministic perfect-matching existence test and return zero if it fails. Thus the remaining case has even \(n\ge2\) and \(Z(G)\ge1\). Set \[ b=1+\frac1n,\qquad K=\min\left\{k\in\mathbb N:k\ge1,\quad n^{n/2}b^{-k}\le\frac\varepsilon{32}\right\}. \tag{56}\] Since \(\log(1+1/n)\ge1/(2n)\), \[ K\le 1+2n\left(\frac n2\log n+\log(32/\varepsilon)\right). \tag{57}\] At stage \(j\in\{0,\ldots,K\}\) the complete logical graph has raw activities \[w^{(j)}_{uv}= \begin{cases}1,&uv\in E,\\ b^{-j},&uv\notin E.\end{cases}\] Write \(Z_j\) for its perfect-matching partition function. Then \[ Z_0=(n-1)!!,\qquad Z_K=(n-1)!!\prod_{j=0}^{K-1}\frac{Z_{j+1}}{Z_j}. \tag{58}\]

For sampling, we additionally maintain positive vertex scales \(s_i\) and use activities \(\lambda_{ij}=w^{(j)}_{ij}s_i s_j\). Every perfect matching has the same vertex factor \(\prod_i s_i\). Consequently \[ Z_\lambda=Z_j\prod_i s_i, \qquad g_\lambda(ij)=\frac{g_{w^{(j)}}(ij)}{s_i s_j}. \tag{59}\] Initialize \(s_i=1/\left\lfloor\sqrt{n-1}\right\rfloor\). At stage zero, \[g_\lambda(ij)=\frac{\left\lfloor\sqrt{n-1}\right\rfloor^2}{n-1}\in[1/4,1],\] so the logical graph is balanced. The update defined below always multiplies each scale by a positive number at most two, whether or not the estimates succeed. Starting from scales at most one and making at most \(K\) updates, we therefore have at every sampling call \[ 0<\lambda_{ij}\le4^K. \tag{60}\] We use the same structural parameters throughout a trial: \[ D_0=10^8(n+1)^4,\qquad p=2K+2,\qquad N=n+4p\binom n2,\qquad D=100N^2D_0. \tag{61}\] For each current \(\lambda\), construct the enlarged graph and marked bag edges of 4.

Observables from one marked sample

Let \(\lambda^+\) denote the activities at the next raw stage, with the current vertex scales retained. If \(P\) is a logical perfect matching, or a logical matching missing one specified pair of vertices, put \[u(P)=b^{-\#\{e\in P:e\notin E\}}.\] Thus \(\operatorname{wt}_{\lambda^+}(P)=u(P)\operatorname{wt}_\lambda(P)\) and \[ \frac12<b^{-n/2}\le u(P)\le1. \tag{62}\] The strict left inequality follows from \((1+1/n)^{n/2}<\mathrm e^{1/2}<2\).

Sample a colored perfect matching from the bag graph and reveal the real, virtual, and probe portions of its edges. Define the following observables, all in \([0,1]\): \[\begin{align*} a&=\mathbf 1_{\{\text{all edges real}\}},\\ d&=\mathbf 1_{\{\text{all edges real}\}}u(P),\\ z_{ij}&=\mathbf 1_{\{\text{probe }ij\text{ and all other edges real}\}}u(P). \end{align*}\] Here \(P\) is decoded only on the indicated event; elsewhere the observable is zero. On the all-real event it is a full logical perfect matching. On the sole-probe event it is a logical perfect matching missing \(i,j\). Write \(I=Z_{\rm bag}/Z'\) for the inflation factor. Under the exact marked bag law, the means satisfy \[ \mathbb Ea=\frac1I,\qquad \mathbb Ed=\frac{Z_{\lambda^+}}{I Z_\lambda},\qquad \mathbb Ez_{ij}=\frac{Z_{\lambda^+}(-ij)}{D I Z_\lambda}. \tag{63}\] The all-real probability in 7 gives the first equality. The other two follow by multiplying the weights in that event, or in the sole-probe identity (18), by \(u(P)\) and summing. In particular, \[ R:=\frac{Z_{j+1}}{Z_j} =\frac{Z_{\lambda^+}}{Z_\lambda} =\frac{\mathbb Ed}{\mathbb Ea},\qquad g_{\lambda^+}(ij)=D\frac{\mathbb Ez_{ij}}{\mathbb Ed}. \tag{64}\] Thus one batch of marked samples will estimate both the next raw partition ratio and the two-hole ratios needed to choose the next vertex scales.

Estimating and restoring balance

Use \[ \alpha=10^{-5},\qquad \sigma=\min\left\{\frac\alpha{100D},\frac\varepsilon{1000K}\right\},\qquad S=\left\lceil 10\sigma^{-2}(n+K+1)\right\rceil,\qquad \xi=\frac1{1000KS}. \tag{65}\] At each stage make \(S\) independent fresh sampler runs, each with marked output law at total variation distance at most \(3\xi\) from the exact marked bag law. 9 implements these runs in bounded bit time. Let bars denote the empirical means of the observables. Set \[ \widehat R= \begin{cases}\bar d/\bar a,&\bar a>0,\\1,&\bar a=0. \end{cases} \tag{66}\] If \(\bar d>0\), form the symmetric array \[ G_{ij}=\max\left\{\alpha, \min\left\{3,D\frac{\bar z_{ij}}{\bar d}+\alpha\right\} \right\} \quad (i\ne j). \tag{67}\] If \(\bar d=0\), set all \(G_{ij}=1\). The following elementary balancing procedure will be useful.

After estimating the current ratio, perform this procedure and replace \(s_i\) by \(s_i v_i\). These scales are used with the next raw stage.

Success probability and output

The form of Hoeffding’s inequality [9] used here is: if \(X_1,\ldots,X_S\) are independent, take values in \([0,1]\), and have a common mean \(\mu\), then, for every \(t>0\), \[ \mathbb P\left[\left|S^{-1}\sum_{i=1}^S X_i-\mu\right|>t\right] \le2\exp(-2St^2). \tag{71}\] The theorem applies separately to each observable; the observables within a single sample need not be independent.

Condition on any history in which all preceding stages succeeded. By 23, the current activities are balanced. For independent ideal samples, a union bound over the at most \(n^2+2\) observables gives failure probability at most \[2(n^2+2)\exp(-2S\sigma^2) \le2(n^2+2)\exp(-20(n+K+1))<\frac1{100K}.\] One elementary verification of the last strict inequality is \(\log(200K(n^2+2))\le6+K+2n<20(n+K+1)\) for \(n\ge2,K\ge1\). The actual fresh runs can be coupled to these ideal samples with disagreement probability at most \(3S\xi=3/(1000K)\). Therefore the conditional probability of failure at the first unsuccessful stage is at most \(13/(1000K)\). Summing over all \(K\) stages proves that a whole trial succeeds with probability at least \(1-13/1000>0.9\).

The trial output is the nonnegative rational \[ Y=(n-1)!!\prod_{j=0}^{K-1}\widehat R_j. \tag{72}\] On success, 23 and \(\prod_i(1-t_i)\ge1-\sum_i t_i\) for \(t_i\in[0,1]\) imply \[(1-\varepsilon/100)Z_K\le Y\le\exp(\varepsilon/100)Z_K.\] Every matching using a nonedge has raw final weight at most \(b^{-K}\). There are at most \((n-1)!!\le n^{n/2}\) perfect matchings on the complete graph, so \[ Z(G)\le Z_K\le Z(G)+n^{n/2}b^{-K} \le(1+\varepsilon/32)Z(G). \tag{73}\] For \(0<\varepsilon<1\), the resulting lower bound exceeds \((1-\varepsilon)Z(G)\), and the upper bound is at most \((1+\varepsilon)Z(G)\). For example, \(\exp(\varepsilon/100)\le1+\varepsilon/50\), which gives \[\exp(\varepsilon/100)(1+\varepsilon/32) \le1+\left(\frac1{50}+\frac1{32}+\frac1{1600}\right)\varepsilon <1+\varepsilon.\]

To amplify, let \(\ell=\left\lceil\log_2(1/\delta)\right\rceil\), computed by rational comparisons with powers of two. Run \(r=10(\ell+1)+1\) independent trials and return their median. This is an odd number of trials. Since each succeeds with probability at least \(0.9\), Hoeffding’s inequality for their success indicators bounds the probability of fewer than half succeeding by \[\exp\bigl(-2(0.4)^2r\bigr) \le\exp\bigl(-3.2(\ell+1)\bigr) <2^{-\ell}\le\delta.\] Whenever a majority succeed, the median belongs to the desired interval. The deterministic zero branch was taken before any trials. This proves the correctness and probability assertions of 1, subject only to the bit implementation established next.

Sizes of the parameters and state spaces

Let \(\mathsf b\) be the total binary input length. We use a standard explicit encoding of a finite graph and of the two rationals. If isolated vertices are specified by a binary vertex count, first return zero when \(n>2|E|\). This test is valid for every nonempty graph admitting a perfect matching and ensures that the nonzero branch has \(n\le2|E|=O(\mathsf b)\). The deterministic existence test [5] therefore has polynomial bit cost on the original graph. It is the only matching-existence oracle used.

[eq:schedule-length] bounds \(K\) polynomially in \(\mathsf b\). It can be found without logarithms by iterating rational powers of \(b\) and comparing the two sides of (56). The parameters \(p,N,D_0,D\) are integers of polynomial numerical size in \(n,K\). The hierarchy has at most \[ H_{\rm tree}=N(2K+1) \tag{74}\] nodes, since its levels are \(0,\ldots,2K\) and every level has at most \(N\) components. The number of colored edges is less than \(N^2\).

Use the ladder parameters \(c,T,L,q,M,\Lambda\) from [eq:energy-constants,eq:ladder-parameters,eq:ladder-update-count,eq:ladder-lambda], and write \[t_{\rm mix}=s(\xi),\] where \(s(\xi)\) is defined in (54). The structural parameters \(N,D,K\) are fixed throughout a trial. Hence the same \(T\) makes the final tier unclamped on every empirical history, by [eq:unconditional-height-bound,eq:ladder-global-range]. It can be found by iterating rational powers of \(c\) until the defining inequality in (42) holds. The bound (44) is polynomial in \(N,K,D\). Ceilings of binary logarithms are found by integer or rational comparisons with powers of two. These bounds and the product-state size \(q(T+1)\) are polynomial in \(\mathsf b\) and \(\varepsilon^{-1}\). Indeed, \[\sigma^{-1}\le C(D+K\varepsilon^{-1}),\qquad S=O\bigl((D^2+K^2\varepsilon^{-2})(n+K+1)\bigr),\qquad \xi^{-1}=1000KS,\] for an absolute constant \(C\). The exponent \(100\) in \(L\) is fixed.

There is no initial-sampling problem at the upper tiers. Pair the logical vertices in their fixed order and expand this matching by the path tilings. All logical activities are positive, so this is an available real perfect matching with assigned colors. Initialize every slot to that matching. The mixing bound applies from this deterministic product state.

Rational bit lengths

We record explicitly why the adaptive scales do not create an arithmetic blowup. Put \[ B_G=1+n\left\lceil\log_2(n+1)\right\rceil+\left\lceil\log_2(S+1)\right\rceil +\left\lceil\log_2(D+1)\right\rceil. \tag{75}\] Absolute multiplicative constants in the bit bounds below are independent of the input.

Every nonzero sample contribution is \((n/(n+1))^a\) for some \(0\le a\le n/2\). These fractions have the common denominator \(Q=(n+1)^{n/2}\). Thus empirical means, their ratios, and the capped entries \(G_{ij}\) have \(O(B_G)\) numerator and denominator bits. This bound does not involve the current scales.

An assignment in (68) chooses one maximizer \(j\) and sets \(v_i=G_{ij}/v_j\). At that moment \(v_j\) is either the initial value two or an expression formed at an earlier update. Following such dependencies strictly decreases the update index, so it produces a chain of at most \(n\) factors. Each \(v_i\) consequently has \(O(nB_G)\) bits. Accumulating at most \(K\) scale updates gives \(O(KnB_G)\) bits per scale. The raw activities have \(O(K\log(n+1))\) bits. These estimates apply even when all observations are unfavorable, because the clipping and branch conventions are unconditional.

Bottleneck entries \(B_{ij}\) and the maxima \(H_i=\max_{j\ne i}B_{ij}\) select from finitely many input activities and one; they do not involve sums over paths. Initialize \(B_{ij}=\max\{1,\lambda_{ij}\}\) and, for each vertex \(k\), update the entries with distinct \(i,j\ne k\) by \[B_{ij}\leftarrow\max\{B_{ij},\min\{B_{ik},B_{kj}\}\}.\] This is the max–min Floyd–Warshall recurrence: after processing \(k\), the entries allow the processed vertices as internal path vertices. It uses \(O(n^3)\) rational comparisons. The path profiles add only powers of two of bit length \(O(K)\). Bag activities are sums of a real activity and a dyadic height divided by \(D\). Finally \(c^{\pm t}\) has \(O(T\log(N+1))\) bits for \(t\le T\).

It follows that every single-edge activity used by the sampler has numerator and denominator lengths at most the integer bound \[ B_W=C_1\bigl(KnB_G+K+\left\lceil\log_2(D+1)\right\rceil +T\left\lceil\log_2(N+1)\right\rceil+1\bigr), \tag{76}\] where \(C_1\) is a sufficiently large absolute positive integer. A product of at most \(N\) such activities has \(O(NB_W)\) bits. In a Metropolis acceptance ratio the unknown normalizing constants cancel. Thus even a cycle exchange at unequal tiers requires only rational products of this size. The real and probe revelation probabilities are likewise rational, respectively the real portion divided by the collapsed activity and \(1/(DW)\) for the probe portion.

By 19, the bottom-tier dynamic program uses counts of \(O(N\log(N+1))\) bits. Its categorical normalizers count the same partial completions and obey the same bound. Label and loop indices use an additional \(O(\log(H_{\rm tree}+1)+\log(N+1))\) bits, since \(b_{\rm tree}\le H_{\rm tree}\) in the constructed hierarchy. The table computation and sampling therefore use polynomially many integer operations with polynomial bit lengths.

Random choices with a fixed number of bits

The exact categorical choices in the unit-tier dynamic program define an ideal transition law. The actual algorithm implements each such choice, along with every proposal and acceptance choice, by the fixed grid in 24. Hence exactness of the dynamic program is a counting identity; sampling from its tables in bounded bit time incurs the error budget derived below.

Every choice in the chain is either uniform on an explicitly indexed finite set, Bernoulli, or proportional to a list of nonnegative integers supplied by the dynamic program. For each used edge, final revelation is a rational categorical choice with at most three outcomes. Subsets and pairings are sampled sequentially, so no exponential list is needed. Uniform choices among the \(M\) scan positions can be implemented by integer arithmetic on the index rather than by enumerating the scan positions.

A bound on the number of random choices in one complete sampler run, including marking its output, is \[ R_{\rm draw}=100(t_{\rm mix}+1)(N+1)^2(K+2). \tag{77}\] To see this, laziness and the scan take two choices per step, and a pair update takes at most eight further choices. By 19, a refresh takes at most \(4(H_{\rm tree}+N+1)\) choices. Allowing both bounds at every step and adding the at most \(N\) choices for final marking gives \[t_{\rm mix}\bigl(10+4(H_{\rm tree}+N+1)\bigr)+N.\] Since \(H_{\rm tree}=N(2K+1)\), this is bounded by (77). The bound therefore includes every categorical choice made by the implemented sampler. The number of options in any one choice is at most \[ R_{\rm opt}=10(M+N^2+1). \tag{78}\] Choose the same bit budget \(k\) for all these choices, with \[ 2^k\ge\frac{2R_{\rm draw}R_{\rm opt}}\xi. \tag{79}\] This integer is found by doubling; \(k\) is polynomially bounded. By 24, the sum of all conditional simulation errors is at most \(\xi\). Successive coupling, while the simulated and exact histories agree, bounds the total variation error of the whole run and its marking by \(\xi\). The exact-kernel mixing error is at most \(\xi\), so the implemented marked sample is within \(2\xi\), and in particular \(3\xi\), of the required law. The number of fair bits used has a deterministic bound \(kR_{\rm draw}\). There is no rejection loop or random convergence test.

Total bit cost

The preceding bounds already give a polynomial-time implementation. To make the accounting explicit, take a common bit bound \[B_* = C\bigl(NB_W+N\log(N+1)+\mathsf b +\log(M+1)+\log(t_{\rm mix}+1) +\log(\xi^{-1}+1)+1\bigr)\] with a sufficiently large absolute constant \(C\) and \(B_W\) from (76). In particular \(NB_W\) dominates \(KB_G\), the bit cost of multiplying the \(K\) empirical ratios. The bound \(B_*\) covers all integers used in weight products, dynamic-program counts, cumulative choices, acceptance comparisons, and output products. All categorical laws just described admit common-denominator integer weights, so computing their cumulative sums respects the same bound. Schoolbook arithmetic and the Euclidean algorithm perform the needed operations within \(O(B_*^3)\) bit time per operation.

A deliberately loose bound of \(O((H_{\rm tree}+1)(N+1)^4)\) arithmetic operations covers graph construction, the bottom-tier tables, and any one transition or refresh; storing the product state additionally takes \(O(q(T+1)N B_*)\) bits. Thus a trial, including all \(K\) stages and \(S\) fresh runs at each stage, has bit cost bounded by a fixed polynomial in \[K,\ S,\ t_{\rm mix},\ H_{\rm tree},\ N,\ q(T+1),\ B_*.\] For example, the sum of \[CKS(t_{\rm mix}+1)(H_{\rm tree}+1)(N+1)^4B_*^3 \quad\text{and}\quad CKSq(T+1)(N+1)B_*^3\] is a valid coarse envelope for arithmetic, repeated initialization, and storage operations. Each displayed parameter is polynomial in \(\mathsf b\) and \(\varepsilon^{-1}\). The \(r=10(\left\lceil\log_2(1/\delta)\right\rceil+1)+1\) repetitions and the rational comparisons used to select their median add a factor polynomial in \(\log\delta^{-1}\) and the input length.

All data structures, loop bounds, and rational comparisons have now been specified by finite procedures. Every auxiliary parameter is bounded by a polynomial of fixed degree, and the bound holds on every random history. Together with 8, this completes the proof of 1.

Constant-fibre reductions

Write \(n=|V|\), \(m=|E|\), and \(d_v=\deg_G(v)\).

A deletion sampler for perfect matchings

The conversion from counts to samples uses the self-reduction framework of Jerrum, Valiant, and Vazirani [12]. We give the deletion procedure explicitly to retain feasibility on every execution and to track its total-variation and bit-time bounds. It uses 1 only through its counting guarantee.