Dimers on Riemann surfaces I: Temperleyan forests

This is the first article in a series of two papers in which we study the Temperleyan dimer model on an arbitrary bounded Riemann surface of finite topolgical type. The end goal of both papers is to prove the convergence of height fluctuations to a universal and conformally invariant scaling limit. In this part we show that the dimer model on the Temperleyan superposition of a graph embedded on the surface and its dual is well posed, provided that we remove an appropriate number of punctures. We further show that the resulting dimer configuration is in bijection with an object which we call Temperleyan forest, whose law is characterised in terms of a certain topological condition. Finally we discuss the relation between height differences and Temperleyan forest, and give a criterion guaranteeing the convergence of the height fluctuations in terms of the Temperleyan forest.


Background and main result
This paper is the first part of a series of two papers in which we show universality and conformal invariance of the Temperleyan dimer model on Riemann surfaces.
Given a finite graph G (with an even number of vertices) and nonnegative weights w = (w e ) e∈E on the edges, a dimer configuration is a perfect matching of the vertices along the edges of G. The dimer model is the random perfect matching obtained by sampling a dimer configuration m on G with probability proportional to the product of the weights of the edges in the matching m: where Z is the partition function. If all the weights are equal to 1, P is simply a uniformly picked perfect matching. Following Thurston, a particularly convenient way to encode the dimer model on a planar bipartite graph is through a notion of height function. Informally, it is a real valued function defined on the faces of the graph, which is uniquely determined by the dimer configuration. This turns the dimer model into a model of discrete random surfaces. A central question in the study of the dimer model is to describe the fluctuations around the mean of the height function. It turns out that the height function of the dimer model is extremely sensitive to boundary conditions (see [26]). In the so-called liquid phase, the Kenyon-Okounkov conjecture predicts that the fluctuations of the height function around its mean are asymptotically given by a Gaussian free field in a certain nontrivial conformal gauge. Despite remarkable progress on this conjecture over the last decade, this remains an outstanding problem in all but a handful of cases: see [4,18,25,30,32] for examples of models where this is already shown, and see [7,34] for an introduction to the Gaussian free field. We emphasise that, given the extreme sensitivity of the model to its boundary conditions, the occurrence of a universal limit such as the Gaussian free field is to some extent very surprising.
The goal of this paper is to describe the fluctuations of the height function for certain graphs embedded on Riemann surfaces of finite topological type, i.e., with finitely many handles and holes. Our assumptions on these graphs ensure that the dimer model is liquid on the entire Riemann surface and the conformal gauge is simply the identity; we call the resulting dimer model Temperleyan, and we say that the sequence of graphs G #δ is a Temperleyan discretisation of the surface. Even though this assumption simplifies the problem, conformal invariance in the scaling limit is by no means obvious; for instance, in the simply connected setting and for the square lattice, this is in fact precisely the content of Kenyon's landmark paper [24].
Before we give a simplified statement of the main theorem, we point out that in the context of a Riemannian surface (M, g), the height function becomes a height one-form h #δ : that is, the gradient ∇h #δ is well-defined but this gradient does not derive itself from a single-valued function on the surface. Instead it can be viewed as the gradient of a (single-valued) function defined on the universal cover of the surface. For concreteness we will denote by h #δ the associated function on the universal cover of the surface. Readers who prefer it can however equivalently think of h #δ as the random one-form on M associated to the dimer configuration. (See Sections 2.4 and 6.1 for precise definitions.) The function h #δ is the central object of this sequence of papers, and we still refer to it as the height function.
Our main theorem, stated informally below, confirms that the dimer model admits a universal scaling limit on Riemann surfaces. This scaling limit is furthermore conformally invariant. Theorem 1.1. Let (M, g) denote a bounded Riemann surface, possibly with a boundary, which is not the sphere or a simply-connected domain. Let G #δ be a sequence of Temperleyan discretisations of M satisfying the invariance principle and a crossing estimate for random walk described in Section 2.3. Let h #δ denote the height function of the dimer model on G #δ . For all smooth compactly supported test functions f , converges as δ → 0 in law and the limit is conformally invariant and independent of the sequence G #δ .
The precise assumptions on the surface M will be given in Section 2.1, the meaning of Temperleyan discretisations and our precise assumptions on the sequence G #δ can be found in Section 2.3. Finally, recall that the smooth test function f is taken on the universal cover M of M , so that the integration takes place on M (equipped e.g. with Lebesgue measure). See Sections 2.4 and 6.1 for details. We refer to Theorem 6.1 for the complete version of the theorem.
Brief outline. Our general approach to study height fluctuations on Riemann surfaces is motivated and inspired by our earlier work [4] connecting the dimer model to imaginary geometry in the simply connected setting. Let us recall the main idea of that work briefly in the original simply connected setting, before explaining how this relates with and differs from the case of Riemann surfaces.
1. The starting point of [4] is Temperley's bijection, which relates the dimer model on a Temperleyan discretisation of a simply connected domain to a pair of dual uniform spanning trees with respectively wired and free boundary conditions. Furthermore, the height function of the dimer model can be determined relatively simply from the associated pair of spanning trees: roughly speaking, the height difference between two points is simply the total winding (sum of turning angles) of the unique path connecting these points in either tree.
2. The second observation is that, as follows from the landmark paper of Lawler, Schramm and Werner [29], paths between points in the wired uniform spanning tree (UST) have a scaling limit which may be described in terms of Schramm-Loewner Evolutions with parameter κ = 2 or SLE 2 for short. 3. Finally, in the continuum, there is a coupling known as imaginary geometry between SLE κ curves and an appropriate multiple of the Gaussian free field. In this coupling, the values of the GFF may informally be identified as the winding of the associated SLE curves (note that this requires careful interpretation, since the Gaussian free field is not defined pointwise and SLE κ paths are not smooth). Nevertheless, this coupling can be thought of as a continuous form of Temperley's bijection in the case κ = 2.
The strategy of [4] is to exploit Temperley's bijection on the one hand, and Imaginary Geometry on the other hand, to show that the winding of paths in a UST is indeed given by the appropriate multiple of the Gaussian free field asymptotically (essentially, one has to show that we can exchange the order of taking limits, or that a diagram commutes). This strategy is very robust; in particular it does not appeal to the more solvable aspects of the dimer model. Indeed, the convergence of loop-erased random walks (which describe the branches of a wired UST at the discrete level) to SLE 2 is known in great generality.
Roughly speaking, we wish to implement a similar approach in the case of Riemann surfaces. However, it should be clear that each part will be substantially different in this new setup.
• To begin with, we observe that Temperley's bijection is local and so may be applied on a Riemann surface to output a pair of random subgraphs that are dual to one another and locally tree-like, and which we call Temperleyan self-dual pair. While in the simply connected case these are of course just dual spanning trees, the nature (i.e. state space) of the Temperleyan self-dual pair on a Riemann surface is a priori not entirely clear, and its law even less so. Our first result will be a characterisation of the Temperleyan self-dual pair. We also describe its law explicitly, which has a simple and concise expression. It turns out that one of the two components, which we call the Temperleyan cycle-rooted spanning forest or Temperleyan CRSF for short, determines the other (and hence the entire dimer configuration) uniquely, and also has a tractable law. Essentially, a Temperleyan CRSF is a random spanning subgraph in which every cycle is noncontractible, and satisfying a certain further topological condition. The introduction and study of Temperleyan CRSFs occupies us throughout both articles in this series.
• A second problem is that the connection between height function and winding is not as straightforward as in the simply connected case; in particular, in addition to the already mentioned fact that only the gradient of the height function is well defined (leading to a multivalued height function), there may not be a path connecting two given points in the Temperleyan CRSF, even when they are neighbours. We thus need to develop a systematic method for computing the height function given the Temperleyan CRSF.
• The next difficulty is to obtain a scaling limit result for our Temperleyan CRSF. Let us immediately note here that there are specific difficulties with defining variants of SLE on surfaces since the simply-connected uniformisation theorem lies at the heart of the standard definition of SLE. In addition, if we were to follow the blueprint of [4] one would also need to establish a version of Imaginary geometry for Riemann surfaces, an endeavour which seems out of reach with current techniques. Fortunately, we are in fact able to entirely bypass this step because it turns out that the (discrete) estimates needed anyway to exchange limits are strong enough to imply convergence regardless of any information in the continuum. Clearly, the downside of this approach that we cannot hope to describe the law of the limit in this manner.
The topological and geometric structures of the surface M raise substantial specific difficulties at each of the steps of this roadmap. In this first paper we develop the relevant discrete/combinatorial arguments corresponding to the first two bullet points above. The second companion paper in this sequence, [5], is devoted to the proof of convergence of the Temperleyan CRSF in the scaling limit. This is the most involved step at the technical level, but can be read essentially independently.

Relation with previous work
One of the first rigourous works on the dimer model on Riemann surfaces is the inspiring paper by Boutillier and de Tilière [9], who were motivated by nonrigourous physics predictions based on the Coulomb gas formalism (see e.g. [16] and the references in [9]). They described the topological part of the fluctuations in the dimer model on the honeycomb lattice on the torus. This was complemented some years later by the fundamental paper of Dubédat [18], in which the full law of the limiting height fluctuations was described for isoradial graphs on the flat torus, and identified as a so-called compactified Gaussian free field. Dubédat also conjectured a greater form of universality beyond the isoradial setting; such a universality statement was subsequently proved by Dubédat and Gheissari [19] but only for the topological part of fluctuations on the torus. Note that in particular our Theorem 1.1 settles this conjecture of Dubédat.
We point out that all techniques in the paper mentioned above, based on Kasteleyn theory (in contrast with our approach), are very specific to the case of the torus and cannot easily be adapted to general Riemann surfaces. The paper of Cimasoni [10] shows that for general surfaces of genus g (and no boundary), the partition function of the dimer model becomes a weighted (in fact, signed) sum of 2 2g determinants of Kasteleyn matrices in which the orientation has been reversed along some of the 2g cycles forming a basis of the first homology group of the surface. The signs in front of the determinants themselves form a nontrivial algebraic topological invariant of the surface which Cimasoni, remarkably, was able to identify in terms of the socalled Arf invariant. (In the case of the torus and the honeycomb lattice, this goes back to the work of Boutillier and de Tilière [9].) This illustrates the difficulty in generalising the Kasteleyn approach to general surfaces. See also [1,11,12,14] for related results.
The above difficulty illustrates the well known fact that the dimer models becomes "less integrable" as the genus of the surface increases, or at any rate integrability becomes harder to exploit. Our approach is not dependent on integrability or solvability and this accounts in large part for the fact that we are able to tackle generic Riemann surfaces (in contrast with methods based on Kasteleyn theory or interacting particle systems). This is not to say that the difficulties caused by the increased complexity of the surface vanish entirely, but they manifest themselves in a different manner, which boil down to the following issue. As we will describe in more details in Section 2, a Temperleyan discretisation of the surface is the graph resulting from the superposition of a graph Γ embedded on the surface and its dual. The resulting graph is (locally) bipartite, with its black vertices being given by the vertices and faces of Γ, and the white vertices being given by the edges of Γ. However, a straightforward calculation using Euler's formula shows that this does not admit a dimer covering unless we remove 2g+b−2 = |χ| white vertices from the graph (where g is the genus and b the number of boundary components, and χ is the Euler characteristic); these can be thought of as punctures in the surface, or monomers. Note that the number of punctures increases with the complexity of the surface, but when the surface is a torus or an annulus no such puncture needs to be removed. Roughly speaking, each puncture or monomer is responsible for an extra topological condition which the Temperleyan forest is required to satisfy; in particular in the case of a torus or an annulus, the Temperleyan forest is not conditioned on any degenerate event and boils down to the simpler notion of cycle-rooted spanning forest studied by Kassel and Kenyon in [23].
The necessity to remove a certain number of punctures from the surface in order to formulate the Temperleyan dimer model has a geometric flavour which is reminiscent of Liouville conformal field theory, in which correlation functions are only well defined with the appropriate number of singularities (known in that context as insertions). The required number of insertions is specified by the so-called Seiberg bounds, see e.g. [15] on the sphere and [20] for Riemann surfaces (see also [7] for an introduction). As is the case here, the number of punctures is entirely specified by the topology of the surface, but in a way that is in some sense opposite to ours: indeed the required number of insertions is an increasing function of χ in Liouville CFT, whereas it is decreasing here.
This difference can be understood informally as follows. To first order, Liouville CFT may be viewed as yielding a random metric which is approximately hyperbolic (this can for instance be made rigourous in the semiclassical limit where the parameter γ → 0, see [28]). However, without adding conical singularities at a certain number of points, it may not be possible to equip the surface with a metric of constant negative Ricci curvature, as follows from the Gauss-Bonnet theorem. By contrast, the Temperleyan dimer model is trying to equip the surface with a flat metric, a fact which will be particularly apparent when we relate the height function to the winding of the Temperleyan forest in Section 4 (note that the winding is computed on the universal cover, where curvature effects are ignored). But on a hyperbolic surface (i.e., when χ < 0), no flat metric exists, again by the Gauss-Bonnet theorem. In summary, Liouville CFT may be viewed as trying to build a hyperbolic metric on a surface which may be parabolic or elliptic, whereas the Temperleyan dimer model tries to impose a flat metric on a surface which may be hyperbolic. In both cases the difficulty is resolved by puncturing the surface at an appropriate number of singularities. (The above parallel between Liouville CFT and dimer model is also to be compared with the two couplings between SLE and Gaussian free fields given respectively by the so-called forward and reverse couplings; the forward coupling describes Liouville conformal field theory whereas the reverse coupling describes imaginary geometry).

Overview of the proof and organisation of the paper
• In Section 2, we recall some basic facts about dimers on surfaces. Sections 2.1 and 2. 2 contain small reminders about Riemann surfaces in two dimensions. We introduce precisely our discrete setup in Section 2.3. As one of the main difference with the classical case is that the height function now becomes multivalued we recall in Section 2.4 the language of height one-forms on graphs which is a classical way to handle such multivalued functions. In Section 2.5, we introduce another way to deal with multivalued functions by lifting them to single-valued functions on the universal cover. This is in practice easier to deal with, and we mainly work with this lift in this article. We also illustrate how an application of the Hodge decomposition theorem allows us to decompose the (multivalued) height function on the surface into a single-valued function (scalar component) and a canonical representative of the "multivalued (or topological) part", the so-called instanton component. In Section 2.6, we recall the notions of windings of curves which were developed in our previous article [4]. • In Section 3, we introduce the generalisation of Temperley's bijection to Riemann surfaces. Section 3.1 contains the bijection itself which, as mentioned above, uses a new type of object which we call a Temperleyan CRSF or Temperleyan forest for short. We also provide in Section 3.2 a simple criterion for a CRSF to be Temerleyan which is the basis of [5]. • In Section 4, we carefully develop the relation between the height differences and the winding of Temperleyan forests; in spite of the curvature of the surface such a connection remains true if one considers a natural embedding of the graph on the universal cover of the surface. This is in the spirit of the work of Kenyon-Propp-Wilson [27], who treated the (planar) simply connected case. In that case of course the bijection is with a uniform spanning tree. Compared to that work there are two additional points that we need to handle. The first one is that the edges of the graph cannot be assumed to be straight lines as in [27]. The second and more significant one is that there are additional terms coming from the fact that the forest is not connected: given points x and y on the universal cover, the height difference h(x) − h(y) (which is unambiguously defined) is essentially given by the intrinsic winding of any path connecting x to y, plus additional discontinuities every time the path jumps between components. The resulting key formula is stated in Theorem 4.10 (see also Lemma 4.6). • In Section 5, we extend our local and global coupling results from [4] to the framework of Riemann surfaces. This coupling is a key ingredient of our approach, which allows us to show that given a finite number of points (z i ) 1≤i≤k on the surface (or, rather, their lifts to the universal cover), the respective geometries of the CRSF in a neighbourhood of these k points decouple and are independent. Clearly, such an independence can only be expected to hold at distances smaller than the distances between the z i 's, but if the points are macroscopically far apart this still leaves a lot of room. The strength of this coupling is that it holds at an essentially optimal scale. More precisely, the independence property holds in neighbourhoods whose size is random and comparable to the distance between the z i 's, with an exponential control on the number of additional finer scales required for the independence to hold. The argument in [4] was based on relatively soft estimates about LERW, and it is not difficult to adapt them to this setting. There is one additional technicality however, which can be summarised as follows. A basic estimate needed for the argument is a polynomial control on the probability that a LERW starting from u comes close to a given point v (Proposition 4.11 in [4]). In [4] we could exploit the fact that if the random walk starting from u comes within e −n of v then it has a positive probability to make a loop around v at all scales e −n+1 , . . . , 1 subsequently, which erases the portion of the walk coming close to v in the loop-erasure. This results in an exponentially small probability that the looperasure comes with distance e −n of v. However, in the setup of Riemann surfaces, the walks also stop when it makes a non-contractible loop. Therefore if a random walk comes within distance e −n of v, there is a possibility that the walk completes a noncontractible loop shortly after, say at distance e −n+m from v, with m n. In this case there are in principle no remaining constraints on larger scales, and we do not get a sufficiently good bound on the probability that the path comes close to v. The argument needed to overcome this technical difficulty is given in Lemma 5.3.
• In Section 6, we finish by proving our main convergence result Theorem 6.1. Sections 6.2 and 6.3 contain some a priori technical estimates on winding of the CRSF branches on the universal cover. Finally in Section 6.4, we finish the proof of Theorem 6.1. This proof shares similarities with the argument in [4], but with the difference that we cannot rely on imaginary geometry to identify the limit. • Finally, appendix A contains some geometric facts about spines in the CRSF (that is, the lift to the universal cover of a cycle of the CRSF). It is shown that such spines, in the hyperbolic case, form simple chords in the unit disc (simple paths joining two boundary points) or a loop touching the boundary at a single point. This basic fact (stated as Lemma 4.9 but proved in the appendix) underlies much of the discussion connecting winding to height differences. The proof relies on some arguments from the classical theory of Riemann surfaces. • M is of finite topological type, meaning that the fundamental group π 1 (M ) is finitely generated. In other words, we assume that the surface has finitely many "handles" and "holes".
• M can be compactified by specifying a boundary ∂M . We denote by M the compactified Riemann surface with the boundary. More precisely, every point is either in the interior and hence has a local chart homeomorphic to C or is on the boundary and has a local chart homeomorphic to the closed upper half planeH. Also there are finitely many such charts which cover the boundary. Note that this condition implies that M has no punctures.
(However for future reference, we note here that we will later introduce punctures on M which will align with the removed vertices of Temperleyan graphs. The resulting punctured manifold will be denoted by M .) • The metric in M extends continuously to ∂M . In other words, recall that the metric inside M can be represented locally in isothermal coordinates as e ρ |dz| 2 for a smooth function ρ and we assume that ρ extends continuously to ∂M .
We say M is nice if M satisfies the above properties. Note that continuous extension of the metric to the boundary means we exclude surfaces such as the hyperbolic plane, and will technically simplify certain topological issues later when we deal with the Schramm topology.
Classification of surfaces Riemann surfaces can be classified into the following classes depending on their conformal type (see e.g. [17] for an account of the classical theory): • Elliptic: this class consists only of the Riemann sphere, i.e., M ≡Ĉ • Parabolic: this class includes the torus, i.e., M ≡ T = C/(Z + τ Z) where (τ > 0)), the cylinder M ≡ C \ {0}, and the complex plane itself (or the Riemann sphere minus a point), M ≡ C.
• Hyperbolic: this class contains everything else. This includes examples such as the twotorus, the annulus, as well as proper simply connected domains in the complex plane, etc.
The proofs in this paper (and its companion [5]) are concerned with the hyperbolic case (subject to the above conditions) as well as the case of the torus so from now on we always assume that M is such a manifold. We note that the case of simply connected domain in C is covered in our previous work [4].

Universal cover
The universal coverM of the Riemann surface M will play an important role in our analysis.
In the cases that we will analyse, due to the above classification, the universal cover can be classified as follows: •M ≡ C in the parabolic case of the torus; •M ≡ D (the unit disc) in all remaining non-simply connected hyperbolic cases.
From now on we assume, unless otherwise explicitly stated, that we are only dealing with these cases. We also recall here the classical Riemann uniformisation theorem: M is conformally equivalent toM /F where F is a discrete subgroup of the Möbius group. In case of the torus, this discrete subgroup is isomorphic to Z 2 and the generators specify translations in the two directions of the torus. In the hyperbolic case, this class of subgroups is much bigger and are known as Fuchsian groups (see e.g. [17] and appendix A). This particular representation of a surface, sometimes called a Fuschian model, will be particularly convenient to describe the scaling limits of Temperleryan or cycle-rooted spanning forests in the universal cover, rather than on the surface itself. This is essential for our purposes since the lack of curvature in the cover allows us to compute winding in a straightforward manner.
Let p :M → M be the covering map associated with the Riemann uniformisation theorem. That is, in the hyperbolic case, we write M = D/F and p : D → M is the canonical projection which is then conformal. In the case of the torus, p : C → M is the standard projection from the plane onto the torus and is also conformal. We recall here that the covering map acts discretely, in the sense that for every z ∈ M , there exists a neighbourhood N containing z so that p is injective in every component of p −1 (N ).
We will require the graphs we are working with to be embedded on M in a way that preserve its topology and geometry. Essentially the topological conditions will be that edges cannot cross, all interior faces must be homeomorphic to discs and the boundary of the manifold must correspond to simple loops on the graph, as outlined in the next section.

Discretization setup
Let Γ be a graph which is properly embedded on a Riemann surface M with b holes and g handles (see Section 2.1 for precise assumptions on the M ). Let Γ † be its planar dual. A face of Γ (resp. Γ † ) is a connected component in M of the complement of the embedded graph Γ (resp. Γ † ). We assume that Γ † is embedded on M so that every vertex of Γ † is in a face of Γ and every edge e of Γ is crossed by a single edge e † of Γ † which joins the two faces incident to e. (Later, when we consider a random walk on Γ we will allow the jump probabilities to be non-uniform and non-reversible, but for now Γ and Γ † are both unweighted, undirected graphs).
We assume that Γ and Γ † are faithful to the topology of the surface in the following sense. We require all faces of Γ and Γ † whose boundary does not intersect ∂M to be homeomorphic to a disc. We also require all faces of Γ † whose boundary intersects ∂M (we call those faces outer faces) to be homeomorphic to an annulus. In other words, we have a boundary vertex of Γ corresponding to each hole (let ∂Γ denote the set of such vertices), and each hole is surrounded by a cycle in Γ † which we call boundary cycle and which we denote by ∂Γ † .
We assign an oriented weight w(e) = w Γ (e) to every oriented edge of Γ. This defines a continuous time Markov chain on Γ where for every edge (x, y) in Γ, the chain jumps from x to y at rate w(x, y), and with a slight abuse of terminology we will from now on refer to this chain as the random walk on Γ. Note that we defined it in continuous time for simplicity but we will in fact always be interested only by the geometry of its paths so its time parametrisation is mostly irrelevant. With a slight abuse of notations, we will often identify the trajectory of the random walk (or other processes such as loop-erased paths) with the continuous path obtained from the set of edges used by the walk, either seen as a path or as a subset of M . Similarly, with a slight abuse of notation, we will use direct image notations to denote sub-paths of the random walk, i.e if s < t are two times, we will denote by X[s, t] the path followed by the random walk between times s and t. Finally we will always stop the random walk whenever it reaches a boundary vertex so we assume that all oriented edges starting at such a vertex have weight 0. We emphasize that there are no weights on the dual graph Γ † and indeed we will never consider the random walk on the dual graph.
Γ and Γ †Ĝ G Figure 1: The graphs Γ, Γ † , and G. In the leftmost picture Γ is pictured in red and Γ † in black. In the rightmost picture note that vertices in B appear as either a black or a red dot. The graph is Temperleyan in the sense that all the B vertices on the boundary are of one type only, namely black (i.e., vertices of Γ).
Remark 2.1. To explain the asymmetry between Γ and Γ † , we point out that Γ can be thought of as a graph that is wired at the boundary of ∂M whereas Γ † has free boundary conditions.
In order to apply locally the bijection between dimers and (a variant of) spanning trees, we need to define the superposition of Γ and Γ † . We introduce a new set of vertices W at each point where some edge e and dual edge e † intersect. Define the graphĜ whose vertices are V = V (Γ) ∪ V (Γ † ) ∪ W and the edges given by joining every vertex in W which is on the edge e with the endpoints of e and e † (in other words, each pair of edges e and e † corresponds to four edges inĜ). We callĜ the superposition graph. Then define G to be the graph where we remove fromĜ the boundary vertices ∂Γ and its corresponding half edges. See Figure 1 for an illustration. Clearly, G is a bipartite graph where we take W to be the set of white vertices and the rest to be black. Also every non-outer face of G is bounded by 4 edges (i.e. every non-outer face is a quadrangle). We can naturally define unoriented weights on G from the ones on Γ. If w is the white vertex in the middle of the edge (x, y) of Γ, we set w G (x, w) = w Γ (x, y) and w G (y, w) = w Γ (y, x) while w G (e) = 1 if e is part of an edge of Γ † .
We will prove in Section 3 that if the manifold M is not a torus or an annulus, one must introduce punctures in the graph G for it to even admit dimer coverings. More precisely we will remove finitely many white vertices from G and call a removed vertex and the edges associated with it a puncture (equivalently, the puncture may be thought of as forcing the vertex to be a monomer instead of in a dimer edge). Note that each puncture creates an octagon G consisting of four white vertices, two black vertices from Γ and two black vertices from Γ † (see Figure 2). We will call the resulting graphs respectively G , Γ , (Γ † ) and we call the manifold obtained by removing the white vertices M . We will see later in Section 3 that the number of punctures we remove must be k = |χ| where χ is the Euler characteristic of M . When this is the case we will call G a Temperleyan graph on M or Temperleyan discretisation of M .
When considering scaling limits, we will of course consider a sequence of graphs Γ #δ , (Γ † ) #δ , G #δ from the setup above. We now introduce our assumptions on such a sequence.
Let p :M → M be a conformal lift to the universal coverM , which is either the unit disc or the whole plane 1 . In the end the choice of this lift does not affect the following assumptions, see Remark 2.3 for a precise discussion about this lift and also about the lift of the punctured  (i) (Bounded density) There exists a constant C independent of δ such that for any x ∈ M , the number of vertices of Γ #δ in the ball {z ∈ M : d M (x, z) < δ} is smaller than C.
(ii) (Good embedding) The edges of the graph are embedded as smooth curves and for every compact set K ⊂M , the intrinsic winding of every edge in the liftΓ #δ intersecting K is bounded by a constant C = C K depending only on K. (Note that this allows edges to wind quite a bit near holes.) (iii) (Invariance principle) As δ → 0, the continuous time random walk {X t } t≥0 onΓ #δ started from a nearest vertex to 0 satisfies: where (B t , t ≥ 0) is a two dimensional standard Brownian motion inM (killed when it leavesM , ifM = D) started from 0, and φ is a nondecreasing, continuous, possibly random function satisfying φ(0) = 0 and φ(∞) = ∞. The above convergence holds in law in Skorokhod topology. We remark that the above condition is equivalent to asserting that simple random walk from some fixed vertex converges to Brownian motion on the Riemann surface itself up to time parametrisation (see e.g. [22]).   Figure 3). The uniform crossing condition is the following. There exist universal constants δ 0 , α 0 > 0 such that for every compact set K ⊂ M , there exists a δ K such that for all δ ∈ (0, δ K ) the following is true. LetK be the lift of K. Let R be a subset of one of the connected components ofK and is a translation of cR or cR where c ≥ δ/δ 0 . Let B 1 (resp. B 2 ) be the same translation of cB 1 , cB 1 (resp. cB 2 , cB 2 ). For all v ∈Γ #δ ∩ B 1 , In what follows, sometimes for a compact set K ⊂M , we will write δ K to mean δ p(K) as defined above.
(v) (punctures) We remove k = |χ| many points from M and call the resulting surface M , where χ is the Euler characteristic of M . We assume that each of the (k many) monomers (white vertices removed fromĜ to obtain G ) converge to a unique puncture as δ → 0.

Remark 2.2.
We also remark here that the invariance principle assumption item iii actually implies something stronger in combination with the other assumptions: for any point x inM , the random walk started from a vertex x #δ nearest to x converges to a Brownian motion started from x up to a time change as above. This is a consequence of the fact that random walk from 0 comes close to x with uniformly positive probability using the crossing estimate and the strong Markov property of Brownian motion.
In case ∂M = ∅, recall that the set of boundary cycles (∂Γ † ) #δ corresponds to the connected components of ∂M . One consequence of the Invariance principle assumption (item iii) is that each boundary cycle converges in the Hausdorff metric (induced by dM ) to the associated component of ∂M .
Sometimes we drop the superscript δ from Γ #δ , (Γ † ) #δ for clarity, when there is no possibility of a confusion. Remark 2.3. By the uniformisation theorem of Riemann surfaces, we know that there exists a conformal map from the Riemann surfaceM /F to M where F is a Fuchsian group which is a discrete subgroup of the group of Möbius transforms onM 2 . In the case of the torus, this subgroup is simply a group of translations isomorphic to Z 2 . In the hyperbolic case F is a subgroup of the group of Möbius transforms of the unit disc D. Such a conformal map is unique up to conformal automorphisms (i.e. Möbius transforms) of the unit disc. In other words, if F, F are two Fuchsian groups such that M is conformally equivalent to both D/F and D/F then there exists a Möbius map φ : We also remark that in item iii, we only require convergence up to time change, and hence this assumption depends only on the conformal type of the metric. Finally, one could be worried about the fact that in item iv, the probability α 0 is uniform over the position and scale of the rectangle despite the distortion between two copies of the same set K ⊂ M . However note that the image of a rectangle by a Möbius transform is made of 4 circular arcs and is crossed by Brownian motion with the same probability as the original rectangle, so our assumption is natural in this sense.
Finally, note also that while we have stated these assumptions on the universal cover of M , we could equivalently stated these assumptions on the universal cover of the punctured surface M instead. The fact that these assumptions are equivalent can be checked using the fact that there is a map from the universal cover of M toM \ p −1 ({x 1 , . . . , x k }) (where x 1 , . . . , x k is the set of punctures) which is locally a conformal bijection.
In summary, we could work with any lift (i.e. any choice of F ) both for the punctured and the non-punctured manifold and we fix a particular choice of this lift and call it p if lifting to M or p if lifting to M for concreteness; in fact, with a small abuse of notation we will use p throughout.

Height function and forms
A flow ω is an antisymmetric function on the oriented edges E of G, i.e., for every oriented edge (u, v), ω(u, v) = −ω(v, u). The total flow out of a vertex v is defined to be w∼v ω(v, w). Similarly, the total flow into a vertex v is defined to be w∼v ω(w, v). A flow f is a closed 1-form if the sum over any oriented contractible cycle is 0: i.e., for any oriented cycle It is clear that if M is simply connected, then a closed 1-form also defines a function on the vertices of G up to a global constant. Suppose G is bipartite. We now associate to any dimer configuration m on G a closed 1-form on E. Let m be a dimer configuration on G, and let e = (w, b) be an oriented edge, where w is a white vertex and b a black vertex. We define the flow ω m by setting ω m ( e) = 1 {e∈m} . Also, ω m is defined in an antisymmetric way: ω m ((b, w)) = −ω m ((w, b)). Note that the total flow out of a white vertex is 1 and that out of a black vertex is -1.
To any flow ω on oriented edges, one can associated a dual flow ω † defined on the oriented edges of the dual graph G † , where if e † crosses the edge e = (w, b) with w on its right and b on its left, then we set ω † (e † ) = ω(e). Note also that if ω is divergence free (i.e., the flow out of every vertex is 0), then ω † is a closed 1-form on E † .
Consider any reference flow ω 0 which has total flow out of white vertex equal to 1 and total flow out of a black vertex equal to −1. Then ω = ω m − ω 0 defines a divergence free flow on E. We call ω † the height 1-form corresponding to m with reference flow ω 0 .
When G is embedded on M in such a way that no cycle in G † is noncontractible, every closed 1-form ω on E † becomes exact: i.e., there exists a function on the faces F (G) of G, h : F (G) → R, so that for any two adjacent faces f, f , Observe further that this function is defined only up to a global constant. The function h is then called the height function of the dimer m, admitting an abuse of terminology.
We recall the following simple but useful observation. A path in G (or G † ) is a sequence of vertices (v 0 , . . . , v n ) (or faces (f 0 , . . . , f n )) of G so that v i is adjacent to v i+1 (or f i is adjacent to f i+1 in G † ) for all 0 ≤ i ≤ n − 1.
We now turn to the definition of height function in the more complicated case when M is no longer assumed to be simply connected. In that case, when we sum the values of the height 1-form (defined above) along any noncontractible cycle, we may get a nonzero value. One can use the Hodge decomposition theorem, to isolate out the part of the height 1-form which is encoded by the topology of the underlying surface. The Hodge decomposition theorem works in great generality, but in the present context, it takes the following simple form. For any function f on the vertices of G we define df to be the closed 1-form defined on E as A harmonic 1-form h is a closed 1-form which is divergence free, so that v∼u h(u, v) = 0. Theorem 2.5 (Hodge decomposition [2,8]). For any closed 1-form ω on G (or G † ), there exist a function f on the vertices of G and a harmonic 1-form h defined on E such that and f is unique up to an additive global constant, and h is unique. Furthermore, h is completely determined by summing ω over a finite set of oriented noncontractible cycles which forms the basis of the first homology group of M .
In this paper, we will analyze (f, h) corresponding to the divergence free flow ω m − ω 0 , where m is dimer configuration chosen from the law (1.1) subject to certain natural conditions and ω 0 is a carefully chosen reference flow (in fact, we will consider the height function see Section 6 for a precise statement). We will call h the instanton component. We remark that changing the reference flow changes (f, h) by a deterministic additive factor, and in particular does not affect the fluctuations of (f, h) around their mean. To be more precise, our main theorem (Theorem 6.1) will be stated in terms of the single-valued function associated to ω m − ω 0 on the universal cover of M . See Remark 6.15 for a statement concerning both the scalar and instanton parts of the height function.

Height function on the universal cover
Throughout the paper, rather than working with the scalar and instanton components of the height 1-form, it will be more convenient to lift the height 1-form ω to the universal cover of M . Since the latter is always simply connected, this allows us to work with actual functions without having to worry about the Hodge decomposition Theorem 2.5. We will then check that the convergence of height function on the universal cover implies convergence of each of the components in the Hodge decomposition.
Our assumptions on the graph G where the dimer model lives will be such thatG = p −1 (G) is a planar graph embedded onM . Moreover, the height 1-form ω on the dual edges of G lifts to a height 1-formω on the dual edges ofG. SinceM is simply connected, and sinceω is a closed one-form on the dual edges ofG (this is a local property, so remains true when we lift toG), we can define a height function h = h(m, G) (up to a global constant) on the dual graphG † . The instanton component h can be related to the height function h on the universal cover by summing up the value ofω along any path in the dual graph ofG corresponding to a noncontractible loop in the dual edges of G. This is easier to explain on an example. Example 2.6. If M is the flat torus T := C/(Z+τ Z) for some complex number τ with (τ ) > 0, then the universal cover is the complex plane C. The universal cover can be thought of as many copies of the fundamental domain (a parallelogram determined by 1 and τ ).
Fix v 0 in the fundamental domain. Then by periodicity of dh, the height function h onG evaluated at a point v = v 0 + m + τ n (where m, n ∈ Z) is given by for some a, b ∈ R which do not depend on either v 0 , m, n. Let us describe what a, b are. Consider the two loops on the torus described by L 1 := (t+1/2τ : t ∈ [0, 1]) and L 2 := (1/2+tτ : t ∈ [0, 1]) in the fundamental domain (i.e., L 1 and L 2 are the two noncontractible loops in the torus which form the basis of the homology group). Then a is the sum of the values of ω along any loop in E † which is homotopic to L 1 , whereas b is the sum of the values of ω along any loop which is homotopic to L 2 . Clearly, the choice of these curves in G † do not matter since the height 1-form is closed. Furthermore, in the Hodge decomposition of Theorem 2.5, the harmonic 1-form h is uniquely determined by the numbers a and b.

Intrinsic and topological winding
The goal of this section is to recall several notions of windings of curves drawn in the plane, which we use in this paper. We refer to [4] for a more detailed exposition. A self-avoiding (or simple) curve in C is an injective continuous map γ : The topological winding of such a curve γ around a point p / ∈ γ[0, T ] is defined as follows. We first write where the function θ(t) : [0, ∞) → R is taken to be continuous, which means that it is unique modulo a global additive constant multiple of 2π. We define the winding for an interval of time (note that this is uniquely defined). Notice that if the curve has a derivative at an endpoint of γ, we can take p to be this endpoint by defining With this definition, winding is additive: 3) The notion of intrinsic winding we describe now, also discussed in [4], is perhaps a more natural definition of windings of curves. This notion is the continuous analogue of the discrete winding of non backtracking paths in Z 2 which can be defined just by the number of anticlockwise turns minus the number of clockwise turns. Notice that we do not need to specify a reference point with respect to which we calculate the winding, hence our name "intrinsic" for this notion.
We call a curve smooth if the map γ is smooth (continuously differentiable). Suppose γ is smooth and for all t, γ (t) = 0. We write The total intrinsic winding is again defined to be lim t→T W int (γ, [0, t]) provided this limit exists. Note that this definition does not depend on the parametrisation of γ (except for the assumption of non-zero derivative). The following topological lemma from [4] connects the intrinsic and topological windings for smooth curves.
Lemma 2.8. Let D be a domain and γ ⊂D a simple smooth curve (or piecewise smooth with smooth endpoints). Let ψ be a conformal map on D and let arg ψ (D) be any realisation of argument on ψ (D). Then 3 Temperley's bijection on Riemann surfaces

Notion of Temperleyan cycle rooted spanning forest; bijection
Let Γ be a graph embedded on a surface with a certain specified set of boundary vertices. Before introducing the notion of Temperleyn forest, we start with the simpler notion of Cycle-Rooted Spanning Forest.
Definition 3.1. A wired oriented cycle rooted spanning forest (which we abbreviate: wired oriented CRSF) of Γ with the specified boundary is an oriented subgraph T of Γ where • Every non-boundary vertex of Γ has exactly one outgoing edge in T . Every boundary vertex has no outgoing edge. (As a result, any cycle of T must have a unique orientation).
• Every cycle of T is noncontractible.
To each CRSF t, we may associate the weight e∈t w(e). We will consider the probability measure on CRSFs such that the probability of t is by definition proportional to its weight.
This is equivalent to the notion of essential CRSF on a graph with wired boundary introduced by Kassel and Kenyon [23]. Ignoring the orientation of T gives an unoriented graph, its connected components will simply be called the connected components of T without any additional precision. Note that if T is a wired oriented CRSF, every connected component of T contains at most one cycle: more precisely, every boundary component must have zero cycles, while every non-boundary component contains exactly one cycle.
We will refer to the set of all noncontractible cycles of a wired oriented CRSF to mean the set of unique cycles corresponding to each (non-boundary) component of the wired oriented CRSF.
Let us come back to the setup of Section 2.3 and recall that we had the graph Γ, its dual Γ † and the superposition graph G all embedded nicely in a Riemann surfaces with g handles and b holes. Compared to the planar setting, there is an immediate topological difficulty, which is that this graph G in general does not admit a dimer cover. Indeed by Euler's formula, if G has v vertices, e edges and f many contractible faces, then v−e+f = χ = 2−2g−b. On the other hand if G admits a dimer cover, one must have e = v + f since e is the number of white vertices which must match v + f (the number of vertices contractible faces of Γ † ). Thus we need to remove k = |χ| many edges (where χ := 2 − 2g − b is the Euler's characteristic) from the superposition graph for it to admit a dimer cover (see Figure 2). These removed edges can be thought of as creating punctures in the surface. Call the graph with these punctures G . (Recall that as per item v, we assume throughout that the punctures are macroscopically far apart in M and converge as δ → 0 to a set of pairwise distinct points in M ). Note that if M is a torus or an annulus then χ = 0 so no punctures need to be removed. See also Ciucu-Krattenthaler [13] and Dubédat [18] for other situations where punctured dimers arise.
We now come to the crucial definition of this paper, in which we introduce the notion of Temperlayan forest. Recall that Temperleyan forests will be shown to be in bijection with a dimer cover on G (and the bijection will also verify as a corollary that G indeed has a dimer cover). For every wired oriented CRSF T of Γ with boundary ∂Γ, one obtains a natural dual free cycle rooted spanning forest T † (abbreviated free oriented CRSF) of (Γ ) † as follows. The vertices of T † are given by the vertices of (Γ ) † (i.e., it spans (Γ ) † ) and an edge e † is present in T † if and only if its dual e is absent in T . Note a priori that T † does not come with an orientation. This is the reason why not every CRSF is associated with a dimer configuration (recall that in the classical simply connected case both the primal and the dual trees come with a canonical orientation and the dimer configuration is obtained from the pair of dual oriented spanning trees in Temperley's bijection by placing a dimer on the "first half" of each oriented edge in both trees). The following adds a condition to the CRSF which allows to assign an unambiguous orientation to the dual T † of T . Definition 3.2. We say that the wired oriented CRSF T is Temperleyan if each connected component of T † contains exactly one cycle.
An example of a wired oriented CRSF T that is not Temperleyan is provided in Figure 4. Temperleyan forests also come with a natural law, which is simply the law of the CRSF conditioned on the event that each component of the dual contains exactly one cycle. We will soon formulate a more explicit criterion for a CRSF T to be Temperleyan (see Section 3.2 and in particular Theorem 3.6). For now, observe that if T is a Temperleyan wired oriented CRSF and T † is its dual, we can assign an orientation to each cycle in each component of T † T T † Figure 4: A non-Temperleyan CRSF (in blue). The surface M is the "pair of pants": a domain of the plane with two holes (in grey in the picture). In this example, T does not contain any cycle, and hence any connected component flows to the boundary of M . Its dual T † must contain a component with two cycles which overlap: the cycles go around each of the two holes, and must touch each other, as otherwise there would have to be a path in T separating them, but this is impossible as such as path would have to connect two distinct boundary points. So T is not Temperleyan. Note that 2g + b = 3 > 2 here. See Lemma 3.5 for a more general argument.
arbitrarily from one of the two possible choices. Then we orient all other edges of T † towards the cycle of that component. We let T be the set of pairs (t, t † ) where t is a Temperleyan wired oriented CRSF, t † is its dual (hence a free CRSF) for which an orientation of its cycles has been specified, which allows us to view t † also as a free, oriented, CRSF such that each vertex has a single outgoing edge attached to it. Note that if (t, t † ) ∈ T then we call t a Temperleyan CRSF and (t, t † ) a Temperleyan pair or Temperleyan self-dual pair. We hope this terminology will not be too confusing: t is the object we will mostly work with, and t determines t † uniquely up to the orientation of its cycles.
We now state the Temperleyan bijection for general surfaces. Recall that we assign oriented weights w e to each edge e in Γ , no weight (or unit weight) to edges of Γ † and that this turns G into a weighted unoriented graph. Indeed, if e = (x, y) is an oriented edge of Γ , let w denote the white vertex in the middle of e. Then we assign to the unoriented edge {x, w} of G the weight w (x,y) and to the edge {w, y} of G the weight w (y,x) .
For every Temperleyan pair (t, t † ) ∈ T , we define the measure where Z Temp is the partition function. Furthermore if (T , T † ) has the law (3.1) then m = ψ((T , T † )) has law (1.1) with unoriented weights on G described above.
We emphasise that it is unclear at this point whether the measure P Temp (and consequently the measure on dimers) is not trivially zero as we have only deduced that the construction of G is necessary, but have not deduced that it is sufficient. We prove this later in Proposition 3.4 using a pants decomposition argument. Proof of Theorem 3.3. We assume that the set of Temperleyan CRSF and the set of dimer covers are both nonempty (this assumption is validated in Proposition 3.4). Given a Temperleyan pair (t, t † ), we obtain a configuration of edges m = ψ((t, t † )) as follows: for every oriented edge e ∈ t (resp. e ∈ t † ), we can write e = e 1 ∪ e 2 where e 1 , e 2 are the first and second halves of e (recall that e is oriented, whether in t or in t † ), and let e 1 ∈ m (see Figure 5). Note that m is a matching on G because every (non boundary) vertex has a unique outgoing edge in either t or t † . Furthermore, since t ∪ t † spans the black vertices of G , the matching is a perfect matching.
Also ψ is injective: if (t 1 , t † 1 ) and (t 2 , t † 2 ) are two distinct elements of T , then there must a black vertex v on G (i.e., a vertex of Γ or (Γ † ) ) such that the unique outgoing edge from v in t 1 or t † 1 is different from the unique outgoing edge from v in t 2 or t † 2 . Hence v will be matched to two distinct white vertices in ψ((t 1 , t † 1 )) and ψ((t 2 , t † 2 )). We now check ψ is onto. Given a matching m of G , we can obtain a pair (t, t † ) by extending the matched edges: for every nonboundary black vertex v of G (i.e., a vertex of Γ or (Γ † ) not on a boundary cycle), let the unique outgoing edge from v be the edge of Γ or (Γ † ) containing the white vertex to which v is matched in m. The fact that neither t nor t † contain contractible cycles is the same as in the standard, planar case: if say t has a contractible cycle C, then since Γ \ C has a dimer cover, an elementary counting argument shows that v − e + f = 0 where v, e, f are the number of vertices, edges and faces of the contractible component Γ \ C. On the other hand, Euler's formula implies v − e + f = 1 (again excluding the outer face). Since every vertex v of t † has a unique outgoing edge, t † must have one cycle per component and so t is Temperleyan. This concludes the proof.

Criterion for a wired CRSF to be Temperleyan
In this Section, we prove the following lemma: Proposition 3.4. The graph G obtained in Section 2.3 has a dimer cover. In particular, the measure P Temp defined in eq. (3.1) is a probability measure.
As a by product of the proof of Proposition 3.4, we derive a simple criterion for a CRSF to be Temperleyan (Theorem 3.6).
First of all note that if we can find a Temperleyan CRSF of Γ then by Theorem 3.3, G has a dimer cover. Next note that if M has the topology of an annulus (with all the nice properties of Section 2.1), then finding a Temperleyan oriented CRSF is straightforward. Indeed, any wired spanning forest in annulus (where both boundaries are wired) is Temperleyan: the dual is a graph containing a single cycle separating the two components touching each boundary. Also notice that for a torus, every oriented CRSF is Temperleyan [19]: essentially an oriented CRSF must contain a cycle (since there is no boundary on a torus) and cutting along this cycle gives us a (bounded) cylinder or equivalently an annulus. The converse is also obviously true when we do not remove edges: Lemma 3.5. Let M be nice with g handles and b boundary components and Γ, Γ † , G be embedded as above (i.e. without punctures). A wired Temperleyan oriented CRSF of Γ exists if and only if M has the topology of either a torus or an annulus. Furthermore in these cases, all oriented CRSF are Temperleyan.
Proof. This is simply a consequence of the extended Tempereley's bijection (Theorem 3.3), because a Temperleyan oriented CRSF and an oriented dual correspond to a dimer configuration (by endowing each cycle with arbitrary orientation and orienting every other edge towards the unique cycle of its component). However, a dimer configuration only exists if and only if χ = 0 or equivalently 2g + b = 2 (as discussed in Section 3.1). This equation has only two feasible solutions: g = 1, b = 0 (i.e. a torus) and g = 0, b = 2 (i.e. an annulus). The last part was already argued above the statement.
Proof of Proposition 3.4. In light of Lemma 3.5, we assume M is neither a torus, nor an annulus. Now recall the pants decomposition [21]: every topological surface of genus g with b boundary components can be decomposed as a finite union of k = 2g + b − 2 pairs of pants. That is, we can find continuous paths on M forming simple noncontractible cycles, so that if we cut open M along these cycles, each component is homeomorphic to a pair of pants. Without loss of generality we can assume that these paths consist of edges from Γ : indeed, any maximal collection of simple closed curves which are disjoint, pairwise non-homotopic as well as non-homotopic to a boundary component is a suitable collection of such paths (Theorem 9.7.1 in [31]). Clearly, the restriction of Γ and (Γ † ) to each such component can be viewed as a pair of dual graphs embedded in a manifold with the topology of a pair of pants where boundary components coming from a cut are described by a single boundary cycle in Γ . Removing these boundary cycles and the attached half-edges results in a pair of dual graphs faithfully embedded on the pair of pants, exactly as described in Section 2.3 (in particular, each boundary component is now associated with a boundary cycle in (Γ † ) and the boundary condition remains 'wired' for Γ).
Furthermore, since k = 2g + b − 2 is also the number of white monomers (punctures) which we remove to get G , we can assume without loss of generality that there is exactly one white vertex removed from each component having the topology of the pair of pants, and that such a vertex does not belong to the boundary of Γ restricted to that component.
We now claim it suffices to prove Proposition 3.4 when M has the topology of a pair of pants. Indeed, for each cycle arising in the pants decomposition, we fix an arbitrary orientation of this cycle. This defines as in Theorem 3.3 a dimer configuration on the edges of the cycle. Recall that the complement of these cycles define faithfully embedded dual graphs in a number of pair of pants. In each such pair of pants P , we superpose a dimer configuration on the graph G restricted to P . Superposing these dimer configurations in each pair of pants and on the separating cycles gives rise to a dimer configuration on the whole graph G : note that there are no conflicts on the separating cycles, since when we removed these cycles to obtain the pair of pants we also removed the half edges attached to them, so these will never be occupied by dimers. Figure 6: An illustration of the proof of Proposition 3.4 for a pair of pants. The skeleton consists of p and q, which decompose the surface into a number of disconnected annuli: three in the first case, two in the second. The dotted edges are the ones removed from G to get G .
We thus now assume that M is a pair of pants with a pair of dual graphs faithfully embedded onto it as described in Section 2.3. Let v 1 and v 2 be the two vertices of Γ which are in the octagon formed because of the removal of white vertex and its attached half-edges. Now consider two disjoint oriented paths p and q in Γ starting from v 1 and v 2 forming a cycle around two boundaries as shown in Figure 6. If we cut along the loop formed by the paths p, q and the octagon, we get three annuli with faithfully embedded dual graphs (note again that removing the cycles formed by the paths p and q and the attached half-edges means that the each boundary component in the resulting annuli is associated to a boundary cycle in (Γ † ) , as in Section 2.3). Thus, by fixing an orientation of p, q and introducing a matching as before, we are back to the annulus case. This we have dealt before in Lemma 3.5, and hence the proof of Proposition 3.4 is complete.
We now deduce from the above an extremely convenient criterion for a wired oriented CRSF to be Temperleyan. Let us supopse M is not the torus or an annulus, whence k = 2g + b − 2 > 0 (we already know that every CRSF is Temperleyan otherwise). Define the branch starting from a primal vertex v of Γ to be the path obtained by going along the unique outgoing edge from each vertex (which necessarily ends when a loop is formed or a boundary is hit). Recall that, at a puncture (i.e., a white monomer), there are exactly two vertices u i , v i of Γ on either side of the puncture. Let B i1 , B i2 be the branches in Γ of u i , v i for 1 ≤ i ≤ k = 2g + b − 2. We call the union of these branches, together with the removed edges e i , the skeleton s of T . That is, Simply put, the skeleton of T is the union of the branches emanating out of the punctures. (Of course, in the case of an annulus, as already noted there are no punctures and so s = ∅, so the condition is automatically satisfied, as we already know.) The criterion for a CRSF to be Temperleyan is thus just to say that the skeleton cuts the surface into disjoint annuli. Figure 6 gives two examples on a surface M (the 'pair of pants') where M \ s consists of topological annuli (three annuli in the first example, and two in the second example). In the companion paper [5] we denote this conditioning event by A M . A Temperleyan forest may thus be viewed as a CRSF conditioned on the topological condition A M , which is asymptotically degenerate (i.e., the probability of A M for a wired CRSF tends to zero as the mesh size δ tends to zero).
Proof. First note that a component can only have the topology of a torus if the manifold is a torus in which case there is nothing to prove by Lemma 3.5.
For the general case, let T be a Temperleyan CRSF and let T † be its dual with a choice of orientation. Note that the vertices (in G ) of s are all matched with each other in the dimer configuration associated to (T , T † ). Therefore in each component of M \ s, all vertices are also matched with each other. By Lemma 3.5, this implies that these components are annuli.
Conversely, note that by definition T † cannot cross s, and so T can be restricted to each component of M \ s to form a wired CRSF. By the other implication in Lemma 3.5, it follows that T is Temperleyan in each such component and so is globally Temperleyan.
The significance of this criterion is as follows. By Theorem 3.6, a Temperleyan forest can be thought of as a wired CRSF conditioned on the event in the statement of that proposition. However, wired CRSF are easier objects to understand owing to the fact that they can be sampled through a version of Wilson's algorithm, as will be recalled in Section 3.3.

Wilson's algorithm to generate wired oriented CRSF
Recall the measure P Temp , i.e., the law of (T , T † ) on T under the measure eq. (3.1). We will not study directly P Temp but rather a measure on T which can be sampled through Wilson's algorithm and is defined as follows: we first sample a wired oriented CRSF of Γ with law and given T = t, we pick T † an oriented dual uniformly among all possibilities of orientation of the dual. Thus P Wils can be viewed also (with a small abuse of notation) as a mesure on T . The two laws look similar but are in fact different due to the fact that any cycle of the dual t † of a Temperleyan oriented CRSF t can be oriented in two possible ways to determine a dual pair (t, t † ) ∈ T . We deduce the following relationship: Lemma 3.7. Let (t, t † ) be a Temperleyan wired oriented CRSF such that t † contains exactly k † noncontractible cycles. Then the Radon-Nikodym derivative satisfies In particular, conditioned on having k † noncontractible cycles for the dual forest, the law P Temp and P Wils coincide.
Let Γ, Γ † be faithfully embedded on a nice Riemann surface M . We now describe Wilson's algorithm to generate a wired (but not necessarily Temperleyan) oriented CRSF on Γ. We prescribe an ordering of the vertices (v 0 , v 1 , . . .) of Γ.
• We start from v 0 and perform a loop-erased random walk until a noncontractible cycle is created or a boundary vertex (i.e., a vertex in ∂Γ) is hit.
• We start from the next vertex in the ordering which is not included in what we sampled so far and start a loop-erased random walk from it. We stop if we create a noncontractible cycle or hit the part of vertices we have sampled before.
There is a natural orientation of the subgraph created since from every non-boundary vertex there is exactly one outgoing edge through which the loop erased walk exits a vertex after visiting it. LetP Wils be the law of the resulting wired oriented CRSF.
Proposition 3.8. We haveP In particular,P Wils generates a wired oriented CRSF of Γ as described by Definition 3.1 with law given by (3.2). Furthermore, conditionally on being Temperleyan,P Wils coincides with the first marginal of P Wils .

Winding and height function
In this section, we explain the connection between winding of the Temperleyan forest and height function. In order to account for the curvature of the surface it will be important to work on the universal cover and the lifts of both the dimer configuation and the Temperleyan CRSF to it. This will also have the advantage that the dimer height one-form (defined properly below) becomes an actual function on the faces of the lift of the Temperleyan graph embedded on the surface. We refer the reader to Sections 2.4-2.6 for relevant definitions and notations. The theory we develop below is similar to [27] but with a few important modifications, related in particular to the fact that the Temperleyan forest is typically not connected. See also [35] for a version on the torus.
Recall thatM = C only for the case of the flat torus in this article. We recall the notation G for a Temperleyan graph embedded in M and also recall Γ , (Γ † ) . Now we need to lift the graph G and define a height function which is consistent with this lift.
At this point, we need to spare a few words related to the removal of white vertices from the graph G to obtain a graph with a dimer cover. This operation can be interpreted as inserting certain discrete version of magnetic operators on the free field (e.g. in the sense of [18]). If we want to interpret the height function as winding, the height function would be additively multivalued where it picks up an additional ±2π winding when it goes around a removed white vertex. This motivates us to introduce a puncture corresponding to each face obtained by removal of a white vertex (cf. Figure 2) and call the new manifold M . Recall that we assume each of the punctures converge to fixed distinct points in the manifold (Section 2.3). We now treat every face with a puncture as an outer face.
We lift G to the universal coverM of M and call itG . Note that the lift of every non-outer face of G is a quadrangle. For every face f ofG , we fix once and for all a diagonal d(f ) joining the two black vertices. We assume without loss of generality that the diagonal is a smooth curve lying completely inside the quadrangle.
We take a dimer configuration of G which corresponds to a dimer configuration onG . We now wish to relate the dimer height function to the winding of a wired oriented CRSF branch adjacent to a path joining the faces, up to certain explicit terms describing the effect of jumping from one component of the CRSF to another. In [27], this connection was established for trees with straight line embeddings in the simply connected case. In [35], the toroidal case was treated, but with only straight line embeddings. In what follows, we need to define the height function properly not necessarily for just straight line embeddings, but for any arbitrary embedding which is smooth except at the vertices. We also need to deal with the fact that the CRSF might have many connected components. The main result of this section, which summarises the desired relationship is stated in Theorem 4.10.

Winding field of embedded trees and choice of reference flow
To describe properly the connection between winding field of trees and height function of dimers, one needs to get past certain technicalities. One technicality with this connection is that the height function is defined on the faces of the graph and the winding should ideally be computed between vertices of the spanning forest. Another issue we have is that now we need to deal with a spanning forest, and hence the winding between vertices belonging to different connected components must be properly defined. Finally, the definitions should ideally be symmetric with respect to the primal and dual trees.
In this subsection, we temporarily forget about dimers and Temperleyan forests, and focus on how to compute the winding field of a deterministic, infinite, one-ended, spanning tree and its dual tree which are embedded smoothly in the complex plane. This will create a simple setup for connecting dimer height function with winding of Temperleyan spanning trees. Hence the rest of this subsection could be read independently of the rest of the paper.
Consider an infinite one-ended tree T embedded in C with the edges embedded smoothly (although having points of non-differentiability at the vertices, which we call corners, are allowed, so that the branches are only piecewise smooth). Recall that the intrinsic winding of any finite branch in this tree is well-defined. Since the infinite tree is one-ended, we can orient the tree towards that unique end. From every vertex x on the tree, we can define an infinite branch γ x from x to infinity. Suppose for now (4.1) Then we can define a winding field {h T (x) : x ∈ T } of the tree to be simply h T (x) := W int (γ x ).
The following elementary lemma expresses the height in a way which later allows us to extend the definition of h T even if (4.1) is not satisfied. If x / ∈ γ y and y / ∈ γ x , notice that γ x and γ y eventually merge and γ y merges either to the right or to the left of γ x (since the tree is embedded in C, this makes sense). Let γ xy be the unique path connecting x and y in T .
Lemma 4.1. In the above setup, if γ y merges with γ x to its right then If γ y merges with γ x to its left, then If y ∈ γ x and y is not a corner (or x ∈ γ(y) and x is not a corner), then Proof. Notice that the last assertion follows simply from additivity of intrinsic winding. Indeed, for example if y ∈ γ x , there is no discontinuity of intrinsic winding at y since y is not a corner.
For the rest, take m ∈ γ x ∩ γ y . Notice that by additivity of winding, . This is clear since we need to do a half turn to move from γ xm to γ my at m and the turn is clockwise if γ y merges to the right of γ x and anticlockwise otherwise.
We remark that the last assertion of the above lemma works even at corners by adding the appropriate angle at the corner to match the winding field difference. In what follows, we avoid using winding field at corners, hence this ambiguity would not bother us. We also remark that such a definition easily extends to a finite tree with the branches oriented towards a fixed vertex in the tree in place of being oriented towards "infinity". Finally the formula for h(x) − h(y) Figure 7: The definition of the height function in terms of winding. described in Lemma 4.1 can be taken to be the definition of the winding field (or rather its gradient) for any tree embedded with piecewise smooth edges, even if (4.1) is not satisfied. This will be the typical situation for our setup.
We will now extend the definition of the winding field of T to the faces of a graph spanned by T . For context, we remind the reader that in Section 3.1, a bijection between Temperleyan CRSF and dimer configurations was established. Also recall that the height function of a dimer configuration is defined on the faces of the graph. Thus it is necessary to do this extension carefully so that the dimer height differences become exactly the same as the winding field, perhaps up to some global topological error term coming from the jumps between various components in the forest.
With this in mind, take an infinite locally finite graph Γ embedded smoothly in C, except perhaps at the vertices where we might have corners and let T be a one ended spanning tree of Γ (we emphasise that we are still considering a spanning tree for the moment and not yet a forest). Therefore the dual spanning tree T † of the dual graph Γ † is also one-ended. Let G be the superposition graph (in C), as introduced in Section 2.3 (the appropriate infinite version of the graph G there). It will be useful to augment the trees T and T † to include diagonals, as follows. We fix a point on the diagonal d(f ) once and for all, which we denote by m(f ) and sometimes refer to as midpoint of the diagonal by an abuse of terminology. For each face, add to T (resp. T † ) the portion of the diagonal d(f ) connecting the point m(f ) to the unique primal (dual) vertex touching d(f ). This way the primal and dual trees meet in each face f at the (smooth) point m(f ) on the diagonal d(f ) of that face. With a small abuse of notation, we will still denote T and T † these augmented trees.
Note that we can orient T , T † towards their unique ends. This allows us to define two winding fields as above (one for T and one for T † ). Having oriented T and T † we can also apply the local operation of Figure 5 as described in Section 3.1 to obtain a dimer matching of G (recall that for the moment T and T † are still a dual pair of spanning trees, not forests). We will first need to define a suitable reference flow on G, which will then allow us to speak of the height function associated to the dimer configuration and then show the relation between this dimer height function and the two above winding fields.
where   Proof. Let w be a white vertex. Notice that in b∼w ω ref (wb), the oriented diagonals form a clockwise loop whose total winding is −2π. Adding +π for each of the surrounding black vertices, and dividing by 2π, we see that the total flow out of w is indeed 1 as desired.
For a black vertex b, the argument is better explained by considering a picture (see Figure 8). Fatten the "star" formed by the half-diagonals incident to b into a star shaped domain. Then notice that the total flow out of b is simply the limit of the total winding, divided by 2π, of the boundary of this domain (again in the clockwise orientation this time), as the domain thins into the "star". Indeed, the −π term in the definition of ω ref in (4.2) counts the half-turn as we move from the left side to the right side of a half-diagonal. Since the total winding of such a curve is −2π, it follows that the total flow out of b is −1, as desired.
We are now ready to relate the three notions of height function defined by the pair of dual spanning trees T , T † . To do so, note that we can extend the definition of the winding field of both T and T † from Lemma 4.1 to the augmented trees.
Remark 4.5. Note that having augmented the trees as explained above, the winding fields in T and T † now play a completely symmetric role.
Proof. Define γ v for the branch of T starting from v and similarly define γ † v . For a face f , we define γ f (respectively γ † f ) to be the branch of the augmented tree T (resp. T † ) starting from m(f ). Fix faces f and f and assume without loss of generality that γ f is to the right of γ f . Note that this means that γ † f is to the left of γ † f . Also note that because γ f f concatenated with γ † f f forms a clockwise loop and there is no jump at m(f ) and m(f ) because by assumption the midpoint m(f ) is a smooth point of d(f ). The first equality easily follows by applying the definition of winding field from Lemma 4.1.
For the last equality let f l , f r be adjacent faces. Let the common (oriented) edge be (bw) with b being the black and w the white vertex and let f r lie to its right. We assume without loss of generality that b is a primal vertex. From the definition of ω ref and recalling the sign convention of the flow defining the height function (Section 2.4) Note that that γ fr and γ f l merge at b. Also note from Temperleyan bijection that if (bw) is occupied by a dimer, then γ b starts by using the (half) edge bw which implies that γ f l lies to the left of γ fr . Otherwise, γ f l lies to the left of γ fr . We conclude using the definition of winding field from Lemma 4.1 and the above equation.
As a step towards extending the correspondence of Proposition 4.4 to forests, we now explain how to read off the height change along a path in the graph which does not necessarily follow the tree branches. Note that although the dimer height function is independent of the chosen path, it is not clear how to see this path independence looking at just the winding field. For an arbitrary path in the graph, there could be several "jumps" over edges in the dual tree and it turns out that these jumps contribute an extra ±π on top of the winding, depending on orientation. It will turn out that in case of forests, there is an analogous contribution for jumping over dual components which separate the primal connected components (and vice-versa). Let P be a self avoiding path in Γ. We can partition this path into segments belonging to T separated by edges not belonging to T (we remind the reader that T is still a tree, not yet a forest). Let us call these segments (P 1 , . . . , P k ) and let (e 1 , . . . , e k−1 ) be the edges in Γ separating them. Note that we allow P i to be a single vertex. Let (x i , y i ) be the starting and ending points of P i . Now we modify P i as follows. Observe that the oriented edge (y i , x i+1 ) has two faces of G to its left. Call the one incident to y i , f − i and the one incident to the midpoint of the diagonals of these faces. We pick a face incident to x 1 and y k and add the segments joining these vertices and the midpoint of their diagonals. Call these endpoints m + 0 , m − k respectively. We also join y i to x i+1 using the diagonal segments of f − i and f + i . Finally we delete the edges e i . This completes the modification of the path P , which we are still going to denote as P (see Figure 9).
We also partition P into segments as Proof. By summing along each segment, the first equalty is simply a consequence of the definition and additivity of h T and Proposition 4.4. The second equality is also an implication of the second equality of Proposition 4.4.
Extension to forests. We now extend the definition of winding field to forests before relating it to a dimer height function. Suppose T is now a spanning forest of Γ where each component is infinite. Let T † be its dual, and note that its components are also necessarily infinite, and it is also a forest. We fix an end of each component of primal and dual, and orient the trees towards that end. Thus we obtain an oriented pair of dual forests. Note again that by the local operation of Temperleyan bijection, we can find a dimer configuration on the superposition graph G corresponding to (T , T † ). Let h dim be the dimer height function corresponding to this dimer configuration with reference flow given by ω ref defined in Definition 4.2.
We augment (T , T † ) in the same way as in the connected case: in each face join the midpoint of the diagonals to the corresponding primal vertices in T and dual vertices in T † .
Definition 4.7. In the above setup, define the winding field of (T , T † ) to be the same as that given by the first equality of the formula in Lemma 4.6. This is well defined up to a global shift by R.
Notice that in the definition of δ i above (the error term we need to add when "jumping" over a dual edge) is defined using the dual tree, and hence there is no issue of connectivity when dealing with spanning forests rather than spanning trees. It is straightforward to see that the above notion of winding field for forests is well-defined as the second equality of Lemma 4.6 also tells us that this quantity is the same as the dimer height function (up to a global shift), and in particular is independent of the path P . The proof of Lemma 4.6 relies on additivity of winding for each edge and does not depend on the fact that (T , T † ) is a tree-dual tree pair and not a forest.
Remark 4.8. Note that there is no assumption about the number of ends of (T , T † ) in the above definition: it makes sense as long as an orientation is fixed towards a unique end in each component. Later we will see that because of toplogical reasons, each component (either of primal or dual) of the relevant (random) lift of the Temperleyan forests we need to deal with has exactly two ends almost surely. This leaves us with a choice of which end to orient for each connected component, and once this choice is fixed, the above definition can be applied to define the winding field, which is going to be equal to the dimer height function (up to a global shift). But since the above definition for general forests with arbitrary ends require no additional extra effort, we decided to keep this definition.

Winding of CRSF and height function
In this section, we show how to extend the connection between winding and height function (established in [27]) from the simply connected setting in the Euclidean plane with straight line embeddings to height one forms on general surfaces. Firstly we work in the universal cover if the surface, so that the CRSF is lifted to a spanning forest for which we have already worked out the connection in the previous section. However there are many ways to choose a path between vertices in different components. Most of the work in this section will actually be to specify a convenient one so that the winding computation is as simple as possible.
Recall Γ, Γ † from Section 2.3. We call the vertices in the lift of (Γ † ) the dual vertices and those in the lift of Γ primal vertices. Recall that the dual of a wired Temperleyan CRSF of Γ is a free CRSF of (Γ † ) . In this setup, the punctures are mapped to the boundary of the disc. Furthermore, each component of both the wired Temperleyan CRSF and its dual contains either a single cycle or contains a boundary component. In the rest of this section we drop the adjective Temperleyan and refer simply to the wired and free CRSF.
Roughly this choice of path is as follows. As a preparatory step, we first observe that the cycles in the CRSF, when lifted, become bi-infinite paths, which we call spines. Next we show in Lemma 4.9 that: • In the hyperbolic case, the spines converge to some point on the boundary of the hyperbolic disc (i.e., on the unit circle), • In the parabolic case of the torus, the spines converge to infinity in some fixed asymptotic direction.
We could not find a soft probabilistic argument for this and had to rely on tools coming from the classical theory of Riemann surfaces (essentially, uniformisation and Fuchsian theory). We will postpone the proof of Lemma 4.9 to the appendix in order to maintain the flow in the paper. Armed with Lemma 4.9, given x and y, we choose the path between them as follows. Suppose we are in the hyperbolic case. We first move from x along the spine of its tree component to infinity and do the same for y. These paths converge to the unit circle at ζ x , ζ y . We join ζ x , ζ y along an arc on the unit circle (there are two choices of this arc, and a canonical choice is possible). This defines a path and we prove in Theorem 4.10 that the height gap between x and y can be written as the winding of this path plus a sum of terms of the form ±π depending on the "jumps" over the primal and dual spines lying in between. Note that this is a non-trivial fact as the path is not completely in the graph, but "goes through infinity". An analogous path is chosen in the torus by choosing points along the corresponding spines in place of ζ x , ζ y and then taking a limit. This choice is particularly useful as we do not need to compute the winding coming from the components lying in between trees containing x and y, which would be hard to keep track of using Wilson's algorithm. When relating the height function to winding, it is convenient to replace each boundary vertex of ∂Γ by a cycle surrounding the corresponding hole in the surface. In this manner, all paths of the wired CRSF that end up in the boundary can also be considered to form loopscall this also a boundary loop in what follows. Hence all components of both the wired and free CRSF can be viewed as a set of paths flowing towards a single noncontractible cycle.
It is clear from the unique path lifting property that each oriented loop in the CRSF (including any "boundary loop" mentioned above, if there is any) corresponds to a bi-infinite simple path in the lift (the one obtained by going along the loop infinitely many times in clockwise and anti-clockwise direction). We call this a spine. Note also that each loop corresponds in the lift to multiple spines (except in the case of the annulus). We now state a useful lemma about the geometry of spines. The proof of Lemma 4.9 uses classical tools from the theory of Riemann surfaces, so we postpone to appendix A in order to not disrupt the rest of the argument too much. Lemma 4.9 allows us to describe the geometry of spines and dual spines. In the case of the torus, it is clear that spines and dual spines form alternating bi-infinite paths, all in the same asymptotic direction. When the universal cover is the disc, each (oriented) spine separates D into two simply connected components which lie to its right and left. Given two spines S and S , we have a well defined region Ω S,S between them which is bounded by S, S and two portions of ∂D (if S and S are loops, then these portions of ∂D are understood as only prime ends associated with the same point). We say that a spine or dual spine S separates S from S if it connects these two portions of ∂D. Since the graphG only has accumulation points on ∂D, it is easy to see that two spines can only be separated by finitely many others. Furthermore, two adjacent dual spines must be separated by a primal spine. Let us call the lifts also (T , T † ) with an abuse of notation. Also augment these trees by joining the midpoint of the diagonals to the respective endpoints, as explained in Section 4.1. Now pick two faces f, f and let S, S be the spines of the components of T containing m(f ), m(f ).
Let Ω S,S be the component between them as above. Draw a simple curve connecting S and S in Ω S,S , and oriented from S to S . For each primal spine σ ⊂ Ω S,S (hence not including S and S ), let ε σ be the algebraic number of times σ is crossed by this curve (with the convention of counting +1 if σ is crossed from its left to its right by the curve and -1 otherwise). If σ is not crossed, define ε σ to be 0. Define δ σ similarly for primal spines. Notice that is a topological term which does not depend on the choice of the curve (where the sum is over all primal and dual spines contained in Ω S,S ). Also, its only dependence on f and f is through the spines S and S . We will need a slightly modified definition of ε S and ε S . For any face f , let γ f be the infinite oriented path in T started from m(f ) and going off to infinity along the unique outgoing oriented edges. Let ζ be the limit point of γ f (which exists due to Lemma 4.9). Let ζ be the limiting endpoint of S which lies in the same connected component of ∂M ∩ Ω S,S as ζ. Note that we are not defining ζ to be the limit point of γ f as it could be the case that ζ and ζ might lie in different boundary components of ∂M ∩ Ω S,S (depending on the orientation of S ). In case S, S are loops through the same boundary point, we want ζ and ζ to be in the same prime end of Ω S,S . Now we define ε S and ε S . Define ε S = +1. Now we have a few cases for ε S .
• If the limit point of γ f lies in the same component of ∂M ∩ Ω S,S , then define ε S = 1.
• If not, the there are two cases. If f lies in Ω S,S , then define ε S = −1. Otherwise, define ε S = +1.
The above choice of ε S , ε S is made in accordance with Lemma 4.6. Indeed imagine a path not going through infinity, but which starts at f moves up along the spine S in the direction of its orientation, joins two vertices in S, S very close to the boundary arc and then moves down along S (either in the direction of the orientation of S or the opposite, depending on the orientation of S itself). Note that for this path, the choice of the ε S , ε S exactly matches with the choice defined for Lemma 4.6. We refer the readers to the second figure of Figure 11 for the case when the limit point of γ f does not lie in the same component of ∂M ∩ Ω S,S and f lies in Ω S,S , in which case ε S = −1.
Let (ζ, ζ ) be the arc joining ζ and ζ (which could be a single point in case ζ = ζ ) in ∂M in the hyperbolic case. Let h dim be the dimer configuration corresponding to the pair (T , T † ) with reference flow given by Definition 4.2. We can now state the following final result relating the dimer height function to the winding of trees in the general setting we consider: Theorem 4.10. In the hyperbolic case, let γ := γ f,f be the curve formed by concatenating γ f , (ζ, ζ ) and the path in T joining ζ to m(f ). Orient γ from m(f ) to m(f ). We have the following deterministic relation: Here the sum σ∈Ω S,S ε σ + δ σ is as in (4.3) but also includes S and S .
Remark 4.11. The same statement holds in the case of the torus, but with a few obvious modifications, where we think of ζ and ζ being the same point at infinity (so starting from f and f both paths go to infinity in the same direction). Then (4.4) is also true but obviously without the winding term between ζ and ζ .
Proof. Consider the winding field h corresponding to (T , T † ) defined as in Lemma 4.6. See also the discussion following that lemma for the definition when T is a forest. Notice that the winding field differences between the diagonal midpoints are exactly the dimer height differences between the corresponding faces, hence it is enough to prove the lemma for the winding field.
Let us consider the hyperbolic case first. Take two vertices x, x in S, S respectively. Observe that we can find a path joining x, x which goes along the primal components, and moves from one component to the other by "jumping" over dual spines (i.e. going along edges whose dual belong to a dual spine). Also one can make sure that the path is minimal, in the sense that every component between those of x, x (i.e. components which separate S, S ) is visited at most once. Call this path (x, x ). Now letting x → ζ and x → ζ , we can also ensure that the portion of Ω S,S bounded by (x, x ), (ζ, ζ ) and the spines S, S contains none of f or f .
Letγ be the path obtained by concatenating (m(f ), x), (x, x ), (x , m(f )). From Lemma 4.6, it is clear that Indeed, notice that since the path is minimal, the ε S , δ S terms are defined so as to match with the definition of ε i , δ i in Lemma 4.6. Furthermore, W int (γ) = W int (γ) because of the above choice of x, x . We finish the proof using Lemma 2.7. For the torus case, the same proof applies by noticing that W ((x, x ), m(f ))+W ((x, x ), m(f )) converges to 0 as x and x go to infinity along S, S in the same asymptotic direction.
Finally, the fact that the term π S (ε S + δ S ) converges follows simply from the fact that the number of noncontractible components is a.s. finite in the limiting CRSF (and hence so is the number of spines which separate S, S into different components).

Convergence of Instanton component
We now prove a lemma about the instanton component of the height one-form which shows that it is a function of only the noncontractible loops and nothing else. Although this may be intuitively clear, it is not so easy to see using the connection between height difference and winding established in Theorem 4.10. Indeed, consider a noncontractible (continuous) closed curve λ in the manifold. To compute the height gap we first lift λ to the cover p −1 (λ) and we then compute the height gap between two copies of the starting point of the loop. A priori this might depend on the winding of the spines corresponding to the start and end points of p −1 (λ). However, it is clear in the case of the torus that the winding terms cancel out as the spines of the starting and ending points are translates of each other. It turns out the same is true in the hyperbolic case but a slightly more involved argument is needed, as the map from one fundamental domain to another is not just a translation but a Möbius map which perturbs the winding. Thus, the proof of the following lemma requires this input from Riemannian geometry and the assumptions on the graph in Section 2.3, so we encourage the reader to read the introduction of appendix A and Section 2.3 before reading this proof.
Take an ordered finite set of continuous simple loops which forms the basis of the first homology group of M , all endowed with a fixed orientation. Let H be the finite set of numbers which denote the height change along these loops. It is well-known (see Theorem 2.5) that h is completely determined by H (at the discrete level). In the following lemma, we write a superscript δ to account for the dependance in δ as in Section 2.3. Furthermore, H #δ also converges and the limit is measurable with respect to the limit of (C #δ , (C † ) #δ ) as δ → 0.
Proof. Using the main result of the companion paper [5, Theorem 1.2], we know that the Temperleyan forest converges in Schramm topology, in particular, the skeleton converges. Therefore the measurability part is clear once we show that h #δ is a function of (C #δ , (C † ) #δ ) only.
Take a set of loops which forms the basis of the homology group of M . It is enough to prove that the height change along one such loop λ is a function of (C, C † ) only.
Let us first consider the hyperbolic case. Lift the loop to the universal cover and compute the winding between its endpoints using Theorem 4.10 (see Figure 11). Note that it is enough to show that the winding term in Theorem 4.10 depends only on (C, C † ). To this end, we borrow the notations from Theorem 4.10. We will actually show that the only contribution of the winding term in Theorem 4.10 comes from the winding of the arc (ζ, ζ ). Clearly, the winding of the interval (ζ, ζ ) is determined by the spines and hence (C, C † ). Using the isomorphism between the Fuchsian group corresponding to M and the fundamental group of M , (see appendix A), we see that there exists a Möbius transform Φ which maps γ f to γ f and which depends only on the loop λ. In the hyperbolic case, observe that a Möbius map extends slightly beyond the unit disc and so it makes sense to talk about Φ (ζ). Using Lemma 2.8, we see that the difference in winding between γ f and γ f is arg Φ (M ) (Φ (ζ)) − arg Φ (M ) (Φ (m(f ))) which is completely determined by (C, C † ).
In the toroidal case, as mentioned before, the mapping Φ is simply a translation (i.e. Φ is a constant), so the lemma is immediate using the same argument as above. This completes the proof.
The next lemma, combined with Lemma 4.12, shows in a precise sense that the instanton component converges in law (see also Remark 6.15).
Lemma 4.13. In the setting of Section 2.3, suppose x, y ∈M and that f, f are two faces closest to x, y (breaking ties arbitrarily). Then σ∈Ω S,S (ε σ + δ σ ) converges in law. In fact, this convergence in law holds jointly for an arbitrary number of pairs of faces approximating given points inM .
Proof. By [5, Theorem 1.2], the Temperleyan forest T #δ converges to a continuum limit T (under the assumptions of this theorem). Consequently the nontrivial primal and dual cycles (C, C † ) also converge (with respect to the Hausdorff metric) to the nontrivial cycles induced by T . Moreover, since the number K of nontrivial cycles has superexponential tail as shown in the companion paper [5], no two cycles are likely to be close to one another by Wilson's algorithm and the crossing assumption. It is easy to deduce that S (ε S + δ S ) converges and the limit is measurable with respect to the limit of (C #δ , (C † ) #δ ). The joint convergence for an arbitrary finite number of pair of faces also follows in the same manner.

Local coupling
In this section we prove a local discrete coupling result which extends ideas of [4] to the setup of Riemann surfaces. Roughly speaking, the goal of such a result is to show that the local geometry in a small neighbourhood of a Temperleyan CRSF is given by that of a uniform spanning tree in the surface or, alternatively (and more usefully), in some reference planar domain. Moreover, locally around a finite number of given points, the configurations can be coupled to independent such USTs.
Recall that the actual Temperleyan CRSF's cannot be completely sampled using the standard Wilson's algorithm. However, due to Theorem 3.6, given the skeleton s emanating from either side of the punctures, the rest of the Temperleyan CRSF can be sampled from Wilson's algorithm. The main result in this section is Proposition 5.9.
The argument follows the same line of arguments as in [4] so let us first recall the strategy there. This consists of two main steps. Consider k points v 1 , . . . , v k in (Γ ) #δ . In the first step, we choose cutsets at a small but macroscopic distance around each of the k points, such that the cutsets separate the points from each other and from the rest of the graph. We reveal all the branches emanating from these cutsets. This leaves k unexplored subgraphs Γ #δ i , one around each point. In this step the key point is to make sure that the Γ #δ i are macroscopic (e.g., contain a ball of a radius roughly of the same order of magnitude as the distance to the cutsets). Clearly, the conditional law of the tree in each Γ #δ i is that of a wired UST. Moreover, these wired USTs are also independent conditionally given the cutset exploration. (Of course, unconditionally there is still some dependence). The second step is then to say that in each Γ #δ i , one can couple a wired UST with a full plane one, which shows, among other things, that the unconditional distribution is close to being independent.
To adapt this strategy to Riemann surfaces, if the Γ #δ i are sufficiently small so that it is simply connected, then the conditional law of the CRSF in each Γ #δ i will also be that of a wired UST, so that the second step can be used directly. For the first step however (cutset exploration), we will need to redo parts of the proof to take into account the possible loops in the CRSF.
Lemma 5.2. Let K 0 ⊂ K 0 ⊂ M be open connected sets and let K ⊂ K 0 , K be compact sets which are the closures of K 0 , K 0 . Also assume K contains a loop which is noncontractible in M . Then there exists a δ K,K > 0 and α K,K > 0 depending only on K, K such that for all δ < δ K,K and all v ∈ Γ #δ such that v ∈ K, simple random walk started from v has probability at least α K,K of forming a noncontractible loop before exiting K .

Cutset exploration
We now describe the construction more precisely. (We point out that in our setup it is natural to define these quantities with respect to Euclidean geometry onM and project to M , because our assumption and in particular the uniform crossing condition are stated with respect to this geometry; whereas the intrinsic geometry of M does not really appear in our setup). We denote by B euc (ṽ, r) the Euclidean ball of radius r aroundṽ and for r > r, let A euc (ṽ, r, r ) = B euc (ṽ, r ) \ B euc (ṽ, r) be the Euclidean annulus of inradius r and outer radius r . Let H i be a set of vertices in p(A euc (ṽ i , r/2, r)) ∩ (Γ ) #δ which disconnects v i from ∂N vi . The cutset exploration is done in two steps a. First sample the skeleton according to the correct law (i.e. according to CRSF branches with the global conditioning that the branches divide the manifold into annuli).
b. Then, reveal all the branches from H i , 1 ≤ i ≤ k using (unconditioned) Wilson's algorithm, resulting in a subgraph T #δ Hi , 1 ≤ i ≤ k. We say a vertex v j has cutset isolation radius 6 −k r at scale r if p(B euc (ṽ j , 6 −k r)) does not contain any vertex from T #δ Hi . ṽ X τ nc Let us define J vj as the the minimum value such that v j has isolation radius 6 −Jv j r and let J = max j J vj . We want to show that J has exponential tail (which results in polynomial tail of the isolation radius). To this end, we rely on the following bound on the distance between a loop-erased walk and a point, which is a version of Proposition 4.11 in [4]. Lemma 5.3. Let γ be a loop-erased random walk starting from a vertex v until it hits (∂Γ ) #δ or a noncontractible loop is created. Let u = v be another vertex and letũ be one of the images under p −1 of u. Let r be small enough such that U :=B euc (ũ, r), is contained in the pre-image of N u which containsũ. Then there exist constants α, c > 0 (depending only on the initial assumptions of the graph) such that for all δ < δ U and all n ∈ (0, log 2 ( Crδ0 δ ) − 1) where δ U , δ 0 are as in the crossing condition (assumption iv), P(γ ∩ p(B euc (ũ, 2 −n r)) = ∅) < cα n .
Proof. This is almost identical to Lemma 4.11 in [4] except we now need to take into account the topology as well. To emphasise the differences with Lemma 4.11 in [4], we recall the strategy there. The idea is that if the loop erased walk on the manifold comes inside p(B euc (ũ, 2 −n r)) then some lift of the random walk must necessarily come within B euc (ũ, 2 −n r) -call this region (exponential) scale n. Furthermore, after the last such visit, this (lift of the) random walk must cross n many annuli without making a loop aroundũ (which we called a "full turn"). This has a probability bounded by e −cn .
In the current situation however, the random walk might form a noncontractible loop before exiting U , and therefore its relevant lift described above does not necessarily have to cross n many annuli (see e.g. Figure 12) after coming within Euclidean distance 2 −n r ofũ before we stop it. Thus we cannot simply apply the argument of the previous paper. We overcome this using the following idea which we first explain informally. The random walk on the manifold may visit p(B euc (ũ, r)) multiple times. Consider one such time when it enters p(B euc (ũ, r)) (not necessarily the first time) and suppose its relevant liftX which also enters B euc (ũ, r) at that time. Let X be the portion of the loop-erasure of X on M at that time such that, if the random walk hits X later on, then a noncontractible loop in M is formed. Assume that X comes inside p(B euc (ũ, 2 −m r)). In order for the whole loop-erasure γ to come inside p(B euc (ũ, 2 −n r)) the following events must take place on one such visit. First, the liftX has to cross the first m scales without performing a full turn (otherwise it would close a noncontractible loop on the manifold and stop); this has a probability bounded by e −cm . ThenX needs to come within Euclidean distance 2 −n r ofũ, but we bound crudely the probability of this event by 1. Finallỹ X needs to come back out after the last visit to B euc (ũ, 2 −n r) in such a way that no full turn occurs between distances 2 −n r and 2 −m r (for the same reason as in the simply connected case). The probability of this event can be bounded by e −c(n−m) . The intersection of these three events gives the right overall upper bound: e −cm e −c(n−m) = e −cn .
Let us fill in the details. Let C i be the Euclidean circle of radius r i := 2 −i r aroundũ for i ≥ 0. We start a random walk X from v (to emphasise again: X is on the manifold M ). Let {τ k } be a sequence of stopping times defined inductively as follows. Set τ 0 = 0. Having defined τ k to be a time when the random walk X crosses or hits some circle p(C i(k) ), define τ k+1 to be the smallest time after τ k when the random walk crosses or hits either p(C i(k)−1 ) or p(C i(k)+1 ); this defines τ k by induction for every k ≥ 1. If i(k) = 0, define τ k+1 to be the smallest time after τ k when the random walk crosses or hits p(C i(k)+1 ). Let k exit be the smallest integer such that in the interval [τ kexit , τ kexit+1 ] either a noncontractible loop is created by the loop-erasure of the random walk or the random walk hits the boundary ∂Γ #δ . Let κ 0 , κ 1 , . . . , κ N be the set of indices k < k exit +1 when X τ k has crossed or hit p(C 0 ). Also note that the portions X[τ κ0+1 , τ κ1 ], X[τ κ1+1 , τ κ2 ], . . . are contained inside p(C 0 ); while the portions X[τ κ0 , τ κ0+1 ], X[τ κ1 , τ κ1+1 ], . . . are contained outside p(C 1 ). Our first claim is that N has exponential tail, i.e., ∃α 1 ∈ (0, 1), such that for all n ≥ 1, P(N > n) < α n 1 .

(5.2)
Indeed, using Lemma 5.2, every time the walk hits p(C 0 ), there is a positive probability independent of the past that the walk creates a noncontractible loop before returning to p(C 1 ). Therefore, N has geometric tail. This proves (5.2). Let S be the set of points {X(τ k )} k≥1 . Note that given S, the pieces of random walk {X[τ k , τ k+1 ]} are independent of each other. We call such pieces elementary pieces of the walk. If i(k) = 0, X[τ k , τ k+1 ] is (given S) a random walk starting from X τ k conditioned to exit the annulus p(A euc (ũ, r i(k)−1 , r i(k)+1 )) at X τ k+1 . If i(k) = 0, X[τ k , τ k+1 ] is a random walk which is conditioned to next hit p(B euc (ũ, r/2)) at X τ k+1 . (Note that this property is lost if we solely condition on {X(τ k )} k≤κ N since the conditioning on N is complicated.) By Lemma 4.7 in [4], conditionally on S, each random walk piece X[τ k , τ k+1 ] has a uniformly positive probability to do a full turn in the annulus p(A euc (ũ, r i(k)−1 , r i(k) )) for δ ≤ δ U and given range of n (where δ U comes from the uniform crossing assumption). Here we define a full turn to be the event that the random walk crosses every curve joining the inner and outer boundary in the specified annulus. Indeed, although Lemma 4.7 in [4] gives the uniform positive probability estimate for crossing in the Euclidean plane, the estimate is valid here by considering the relevant lift of X[τ k , τ k+1 ] which is insideB euc (ũ, r) = U and applying the uniform crossing estimate in the universal cover. Let us also point out that δ is chosen small enough so that the uniform crossing estimate is valid inside U .
We now define certain low probability events E 1 , . . . , E n such that one of them must take place if the loop-erasure of X is to enter B n := p(B euc (ũ, 2 −n r)). Let us condition on S. We say the event E j occurs if: such that X τ k crosses (or hits) p(C n ).
(ii) Let X j be the portion of the loop-erasure of X[0, τ κj−1+1 ] such that if the walk X[τ κj−1+1 , τ κj ] intersects X j , a noncontractible loop would be created. Let m be the maximum index such that X intersects the circle p(C m ). Then for the event E j to hold we also require, in addition to the previous point, that: for any κ j−1 + 1 ≤ k ≤ λ j if i(k) < m, then the walk X[τ k , τ k+1 ] does not perform a full turn in p(A euc (ũ, r i(k)−1 , r i(k) )). If no such m exists (e.g., if X j is empty), we do not require anything further.
(iii) Let j be the last k before κ j such that X τ k crosses (or hits) p(C n ). Let j be the first time after j that the walk intersects p(C m+1 ). Then in addition to the previous two points, for E j to hold we require that the walk X[τ k , τ k+1 ] does not perform a full turn for any From the discussion above, it is clear that we have the following lemma.
Lemma 5.4. We have Thus all we need to show is that the event on the right hand side of Lemma 5.4 has exponential tail bound. Now we claim that P(E j |S) ≤ e −c n and so P(E j ) ≤ e −c n . This is justified using the uniform positive probability of the walk X performing a full turn, even conditionally given S. Indeed, conditioned on S, we have that the event (ii) in the definition of E j has probability at most e −c m , and conditioned on event in (ii), the event in (iii) has probability at most e −c (n−m) since conditionally given S, the random walk portions X[τ k , τ k+1 ] are independent. Thus the overall probability given S is at most e −cn with c = c ∧c . All in all, using eq. (5.2), Lemma 5.4 and a union bound, we obtain thereby concluding the proof.
Lemma 5.5. Let s denote the skeleton. Let u be a vertex which is not one of the punctures and letũ be one of the images under p −1 of u. Let s be small enough such that U :=B euc (ũ, s), is contained in the pre-image of N u which containsũ. Then there exist constants α, c > 0 (depending only on the initial assumptions of the graph) such that for all δ < δ U and all n ∈ (0, log 2 ( Csδ0 δ ) − 1) where δ U , δ 0 are as in the crossing condition (assumption iv), P(s ∩ p(B euc (ũ, 2 −n s)) = ∅) < cα n .
Proof. This follows from [5, Proposition 6.5] and Lemma 5.3. Let us explain why this Proposition implies the desired statement in the case k = 1. For this we need to introduce a few notations used in [5]. Letλ M denote the law s. Let η r = (η 1 r , η 2 r ) denote the portion of the branches of the skeleton until they exit a ball of radius r around the puncture, for a very small choice of r. For each such path η i r (with i = 1, 2) run independent simple random walks from their tips, where each walk is conditioned to exit the ball of radius r > r before hitting η i r , for another small choice of r . After hitting r the walks are allowed to run unconditionally until they either hit η r or (∂Γ ) #δ or their loop erasure creates a noncontractible loop. Letμ r,r denote the joint law of the loop erasures of the resulting pair of walks. For any event E which is measurable with respect to s \ η r , [5, Proposition 6.4] asserts that where the supremum is over all possible η r/2 (resp. η r/4 ) It is easy to see using the Markov property of the independent random walks used to samplẽ µ r/2,r that after hitting r, the result of Lemma 5.3 kicks in, yielding the required exponential bound forμ r/2,r (E). The same argument obviously applies toμ r/4,r/2 (E), yielding the desired statement.
The application of [5, Proposition 6.5] for k > 1 is analogous. We skip the details here.
We now state the result showing an exponential tail of J, which is a combination of Lemmas 5.3 and 5.5 and Schramm's lemma, and is identical to the proof of Lemma 4.20 in [4] (see also the proof of Theorem 4.21 in [4]); as it is identical we skip the proof here.
Lemma 5.6. There exist constants c, c > 0 such that the following holds. LetD be a compact set containing B(ṽ i , r) for 1 ≤ i ≤ k where r is as in (5.1). Then for all δ ∈ (0, δD) and for all m ∈ (0, log 6 (δ 0 r/δ) − 1), Finally, we state a lemma which says that with exponentially high probability, a branch of the CRSF, after entering an exponential scale t does not backtrack to a smaller scale. Such a lemma for SLE curves can be found in [33] (see also Lemma 3.4 in [4] in the simply connected case).
Lemma 5.7. Fix u, U, δ U as in Lemma 5.3. Suppose γ is the loop-erasure of a simple random walk in Γ #δ started from vertex u until it hits the boundary or creates a noncontractible loop. Supposeγ is the lift of γ started fromũ (parametrised fromũ to its endpoint). There exist constants c, c > 0 (depending only on the initial assumptions of the graph) such that for all δ < δ U and all n ∈ (0, log 2 ( Crδ0 δ ) − 1) Proof. Let E be the event thatγ enters B euc (ũ, r2 −n ) after exiting B euc (ũ, r2 −n/2 ). Let r be as in Lemma 5.3 and assume n is even without loss of generality. The argument for this is very similar to Lemma 5.3 and in fact simpler, so we content ourselves with a sketch. Let r, C i , B i be as in the proof of Lemma 5.3. We look at the liftX of the simple random walk X started at u and we stop when either X hits the boundary or a noncontractible loop is created. Let τ k be the set of stopping times defined as in the proof of Lemma 5.3 but for the liftX instead of X and let S be the set {X(τ k )} k≥1 . Observe that lift of the loop erasure of X is the loop erasure ofX since erasing contractible loops commutes with lifting to the universal cover. If the loop erasure ofX has to backtrack to scale n, the random walk has to necessarily enter scale n after leaving scale n/2. Let J be the largest j such thatX τj crosses or hits C n after leaving C n/2 and let I be the largest index i smaller than J whenX τi crosses or hits C n/2 . Conditioned on S, if E occurs, thenX enters B n at least once after leaving B n/2 , and none of the elementary pieces of the walk between τ I and τ J can perform a full turn. But again, conditioned on S there is a uniformly positive probability to do a full turn for each elementary piece. Since there are at least n/4 such elementary pieces contained in [τ I , τ J ], we conclude the proof of the lemma applying the upper bound on the full-turn estimate on each elementary piece.

Full coupling
The results of the previous section covered the first step in the proof of the coupling. As we mentioned above, the second step is identical to the simply connected case so we only recall the main statements. We first recall the result which we will need from [4]. We use the notations and assumptions already in force. LetD ⊂D ⊂M be two simply connected compact domains and fixṽ ∈D #δ . Let TD , TD denote the wired UST respectively in (D ) #δ andD #δ . Let rṽ denote R(ṽ,Ñṽ) as in eq. (5.1) (the largest Euclidean radius so that p is injective).
Lemma 5.8 (Theorem 4.21, [4]). There exists c, c > 0 such that the following holds. Fix v,D,D as above. There exists a coupling between TD and TD and a random variable R > 0 such that Furthermore, for all δ < δD and for all m ∈ (0, log 6 (δ 0 rṽ/δ) − 1), if we write R = 6 −Kv rṽ ,D where rṽ ,D is the minimum of rṽ and the Euclidean distance betweenṽ and ∂D, then We now put together the cutset exploration with the above one-point coupling from the simply connected case as follows. Recall the notations from Section 5.1. Let T #δ Hi be the branches revealed in the cutset exploration around each v i . Let T #δ containing v i . Observe that given T #δ H , the law of T #δ ∩ Γ #δ i are independent wired UST in Γ #δ i (with the natural boundary), by the generalised Wilson algorithm. Applying the coupling of Lemma 5.8 in the lift of Γ #δ i toM and some fixed compactD inM containing some liftsṽ 1 , . . . ,ṽ k of points v 1 , . . . , v k , we obtain a coupling of the oriented CRSF in Γ #δ and independent wired USTs (TD i ) 1≤i≤k inD. Furthermore, the r in (5.1) and rṽ in Lemma 5.8 differs by a constant multiplicative factor which depends only on the choice of lifts of the vertices. Furthermore, we may (and will) choose the (TD i ) 1≤i≤k to be independent of T #δ H . In fact, note that for each 1 ≤ i ≤ k, TD i may be chosen independent of the restriction of T #δ to ∪ j =i p(B(ṽ j , r/2)).
Proposition 5.9. The above coupling has the following properties: there exists random variables Furthermore if we write R i = 6 −Iv i r i where r i is the minimum of r as in (5.1) and the distance betweenṽ i and ∂D, then for all δ ≤ δD, and for all 1 ≤ i ≤ k, for all n ∈ (0, log 6 (δ 0 r i /δ) − 1), for some constants c, c > 0 (depending only on the initial assumptions on the graph). In particular, P(I vi ≥ n) ≤ ce −c n ∨ δ c . The set of noncontractible loops of T #δ is measurable with respect to T #δ H . In particular, (TD i ) 1≤i≤k are also independent of the set of oriented noncontractible loops in the oriented CRSF.
Observe that when I vi is very big or r i is very small, it is possible that the ball B(v i , R i ) is reduced to a point so the statement (5.4) is trivial (that is, (5.4) holds for a single vertex). This happens with a probability which is at most δ c for some c .
Proof. This follows immediately from Lemma 5.6 and Lemma 5.8 once we observe that I vi is within O(1) of the sum of J vi in Lemma 5.6 and K vi in Lemma 5.8, both of which have exponential tails.
The proof of the final assertion is a topological fact. Indeed, conditioned on T #δ H , it is a deterministic fact that none of the oriented noncontractible loops pass through the portion of the CRSF in Γ #δ i (which is a wired UST). Indeed, otherwise we would have a path joining two points of the wired boundary of a wired UST which is impossible. Thus, the wired UST in each Γ #δ i is conditionally independent of the noncontractible loops and hence so are TD i s.

Convergence of height function and forms
In this section, we precisely state our main result (this is Theorem 6.1) and then prove it. Recalling the sketch from Section 1.3, we see that what remains to be done is to go from the convergence of the Temperleyan CRSF to the convergence of its winding field using the coupling of Section 5. The global idea is the same as in [4] but some of the estimates on winding of LERW have to be redone. The primary difficulty is that spines are made of multiple copies of LERW, and hence Wilson's algorithm cannot be used to estimate the winding of all the copies. These are dealt with respectively in Sections 6.2 and 6.3. The proof of the main result is finally concluded in Section 6.4.

Precise statement of the result
We first provide a precise statement of what we prove. From now on, we work with the manifold M which is obtained by removing 2g + b − 2 points from its interior. Recall that the white vertices removed from G #δ to obtain a Temperleyan graph in Theorem 3.6 will converge to these points. Denote byM the universal cover of M ; in the sequel "lift" refers to lift onM . Let h #δ be the height function defined as in Section 4 sampled using eq. (1.1). Recall that h #δ is a function from the the dual of (G ) #δ to R, that is defined up to a global additive constant and that a choice of reference flow is used in the definition that depends only on the embedding. Recall the probability measures P Wils and P Temp from Section 3.3.
First, we extend h #δ to a function h #δ ext :M → R by defining h #δ ext (x) to be the value of h #δ on the face containing x. Recall the graphs (Γ #δ ) , (Γ †,#δ ) from Section 2.3. Also recall that X denotes X − E(X). Let T #δ be the Temperleyan forest sampled using P Temp . converges jointly in law as δ → 0. The first coordinate also converges in the sense of all moments. Furthermore, the limit of the first coordinate is measurable with respect to the limit T of T #δ , is universal (in the sense that it does not depend on the graph sequence (G ) #δ ), and is conformally invariant.
To clarify, convergence in the sense of all moments means that for all i, denoting X #δ := f (x)h #δ ext (x) dµ(x), E(|X #δ | i ) converges as δ → 0. Notice also that since f dµ = 0, the fact thath #δ ext (x) is defined only up to a global additive constant is irrelevant. Let us explain briefly what is meant by conformal invariance in our setting. Suppose |χ| = k and (M, x 1 , x 2 , . . . , x k ) and (N, y 1 , . . . , y k ) are conformally equivalent in the sense that there is a map φ : M → N so that φ(x i ) = y i for all 1 ≤ i ≤ k and φ is a conformal bijection between the two Riemann surfaces. Let h M (resp. h N ) denote the limit field obtained in the setup of Theorem 6.1 with the punctures being x i and y i for 1 ≤ i ≤ k. Then h N = h M • φ in the sense that for any compactly supported test function f on the universal cover as in Theorem 6.1, Remark 6.2. Before we start giving the details, we point out that the theorem follows if we can establish the convergence of the Temperleyan CRSF (in the Schramm sense) to a limit in which all curves are a.s. simple, together with a few additional ingredients. In addition to the assumptions on the discretisation of the surface described in Section 2.3 and in particular the crossing assumption, these are the following: for each fixed vertex v ∈ (Γ ) #δ , let γ v the branch of the Temperleyan forest starting from v. Then we require that γ v satisfies Lemma 5.3 (the branch does not come close to a given point) as well as the moments on winding of Lemma 6.8 and Lemma 6.9. Note that the convergence to Brownian motion in item iii is not really needed here; instead, convergence of random walk to a process which is absolutely continuous with respect to Brownian motion would be sufficient. Stated this way, we note that Theorem 6.1 is novel even in the simply connected case. For instance, this result is applied in [3] and in the forthcoming [6].

Some a priori tail estimates on winding
The goal of this section is to obtain tail estimates on the winding of a branch of a Temperleyan forest. This will be achieved in Lemma 6.6 which is the main result of this section. Conceptually the arguments are similar to Section 4 of [4]. However, as in Section 5, there are additional difficulties linked with the fact that Wilson's algorithm can stop because of a noncontractible loop being formed. Furthermore we ultimately want to consider the winding of entire spines, yet only one copy of the lift of a loop is directly connected to a loop erased random walk path. We first treat the part directly obtained by loop-erasure in this section and defer estimates on the copies for the next. Throughout this section, we deal with a CRSF sampled from P Wils (cf. Section 3.3). First we sample the skeleton (cf. Theorem 3.6) which we denote by s. Lets be the lift of s toM . Recall that s decomposes M into finite number of disjoint annuli and conditionally on s, in each of them we can sample the rest of the Temperleyan CRSF using Wilson's algorithm.
We now work conditionally on s. Let γ v be the branch of the CRSF started from v (which can thus be sampled using Wilson's algorithm). For any vertexṽ ∈ p −1 (v), letγṽ be the lift of γ v starting fromṽ and up until the time when γ v closes a noncontractible loop or hits the boundary ∂Γ #δ ∪ s. Let τ nc be the stopping time when the simple random walk generating γṽ creates a noncontractible loop or hits the boundary. Recalling the definition of spines from Section 4, note thatγṽ will include a part of a spine as soon as the random walk does not hit ∂Γ #δ ∪ s when it stops (see Figure 12).
Let us orientγṽ starting fromṽ and going away from it in some continuous manner in [0, 1]. For t > 0, ifγṽ exits B(ṽ, e −t−1 ), let t 1 be the first timeγṽ exits B(ṽ, e −t−1 ) and let t 2 be the last time it exits B(ṽ, e −t ). In this case, ifγṽ ends in B(ṽ, e −t ) we set t 2 = 1.
We first state a simple deterministic lemma connecting the winding of curves avoiding each other. Lemma 6.3. Consider the annulus A = {z ∈ C : r < |z − x| < R} where r > 0 and R ∈ (r, ∞] and let γ 0 be a simple curve in A connecting the outer and inner boundaries of A, assumed to be parametrised in [0, 1]. Then for any simple curve γ in A \ γ 0 , we have For R = ∞, the annulus is interpreted as the complement of a ball.
Note that in the above lemma we do not need any assumption on the regularity of the curves beyond continuity.
Proof. Join each endpoint of γ to its nearest points in γ 0 using a geodesic in A (with its intrinsic metric inherited from the Euclidean plane). Assume these points are γ 0 [s 1 ] and γ 0 [s 2 ]. Consider the loop L obtained by γ, these two geodesics and γ 0 [s 1 , s 2 ]. Note that L might be non-simple, but any loop which is created by the union of γ and the two geodesics lie in the annulus A slitted by γ 0 which is simply connected and do not contain x. Hence all these loops contribute winding 0 around x. Furthermore, it is easy to see that erasing all these loops (say in a chronological manner started from one of the tips of γ) creates a simple loop. Since the winding of L after erasing these loops is 0 or ±2π, the total winding of L around x is also either 0 or ±2π. It is easy to observe that the geodesic joining any two points in the annulus lie in the smaller half annulus containing the two points, and hence has winding at most π. Using this observation, the proof of the lemma is complete.
Assume without loss of generality and to simplify notations that t is an integer. Let rṽ be such that B(ṽ, 10rṽ) does not intersects and p is injective in B(ṽ, rṽ) (this is a slight modification of the previous definition of rṽ). Let e −t0 = rṽ. We consider concentric circles C j of radius {e −t0−j } 0≤j≤t+2 (with B j the ball inside it) aroundṽ. Take a random walk X in (Γ ) #δ starting fromṽ stopped when it hits (∂Γ ) #δ ∪ s. In case (∂Γ ) #δ ∪ s = ∅ (i.e., in the case of the torus), the random walk continues forever. LetX be the lift of this walk starting fromṽ. Now let {τ k } k≥0 be the set of stopping times as described in the proof of Lemma 5.3 for the random walkX. That is, if we hit or cross the circle C i(k) at time τ k , we wait until we hit or cross C i(k)±1 at time τ k+1 . If i(k) = 0 (or t + 2), we wait until we hit or cross C 1 (or C t+1 ). Let B i denote the disc inside C i .
Let τ i1 , τ i2 , . . . be the successive times in the sequence (τ k ) k≥1 when the walk hits C 0 . Note that the interval [τ ij , τ ij +1 ) is spent completely outside B 1 . We now claim that the random walk cannot wind too much outside B 1 . Lemma 6.4. There exist c, c > 0 such that for all δ < δ p(B0) , n ≥ 1, j ≥ 1 and u ∈ B 1 such that P(X τi j +1 = u) > 0, Here the supremum is over all continuous paths obtained by erasing portions ofX[τ ij , τ ij +1 ].
Proof. The proof of Lemma 6.4 is very similar to Lemma 4.8 in [4]. But it needs an input from Riemannian geometry to control the winding of the spines near the boundary of the universal cover. We postpone this proof to appendix A. Lemma 6.4 takes care of the winding of the excursions outside B 1 . For excursions inside, most of the technical work was done in [4]. LetỸ j be the loop-erasure ofX[0, j]. For any j, we parametriseỸ j in some continuous way away fromṽ. Fix m ≥ 1 and let t 1 be the first timẽ Y τi m exits B(ṽ, e −t−1 ) and t 2 be the last timeỸ τi m exits B(ṽ, e −t ).
Lemma 6.5. There exist C, c > 0 such that for all δ < δ p(B0) , t ∈ (t 0 , log(Cδ 0 /δ)), and n ≥ 1, m ≥ 1, We emphasise that in the above lemma, m is non-random as it is going to be important in what follows.
Proof. The proof of this lemma is identical to that of [4], Lemma 4.13 and the only difference is in the setup, which we point out. In [4], we were working with a simply connected domain and we waited until the random walk exited it. The argument proceeds by conditioning on the positions S := {X τ k } k≥1 and arguing that in each interval of walk between successive points in S, the winding of any continuous subpath has exponential tail. Then we proved that we only need to look at a random number of intervals which itself has exponential tail. In the present setup, we can condition on S so that τ im < ∞ is satisfied. The exponential tail of winding inside any inner annulus follows from Lemma 4.7 in [4] (this is exactly the same as the simply connected case) and for the outer annulus, we use Lemma 6.4. The number of relevant intervals to consider also has exponential tail following verbatim the proof of Lemma 4.15 in [4].
With these lemmas we can now state and prove the main result of this section, which controls the (topological) winding in an annulus. Lemma 6.6. There exist constants C, C , c > 0 so that for all n ≥ 1, for allṽ, for all δ < δB (ṽ,rṽ) and for all 0 < t < log(C rṽδ 0 /δ), where the maximum is taken over all continuous segments fromγṽ[t 1 , t 2 ]. Also P max A few remarks are in order. Firstly, the stretched exponential tail is an artefact of the proof and we believe that an exponential tail bound can be proved with more care. Secondly, note that the intersection in the argument of the max above is used so that we don't need to look into what happens very close toṽ. This is a technicality which simplifies the proof, but later it is not going to matter. In the end, we can decompose the whole pathγṽ into a disjoint union over t of [t 1 , t 2 ] where t 1 is last exit of B(ṽ, e −t−1 ) and t 2 is the last exit of B(ṽ, e −t ) which will accomplish the desired moment bound of truncated winding by exploiting the moment bound for each of these segments. We also emphasise that, although the noncontractible loop in the component containing v could be a branch from the skeleton,γṽ itself does not contain any portion of the skeleton. This is simply because we have sampled the skeleton first before defining γ v .
Proof of Lemma 6.6. The main idea is to use a union bound for the winding of each excursion between annuli. One difficulty with working withγṽ is that if we condition on where or when a noncontractible loop is formed then we break the independence between the pieces of random walk. Indeed, a walk conditioned on not forming a noncontractible loop will avoid certain portions of its previous trajectory, which a priori might bias it to have a lot of winding. Recall that at the point when a non-contractible loop is formed, the walk could be very close to the starting point. Thus the idea is to run the random walk until the boundary (∂Γ ) #δ is hit (which could be potentially longer than dictated by Wilson's algorithm) and control the winding of its loop erasure uniformly over all scales using Lemmas 6.4 and 6.5. After this, we need to separately compare the winding made by the last excursion into the ball to the previous estimate.
Recall that t 0 satisfies e −t0 =rṽ. Note that (6.2) follows easily from Lemma 6.4. Indeed, let us call τ nc the first time where a noncontractible loop is created and τ ∂ be the first hitting of the boundary. Fix N = max{j : τ ij < τ nc ∧ τ ∂ }. Recall that N has exponential tail since every time we hit C 0 , conditioned on what happened before, we have a positive probability to create a noncontractible loop or hit the boundary before hitting C 1 , by Lemma 5.2. Thus we can work on the event N ≤ n at a cost which is exponentially small in n. Now Lemma 6.4 entails that on each of the n pieces, the winding has uniform exponential tail. This completes the proof of (6.2). We emphasise here that to bound portions of the loop erasure outside radius e −1 , we can run the random walk until it hits the boundary (∂Γ ) #δ and not worry about the case when the walk stops deep inside the ball, as the loop erasure outside radius e −1 is already taken care of. Now we turn to (6.1). As before, let S := {X τ k } k≥1 . The idea is to compare the winding of γṽ to the winding ofỸ τi m because there the conditioning on S preserves the independence and Lemma 6.5 controls the winding. The main work will be to do a case-by-case analysis of what can be erased and created between τ im and the creation of a noncontractible loop.
Let N be as above and note that using the same idea of exponential tail of N and exponential tail of winding up to a fixed number of hits of the outermost circle (Lemma 6.5), we get P sup Let λ 1 and λ 2 be respectively the first timeỸ τi N exits B(ṽ, e −t−1 ) and the last timeỸ τi N exits B(ṽ, e −t ). To ease notations, from now on, all maximums in this proof come with the additional condition that the paths stay outside B(ṽ, e −t−1 ) without writing it explicitly. We will also assume without loss of generality that the parametrisations ofγṽ andỸ τi N are identical up to the first point their traces differ.
If a noncontractible cycle is created or the boundary is hit in the interval [τ i N , τ (i N +1) ] then, sinceX[τ i N , τ (i N +1) ] does not intersect B(ṽ, e −t0−1 ), only pieces ofỸ τi N [0, λ 2 ] can be erased and nothing can be added in the time interval [τ i N , τ (i N +1) ] to getγṽ[0, t 2 ]. In particular we see that in this case and we are done using (6.3).
The only other possibility is that a noncontractible cycle is created between τ (i N +1) and τ i (N +1) (otherwise that would contradict the maximality of N ). Since p is injective in B(ṽ, e −t0 ), this occurs if and only if the walk hits a copy inside B(ṽ, e −t0 ) of a portion of the walk that is further away fromṽ, as illustrated schematically in Figure 13. From now on we assume we are on this event.
First we make a topological observation. Let β be the first exit time of B(ṽ, e −t0 ) byỸ τi N . We claim thatX[τ i N +1 , τ nc ] cannot hitỸ τi N [0, β]. Indeed, suppose by contradiction that it does so at some time T ∈ [τ i N +1 , τ nc ]. Then, after erasure, we are left with a path completely contained in B(ṽ, e −t0 ) where p is injective. No noncontractible loop can then be created before τ i N +1 , which contradicts the maximality of N as explained above.
Let k nc be the index during which the noncontractible loop happens, i.e, the unique index k such that τ nc ∈ [τ k , τ k+1 ]. By the topological claim in the previous paragraph, no full turn can occur in any of the intervals [τ k , τ k+1 ] for i N ≤ k < k nc . However, by Corollary 4.5 in [4], a full turn can occur in any interval [τ k , τ k+1 ] independently with uniformly positive probability given S. Hence, there exist constants c, C depending only on the crossing estimate such that P(k nc − i N > n) ≤ P(k nc − i N > n, N ≤ n) + P(N > n) ≤ Ce −cn .
Combining the above with Lemma 4.7 in [4] (which bounds the winding of the random walk during an interval of the type [τ i , τ i+1 ]), we obtain the following stretched exponential tail bound: for some C, c > 0 depending only on the crossing estimate and where, as before, Y is any continuous portion obtained from X[τ i N , τ nc ] which preserves the order of the random walk path. In particular, this gives a good control of the winding of the piece ofγṽ added to Y τi N . However, we are still not done as illustrated in Figure 13: the times t 1 , t 2 (which were first entry, last exit times forγṽ) may be different from λ 1 , λ 2 (which were the first entry, last exit ofỸ τi N ), so we need additional arguments. Let be the last time inỸ τi N which is not erased. Note that with this notation we showed above that λ E ≥ β. This implies immediately that t 1 = λ 1 . Now we need to consider several cases depending on where λ E is in relation to λ 2 .
• If λ E ≥ λ 2 andγṽ(λ E , t 2 ] (the loop erasure after λ E ) does not enter B(ṽ, e −t ) theñ and there is nothing more to prove.
• If λ E ≥ λ 2 andγṽ(λ E , t 2 ] enters B(ṽ, e −t ) (this is the case pictured in Figure 13), then we decomposeγṽ as the union ofỸ τi N [λ 1 , . The winding of any continuous portion of the first part has exponential tail by (6.3). The winding of any continuous portion of the last part has stretched exponential tail since it is bounded by ∆ as in the first case. So we only need to take care of the middle part.
To that end, we decomposeỸ τi N [λ 2 , λ E ] into excursions inside and outside B 0 as follows: where β 0 = β and for k ≥ 1, β 2k is the first exit of B(ṽ, e −t0 ) after β 2k−1 and β 2k−1 is the first entry into B(ṽ, e −t0−1 ) after β 2k−2 . Note that any portion ofỸ τi N [β 2k , β 2k+1 ] is outside B(ṽ, e −t0−1 ) so its winding has exponential tail using Lemma 6.4. Note also that Y τi N [β 2k+1 , β 2k+2 ] lies inside the annulus bounded between C t and C 0 and never intersects the curveX[τ i N , τ nc ] by maximality of λ E . Let Y be the loop erasure of the portion ofX from τ i N until its first hit of C t . Hence using Lemma 6.3 Note that the right hand side has stretched exponential tail by (6.4). Thus each term in the decomposition (6.5) has stretched exponential tail and clearly the number of terms is bounded by N which itself has exponential tail. Combining these two, we can conclude that the absolute value of winding of any continuous portion ofỸ τi N [λ 2 , λ E ] has stretched exponential tail (with exponent 1/3), thereby completing the proof of this case.
This concludes the proof of Lemma 6.6.
6.3 From partial path to the full spine.
As shown in Theorem 4.10, the pathγṽ defined above is only a part of what is needed to compute the height function. Recall the notion of spine from Section 4.2 which is just a bi-infinite path in the cover. As observed before, when we follow the outgoing edges in the Temperleyan forest from any givenṽ, we always end up in a unique spine (since boundary loops have been added around every hole).
Recall from Theorem 4.10 that we are interested in the winding of the path starting fromṽ, and then moving along the spine to infinity or the boundary of the disc. Observe that the initial portion of this path isγṽ ( Figure 12) and then it moves along copies of the noncontractible loop in the component ofṽ in the Temperleyan CRSF. In this section, we will call this path pṽ (i.e., γṽ followed by a semi-infinite piece of the spine). Observe that this notation is in contrast to the notation used in Theorem 4.10, where the path pṽ was called γ f , whereas here we emphasise thatγṽ denotes what is sampled from Wilson's algorithm given the skeleton. When we apply Wilson's algorithm to sampleγṽ, given the skeleton, we may form a new noncontractible loop, in which case we have discovered a portion of a spine (the unique spine attached toṽ). The winding of this portion is controlled by Lemma 6.6. However, we also need to control the winding of the rest of pṽ. The purpose of this section is precisely to achieve this (done in Lemma 6.7). We also need to control the winding of the other semi-infinite piece of the spine attached toṽ, which is done in Lemma 6.8. Finally, if we want to control the height gap betweenũ andṽ, we need to control the winding of pũ aroundṽ as well, which is done in Lemma 6.9.
Let us parametrise pṽ in [0, ∞) in any continuous manner away fromṽ. Let t 0 be as in Section 6.2, i.e., e −t0 =rṽ. For t ≥ t 0 , let t 1 be the first exit time of pṽ from B(ṽ, e −t−1 ) and let t 2 be the last exit time of the same from B(ṽ, e −t ). Note again that here, it is possible that the unique spine attached toṽ could be the same as the spine corresponding to a skeleton branch. Lemma 6.7. For all k ≥ 1, there exists a constant m > 0 so that for all δ < δB (ṽ,e −t 0 ) and for all t 0 ≤ t < log(C δ 0 /δ), where the supremum is taken over all continuous segments. Also Proof. In this proof, we write pṽ = p to lighten notation. Let us first consider the case when p does not contain a portion of the spine corresponding to the skeleton. We parametriseγṽ in [0, 1] as before and assume that the parametrisation ofγṽ and p is the same until they start to differ. We dropṽ and writeγṽ = γ throughout the rest of this proof for notational clarity. To differentiate the first and last exit times for these two parametrisations we write t 1 (p) (resp. t 2 (p)) for the first exit of B(ṽ, e −t−1 ) (resp. last exit of B(ṽ, e −t )) of p. We also write t 1 = t 1 (γ) and t 2 = t 2 (γ) for clarity.
Note that it is always the case that t 1 (γ) = t 1 (p) and that t 2 (γ) ≤ t 2 (p). On the event t 2 (p) = t 2 (γ), the maximum in this lemma is over the same set as the one in Lemma 6.6 so there is nothing to prove. We focus now on the case where t 2 (p) > t 2 (γ), meaning that a part of the spine comes back close to x on the coverM after the time that a noncontractible loop was created.
We start with the case of the torus. Let us denote by S the spine attached toṽ (meaning the full bi-infinite path). Since S is a periodic path, it has a well-defined direction d (say represented by a unit vector in R 2 ) and separates the plane in two sets right and left of S. Let us assume without loss of generality that the direction d is horizontal and thatṽ is below S. Let τ be a time at which the vertical coordinate of γ reaches its maximum (for topological reasons, this has to occur on the spine S) and let us define γ by appending to γ[0, τ ] a vertical segment going up to infinity. Now it is easy to see that the winding of a straight segment around a point is bounded by π, so max By Lemma 6.6 (theγṽ there is γ here) applied to all scales up to scale t, the right hand side is a sum of O(t) variables with uniform stretched exponential tails. Let us now bound P(t 2 (p) > t 2 (γ); p ∩ s = ∅). Note that on the event t 2 (p) > t 2 (γ),ṽ has to be at distance less that e −t of its spine. Using the fact that we can sample points in any order in Wilson's algorithm, we can bound that probability by first doing a cutset exploration around v at scale rṽ. After this cutset exploration, the probability that any of the branches sampled intersects B(ṽ, e −t ) is at most Ce −c(t−t0) for some constants C, c (see Lemma 5.6). At the same time, after the cutset exploration aroundṽ the spine attached toṽ is necessarily fully sampled (see e.g., Proposition 5.9, final assertion). Hence Thus overall, either max Y⊆p[t1,t2]∩B(ṽ,e −t−1 ) c |W (Y,ṽ)| is the same as the variable in Lemma 6.6, or with probability at most Ce −c(t−t0) , it is at most a sum of O(t) variables with uniform stretched exponential tail. The moment bound now easily follows.
For the case where the universal cover is the unit disc, a similar argument can be done provided we can construct a path γ which is nice and avoids the spine. Let us also introduce τ S as the first time γ hits the spine and τ nc as the ending time of γ. As recalled in more details in appendix A, it is a well known fact from Riemann surfaces that there exists a Möbius transform φ = φ γ,ṽ such that where φ (n) = φ • . . . • φ and the union is disjoint. In other words, we keep applying the same Möbius transform to obtain all the copies. Furthermore φ can be written as φ = Φ −1 •µ•Φ where Φ is a Möbius map from the unit disc to the upper half plane and µ is either the translation by ±1 or a multiplication by λ > 1.
If it is a translation we use the same argument as before. Otherwise, it is a scaling by λ > 1 in which case φ(S) necessarily converges to infinity. Furthermore, let τ be the time at which (Φ(γ)) reaches its minimum, we define γ by appending a straight vertical segment from Φ(γ(τ )) to R (note that this segment may not intersect Φ(γ) nor the subsequent scalings of the image of the portion of the spine Φ(γ([τ S , τ nc ])) and then mapping Φ(γ[0, τ ] ∪ ) back to the disc by Φ −1 . By construction, it is trivial to check that Φ(γ ) does not intersects Φ(p(τ, ∞)) so γ does not intersect p(τ, ∞). On the other hand, Φ −1 is the image of a line segment by a Möbius transform so it is a circular arc and therefore its winding around any point is bounded by 2π. We can then conclude using exactly the same reasoning as in the torus case, with only the constant π replaced by 2π on the right hand side of (6.6) and hence obtain the analogue of (6.7).
Finally, let us consider the case when p contains a portion of the spine corresponding to the skeleton. Note that the whole argument was geometric and the only place where we used the probability was in (6.8) which is provided by Lemma 5.5. Lemma 6.8. For all k ≥ 1, there exist a constant C > 0 such that the following holds. Fix a compact set K ⊂M . Suppose we are in the setup of Lemma 6.7. Let p = pṽ be as above andp be the curve which starts atṽ, hits the spine, and then goes to infinity in the direction opposite to the orientation of the spine. Orient and parametrise both p,p fromṽ to ∞ so that t 2 is the last time they exit B(ṽ, e −t ). Then for all δ < δ K , t ∈ [t 0 , ∞],ṽ ∈ K, We emphasize here that the spine might be identical to a spine of a skeleton branch.
Proof. First we tackle the case of p. If t < log(C δ 0 /δ), we need to add O(t) many variables given by Lemma 6.7 and hence we are done by Minkowski's inequality. If t > log(C δ 0 /δ), first we apply the previous bound up to t = log(C δ 0 /δ). For the remainder, we can simply bound the winding by the volume of the ball of radius C δ/δ 0 aroundṽ which is O(1) by our assumptions on the graph (see assumption (i) in Section 2.3). Indeed, this quantity is bounded by the number of times p crosses a straight line joiningṽ to the boundary of the ball, which is simply bounded by the volume of the ball. Forp, we can use Lemma 6.3 and the bound for p and the proof is complete.
If |ũ −ṽ| ≥ e −T , we take t = − log |ũ −ṽ|. Then using Lemma 6.3, we can write where k 2 is the last exit from the B(ṽ, e −k ) by pũ. Using Lemma 5.3 or Lemma 5.5 whichever is appropriate, we get an exponential bound on the expectation of the indicator event above (notice for scales k > log(C δ 0 /δ), the probability actually becomes 0 so this is a finite sum). Therefore, we can again use Lemma 6.8 and Cauchy-Schwarz to conclude.

Convergence of height function
In this section we prove Theorem 6.1. Let us describe informally the general structure of the proof. Most of the work will be in the universal cover to obtain the convergence of expressions of the form E i (h(z i ) − h(w i )) for distinct pointsz i andw i (Lemma 6.11). By integrating this expression, we will then obtain the convergence of the height function, seen as a function on the universal cover.
As mentioned in the introduction, for the first part, the idea is to introduce a regularised height function h #δ t which is a continuous function of the CRSF so that the convergence of h #δ t as δ → 0 is immediate. This leaves us with two issues: the comparison between h #δ and h #δ t for fixed δ and the convergence of h t as t → ∞ in the limit. Both questions are actually solved simultaneously by the estimate of Lemma 6.11 which compares h #δ and h #δ t with an error term that becomes small both in δ and t independently.
Setup and notations. Recall that our goal is to prove that M h #δ ext (z)f (z)dµ(z) (6.9) converges in law and also in the sense of moments as δ → 0. (Recall thatX denote the centered random variableX = X − E(X) whenever this expectation is well defined). Implicitly, h #δ is sampled according to the dimer law (1.1). In other words, by Theorem 3.3, we need to first sample a Temperleyan pair (T , T † ) and apply the bijection ψ in that theorem.
However, it turns out to be more convenient in the proof to work with a pair (T , T † ) sampled from P Wils . To this pair (T , T † ) we can apply the bijection ψ from Theorem 3.3 and obtain a dimer configuration m whose centered height function can be studied. We will first prove (6.9) for P Wils and then later explain how this implies the same for P Temp .
Recall that since M f (z)dµ(z) = 0, the expression in (6.9) is in fact well-defined. We wish to compute the moments of this integral but only in terms of height differences since the field is a priori defined only up to constant. To do this, we use the following trick. Note that since f dµ = 0, f + dµ = f − dµ =: Z(f ) where f + = max{f, 0} and f − = − min{f, 0}. Now we can write . Pick k distinct pairs of points (z 1 , w 1 ), . . . , (z k , w k ) ∈K and let (f (z 1 ), f (w 1 )), . . . , (f (z k ), f (w k )) be the faces containing them. Let z #δ i , w #δ i be the midpoint of the diagonals of f (z i ), f (w i ). Now recall Theorem 4.10 which relates the dimer height difference between two faces f and f to the winding of a specific path γ f,f connecting m(f ) and m(f ) and additional terms of the form ±π associated with jumping over components of the CRSF. Let γ #δ i be the path γ f,f when f = f (z i ) and f = f (w i ), and orient it from z #δ i to w #δ i . Then with these notations, Theorem 4.10 says that where Ψ #δ i is the π S (ε S + δ S ) term in Theorem 4.10. We drop the superscript δ from now on for clarity. Thus from now on we focus on proving convergence of Clearly, (6.13) can be expanded as a sum of 2 k many terms of the form: where S is a subset of vertices of X with distinct indices and γ x is γ i for some i such that x = z i or w i (of course, the products in the expansion of (6.14) has further restrictions, but we ignore that for clarity). We are interested in estimating and proving convergence of (6.14). To that end we employ an idea similar to that in [4]: we truncate the CRSF branch at a macroscopic scale and deal with the truncated macroscopic winding and the remaining microscopic winding part separately.
We now fill in the details. Parametrise γ = γ x . Define (γ x (t)) t≥1 to be γ x at the first entry time into the ball B(x, e −t ) (Note that at the moment, B(x, e −t ) might overlap for different x ∈ X .) We emphasise here that we parametrise γ in the opposite direction compared to Sections 6.2 and 6.3 to be consistent with [4]. With an abuse of notation, we will denote by γ x (−∞, t] the whole path from the opposite end of γ up until γ x (t). Define the regularised term and the error term as When we want to emphasise the role of z, w, we write γ zw , Ψ zw , F z (γ zw ) in place of the above. We start with a general moment bound for the truncated winding.
Lemma 6.10. For all m, there exists constants c = c(K, m), α = α(m) such that for all z, w, δ < δK, t ≥ 0, Proof. The proof is immediate from Lemmas 6.8 and 6.9 and the fact that Ψ zw is bounded by a constant c(K) times the number of noncontractible loops in the CRSF. Indeed, this is clear from the fact that for a fixed compact setK ⊂M there exists a constant c(K) such that any curve P ⊂K will cross at most c(K) copies of a given spine (or lift of a loop). Since the number of noncontractible loops has superexponential tail ([5, Theorem 3.1]), the moment bound of Ψ zw is immediate.
By compactness, choose rK small enough so that p is injective in B(z, rK) for all z ∈K. Now let r X be defined as in (5.1) but for the set of vertices p(X ) (which are all distinct by assumption on X ). We observe that r X ≥ c(K) min x =y∈X |x − y| for some constant c(K). We set ∆ = ∆(X ) = 1 10 (r X ∧ rK) (6.17) Lemma 6.11. There exist constants c, c > 0 such that for all m, m ≥ 1 there exists α > 0 such that the following holds. Let S ⊂ X , S ⊂ X be disjoint. Also assume that |S| = m, |S | = m ≥ 0. Let ∆ be as in (6.17). Then for all δ < δK and t ≥ log(rK/∆) Proof. For simplicity, we first assume S is empty. Recall that we can sample a CRSF from P Wils by first sampling the branches s and then the rest by Wilson's algorithm. Note that Lemma 5.5 tells us that the isolation radius corresponding of each vertex in S with respect to the skeleton has polynomial tail. Since after sampling the skeleton, the rest of the branches are sampled simply by Wilson's algorithm, we can conclude that the application of Proposition 5.9 is valid with ∆ in place of r i . We perform the coupling in Section 5.2 (in particular Proposition 5.9) with points p(X ) and their lift X , and a compact domain D ⊂M containing all points in X so that the minimal Euclidean distance between any point in X and ∂D is at least r X . Note that we can choose r X for the r in (5.1) there. Call (T D x ) x∈X the resulting independent UST in D. For x ∈ X let γ D x be the branch in the UST T D x joining x to ∂D. Let γ D x (t) be parametrised so that γ D x enters B(x, e −t rK) for the first time (going from the outside to x) at time t.
Let R x be the isolation radius of x in the application of Proposition 5.9 and write R x = e −Ix rK. Let I = max x∈X I x . (Note that for notational clarity, I here is shifted by log(rK/∆) from that in Section 5.2.) We now decompose Therefore, we need to deal with expectation of products of α x , ξ x for different indices x. Let C be the sigma algebra generated by the cutset exploration. We will first compute the conditional expectation and then take the overall expectation. Note that T D x are independent for different x and also independent of C and hence any term involving α x is 0. Thus we only have to deal with terms involving ξ x , x ∈ S.
Let Λ be the last time γ x enters B(x, e −I rK). In the coupling, γ x , γ D x agree inside B(x, e −I rK) so in particular γ x (Λ, ∞) = γ D x (Λ, ∞), so Let G be the event that ξ x without the bar (meaning without the expectation terms) is 0.
Observe that if G does not occur then one of the following events happen. Either which has probability at most e −c(t−log(rK /∆)) ∨ δ c ≤ (e −c(t−log(rK /∆)) + δ c ) (Proposition 5.9). Otherwise, γ x has to exit B(x, e −I rK) after hitting e −t rK. This also has probability at most e −c(t−log(rK /∆)) (Lemma 5.7). Now we bound the moments of ξ x . Notice that we can write where we first used Cauchy-Schwarz and then used Lemma 6.8 to bound the moment and Proposition 5.9. Thus overall E Wils (|ξ x | k ) ≤ |1 + t| k (e −c(t−log(rK /∆)) + δ c ).
We have a product of at most m terms and so using Hölder's inequality, taking α = m works. Finally, if S is non-empty, the proof is exactly the same as the vertices in S are distinct from S and hence local independence from the coupling still holds. We then use Lemma 6.10 to conclude using Cauchy-Schwarz. Details are left to the reader. Corollary 6.12. Let S ⊂ X containing vertices with distinct indices such that |S| = m. There exists a constant c = c(K, m), α = α(m) such that for all δ < δK, with t = 2 log(rK/∆). Then this is a straightforward application of Lemmas 6.10 and 6.11. Now we prove a convergence of the height function integrated against f in the sense of moments (still for the Wilson law P Wils ). Recall that X δ converges in the sense of moments as δ → 0 if for all k, E(X k δ ) converges as δ → 0.
Lemma 6.13. K h #δ ext (z)f (z)dµ(z) converges in the sense of moments under the law P Wils . Furthermore, the limit does not depend on the sequence (G ) #δ .
Proof. Using eqs. (6.11), (6.13) and (6.14), we need to prove the convergence as δ → 0 of K2k E Wils ( x∈SF x (γ #δ x ))f We use Fubini to bring the expectation inside the integral in the above display. We first observe that integrating (6.20) over all sets of vertices S so that p(x) = p(y) for all x = y ∈ S is enough (recall p :M → M is the covering map). Indeed, using Lemma 6.10 and the fact that f ∞ < ∞, we see that the integrand can be bounded by O(log α (1/δ)). Since the volume of {S ⊂K 2k : p(x) = p(y) for some x, y ∈ S} is O(δ) (indeed the number of preimages of any vertex inK is bounded by a constant depending only onK), we see that the integral over this set is O(δ log α (1/δ)).
Thus we now concentrate on the integral (6.20) over A(K) = {sets of vertices S so that p(x) = p(y) for all x = y ∈ S} (6.21) We now use Corollary 6.12 and dominated convergence theorem. Since f ∞ < ∞ and log m (∆) is integrable for any m > 0, we need to prove convergence of the expectation inside the integral above. Now we claim that the regularised part E Wils ( x∈SF x (γ #δ x , t)) (6.22) converges as δ → 0. This follows from our assumption of a.s. convergence of Temperleyan CRSF and because the term inside the expectation is a.s. continuous function of the Temperleyan CRSF. The same can be said about the skeleton using the absolute continuity statement given by [5, Proposition 6.5] (see also (5.3)). Indeed, this follows from a.s. continuity properties of SLE 2 and the fact that the CRSF branches are made a.s. from finitely many chunks of SLE 2 , in particular for a fixed t, the SLE curve is a.s. not a tangent to the boundary of circle of radius e −t (we refer to [5, Section 3]) for details. Lemma 6.10 tells us that the random variable in (6.22) is uniformly integrable which completes the proof of convergence of the regularised term (6.22). Now given a fixed S, choose t ≥ log(rK/∆). Now we claim that that the error term satisfies | E Wils ( x∈SF x (γ x )) − E Wils ( x∈SF x (γ x , t))| < c|1 + t| α (e −c(t−log(rK /∆)) + δ c ) (6.23) Indeed, (6.23) follows by writingF x (γ x ) =F x (γ x , t) +ē x (t) and then expanding and using the bounds in Lemmas 6.10 and 6.11. To finish the proof of the lemma, fix ε > 0. First choose t large enough and then δ small enough so that the right hand side of (6.23) is less than ε. Next choose a smaller δ if needed so that for all δ < δ, via the convergence of the regularised term. This completes the proof as ε is arbitrary. Observe that (6.23) completes the proof that the limit does not depend on the sequence ((G ) #δ ) δ>0 since the main term (6.22) is measurable with respect to the limiting continuum Temperleyan CRSF which is universal by [5,Theorem 1.1].
It remains to prove convergence in law (still for P Wils for now). For this, we need to alter the definition of truncation slightly. Given t > 0, take a cover {B euc (x, rKe −10t ) : x ∈K} ofK and take a finite subcover. Let c(z), c(w) denote the center of one of the balls (chosen arbitrarily) in the finite subcover in which z, w belong. Define F x (t), e x (t) to be the same as F x (t) and e x (t) but we cut off the first time γ x enters B(c(x), e −t ) (so compared to the above, the centre of the cutoff ball is shifted by an amount which is at most e −10t ). However, in this definition, we still measure winding around x not c(x). Lemma 6.14. The statements of Lemmas 6.10 and 6.11 still hold if we replace F, e by F, e everywhere.
Proof. For Lemma 6.10, the proof identically follows from Lemmas 6.8 and 6.9 by noticing that the supremum over all continuous subpaths contain the portion of the branch until it first hits the shifted ball. Notice also that we shift only by an additive term in the exponential scale which is O(t).
For Lemma 6.11, the proof is also identical, in particular we still consider the coupling around the points in S. Because we still shift only by O(t) in the exponential scale, the proof readily follows.
We can now complete the proof of Theorem 6.1.
Proof of Theorem 6.1. Again, we first prove this under P Wils before explaining how to extend this and concluding the proof of Theorem 6.1 to P Temp . We write r = rK and pick a t > 0 to be taken large later. Recall from eqs. (6.10) and (6.12), Here Z(f ) is a deterministic constant so we can assume it to be 1 from now on without loss of generality. Firstly we recall that we only need to integrate over the set A(K) as in (6.21) as the integral over the remaining part is O(δ log c (1/δ)). Let us denote by Y = Y (δ) the above integral over A(K). Introduce Note that X(δ, t) is a sum of regularised winding of finitely many branches and hence converges in law (as δ → 0 and t is fixed) by our assumption. Now we show that for all t > 0 and δ < C δ 0 e −t , E(Y (δ) − X(δ, t)) 2 ≤ c(1 + t) α e −c t .
We expand the above square to get an integral over (z, w, z , w ) ∈ A(K) 2 . We can again without loss of generality restrict to the the set A 2 (K) such that the projection p of all four points (z, w, z , w ) maps to pairwise distinct points on M . In that case let ∆ = ∆(z, w, z , w ) be as in (6.17). We now argue that this integral over the set ∆ ≤ e −t/11 is exponentially small in t. Indeed, we are integrating over a set which has exponentially small measure with respect to µ 4 and furthermore the moments are integrable using Corollary 6.12 and Lemma 6.14.
For the rest of the integral, we write (W (γ #δ zw , z) +W (γ #δ zw , w) +Ψ zw ) =F zw (t) +F wz (t) +ē zw (t) +ē wz (t) and we expand again the product inside the integral to get products of terms of the form e zw (t) +ē wz (t);F zw (t) +F wz (t) −F c(z)c(w) (t) −F c(w)c(z) (t) Note here that the indicator over |c(z) − c(w)| > e −t/10 is included in ∆ > e −t/11 so we can get rid of the indicator. Products containing at least one e are small because of Lemma 6.14. Also note that on ∆ > e −t/11 , |z − w| > c(K)e −t/11 . Thus we need to show E(|F zw (t) −F c(w)c(z) (t)| 2 1 |z−w|>c(K)e −t/11 ) ≤ c(1 + t) α e −c t .
It is well-known (see Theorem 2.5) that the instanton component of the above one-form is completely determined by H #δ . Then we have: f (x)h #δ ext (x) dµ(x), H #δ , T #δ converges jointly in law as δ → 0. The first two coordinates also jointly converge in the sense of all moments. Furthermore, lim δ→0 f (x)h #δ ext (x) dµ(x), H #δ is measurable with respect to the limit T of T #δ . We do not go into the details of these claims for the sake of brevity.

A Geometry of spines
In this section we prove Lemma 4.9 and Lemma 6.4. But before getting into the proofs, we remind the readers certain basic facts from the classical theory of Riemann surfaces. By the uniformisation theorem of Riemann surfaces, recall that there exists a conformal map from the Riemann surface D/Γ to M where Γ is a Fuchsian group which is a discrete subgroup of the group of Möbius transforms. Furthermore, such a conformal map is unique up to conformal automorphisms (i.e. Möbius transforms) of the unit disc. In other words, if Γ, Γ are two Fuchsian groups such that M is conformally equivalent to both D/Γ and D/Γ then there exists a Möbius map φ : D → D such that Γ = φ −1 • Γ • φ. Since we have fixed a canonical lift, we have defined Γ uniquely.
It is further known that Γ is isomorphic to the fundamental group π 1 (M , x) (topologically, this is also known as the group of Deck transformations). This connection is described as follows: choose a base point in the manifold x 0 and a particular lift ofx 0 , then any simple loop based in x 0 can be lifted to a simple curve˜ in D starting fromx 0 , with the endpointx 1 depending only on the homology class of the loop. Then to we associate the map φ ,x0 ∈ Γ that sendsx 0 tox 1 . Note then that φ ,x0 (˜ ) is a curve which projects via p to the same curve in M since M is conformally equivalent to D/Γ (and in particular is homeomorphic) and it does not intersect since is a simple loop. Furthermore, since the endpoint of φ ,x0 (˜ ) is φ ,x0 • φ ,x0 (x 0 ), by the unique path lifting property, the curve˜ ∪ φ ,x0 (˜ ) is the unique lift obtained by going around twice in the same direction. Iterating, we obtain that the infinite path obtained by going around in the same direction is given by ∪ ∞ n=0 φ (n) ,x0 (˜ ) where φ (n) ,x0 is the n-fold composition of φ ,x0 and that the union is a disjoint union.
We can actually say more using the classification of Möbius maps according to their trace. Recall that Möbius maps preserving the unit disc have the form φ(z) = e iθ z − ā az − 1 ; |a| < 1; θ ∈ [0, 2π) and can be classified (up to conjugation with Möbius transforms) depending upon the behaviour of the trace: • If |e iθ − 1| 2 = 4(|a| 2 − 1) 2 , then φ(z) is conjugate to either z + 1 or z − 1 seen as maps from H to H.
Thus there is a Möbius map Φ from D to D or H such that φ = Φ −1 • µ φ • Φ where µ φ is either a translation by ±1 or a scaling by λ (in case Φ maps D to H) or a rotation (in case Φ maps D to D). We argue for any loop , µ φ ,x 0 cannot be a rotation. Indeed, if a map φ ,x0 was conjugate to a rotation, then its iterates would be either periodic, or a dense set (on the image of a circle). Being periodic is forbidden because of path lifting property and being dense is forbidden because Γ is discrete. From the two remaining cases, we can complete the proof of Lemma 4.9.
Proof of Lemma 4.9. Let S be a spine, we can choose a point x 0 on p(S) and a liftx 0 on S. Then we see that the previous general theory applies so we can find an open path˜ (a certain lift of the path going once around a noncontractible loop) and a map φ p(S),x0 such that We now prove Lemma 6.4 which provides bounds on the macroscopic winding of the spines. We repeat the statement of the Lemma here for convenience.
Here the supremum is over all continuous paths obtained by erasing portions ofX[τ ij , τ ij +1 ].
Proof. We are going to borrow the notations from Section 6.2. Take a noncontractible simple loop in M through v and find a continuous path #δ in G #δ which approximates this loop in the sense that it stays within distance cδ/δ 0 from (this is guaranteed to exist by the uniform crossing estimate for small enough δ depending on ). Notice that the lift starting fromṽ of a curve which winds infinitely many times in the clockwise (resp. anti-clockwise) direction of #δ defines an infinite path˜ + (resp.˜ − ) which converge to the boundary in the hyperbolic case (Lemma 4.9) or goes to infinity in an asymptotic direction in case of the torus. Let˜ ± denote the portion of˜ ± from the last exit of B 1 to infinity (given an arbitrary parametrisation starting fromṽ). Notice that (˜ − ∪˜ + ) divides the annulusM \ B 1 into two simply connected domains. By compactness we can find a constant C such that the winding of˜ + aroundṽ is bounded by C, uniformly over all pointsṽ ∈K for a suitable choice of loop . Let τ +,1 be the first hitting time of˜ + in the interval [τ ij , τ ij +1 ] and by induction define τ −,i to be the first hitting time of − after τ +,i and define τ +,i the hitting time of˜ + after τ −,i−1 . Let I + the number of τ +,i before τ ij +1 , i.e I + = |{i|τ +,i ≤ τ ij +1 }|. We now use the deterministic bound sup Y⊂X[τi j ,τi j +1] |W (Y,ṽ)| ≤ (C + 2π)(I + + 1).
and observe that conditionally on either {X τi j +1 = u} or {X τi j +1 ∈ ∂(G ) #δ ∪ s}, I + has an exponential tail. The proof of this fact is essentially the same as the second item of Lemma 4.8 in [4]. Indeed we need to show that once the random walk intersects˜ + outside B 1 , there is a positive probability (uniform in δ) for the walk to create a noncontractible loop without intersecting˜ − (as in the proof of Lemma 5.2). This is intuitively clear, but a complete proof of this needs an input from Riemannian geometry. The issue at hand is that too close to the boundary, the uniform crossing ceases to hold uniformly in δ. We control the winding of this portion of this walk as follows. Consider a compact set A ⊂ M containing which forms an approximation of , in the sense that topologically K is an annulus with being noncontractible in K. Recall from the proof of Lemma 5.2 that the simple random walk has a uniform positive probability to create a noncontractible loop inside K by winding around exactly once and we stop the simple random walk if this happens. But in this process, we can assume without loss of generality on that the lift of the walk stays inside at most four consecutive copies of K (where these four copies are mapped to each other by Möbius transforms). So if the walk is on˜ + at a copy ofÃ which is more than four copies away fromṽ and˜ − , the lift of that walk cannot intersect − in this process. On the other hand, applying four copies of the corresponding Möbius map to an arbitraryṽ ∈K yields by compactness ofK a slightly bigger compactK . Applying the uniform crossing property for the given choice of δ toK , we now simply apply the argument of Lemma 4.8 in [4].