The Concept of Spin Ice Graphs and a Field Theory for their Topological Monopoles and Charges

Spin ices can now be fabricated in a variety of geometries which control their collective behavior and exotic properties. Therefore, their proper framework is graph theory. We relate spin ice notions such as ice rule, ice manifold, Coulomb phases, charges and monopoles, to graph-theoretical notions, such as balance, in/out-degrees, and Eulerianicity. We then propose a field-theoretical treatment in which topological charges and monopoles are the degrees of freedom while the binary spins are subsumed into entropic interaction among charges. We show that for a spin ice on a graph in a Gaussian approximation, the kernel of the entropic interaction is the inverse of the graph Laplacian, and we compute screening functions from the graph spectra as Green operators for the screened Poisson problem on a graph. We then apply the treatment on star graphs, tournaments, cycles, and regular spin ice in different dimensions.


I. INTRODUCTION
Since the Bernal-Fowler ice rule 1 was invoked by Pauling 2 to explain the zero-point entropy of water ice 3,4 the concept has come to describe a variety of other materials, such as pyroclore rare-earth spin ice antiferromagnets [5][6][7] , artificial magnetic spin ice antiferromagnets 8-12 , or artificial particle-based ices 13 . Often 5-8,14-20 , but not always 9,21,22 , frustration impedes ordering in these materials, and leads to degenerate states of constrained disorder, or ice manifolds, of interesting topological properties.
Indeed, the ice rule is a topological concept related to the local minimization of a topological charge 16,23 . As such, it has a wide applicability, and artificial spin ice materials are being designed 24 for a variety of emergent behaviors not necessarily found in natural magnets 16,17,[24][25][26][27][28][29][30] . Violations of the ice manifold are topological excitations 31,32 that, depending on the geometry and local degeneracy of the system, can be deconfined. Further, in magnetic materials these topological charges are also magnetic charges, and if deconfinable, magnetic monopoles [31][32][33][34][35][36] . As such, they interact via a Coulomb law 32 , are sources and sinks of the H field, can pin superconductive vortices in spin ice/superconductors heterostructures 37 , and could perhaps exert a localized and mobile magnetic proximity effect 38 in heterostructures that interface two-dimensional spin ices to transition metal dichalcogenides or Dirac materials. Finally, in artificial realizations, these topological objects can also be read and written 39,40 , and might therefore function as binary, mobile information carriers for spintronics for for neuromorphic computation 41,42 .
In this work, first we extend the concept of spin ice on general graphs, and relate the two fields and their jargons, because many concepts of spin ice physics have been investigated under different names and with different aims in graph theory.
Then we propose here a unifying framework for degenerate spin ice on a general graph. We treat these charge excitations as degrees of freedom, subsuming the underlying spin vacuum into entropic interactions among charges. This work has many motivations.
Thirdly, spin ice can now be realized in the quantum dots of a quantum annealer 81 on a wide variety of graphs.
Finally, many common notions in spin ice physics posses a direct equivalent in graph theory. For instance: spin ice charges are the degree-excess between indegrees and outdegrees; the degeneracy of spin ice on a complete graph is simply the number of its regular tournaments; an ice rule configuration realizes a graph that can be unraveled in an Eulerian path, etcetera. A wealth of work has been done in graph theory, some of which could be easily translated to spin ice physics. This work represent an invitation extended both to graph theorist, to contribute to the field of spin ice physics, and to physicists working in spin ice, to broaden their mathematical approaches to these systems.

A. Spin Ice
Ice comes in about eighteen crystalline forms 82 , all of them involving oxygen atoms residing at the center of tetrahedra, sharing four hydrogen atoms with four nearest neighbor oxygen atoms (Fig. 1). Two of such hydrogens are covalently bonded to the oxygen of their molecule, and two realize hydrogen bonds with oxygens of different molecules. Thus, two are "in", two are "out" of the tetrahedron, and this is the so called ice-rule introduced by Bernal and Fowler 1 . Each tetrahedron therefore realizes 6 admissible configurations out of the 2 4 = 16 ideally possible, and the collective degeneracy of the ice grows exponentially in the number of tetrahedra N as W N . This leads to a non-zero residual entropy per tetrahedron k B ln W . Pauling famously estimated W = 3/2, remarkably close to both the experimental 3,4 and the numerically obtained value (W = 1.50685 ± 0.00015 83 ).
One can associate to the ice rule a spin model (Fig 1c) where binary spins are assigned to the bonds between molecules pointing toward the proton. Then the ice rule dictates that two spins point in, two point out, as is the case of magnetic spin ices. In rare earth titanates such as Ho 2 Ti 2 O 7 and Dy 2 Ti 2 O 7 , the magnetic cations Ho 3 + and Dy 3 + carry a very large magnetic moment, µ ∼ 10µ B . At low temperature they can be considered binary, classical Ising spins constrained to point along the directions of the lattice bonds which form a pyrochlore lattice, and are thus expected to interact as magnetic dipoles obeying the ice rule 5,84 .
While this ice manifold represents an interesting manifestation of a Coulomb phase which has been amply studied 57,59-61 and which can be considered a prototype for classical topological order 62,63,85 , the exotic behaviors of spin ices proceed not so much from said topological structure, but rather from how it is broken, e.g. via fractionalization into monopoles 31,32,63 , and how much of it is instead retained, e.g. via spin fragmentation [86][87][88] .

B. Spins on a Graph
We consider the most general case of a spin ice on a connected, undirected, simple graph 89 G. The phase space for G is the set of all directed graphs, o digraphs, that can be built on G by assigning an orientation to its edges.  1. a) The crystalline structure of ice shows proton disorder in the allocation of hydrogen atoms among oxygen centered tetrahedra (image from ref 82 ). b) In water ice oxygen atoms sit at the center of tetrahedra, connected to each other by a hydrogen atom. Two of such protons are close (covalently bonded) to the oxygen at the center, two are further away. c) One might replace this picture with spins pointing in or out depending on whether the proton is close or far away. Then two spins point in, two point out. This corresponds to the disposition of magnetic moments on pyrochlore spin ices, rare earth titanates whose magnetic ensemble does not order at low temperature, because of frustration, and, much like water ice, provides non-zero low temperature entropy density. Consider an undirected, simple graph G of a number N l of edges labeled by l, connecting a number N v of vertices labeled by v and of various degree of coordination z v . We call {vv } an edge l among vertices v, v . For such graph, the adjacency matrix 89 is the matrix A vv , such that A vv = 1 if v, v are connected and A vv = 0 otherwise. It contains all the information of the graph. Obviously, A vv is symmetric and We can define binary variables or Ising spins S l on each edge l, as in Fig. 2, via an antisymmetric matrix S vv such that S vv = 0 if v and v do not share an edge, S vv = 1 if they do and the spin points toward v , and S vv = −1 if they do and the spin points toward v.
In the language of graph theory, A vv defines an undirected graph, while an Ising spin structure S vv defines a digraph. The S matrix is related to the non-symmetric adjacency matrix of the corresponding digraph, A dir vv , whose elements have value one if and only if vv are connected by an edge pointing toward v (on simple graphs). Then S is the anti- . For an undirected graph G, we call its phase space P, the set of the digraphs that can be specified on it. Clearly, the cardinality of the phase space is |P| = 2 N l .

Ice Manifolds
An ice manifold is a proper subset of the phase space P that minimizes the topological charge, defined as follows.
Given S vv , its topological charge distribution is the vector Q v defined for each vertex v as the difference between the edges pointing in and out of v, or (1) In graph theory, Q v is thus the difference between indegrees and outdegrees of the directed graph that corresponds to a particular spin configuration on G. Q v can have the values from which we immediately have A digraph S vv (i.e. a spin configuration) for which |Q v | is minimal on each vertex, i.e. has zero charge on all vertices of even coordination and ±1 charge on all vertices of odd coordination, is said to obey the ice-rule.
Then, given a graph G, we call the ice manifold of G the subsets I ⊂ P of its spin phase space P made of all directed graphs that obey the ice rule. Then, |I| is the ice-degeneracy of the graph and s = N −1 l ln |I| is its Pauling entropy per spin. In the language of graph theory, an ice rule configuration corresponds to a balanced or quasi-balanced directed graph. By a famous theorem, if all vertices have even degrees (i.e. even coordination) a balanced digraph is Eulerian. That means that an Eulerian path exists (a walker can follow the arrows and walk the entire graph without passing the same edge twice 90,91 ). Moreover, if the graph is complete and has odd vertices, a spin configuration is called a tournament because it describes winners and losers in a roundrobin tournament 92,93 . Thus, an ice rule configuration on a complete graph describes a regular tournament, one where all players score the same. Therefore, the ice-degeneracy of a complete graph of odd number of vertices counts all the possible outcomes of a roundrobin tournament where every player beats the same number of other players (which is still an open problem).
We see therefore that spin ice physics is a case of the theory of balanced or Eulerian graphs. The latter is greatly consequential, with applications in biology 94 , computer science 95 , logistics 96 , and sociology 97 . In complex networks, digraphs can describe the sharing of a token among agents, and deviation from the ice-rule can quantify the agent satisfaction or frustration, and give measures of fairness 98 .
Depending on the topology of the graph it is not obvious if and when the cardinality of the ice manifold scales exponentially with the size of the graph, thus leading to a non-zero Pauling entropy 2 . Indeed, the one dimensional ferromagnetic Ising model can be mapped into a spin ice on a path graph, or a cycle, and its ice manifold has cardinality two, regardless of the number of vertices. We will show elsewhere that the path graph is the only such case.

Coulomb Phases
The concept of a Coulomb phase appeared first in gauge field theories 56 and was then introduced in the theory of pyrochlore spin ice 57,59,61 . In simple terms, it corresponds to a disordered spin texture that can be coarse-grained to a solenoidal magnetization field. It can be considered a case of classical topological order 62,63,85 where no order parameter exists but instead the disordered states are labeled by gauge fields whose solenoidal nature expresses the constraint over the disorder. Excitations are then violations of the solenoidal requirement.
We can generalize this notion to a graph spin ice by expressing it without coarse graining or gauge fields. We say that two spin assignations S, S ∈ P are charge equivalent if and only if the difference (S vv − S vv )/2 has zero topological charge on every vertex. Chargeequivalence is an equivalence relation and thus induces a partition on the phase space P. We call each class of equivalence in that partition a Coulomb class. Each Coulomb class is labeled by a distribution of charge Q v .
A trivial example: a graph made of only two vertices connected by one edge has a spin phase space of cardinality 2, corresponding to the two orientations of its only spin. Its ice manifold coincides with the spin phase space, and it contains 2 Coulomb classes, each of cardinality 1. For a less trivial example, the ice manifold of pyrochlore ice is a Coulomb class corresponding to charge Q v = 0 everywhere. The same is true by definition for every graph that has only even coordinations.
Clearly, not all ice manifolds are Coulomb classes. For instance the ice manifold of the honeycomb graph, so-called kagome ice, is not: for any spin configuration there is always at least one spin, and in fact an extensive number of spins, that can be flipped individually, thus changing charge configurations (and Coulomb class) without violating the ice rule. And yet, the kagome ice manifold can be partitioned into Coulomb classes. Crucially, its Ice II phase [18][19][20]53,99 corresponds to two Coulomb classes, each of charge alternating in sign ±1 on nearby vertices. Indeed, the spontaneous symmetry breaking between the two Coulomb classes drives a second order transition to charge ordering (of the Ising universality class).
A Coulomb class becomes important if it represents a low energy state. In such case it defines a Coulomb phase. Thus, the ice manifold of spin ices of even coordination, including therefore pyrochlores, is a Coulomb class and also a Coulomb phase. In kagome ice, each of the two Coulomb classes of the Ice-II phase are Coulomb phases, which might explain why it is so hard to observe it experimentally 20,25,100,101 , whereas the Ice-I phase is found easily 8,14,102 .
Indeed, a Coulomb class imposes topological constraints on the kinetics within the class: it prohibits single spin flips and requires collaborative flips (corresponding to Euler trails) that might be extremely unlikely in a realistic dynamics. Thus, in a Coulomb phase, all kinetics must happen above that phase, by breaking the topology of the Coulomb class, that is by changing its defining distribution of charges. In pyrochlore ices or in degenerate artificial square ices, where the entire ice manifold is a Coulomb phase, the breakage of topological protection consists in the appearance of monopoles, that is ±2 charges. In such systems all single spin kinetics is then monopole kinetics. In the Ice-II phase of kagome ice, also a Coulomb phase as noted by McDonald et al. 100 , kinetics consists in the breaking of the ±1 charge alternation.
Because properties of topological order pertain to Coulomb classes, which are defined by charges, we develop a formalism in which charges are the relevant degrees of freedom. We reach an effective free energy for charge distribution, which contains the entropic effect of the corresponding Coulomb class.

III. FIELD THEORY OF EMERGENT CHARGES
We need a Hamiltonian whose ground state is the ice manifold, and thus allows us to deal with breakages of the topological protection in terms of energy costs. The following is minimal where > 0 is an energy. Note that H more generally measures how "unbalanced" a digraph is. Note also that for every digraph v Q 2 v is larger than or equal to the number of vertices with odd coordination.
The corresponding partition function reads where β = 1/T , H vv is an antisymmetric matrix modeling an external "magnetic" field defined on the edges with respect to the orientation l = vv , so that H vv = 0 if and only if v and v share an edge, and H has dimensions of an energy (the Zeeman energy), and V v is an external potential acting on the charges. Clearly then are the one-and two-point correlation functions for spins and charges.
The Concept of Spin Ice Graphs and Field Theory for their Topological Monopoles and Charges 7 To produce a field theory we need to remove the discrete variables S l . We do so by a common trick. We insert in the sum of (5) the tautological expression 1 = and then we sum over the spins. We introduce where vv are edges (counted once, so that if vv is counted v v is not). From it we obtain where the matrix ∇φ is defined on the edges only as and is called the gradient matrix of φ. Note that if the graph is embedded in a linear space and vv is the vector pointing from v to v then (5) can now be rewritten as The averages and correlations of q are the same as for Q and therefore, for instance q v = Q v and also We have thus replaced the binary spin variables S l defined on the edges with a continuum field q v defined on the discrete set of vertices, but we have gained a termΩ[q]. It represents a generalized degeneracy or density of states for the charge distribution q v , emergent from the many possible underlying spin ensembles compatible with q v . It also constraints q v to proper, discretized values.
We can call S = lnΩ[q] the generalized entropy for the charge distribution q. Then, the effective free energy at zero loop 103 for q v is given by the quadratic part of H[q] − T S[q] as we shall see. Also, in absence of a field, from (12), (6) and the parity of the functions involved, q v = 0 for every v and S l = 0 for every l, as one would indeed expect from trivial considerations on the model. To see the same formalism from a different angle we can rewrite the partition function as with and Equations (14), (15) look familiar in the language of quantum field theory. They correspond to a charge field q v , for which iT φ v acts as a bosonic field (of "Lagrangian" F[φ]) mediating an interaction between charges. Again, the interaction is not real but comes from the underlying binary ensemble from which the charge field is an emergent observable. Both pictures come useful in different scenarios, as we shall see.
Integrating (14), over dq v we obtain which shows that in general iφ is real. Indeed from Eqs. (6,17) we have To further clarify the intuitive meaning of iφ as an emergent, entropic field translating the effect of the spin correlations, from Eq. (17) we can immediately deduce It is then natural to introduce an entropic "magnetic field" B e vv = i∇ vv φ, and then when S vv is small (which requires H vv small). Note than in general, from Eq. (18), the entropic magnetic field is related to the gradient of the charges, or where E is the gradient of V . It is time to confess that V v is superabundant, though useful. Indeed, nothing changes by incorporating V v into H vv by replacing in our equations Similarly, if H can be divided into a gradient plus a term irreducible to a gradient, or H vv = H vv + ∇V vv , the equations above are still valid with the substitution V → V + V , H → H . In the following we will assume V = 0.
Finally, taking the divergence of Eq. (19), we obtain In the linearized limit in which Eq. (20) is valid, we obtain q = −βµ q ext + iLφ.
In the previous equation the external charge is defined as µq ext = div[H], where µ is an energy, and the Laplace operatorL is defined via the Laplacian matrix orL =D −Â where D vv = z v δ vv is the degree matrix and z v is the degree or coordination of the vertex v. The Laplacian matrix is the generalization on a graph of the discretized Laplacian operator on a lattice. The reader can easily verify that on a square lattice of edge length a when one takes the usual continuum limit for a → 0 one findsL → −a 2 ∇ 2 . Also, the reader can verify that for a generic ζ v defined on the nodes v as one would expect as the generalization of the notorious ∇· ∇ = ∇ 2 , valid in linear spaces, and which we have used to deduce Eq. (23). Thus, Eq. (23) tells us that charges are the sources and sinks of the entropic potential, and represent a generalization of the at least in a linear approximation. There, the entropic field iφ obeys a generalized Poisson equation and is therefore the Coulomb potential of the charges on a graph. We shall see in the last section that indeed when the graph can be properly embedded in a linear space, the entropic interaction among charges is indeed the standard Coulomb interaction in the proper dimension of the space.

IV. HIGH T APPROXIMATION
By taking the high T , high H limit but keeping H/T finite, the Gaussians in (17) tend to delta functions in φ and we obtain Unsurprisingly, the above is the standard "paramagnetic" partition function for an uncorrelated system. It leads, via Eqs. (6) to the familiar magnetization law for a paramagnet S l = tanh(βH l ).
When H = 0, we get from (26) the correct entropy per spin at high temperature, or s = ln 2. Thus under assumption of disorder, the high T expansion corresponds to an expansion in the small entropic field φ. In doing so we lose the constraint to discrete values of q v imposed by Ω[φ], which is however not relevant at high T .

A. Free Energy and Entropic Interactions
We can expand Eq. (16) at lowest order and find and thus write F in Eq. (15) at the second order as where matrices are expressed as operators acting on q, φ, which are vectors of dimension N v , q 2 = q · q = v q 2 v , and |H| 2 = vv A vv H 2 vv /2. Note that Eq. (29) is correct only because the divergence of H is hereby defined discretely with respect to vertices by Eq. (2). One might at first suspect that it is incorrect, and immagine that, e.g. on a square lattice a continuous solenoidal H would not couple to the spin ice. However, a moment thought shows that even for an uniform H, its divergence, as defined on vertices in Eq. (2), is necessarily non-zero at the boundaries of the lattice even for a uniform field. In other words, the definition of divergence on a graph already accounts for boundary conditions (inclusive of "internal" boundaries 46 ).
By integrating over q v we can express the partition function in Eq. (14) in the high T approximation as coming from a free energy in the entropic field only, or In a similar way, to find an effective free energy for the charges, using Eqs (14), (29), we can integrate instead over the entropic fields φ. However, to do so, we must first consider the spectrum 104 ofL.
The following is well known from the spectral theory of graphs.L is symmetric and thus has real eigenvalues {γ(k) 2 } k=0...kmax with k max ≤ N v − 1, and corresponding N v eigenvectors ψ α v (k) (where α counts the eigenvalue degeneracy). It is immediate to verify that γ 2 (0) = 0 for the uniform eigenvector ψ v (0) = 1/ √ N v . In a simple and connected graph all other eigenvalues are strictly positive.
We can go to the new basis, definingq α (k) = ψ α * (k) · q andφ α (k) = ψ α * (k) · φ. Then in Eqs (14), (29), the integration over dφ(0) merely returns a δ(q(0)), which in "real space" corresponds to δ ( v q v ). This ensures that in the new free energy we sum only over charge configurations of zero net charge-as it should be, since a system of dipoles is charge neutral. All other charge modes have zero net charge. Indeed for any eigenvector of L except the one of zero eigenvalue it is true that v ψ v (k) = 0. This follows immediately From Gaussian integration of F 2 in Eq. (29) in the space orthogonal to ψ(0) (whereL can be inverted) we obtain, in absence of field H, the effective free energy for q in the form which can be interpreted both in real or spectral space. The first term is the usual energy cost to produce charges. The second term tells us that, at quadratic order, the effect of the underlying spin manifold can be subsumed into a pairwise, entropic interaction that corresponds to T L −1 , i.e. the Green operator of the Laplacian.
In regular lattices embedded in a linear space (see below), Eq. (60) implies that in three dimensions (3D) charges interact entropically via a 1/x law, as indeed found numerically. In two dimensions (2D) one expects instead a logarithmic interaction. In both case we have Coulomb potentials in the proper dimension, as anticipated in our discussion of Eq. (23).
This, however, also implies a mismatch in systems of reduced dimensionality, for instance in square ice, whose entropic interaction is the Green function of the 2D laplacian while the real monopole interaction (not considered in this work) is the Green function of the 3D laplacian. Such mismatch leads to lack of screening, as we show elsewhere.

B. Charge Correlations
We defineŴ q as the inverse of the kernel of the free energy for the charge of Eq. (60), or and then, from equipartition we have We note thatŴ q can be written in various ways, includinĝ where we call the correlation length at high temperature, as it was already appreciated in more specific systems, via other means 60,61 (see next section). In Eq. (34) we have introduced the "Green function" (actually, a matrix) of the screened Poisson equation on a graph, or The Concept of Spin Ice Graphs and Field Theory for their Topological Monopoles and Charges 11 which, as we show in the next subsection, controls the screening from an external charge. Note that, from Eq. (30)Ĝ is also the correlation for the entropic field φ, or Finally, going to the spectrum of L we find In the infinite temperature limit the correlation length becomes zero, and from Eq. (38) we obtain From Eq. (39) we have which indeed corresponds to the average square charge of uncorrelated vertices, or q 2 uncorr defined as the square charge obtained from counting arguments. Indeed, considering vertex multiplicities only, and computing the average square charge for a vertex of degree z with each charge 2n − z weighted merely by its vertex multiplicity z n , one obtains Note also that Eq. (39) implies zero correlations among vertices that are not nearest neighbors, but a correlation of −1 among neighboring vertices: since nearby vertices share a spin, they have in average opposite charges. How is the infinite temperature limit approached? From Eq. (38) we obtain, at least formally, the series Note now that L n can be written as sums of products of nÂ andD matrices (e.g. ADDAADAD . . . ). A moment's thought should convince that if v, v are separated by more edges than the number of A matrices in such product, then the vv element of the product is zero. An obvious notion of distance between two vertices on a graph is given by the number of edges in the shortest path (aka graph geodesic) connecting them 105 . It follows that if v and v are at a distance d vv > 1 the first nonzero term of the series in Eq. (42) is given by Interestingly, A k vv is known to be the number of walks of length k between the two vertices v and v . Thus, the coefficient A d vv vv in Eq. (43) is the number of walks between the two vertices v and v of length equal to the distance d vv . It is therefore the number of geodesics connecting the two vertices. Note that this considerations are conditional to the possibility of the expansion of Eq. (42), which might not hold in the thermodynamic limit of certain systems.
By construction our approximation doesn't work for T → 0 where the fluctuations of the entropic field diverge. One could speculate that by including perturbative terms functional forms would not change, except for replacing ξ(T ) → ξ r (T ) where ξ r (T ) is the real correlation length at low T . For instance in pyrochlore ice, correlations become screened-algebraic at low T , suggesting that a quadratic free energy might work with proper renormalizations of the parameters.
In the limit ξ → ∞ From Eq. (38) we find Then, for graphs of even coordination we obtain Note that the equation above corresponds to the Debye screening length for a Coulomb potential whose coupling constant is proportional to T , which is indeed the case of our entropic field. At high temperature, the situation is much different, with ξ 2 = /T . The Debye-Hückel approach applies to strong electrolytes that are fully dissociated, and the disorder brought by higher temperature prevents charges from screening, thus increasing the screening radius. In spin ice, instead, charges carry an energy cost and at higher temperature there are more charges available for the screening. Furthermore, Eq. (45) corresponds in pyrochlore spin ice to the experimentally found 106 exponentially divergent behavior of the correlation length in the proximity of the ice manifold. And because for v = v the correlations tend to ∝L −1 , we call them Coulomb in this limit.

C. Entropic Screening of External Charges
Since we have no real interaction among charges, all screening is entropic. One can consider two cases: screening of an external charge and screening of a pinned charge.
One can define external charges as sources and sinks of the H field, which interact with emergent spin ice charges via the coupling of the H field to the spins S. To understand this formally, consider the term vv S vv H vv in Eq. (5). Imagine that we can write a Helmholtz decomposition on the graph, so that H can be represented as H vv = ∇ vv Ψ + H ⊥ vv where Ψ v is a field and the second term has no divergence. Then we have that is, the potential responsible for the divergence-full part of the magnetic field couples to the emergent spin ice charges, thus generalizing on a graph the notion of magnetic fragmentation in spin ice 86 . Let us then call q ext = µdiv[H] the "external charge" (here µ has dimension of an energy) and assume v q ext,v = 0. Then from Eq. (29), integrating over φ one has Thus, an external charge is screened by the screened Green function of the Laplacian. Note that becauseLΨ = −q ext we can also write From Eq. (47) we obtain the two correlation limits The Concept of Spin Ice Graphs and Field Theory for their Topological Monopoles and Charges 13

D. Entropic Screening of Pinned Charges
The case of a pinned charge is obtained by summing the partition function only over spin configurations corresponding to a fixed charge q pin on a vertexv. This corresponds to inserting a δ(qv − q pin ) in the functional integral. We leave to the reader the simple calculation, whose result is and correctly yields qv = q pin . Note also, from Eq. (34), The difference between the two screenings should not surprise. In the first case an external field interacts locally [see Eq. (46)] with the spin ice, inducing local effects that then propagate via the charge-charge correlation, while the second case is due to correlations of the free charge with the pinned one.

A. Star Graph
A star graph is made of n nodes v i each connected only a central node w. Thus, N l = n, N v = n + 1. The number of ice rule configurations are W = n n/2 when n is even and W = 2 n (n+1)/2 when n is odd. This leads to a Pauling entropy per spin s = ln(2) in the limit of large n. The fact that the Pauling entropy coincides with the entropy at infinite temperature is a feature of spin ices with infinite spins per vertex. In particular, to have a finite energy per spin we need to scale the energy coupling as = J/n.
To find the charge correlations, note immediately that by symmetry it must be q v = 0, q w = 0, and q 2 v = 1 for every v. The only possible correlations are q 2 v = 1, q v q w , q 2 w and q v q v . Using tricks such as q v q w = n −1 v q v q w = n −1 v q v q w and charge neutrality one readily finds and we need to find only q 2 w . Note that Z = exp(βJ/2)Z w , where the latter is given by in the limit of large is the usual binomial entropy]. In the limit of large n we can perform the quadratic expansion around the minimum of f (which is at ρ = 0) obtaining f (ρ) ∼ (1/2)(βJ + 1)ρ 2 , and thus ρ 2 = n −1 (βJ + 1) −1 . Then from q w = nρ and from the Eqs. (52) we have Note that q 2 w → z w at large T and q 2 w , q v q w → 0 at low T , as expected. Note also that the Eq. (43) is valid: from Eq. (54), q v q v = −ξ 2 + O(ξ 4 ) and indeed the number of geodesics between v and v is one.
What are the correlations deduced from the field theory at high T ? We leave to the reader the simple task of computing the matrixŴ q and verifying that W q ww = n 2 n + (n + 1)βJ W q vw = − n n + (n + 1)βJ The equations above correctly reduce to the Eqs. (54) in the limit of large n in which the latter were derived. Therefore, the field theory is exact at any temperature for the star graph in the thermodynamic limit. That should not surprise. We know that field theories become exact for Ising models of binary spins when each spin interact with each other in the same way, that is when the graph connecting interacting spins is a complete graph, which is the case for the star graph.
Finally, pinning a charge Q p on w elicits, from Eq. (50) a screening charge q v = −Q p /n, which can also be obviously deduced from charge cancellation. More interesting is when the charge is pinned on v. Then it must be Q p = ±1, and the screening is: q w = −Q p (1 + βJ) −1 , q v =v = −(Q p /n)βJ/(1 + βJ), going to zero in the thermodynamic limit.

B. Complete Graph
In a complete graph K n all nodes are connected to all others. We have N v = n, N l = n(n − 1)/2. As already explained, a spin configuration configuration is a digraph called a tournament, and there are |P| = 2 n(n−1)/2 tournaments. When n is odd, vertex coordination is even and an ice-rule obeying tournament is called a regular tournament. Thus |I| is the number of regular tournaments with n players, which is still an open problem. There are, however, asymptotic formulas for large n. From McKay's formula 107 we find that the Pauling entropy per spin in the thermodynamic limit is s = ln(2). As in the star graph is coincides to the entropy at infinite temperature.
By the symmetry of the problem there are only two charge correlations, q v q v (when v = v ) and q 2 v , and they are related: . Note that, unlike in a star graph, not all spin interact with other spins. However, its partition function can be factorized in terms corresponding to star graphs in the limit of large n. Indeed, consider K n−1 , then K n can be former by adding a vertex v n and n − 1 edges connecting it to the v 1 . . . v n−1 vertices of K n−1 . Its partition function is therefore where j k = ln e 2β k j=1 qj n−1 , and the suffix indicates that the average is performed on the K n−1 graph. Notice than, again, we must assume = J/n in order to have an extensive energy. Then, with the same transformation as for the star graph, β (2k−n) 2 = nβJρ scales as n, whereas j k is subextensive in n and can be neglected. Thus, for large n, Z n = Z n−1 Z w , with Z w given by Eq. (53). We obtain therefore, as before and from q v q v = − q 2 v /(n − 1) we have Note in particular that while the entropy is extensive in the number of spins, energy scales linearly in the number of vertices. Indeed E = v q 2 v /2 = (n/2)JT /(T + J): the energy per vertex is finite and non zero (except for T = 0) while the energy per spin is always zero.
From Eq. (50) we see that by pinning a charge Q p in a vertex, the charge elicited by entropic screening in all other vertices is q = −Q p /n, which is obvious since there must be charge cancellation. The elicited charge is zero in the thermodynamic limit if Q p is finite, and non-zero if Q p instead scales linearly in n.
Finally, the adiacency matrix for the complete graph K n isÂ =Ĵ −1, and thus the Laplacian isL = (n − 1)1 −Ĵ, whereĴ is the matrix of ones. From that it is immediate to findŴ q as which corresponds to the correlations of Eqs. (57,58) for n large. Again, we see that the field theory is exact at any temperature.

C. Spin Ice on a Lattice
We will treat elsewhere the cases of spin ice on a lattice in full. Here we report general considerations that are common to spin ices that can be properly embedded in a linear space of dimension d, such that in the thermodynamic limit it can be homogenized to a continuum (v → x), and distances are measured in units of the lattice constant such that Eq. (60) becomes Because the Laplacian operator of the graph coincides with that of the linear space, its eigenvectors are plane waves and at small momentum (or large distances) and the eigenvalues are γ(k) 2 k 2 (we measure space in units of the lattice constant). Then we have for the entropic interaction at large distances which are all Coulomb potentials in their proper dimension, as already anticipated in section III. The third line coincides with the numerically verified entropic interaction for pyrochlore spin ice 64,108 .
From Eq. (38) we obtain |q(k)| 2 = k 2 1 + ξ 2 k 2 = ξ −2 − ξ −4 from which we have, via Fourier transform, valid at large distances. The first line corresponds to the exponential form of the screening in a 1D Ising system, where domain walls correspond to charges. We note that ξ thus known exactly as ξ 1D-Ising = −1/ ln tanh(β ) and does not follow the ∼ /T at high T . Though one can verify that at low T we have ξ 1D-Ising ∼ 1/ |q| which of course is different from ξ 2 ∼ 1/ q 2 . In general, Eq. (35) underestimates the correlation in the one-dimensional case because one-dimensional spin ice has an ordered ice manifold. As we show in future work, it enjoys a unique position in such regard.
The second line in Eq. (63) was recently experimentally verified in square spin ice realized on a quantum annealer by pinning a charge and relaxing the system 81 .

VI. CONCLUSIONS
We have proposed the concept of spin ice on a general graph, and related it to well known graph-theoretical concepts of balance and Eulerian paths. We have then developed a field theoretical framework for its excitations, monopoles or charges. We have obtained a series of results that are independent of geometry. The partition function of general Graph Spin Ice can be reformulated exactly as a functional integral over the charge distribution and its entropic field, the latter subsuming the effect of the underlying spin ensemble on emergent charges [Eqs (11)- (16)]. The high T behavior is described by a quadratic free energy of the average charges and contains informations on the graph via the graph Laplacian [Eq. (60)].
In absence of charge interaction and external fields and in the limit of high T , the entropic interaction among charges corresponds to the Green operator of the graph Laplacian. Thus, correlations correspond to the screened Green operator of the graph [Eqs (33)- (38)]. The correlation length is ξ 2 = /T at high T and ξ 2 = 1/ q 2 at low T , results already appreciated in various special spin ice systems 60,106 .
Beside generalizing condensed matter notions to graphs where they can provide intuitive insight or inspiration for broader problems in complex networks, we have shown that many of the properties of spin ice systems follow directly from the graph structure, which is essentially topological.

VII. ACKNOWLEDGEMENTS
We thank Beatrice Nisoli for proofreading. This work was carried out under the auspices of the U.S. DoE through the Los Alamos National Laboratory, operated by Triad National Security, LLC (Contract No. 892333218NCA000001).