Functional Relations on Anisotropic Potts Models: from Biggs Formula to the Tetrahedron Equation

We explore several types of functional relations on the family of multivariate Tutte polynomials: the Biggs formula and the star-triangle (Y −∆) transformation at the critical point n = 2. We deduce the theorem of Matiyasevich and its inverse from the Biggs formula, and we apply this relation to construct the recursion on the parameter n. We provide two different proofs of the Zamolodchikov tetrahedron equation satisfied by the star-triangle transformation in the case of n = 2 multivariate Tutte polynomial, we extend the latter to the case of valency 2 points and show that the Biggs formula and the startriangle transformation commute.


Introduction
The theory of polynomial invariants of graphs in its current state uses many methods and tools of integrable statistical mechanics. This phenomenon demonstrates the inherent intrusion of mathematical physics methods into topology and combinatorics. In this paper, the main subject of research is functional relations in the family of polynomial invariants for framed graphs, in particular for multivariate Tutte polynomials [22], their specializations for Potts models, multivariate chromatic and flow polynomials.
The flow generating function is closely related to the problems of electrical networks on a graph over a finite field. Each flow defines a discrete harmonic function, and non-zero flows can be interpreted as harmonic functions with a completely non-zero gradient. Specifically, we discuss the full flow polynomial which is a linearization of the flow polynomial and, in particular, corresponds to the point of the compactification of the parameter space for the Biggs model.
One of the central tools of the paper is the Biggs formula (Lemma 2.13), which connects n-Potts models for different parameter values as a convolution with some weight over all edge subgraphs. In particular, we offer a new proof of the theorem of Matiyasevich 2.18 about the connection of a flow and chromatic polynomial, as a special case of the Biggs formula. This interpretation allows us to construct an inverse statement of the theorem of Matiyasevich. Moreover, using the connection between the flow and the complete flow polynomial, we obtain a shift of parameters in the Potts models (Theorem 2.25).
The fundamental type of correspondences on the space of the aforementioned invariants is the star-triangle type relations (also known as "wye-delta" relations) and the associated deletioncontraction relations. In a sense, the kinship between these relations is analogous to the role of the tetrahedron equation in the local Yang-Baxter equation. Despite the fact that the invariance of the Ising model with respect to the star-triangle transformation is very well known [2], we have not found in the literature a full proof of the fact that the action of this transformation on the weights of an anisotropic system is a solution of the tetrahedron equation that corresponds to the orthogonal solution of the local Yang-Baxter equation: Theorem 4.1 (parts of this statement were mentioned in [15,17,21]). We offer here two new proofs of this fact. We find them instructive due to their anticipated relation to the theory of positive orthogonal grassmannians [14].
The identification of the Potts model and the multivariate Tutte polynomial allows us to assert the existence of a critical point for the parameter n in the family of Tutte polynomials. Namely, for n = 2, this model has a groupoid symmetry generated by a family of transformations defined by the trigonometric solution of the Zamolodchikov tetrahedron equation. We extend the star-triangle transformation for the graphs of lower valency in Section 5. In this way, we obtain a 14-term correspondence. This extension commute with the Biggs formula. We should mention the relation of this subject with the theory of cluster algebras. We suppose that the multivariate Tutte polynomial on standard graphs at the critical point n = 2 corresponds to the orthogonal version of the Lusztig variety [4] in the case of the unipotent group and the electrical variety [12] for the symplectic group.

Organization of the paper
In Section 2, we concentrate our attention on the Biggs formalism in the Ising and Potts type models. We define the main recurrence relations and also identify the Tutte polynomial with the Potts model. Then, we apply the Biggs formula to the proof of theorem of Matiyasevich and propose its inverse version. We examine in details the recursion of the Potts model with respect to the parameter n.
In Section 3, we show that, if n = 2, then the Potts model is invariant with respect to the star-triangle transformation given by the orthogonal solution for the local Yang-Baxter equation and the corresponding solution for the Zamolodchikov tetrahedron equation. In Section 4 we provide two different proofs for this fact. Both of them are interesting in the context of cluster variables on the space of Ising models. The first proof operates with the space of boundary measurement matrices and the second with the matrix of boundary partition function.
In Section 5, we show that the Biggs formula considered as a correspondence on the set of multivariate Tutte polynomials commutes with the star-triangle transformation.

Biggs interaction models 2.1 n-Potts models and Tutte polynomial
We define the anisotropic Biggs model (interaction model) on an undirected graph G with the set of edges E and the set of vertices V (a graph can have multiple edges and loops) as follows: a state σ is a map σ : V → R, where R is a commutative ring with the unit, the weight of the state σ is defined by the formula where δ(e) = σ(v)−σ(w); the edge e connects the vertices v and w, the functions i e : R → C are even: ∀b ∈ R : i e (b) = i e (−b), the partition function Z(G) of a model is the following sum where the summation is taken over all possible states σ.
Let us consider the most simple Biggs interaction models: Definition 2.1. If R ∼ = Z n and functions i e given as i e (0) = α e , i e (a) = β e , ∀a = 0 ∈ R, we call such model the anisotropic n-Potts model with the set of parameters α e and β e and we denote it by M (G; i e ).
In addition, if the maps i e = i do not depend on edges, then we call such model the isotropic n-Potts model (or just n-Potts model) with parameters α and β. We denote it by M (G; i), also we use the notation M (G; α, β).
this model can be identified with the classic isotropic Ising model [2]. Therefore we will call any anisotropic or isotropic 2-Potts model just Ising model. Definition 2.3. Consider an anisotropic n-Potts model M (G; i e ), we denote its partition function as Z n (G). In addition, if n = 2, we omit index 2 and write just Z(G).
Remark 2.4. For the empty graph, we define the partition function of any n-Potts model to be equal to 1, and for a disjoint set of m points to be equal to n m . Now we will consider the combinatorial properties of the anisotropic n-Potts models (compare with [3, Theorem 3.2]): Theorem 2.5. Consider an anisotropic n-Potts model M (G; i) and its partition function Z n (G).
Let graph G be the disjoint union of graphs G 1 and G 2 , then Let graph G be the joining of graphs G 1 and G 2 by the vertex v, then Consider a graph G and its edge e, where e is neither a bridge nor a loop. Consider the graph G/e obtained by contraction of e, and the graph G\e obtained by deletion of e. Then the following formula holds Z n (G) = (α e − β e )Z n (G/e) + β e Z n (G\e).

Proof .
1. The statement directly follows from Definition 2.1.
2. Let us rewrite the partition function Z n (G): , therefore for any i = j we have the following identity Hence we obtain Let us introduce the partial partition functions X k := then we could rewrite 3. Let the edge e is neither a bridge nor a loop and denote by X the income in the partition function of all states such that the values of the ends of e coincide and by Y another part of the partition function (that of the distinct values of the ends of e), then and we obtain the statement. Now let us recall the definition of the Tutte polynomial of a graph G.
Definition 2.6. Let us define the Tutte polynomial T G (x, y) by the deletion-contraction recurrence relation: 1. If an edge e is neither a bridge nor a loop, then T G (x, y) = T G\e (x, y) + T G/e (x, y).
2. If the graph G consists of i bridges and j loops, then T G (x, y) = x i y j .
Theorem 2.7 ( [9,23]). Let F (G) be a function of a graph G satisfuing the following conditions: F (G) = aF (G\e) + bF (G/e), if an edge e is not a bridge neither a loop. Then where K 2 is a complete graph on two vertices, L is a loop, r(G) = v(G) − k(G) is a rank of G and c(G) = e(G) − r(G) is a corank. Here and below e(G) is the number of edges in the graph G.
Now we are ready to connect the partition function Z n (G) of the isotropic n-Potts model M (G; α, β) with the Tutte polynomial T G (x, y) of the same graph G using a well-known trick (for instance see [3]). Let us consider the weighted partition function where k(G) is the number of connected components in the graph G. It is easy to verify that the weighted partition function Zn(G) n k(G) satisfies Theorem 2.7, therefore the following theorem holds: Theorem 2.8 (Theorem 3.2 [3]). The partition function Z n (G) of the n-Potts model M (G; α, β) coincides with the Tutte polynomial of a graph G up to a multiplicative factor Example 2.9 (the bad coloring polynomial [9]). Consider a graph G and all possible colorings of V (G) in n colors. Define the bad coloring polynomial as here b j (G, n) is the number of colorings such that each of them has exactly j bad edges (we call an edge "bad" if its ends have the same colors). So, easy to see that B G (n, t) = Z n (G), here Z n (G) is the partition function of the n-Potts model M (G; t, 1). Hence, using the Theorem 2.8 we immediately obtain B G (n, t + 1) = n k(G) t r(G) T G t + n t , t + 1 .

n-Potts models and the theorem of Matiyasevich
The connection between the n-Potts models and Tutte polynomials allows us to give a simple proof of the theorem of Matiyasevich about the chromatic and flow polynomials, but at first we introduce a few definitions. Definition 2.11. A graph A is called an edge induced subgraph ( Figure 2) of a graph G, if A is induced by a subset of the set E(G). Every edge induced subgraph A of a graph G could be completed to the spanning subgraph A by adding all the vertices of G which is not contained in the subgraph A.
Definition 2.12. For a n-Potts model M (G; i) we introduce the normalized partition function as follows The edge induced subgraph The spanning subgraph We start with the following lemma, which is a generalization of the high temperature formula for the Ising model: Lemma 2.13 (Biggs formula [5]). Let us consider two n-Potts models M 1 (G; i 1 ) with parameters α 1 , β 1 and M 2 (G; i 2 ) with parameters α 2 , β 2 . Then the normalized partition function Z 1 n (G) of the first model could be expressed in terms of the normalized partition functions of the models of all edge induced subgraphs of the second model: where p = α 1 −β 1 α 2 −β 2 , and q = α 2 β 1 −α 1 β 2 we assume that Z i n (∅) = 1 . Proof . Let us notice that i 1 = p · i 2 + q, therefore i 2 (δ(e)).
In order to complete the proof we consider the following term for a fixed A: We consider further the chromatic and flow polynomials, first of all remain some well-known definitions.
is the set of edges each of them is directed to (respectively from) v.
Next, we formulate one of the classic results of graph theory which can be found for instance in [9]: Theorem 2.17. The number of proper colorings of a graph G in n colors is the following polynomial (called chromatic polynomial) in the variable n: The number of nowhere-zero n-flows of a graph G is independent on the choice of orientation and is obtained by the following polynomial (called flow polynomial) in the variable n: Now we are ready to formulate and prove the theorem of Matiyasevich: [20]). Let us consider a graph G, its chromatic polynomial χ G and its flow polynomial C G , then where the summation goes through all spanning subgraphs A.
Proof . Let us consider two n-Potts models with the special parameters: the model M 1 (G; i 1 ) with the parameters α 1 = 0, β 1 = 1 and the model M 2 (G; i 2 ) with the parameters α 2 = 1 − n, β 2 = 1. By Theorem 2.8 we could express the partition function of the first model in terms of the chromatic polynomial So we have Analogously, we express the partition function of the second model in terms of the flow polynomial Then by Lemma 2.13 after the substitutions (2.2) and (2.3) we obtain where the summation goes through all edge induced subgraphs A .
We finish the proof by noticing that the edge induced subgraph differs from the spanning subgraph by the set of isolated vertices. Therefore we can complete each edge induced subgraph to its corresponding spanning subgraph and then replace the summation over all edge induced subgraph by the summation over all spanning subgraph, because the value of the each flow polynomial C A remains the same and finally we obtain Note that we could produce series of statements that look like Theorem 2.18: Theorem 2.19. Let us consider a graph G, then we can obtain the following formulas
The specification of the first case: M 1 (G; 1 − n, 1) and the same M 2 (G; 0, 1) for the proof of the formula (2.5).
Now it is left to repeat step by step the proof of Theorem 2.18 for these two models.
Remark 2.20. We notice that the formula (2.5) naturally can be considered as "inversion" of Theorem 2.18.

Shifting the order in the Potts models
Biggs Lemma 2.13 allows us to relate the values of the partition functions of the n-Potts models with fixed n, but different values of parameters α and β. The goal of the current subsection is to present a method for connecting partition functions of the n-Potts models for different n.
We will call it shifting order formulas.
The first method is based on the multiplicativity property of the complete flow polynomial.
Definition 2.21. Let G be a graph with the edge set E(G) and the vertex set V (G), let us chose a fixed edge orientation on G. Then, a function f : E → Z n is called an n-flow if the following condition holds ∀v ∈ V (G) : here again M + (v) (respectively M − (v)) is the set of edges each of them is directed to (respectively from) v.
Let us formulate a few well known results concerning a flow polynomial and a number of all n-flows. The proofs could be found for example in [22]. Proposition 2.22. Denote the number of all n-flows on a graph G by F C G (n), then F C G (n) is independent of the choice of an orientation and the following identity holds where the summation goes through all spanning subgraphs A of the graph G.
Proposition 2.23. The number of all n-flows on a graph is the following polynomial (called complete flow polynomial) where e(G), v(G), k(G) are numbers of edges, vertices and connected components in the graph G correspondingly.
Proposition 2.24. The flow polynomial C G (n) of a graph G could be expressed in terms of the complete flow polynomials of its spanning subgraphs by the following identity: The complete flow polynomial F C G (n) is a multiplicative invariant: F C G (n 1 n 2 ) = F C G (n 1 ) × F C G (n 2 ), therefore we are ready to formulate the following theorem: Theorem 2.25. The partition function Z n 1 n 2 (G) of the n 1 n 2 -Potts model M (G; α 1 , β 1 ) could be expressed in terms of the partition functions Z n 1 (A) and Z n 2 (A) of the n 1 -Potts model M 1 (A; α 1 , β 1 ) and n 2 -Potts model M 2 (A; α 1 , β 1 ) of all spanning subgraphs A of the graph G correspondingly.
Proof . Indeed, by Theorem 2.8 and the formula (2.6) we have notice that we used for the second resummations the following simple observation: if X is a spanning subgraph of Y , which is a spanning subgraph of graph Z, so X is a spanning subgraph of a graph Z. We omit this remark below. The Proposition 2.22 implies Finally, with the help of formula (2.7) and Theorem 2.8 we obtain where η A A are some constants, appeared after the resummations.
Remark 2.26 (convolution formula [16]). It seems extremely interesting and fruitful to compare Lemma 2.13 and Theorem 2.25 with the convolution formula here the summation is over all possible subsets of E(G), here G|A is a graph obtained by the restriction of G on the edge subset A and G/A is a graph obtained from G by the contraction of all edges from A (see [9] for more details).
Our second method is based on the Tutte identity for the chromatic polynomial: . Consider a graph G with the set of edge V (G) then the following formula holds Using this fact we can formulate the following theorem: Theorem 2.28. The partition function Z n 1 +n 2 (G) of the n 1 + n 2 -Potts model M (G; α 1 , β 1 ) could be expressed in terms of the partition functions Z n 1 (A) and Z n 2 (A) of the n 1 -Potts model M 1 (A; α 1 , β 1 ) and n 2 -Potts model M 2 (A; α 1 , β 1 ) of all spanning subgraphs A of the graph G correspondingly.
Proof . The proof is very similar to the proof of Theorem 2.25. Again, from Theorem 2.8 and the formula (2.4) we have Let us complete each subgraph A 1 (each A 2 ) to the corresponding spanning subgraph of G by adding isolating vertices where η A A are again some constants, appeared after the resummations.
3 Star-triangle equation for Ising and Potts models

General properties
Let us rewrite the partition function of the anisotropic n-Potts models in the so called Fortuin-Kasteleyn representation: Proposition 3.1 (compare with the formula (2.7) from [22]). Consider the anisotropic n-Potts model M (G; i e ), then its partition function could be expressed as follows

1)
where δ(σ e ) is a value of standard Kronecker delta function of the values of σ on the boundary vertices of the edge e and t e = αe βe is a reduced weight of the edge e. Proof . Indeed, it is easy to see that if σ(v) = σ(w): i e (δ(e)) = i e (σ(v) − σ(w)) = α e = β e + (α e − β e )δ(σ e ) = β e + (α e − β e )δ(σ(v), σ(w)), Also, we introduce the boundary partition function of the n-Potts models: Definition 3.2. Let G be a graph (with possible loops and multiple edges) with the set of vertices V , the set of edges E and the boundary subset S ⊆ V of enumerated vertices: S = {v 1 , v 2 , . . . , v k }. The boundary partition function on G is defined by the following expression where A = {a 1 , a 2 , . . . , a k }, ∀i : a i ∈ Z n is the set of fixed values, and the summation is over such states σ A that σ A (v i ) = a i . Remark 3.3. If n = 2, we will omit the index 2 and will write just Z S(A) (G).
The next Lemma connects boundary and ordinary partition functions: Lemma 3.4. Consider two graphs G 1 = (V 1 , E 1 ) and G 2 = (V 2 , E 2 ) with the only common vertices in the boundary subset S = {v 1 , v 2 , . . . , v n } in V 1 and V 2 . We can glue these graphs and obtain the third graph G Then, the following identity holds where the summation is over all possible sets A.

The case n = 2
In this subsection we consider the case n = 2. Our first goal is to find such conditions that the partition function (3.1) is invariant under the star-triangle transformation which changes the subgraph Ω to the subgraph Ω . We derive these conditions with the use of the boundary partition functions: consider a graph G with the subgraph Ω, then using Lemma 3.4 for graphs Ω and G − Ω we obtain the following identity (Figure 4). After the star-triangle transformation, we obtain a graph G with the following partition function We write them down in detail. Let us note that these conditions do not depend on the states of the vertices (see Figure 5), but depend on the number and the positions of the vertices with equal states. Therefore, we have the following possibilities: two states in the triangle are the same, then the central vertex either has the same state or has the different state, then α 1 β 2 β 3 + β 1 α 2 α 3 → α 1 β 2 β 3 and two more maps after permuting indexes, all states are the same, then α 1 α 2 α 3 + β 1 β 2 β 3 → α 1 α 2 α 3 .
In this way we obtain the following equations After the substitution This set of equations defines a correspondence which preserves the Ising model partition function if we mutate the graph G to G . Let us denote the product β 1 β 2 β 3 by β and the product β 1 β 2 β 3 by β . Then, we obtain the following map from the (t, β)-variables to the (t , β ) variables, we will call it F , F (t 1 , t 2 , t 3 , β) = (t 1 , t 2 , t 3 , β ) : Remark 3.6. Formally speaking to define a map on the space of edge weight adopted to the star-triangle transformation we have to resolve the map F somehow for the parameters β i . For example one can take the following one

Actually, the choice of a resolution is not important in what follows.
Remark 3.7. We choose the positive branch of the root function for real positive values of variables t i for purposes emphasized further. This is relevant to the almost positive version of the orthogonal grassmanian. See [11] for more details about the connection of the Ising model and positive orthogonal grassmanian.

The case n = 2
Let us demonstrate how the method, described above, works for the star-triangle transformation in the case n ≥ 3. Using the same ideas as in the previous subsection we could obtain the following conditions In general, the system (3.5) does not have a solution and the star-triangle transformation is not possible. But, if t i satisfy the special condition, partition function of the n-Potts model is still invariant under the star-triangle transformation.
Proposition 3.8. The system (3.5) together with equation has a solution in terms of prime variables.
Proof . Using the first three and the last equations of (3.5) we immediately obtain the expressions for t i and Substitute these expressions into the fourth equation of (3.5) and obtain the equation By the straightforward computation, we retrieve that the identity (3.6) is the consequence from the equation (3.7). Below we present two nontrivial specialization of the partition function of n-Potts model which are agreed with the system (3.5), (3.6).
Example 3.10. Consider a graph G and equip each e ∈ E(G) with sign + or −. Let us consider the n-Potts model M k (G, α e , β e ) with following parameters: for all e ∈ E equipped with + the parameters α e , β e equal α e = A + = −t − 3 4 , β e = B + = t and n = t + 1 t + 2 (we suppose that parameter t is chosen such that n ∈ N).
Let the graph G has a triangle subgraph, the edges of which have signs +, −, −. The reduced weights of edges are It is easy to see that these t i satisfy the equation (3.6).
We notice that the signed graph G could be considered as the signed Tait graph for a diagram D(K) of a knot K ( [24], the chapter "Knot invariants from edge-interaction models"). Moreover, the value of the Jones polynomial of the knot K at the point n is closely related with the partition function of the n-Potts model M k (G, α e , β e ) (see [24, equation (7.17)]). Thus, the identification of the third Reidemeister move of the diagram D(K) with the star-triangle transformation of the signed graph G is agreed with the star-triangle transformation defined by the system (3.5), (3.6) for the n-Potts model M k (G, α e , β e ).
Our second example is about the models of bond percolation. Firstly, we briefly give their definitions: Definition 3.11 (bond percolation [13]). Consider a graph G. An edge e ∈ E(G) is considered to be open with probability p e or closed with probability 1 − p e . We suppose that all edges might be closed or open independently. One is interested in probabilistic properties of cluster formation (i.e. maximal connected sets of closed edges of the graph G).
Example 3.12. The bond percolation models could be considered as a limit n → 1 of the n-Potts models at the level of the boundary partition functions [7]. This identification corresponds to the specialization of the system (3.5), (3.6) by n → 1.
Substitute t i = 1 p i , t i = 1 p i and n = 1 in (3.5) and (3.6), then After simplifications we obtain the condition for the star-triangle transformation of the bond percolation models (for instance, see [13])

Tetrahedron equation
The tetrahedron equation firstly was considered by A. Zamolodchikov [25] who has constructed its solution in S-form. We consider the following form of the equation where T ijk is an operator acting nontrivially in the tensor product of three vector spaces V i , V j , V k , indexed by i, j and k. Tetrahedron equation is the higher order analog of the Yang-Baxter equation. Both equations are examples of n-simplex equations [18] and play an important role in hypercube combinatorics and higher Bruhat orders. For the complete introduction to the topic see for example [21]. In this section we present two proofs of the main theorem of the paper: These two proofs have a lot of common points and ideas, but have the crucial differences in the last stages. It is interesting to compare proofs for the purpose of combining arguments of boundary partition functions and the technique of correlation functions in the Ising-Potts models.
At first in the next subsection we prove that the change of variables (3.4) corresponds to the variables transform in the trigonometric solution of a local Yang-Baxter equation.

Local Yang-Baxter equation
Let us recall that the following change of variables (t 1 , t 2 , t 3 ) → (t 1 , t 2 , t 3 ) provides an invariance of the Ising model (3.1) under the star-triangle transformation (3.4): (4.2) Following [17], we construct orthogonal hyperbolic 3 × 3 matrices R ij which solve the local Yang-Baxter equation where S(t) is the following involution On the left hand side of (4.3) we have i sinh(log(S(t 2 ))) 0 cosh(log(S(t 2 ))) 0 1 0 cosh(log(S(t 2 ))) 0 −i sinh(log(S(t 2 )))   , (4.6) Proof . It can be proved by a straightforward computation. For example let us write down the result of the product on the left hand side At the first glance the product on the right hand side looks much more cumbersome, but occasionally all terms are simplified and the matrix on the right hand side coincides with (4.8).

Tetrahedron equation, first proof
Let us encode the tetrahedron equation by the Figure 7. The standard graph encodes R-matrices in the following way: in each inner vertex numbered by k, which is the intersection of strands i and j we put the matrix R ij (t k ) which is the 2 × 2 matrix in the 4-dimensional space with basis vectors indexed by a, b, c, d. For instance, R ac (t 5 ) =     i sinh(log(t 5 )) 0 cosh(log(t 5 )) 0 0 1 0 0 cosh(log(t 5 )) 0 −i sinh(log(t 5 )) 0 Let us orient each strand from the left to the right and multiply R-matrices in order of the orientation, for instance for the Figure 7 we have the following product of R-matrices R cd (t 1 )R bd (S(t 2 ))R bc (t 3 )R ad (t 4 )R ac (S(t 5 ))R ab (t 6 ). (4.9) We note that the orientation defines the product (4.9) uniquely. Then let us apply four local Yang-Baxter equations consequently to the inner triangles with vertices numbered (1, 2, 3), (1,4,5), (2,4,6) and (3,5,6) as on the Figure 8. As a result, we have one and the same standard graph as on the Figure 7  At the same time we could apply local Yang-Baxter equations in the opposite direction: firstly to the triangle (3, 5, 6), then (2,4,6), (1,4,5) and (1,2,3). Eventually in this case we will have again the same standard graph.
By Theorem 4.2 the product (4.9) preserves by each local Yang-Baxter equation encoded on the Figure 8. As a result of two sequences of four Local Yang-Baxter equations we obtain an equality of two products of six 4 × 4 R-matrices where the parameters u i , v j , i, j = 1, . . . , 6 depend on the initial variables t i , on the mapping (3.4) and on the involution (4.4).
Let us consider this equation element-wise, and note that we could uniquely express parameters in the right hand side in terms of the parameters on the left hand side , Here U is a matrix on the left hand side and a, b are some invertible functions, come from (4.5), (4.6), (4.7). So we could uniquely determine t 4 from the first equation, t 2 and t 5 from the second and the third, then t 1 and t 6 , and finally t 3 from the element U 2,3 . Let us note that this algebraic proof could be formulated in terms of the paths on the standard graph (Figure 7) with orientation. So the equation (4.10) provides coincidence of the parameters in the vertices given by the two sides of the tetrahedron equation. This finishes the proof.

Involution lemma
Let us consider the map F (t 1 , t 2 , t 3 ) = (t 1 , t 2 , t 3 ), where prime variables defined by (3.4). First of all, we formulate one technical lemma: Lemma 4.3. The following identity holds for all t 1 , t 2 , t 3 : We present the proof in Appendix A.

Towards the tetrahedron equation
Let us consider any graph G with a subgraph Γ 1 which coincides with the leftmost graph on the Figure 9. We can transform the graph G to the graph G with the subgraph Γ 2 which coincides with the rightmost graph on the Figure 9. We could make this mutation by two different chains of star-triangle transformations: F −1 356 F 246 F −1 145 F 123 and F −1 123 F 145 F −1 246 F 356 . Both are figured out on the Figure 8. This observation turns us to the following hypothesis This equality is equivalent to the Zamolodchikov equation (4.12) where Φ ijk = S i S k F ijk S j . Indeed, using Lemma 4.3 and the simple observation that S l F ijk = F ijk S l , l = {i, j, k} we can write down Conjugating both sides of (4.11) by S 2 S 5 we obtain the Zamolodchikov tetrahedron equation.  (1 + (t e − 1)δ(σ e )).

Solution for the tetrahedron equation
The proof of the second part of the statement is straightforward. It follows from the Definition 3.2 and the condition (3.2).
We will prove the Zamolodchikov equation in its equivalent form Let us notice that due to the local nature of the star-triangle transformation and the convolution property of the boundary partition function we have a choice to take some suitable graph to Figure 10. The graphs Γ 1 and Γ 2 .
prove equation (4.13). So let us take the graph Γ 1 from the Figure 10 with the following choice of boundary set S 0 : We will prove that the values of the second-type invariant functions which are preserved by both sides of the equation (4.13) allows us to uniquely reconstruct weights of all edges. Explain this idea in detail, let us consider the left hand side of the equation (4.13) and the map F 123 , then for any A = {a 1 , . . . , a 4 } the following identity holds Figure 11) where Γ is obtained from Γ 1 by the star-triangle transformation (see Figure 11). And therefore we deduce that Repeating these arguments for the remaining maps F ijk from the left hand side of (4.13) we obtain (4.14) In the same fashion, if we consider the right hand side of the equation (4.13), we similarly obtain that Hence, in order to prove the equation (4.13) it is sufficient to prove that we can reconstruct the parameters t i , i = 1, . . . , 6 from the values Z S 0 (A) (Γ 2 )/Z(Γ 2 ) for different values of A in a unique way. We understand the identity (4.14) as a system of 2 4 linear equations with unknowns Z S 0 (A) (Γ 2 ) of the following type F 246 Figure 11. The left hand side of (4.13) . Figure 12. The graph Γ 2 .
which is equivalent to The rank of the system is equal to 15. Indeed, the rank is ≥ 15 and we know that there is a nontrivial solution coming from the boundary partition functions for the graph Γ 2 .
Hence any solution has the form a 2 , a 3 , a 4 ), (4.15) where C is some constant and a i are the states. Now we will prove that the parameters t 1 , . . . , t 6 are reconstructed uniquely from the equation (4.15). Let us introduce some auxiliary variables and rewrite the partition function in the following way: we have 16 states of boundary vertices S 0 = {v 1 , v 2 , v 3 , v 4 }. Each expression Z S 0 (A) (Γ 2 ) is a sum of four terms corresponding to the states of internal vertices v 5 and v 6 . We consider in details the case S 0 = {0, 0, 0, 0}. Let us denote the weights of the states of the square {v 1 , v 6 , v 5 , v 4 } by v, z, y and x ( Figure 13). Then, we obtain the following equations where B = β 3 β 4 β 5 β 6 and v + t 1 t 2 z + t 2 y + t 1 x = C B 1 α 0 (0, 0, 0, 0), B 1 = β 1 β 2 .  Similarly we obtain seven more equations (we omit eighth equation with v 1 = 1 due to the symmetry of the model with respect to the total involution of spins) By straightforward calculation we retrieve Using the auxiliary variables it is easy to see that In the same fashion we obtain v 1 Then we can straightforwardly deduce expressions for the variables t 1 and t 2 (equations (B.1), (B.2) in Appendix B), then obtain expressions for the auxiliary variables from equations (4.17) and finally obtain variables t 3 , t 4 , t 5 , t 6 from the equations (4.16): This completes the proof of the Zamolodchikov equation due to the fact that there is a unique way to choose positive weights for the edges of the model to provide the expected values of boundary partitions function for the graph Γ 2 .

Star-triangle transformation, Biggs formula and conclusion
The main results of the paper represent the functional relations on the space of multivariate Tutte polynomials. This problem is a step of the program of investigation of the framed graph structures and the related statistical models. We examined in details the Biggs formula and applied it to the multivariate case. We also provided a new proof of the theorem of Matiyasevich as a partial case of such formula. The second principal result is the reveal of the tetrahedral symmetry of the multivariate Tutte polynomial at the point n = 2. Therefore, we have a connection between the multivariate Tutte polynomial, functions on Lustig cluster manifolds [4] and its electrical analogues [12,19]. We would like to interpret this property as the critical point of the model described by the multivariate Tutte polynomial, and the tetrahedral symmetry as a longstanding analog of the conformal symmetry of the Ising model at the critical point [8].
Both correspondences are related by the following observation. Let G and G be two graphs related by the star-triangle transformation. Let us consider the case when the partition function is invariant with respect to this transformation (n = 2 or n > 2 and the system (3.5), (3.6) holds) Z n (G ) = Z n (G).
On the other hand the star-triangle transformation provides a groupoid symmetry on a wide class of objects, in our case on the space of Ising models. The Biggs formula allows us to extend this action to the points of valency 1 and 2. And we can obtain the 14-term relation (Theorem 5.1) by comparing the right-hand sides of Biggs formulas for G and G .
Let us explain this idea in details, consider two pairs of n-Potts models: M 1 (G, i 1 e ) and M 1 (G , i 1 e ), M 2 (G, i 2 e ) and M 2 (G , i 2 e ) (for simplicity we denote these models M 1 (G), M 2 (G), M 1 (G ), M 2 (G ) correspondingly). After multiplying both parts of the formula (2.1) for M 1 (G) and M 2 (G) by n v(G) (by n v(G ) for M 1 (G ) and M 2 (G )) we obtain  Figure 14. The 14-term relation.
Proof . We will prove this theorem using induction on ex(G) := e(G) − 3. We show that the base of the induction k = 0 is trivial. Hence, let us consider the n-Potts models M 2 (G) and M 2 (G ) and the special models M 1 (G) and M 1 (G ) such that p 1 p 2 p 3 = p 1 p 2 p 3 . Then, we write the formulas (5.2) and (5.3), after the comparison for each terms using the reasoning above we immediately obtain the result in the case ex(G) = 0.
Then, make the step of induction. Again, let us write down the formulas (5.2) and (5.3): where each A 1 contains the full triangle, each A 2 contains only a part of the triangle and such that ex(A 2 ) = ex(G), and by S 1 we denoted the left hand side of (5.4), where each A 1 contains the full star, each A 2 contains only a part of the star and such that ex(A 2 ) = ex(G ), and by S 2 we denoted the right hand side of (5.4). Then the induction assumption ends the proof.
We consider these results in the context of numerous generalizations, both for other models of statistical physics, and in a purely mathematical direction. In particular, we are interested in applying this technique to the Potts model in the presence of an external magnetic field [10], including an inhomogeneous one. In addition, we are going to develop these methods in a more general algebraic sense, in particular in a non-commutative situation. Partial results of this activity have already been obtained in [6].

A The proof of Lemma 4.3
Proof . We start by reformulating this statement in terms of equivalent rational identities. Let us introduce the x-variables by the following formula