Kazakov-Migdal model on the Graph and Ihara Zeta Function

We propose the Kazakov-Migdal model on graphs and show that, when the parameters of this model are appropriately tuned, the partition function is represented by the unitary matrix integral of an extended Ihara zeta function, which has a series expansion by all non-collapsing Wilson loops with their lengths as weights. The partition function of the model is expressed in two different ways according to the order of integration. A specific unitary matrix integral can be performed at any finite $N$ thanks to this duality. We exactly evaluate the partition function of the parameter-tuned Kazakov-Migdal model on an arbitrary graph in the large $N$ limit and show that it is expressed by the infinite product of the Ihara zeta functions of the graph.


Introduction
The large N limit in gauge theory has so far led to many insights not only in gauge theory but also in string theory and has resulted in AdS/CFT correspondence.In this connection, for various types of matrix models, analysis in the large N approximation has yielded important results in a wide range of physics and mathematics.The matrix models have appeared from various points of view such as toy models of QCD, large N reduction of gauge theories, the exact partition functions of supersymmetric gauge theories, and more.
Kazakov and Migdal proposed a lattice gauge model, which induces the effective gauge field action including multiple Wilson loops [1].The Kazakov-Migdal (KM) model consists of a scalar field in the adjoint representation on sites and gauge field (unitary matrix) on links of the lattice.The action is inspired by the well known Harish-Chandra-Itzykson-Zuber (HCIZ) integral over a unitary group with two Hermite matrices [2,3].In the large N limit, the model is studied by using contemporary matrix model techniques.
Although the original KM model is defined on a square lattice in various dimensions, we can extend it to a model on a generic graph, which describes a kind of discretization of space-time as we show in this paper.By adjusting the coupling constants in the extended KM model on the graph, the partition function of the model reduces to a generating function of the "non-collapsing" Wilson loops corresponding to reduced cycles in graph theory, thanks to the important results of the Ihara zeta function [4][5][6] (see also [7] for a recent review).
The Ihara zeta function associated with a finite graph is an analog of the Riemann zeta function, which plays very important role in the number theory.The Ihara zeta function and its related theorems are very beautiful and meaningful and have been applied mainly to mathematics, but have not been applied much in high energy physics, except for applications to quiver gauge theory [8][9][10] and discretized (supersymmetric) gauge theory [11][12][13][14][15][16][17].
In this paper, we give an essential relationship between the Kazakov-Migdal model on the graph and the Ihara zeta function.We show that, by appropriately tuning the parameters of the model, the partition function is expressed as the unitary matrix integral of an extended Ihara zeta function, which has a series expansion by all non-collapsing Wilson loops with their lengths as weights.We will show that, in the large N limit, the unitary matrices (and also the scalar fields) can be integrated out completely and the partition function is given by an infinite product of the Ihara zeta function.This means that we can exactly solve the KM model on the graph in the large N limit.
The organization of this paper is as follows: In the next section, we give a brief review on the Ihara zeta function for unfamiliar readers.It contains definitions of various terminologies of graph theory.For later discussions, we also extend the Ihara zeta function to include matrices on the edges.In Sect.3, we define the KM model on a generic graph.We show that the partition function of the graph KM model is expressed by the integral of the extended Ihara zeta function over the unitary matrices, after integrating the scalar field first.This model has a different picture depending on the order of integration.If we integrate the unitary matrices first by using the HCIZ integral formula, we obtain a dual integral expression over the eigenvalues of the Hermitian scalar fields.As a byproduct, we can obtain the exact value of the integral dU exp ∞ n=1 q n n | Tr U n | 2 at finite N by using the dual description in cycle graphs.In Sect.4, we perform the integral of the extended Ihara zeta function over the unitary matrices.For finite N , this integral is generally difficult except for the cycle graphs.However, in the large N limit, due to the cluster decomposing nature of the expectation value of the Wilson loops, the integrals can be performed exactly.Our main result is that the partition function of our model is given by an infinite product of the Ihara zeta function.The last section is devoted to the conclusion and discussion.In Appendix A, we show some examples of the Ihara zeta function.In Appendix B, we give a proof of identities associated with symmetric polynomials used in this paper.In Appendix C, we directly carry out the unitary matrix integral appearing in Sect. 3 for small N using Cauchy's integral formula.

Graph Theory and the Ihara Zeta Function
In this section, we give definitions of some objects appeared in graph theory and a brief review of the Ihara zeta function.We then propose an extension of the Ihara zeta function using matrices.

Ihara zeta function
A graph G consists of vertices connecting with edges.We denote the set of vertices and edges as V and E, respectively.The numbers of the vertices and edges (the numbers of the elements in V and E) are also denoted by n V = |V | and n E = |E|.We here assume the graph is connected and directed.So each edge has the direction.A directed edge e = u, v connecting the vertices u, v ∈ V has the "source" (starting point) u = s(e) and "target" (end point) v = t(e).e −1 = v, u is called the inverse of e.We denote the set of the edges with their inverses as In other words, E D is the direct sum of the set of edges E and the set of the inverse of edges E −1 ; The Ihara zeta function corresponding to a graph G is defined as [4,6] where q ∈ C is assumed to have sufficiently small |q|, and [C] runs the equivalence classes of the primitive cycles on the graph G.
Noting that, if C is a primitive cycle, C −1 is also a primitive cycle with the same length as C, a set of the equivalence classes of the primitive cycles of length can be decomposed 1 into Π + Π − , where Π − is the set of the inverse of the elements in Π + such that Π − ≡ C −1 |C ∈ Π + .Since these sets Π + and Π − have the same number of the elements by definition; |Π + | = |Π − |, the Ihara zeta function (2.1) can also be rewritten as In the following, we call the element of Π + as the equivalence class of chiral primitive cycles of length .
For later use, we give another expression of the Ihara zeta function.By using (1−x) Since any reduced cycle is expressed as a positive power of a primitive cycle C and [C n ] (n ∈ N) is consist of (C) different reduced cycles, 2 |Π + | is the number of the reduced cycles of length .Therefore we can regard the Ihara zeta function (2.1) as the generating function of the number of the reduced cycles: where N k is the number of the reduced cycles of length k.

Ihara zeta function via the edge adjacency matrix
It is remarkable that the Ihara zeta function is expressed as an inverse of a polynomial of q.To show it, we define a 2n E × 2n E matrix W called the edge adjacency matrix of G as [18,19] where e, e ∈ E D .
It is obvious from the construction that the element (W k ) ee is the number of the paths of length k whose source and target are s(e) and s(e ), respectively, but the last edge is not e −1 .In particular, the diagonal element of a power of W , (W k ) ee , is the number of the reduced cycles of length n which start from e. Therefore, we see Comparing this expression to (2.3) we obtain [18][19][20] ζ which is the inverse of a polynomial of q as announced.

Ihara's theorem
We define an n V × n V matrix A called the adjacency matrix of G as where v, v ∈ V , which is the number of v, v in E D .We define the degree of the vertex v as which is the summation of the number of edges which start from or end on the vertex v.We then define the degree matrix D of the graph G as In [4], it is shown that the Ihara zeta function is expressed by a rational function, This result is called Ihara's theorem.Some examples of Ihara's theorem are given in Appendix A. Instead of proving this theorem, we will extend the Ihara zeta function and prove a similar theorem with respect to the extended function.

Extended Ihara zeta function
Suppose X e (e ∈ E) are invertible matrices of size K living on each edge.Using these matrices, we extend the Ihara zeta function as where W X is an extension of the edge adjacency matrix (2.4) whose elements are matrices of size K: (2.9) Note that (2.8) reduces to the weighted Ihara zeta function discussed in [21] when K = 1 (weighted by complex scalar numbers).From the construction, the diagonal element of n's power of W n U , (W n U ) ee , is a sum of the product of X e along a reduced cycle C of length n which is denoted by P C (X).For example, if C = e 1 e −1 2 e 3 , P C (X) = X e1 X −1 e2 X e3 .Therefore, we can rewrite the extended Ihara zeta function (2.8) as where R k is the set of reduced cycles of length k.
In the following, we show that the extended Ihara zeta function (2.8) can be expressed as where A X is an extension of the adjacency matrix (2.6) defined by whose elements are matrices of size K as well as W X .Since this relation holds for K = 1 and X = 1, the original Ihara's theorem is a corollary.The following proof is parallel to the one given in [20] for the original Ihara zeta function: Firstly, we define matrices S X and T X of size n E × n V whose elements are matrices of size K; We also define S ≡ (S X )| X=1 and T ≡ (T X )| X=1 .Then we can easily show Using these matrices, we define After a straightforward calculation, we obtain where and then conclude (2.11).

Graph Kazakov-Migdal model
In this section, we propose a matrix integral that adds up all the Wilson loops on the graph with their lengths as weights.This integral is obtained by appropriately adjusting the parameters of the KM model on the graph.As an application, a certain integral over a unitary matrix of finite size is evaluated exactly through elementary linear algebra.

Kazakov-Migdal model on the graph
The KM model [1] is a model defined on the D-dimensional square lattice: where Φ(x) is a scalar field on the site x and U µ (x) is the link variable connecting from x to x + µ.By integrating out the scalar fields, the effective theory becomes a summation of all possible Wilson loops.Therefore this model is expected to induce QCD in an arbitrary dimension.
In this paper, we would like to generalize the square lattice on which the KM model is defined to an arbitrary discretized space-time, which is described by a graph in general.Let us consider a directed graph G and suppose that a Hermitian matrix Φ v and a unitary matrix U e ∈ U (N ) live on each vertex v ∈ V and each edge e ∈ E, respectively.Then we can naturally generalize the action of the KM model (3.1) to where β, m 2 v and q are parameters of the model.We call this model as the graph KM (gKM) model in the following.Note that we can recover the original KM model by choosing the square lattice as the graph structure and setting the parameters as β = N , m 2 v = 2m 2 0 (common for all the vertices) and q = 1.Apart from the induced QCD, this model can also be regarded as a model of a scalar field coupled to the gauge field on a discretized space-time (graph).To see this, we define the matrices, where S U and T U are defined by (2.13) with X e = U e ⊗ U † e , which are regarded as q-deformed covariant incidence matrices associated with forward and backward difference operator on the graph, respectively.Indeed, if we take q = 1 and N = 1 (U e = 1), L q U is nothing but the incidence matrix, and if we take q = 1 only L q U reduces to the covariant incidence matrix discussed in [17].In [17], it is also shown that the covariant incidence matrix plays a role of the first order difference operator and closely related to the Dirac operator on the graph, which reproduces the index theorem on the graph.
Using the q-deformed covariant incidence matrices, the action (3.2) can be written by where ∆ q U without U e is called q-deformed graph Laplacian, which becomes a square of the incidence matrix in the q → 1 limit.The graph Laplacian is a second order difference operator, which is a discrete analog to the Laplace differential operator on the continuum space-time.As discussed in [1], the action (3.4) would correspond to the gauged scalar field theory without an explicit Yang-Mills term on a curved space-time; where ∆ U ≡ D µ D µ is the gauge covariant Laplace operator and R(x) is a scalar function depending on the position x, in a suitable continuum and q → 1 limit.
The advantage of being able to define a theory on an arbitrary graph is not only the ability to discretize a curved space-time.Even in considering the flat space-time, the advantage of having a various kind of discretizations is significant.For example, the flat space-time in two dimensions can be discretized by graphs with a periodic structure as shown in Fig. 1.Of course, the original KM model is included in these, but this extension allows the theory to have greater discrete symmetry.This will be also useful when considering applications to condensed matter physics in particular.

Parameter tuning
Now let us discuss the relationship between the graph Kazakon-Migdal model and the Ihara zeta function on the graph.First of all, since the integrand is Gaussian with respect to Φ v , we can integrate over Φ v of the partition function then it reduces explicitly to where M is a diagonal matrix given by M = m 2 v δ vv ⊗ 1 N 2 and A U is the extended adjacency matrix given by (2.12) with X e = U e ⊗ U † e .Here, motivated from the expression (2.11), if we set then the partition function (3.8) can be expressed through the extended Ihara zeta function; We can eliminate the q-dependence in the overall factor by further setting with a q-independent parameter α.Note that if we use the q-deformed graph Laplacian, the modified mass in (3.5) becomes which represents massless field when deg v = 2 for all v or q = 1.
Since the extended Ihara zeta function can be expressed as (2.10), we see that the partition function of the gKM model (3.7) with the parametrization (3.9) and (3.10) becomes Since Tr P C (U ) is the Wilson loop along a reduced cycle C, this matrix integral is regarded as the generating functional of the integral of the Wilson loops on the graph of constant length 2 .We can of course express Z(q) by using the extended Ihara zeta function; Note that the specific parameter settings (3.9) and (3.10) are essential for the partition function of the gKM model to count the Wilson loops without over-or under-counting.As an example, let us consider the case when β = const.and m 2 v is a constant independent of v, as in the original KM model.In this case, the partition function becomes 2 This kind of the unitary matrix integral also appears in the context of a generalization of the two-dimensional lattice gauge theory [22,23] or the thermodynamics of (supersymmetric) gauge theory [24][25][26].
The point is that Tr A k U includes not only the Wilson loops of length n but also some extra constant coming from the backtrackings.A typical example is the triangle graph whose (extended) adjacency matrix is given by Although it is easy to see Tr A 2 U = 6N 2 , it is obvious that the triangle graph does not have any Wilson loop of length 2; it does not count Wilson loops but rather counts the backtrackings.Of course, if we impose a general q-dependence on m 2 v or β, the relationship between the power of q and the length of the Wilson loop would be more complicated.This suggests that the parameter settings (3.9) and (3.10) are crucial for the partition function of the gKM model (3.7) counts the Wilson loops correctly.

Duality between unitary and Hermitian matrix integrals
We have obtained the expression (3.8) (or (3.12) after the parameter settings) by integrating the Hermitian matrices Φ v first.Instead, we can express the partition function (3.7) by integrating the unitary matrices U e first [1] by using the HCIZ integral formula [27], where U is a unitary matrix of size N , dU is the Haar measure of U (N ) normalized as dU = 1, We apply the formula (3.13) to the partition function of the gKM model (3.7).Since the integrand includes only the eigenvalues of the Hermitian matrices Φ v , we also use a mapping from the matrix integral to the integral over the eigenvalues (see e.g.[28]), for a Hermitian matrix Φ with eigenvalues We then obtain e qφ s(e),i φ t(e),j , ) are the eigenvalues of Φ v rescaled so that the parameter β is included only in the overall factor.In particular, after setting (3.9) and (3.10), the partition becomes as (3.12) and Z G (q) can be expressed in two ways as e qφ s(e),i φ t(e),j .(3.15) In the next section, we evaluate this integral exactly in the large N limit.Before that, however, we consider special cases when we can evaluate this integral at finite N .

The dual integral in cycle graphs
The integrals (3.We first evaluate the integral (3.14).Since the cycle graph C n has only one chiral primitive cycle e 1 • • • e n as a representative of the equivalence class, the length of a general reduced cycle is a multiple of n, and there are 2n reduced cycles of length nk in total, that is, k-power of n chiral primitive cycles and their inverses.Therefore, the integral (3.14) becomes where the integration over the n unitary matrices has been reduced to an integration over a single unitary matrix U in the last line by using the left-right invariance of the Haar measure dU and dU = 1.Note that the similar type of integral, which has more generic coefficients in the exponent, is evaluated by the Frobenius formula in [26].It is remarkable that we will obtain the more explicit result despite finite N thanks to the dual description.We can determine the function Z C1 (q) at finite N by evaluating the integral (3.15) for n = 1 (1-gon); where M σ is a matrix of size N .
The matrix with the i-th and j-th rows and columns of M σ simultaneously replaced is M τij στij , where τ ij ∈ S N is the transposition of i and j, which does not change the determinant.Therefore, assuming that σ is an element of the equivalent class labeled by the partition λ = (λ 1 , • • • , λ r ), we can transform M σ to the following standard form without changing the determinant by acting an appropriate permutation: where K m (m ∈ N) are matrices with size m given by Since it is straightforward to show we see Since the summation ) converges only for |q| < 1, (3.17) is evaluated as where sgn (λ) = r i=1 (−1) λi−1 and n λ is the number of conjugates of a permutation σ ∈ S N labeled by the partition (Young diagram) λ.The partition λ is also given by a set of the length l i and the multiplicity m i of the rows of the Young diagram as λ = (l m1 1 • • • l mr r ).Using this notation, if we define n λ is expressed by N !z λ .Thus, applying the formula proved in Appendix B; we finally obtain for any finite N and |q| < 1.
We make a comment on the region of q.Although the expression (3.16) is well defined only in the region |q| < 1, the dual expression (3.17) is defined also in |q| > 1.Using the formula proved again in Appendix B, we can evaluate (3.17) in the region |q| > 1 as (3.23) In Appendix C, we evaluate the unitary matrix integral (3.16) directly by using Cauchy's integral formula and reproduce the result (3.21) for |q| < 1. Interestingly enough, in the same appendix, we also evaluate the integral for |q| > 1 ignoring the convergence, which reproduces (3.23).We thus expect that the unitary matrix integral (3.16) can be analytically connected to the region |q| > 1.We will also discuss this point in the last section.
4 The matrix integral in the large N limit In this section, we evaluate Z G (q) in the large N limit.In preparation for this, we start with showing some properties satisfied by unitary matrix integrals at large N .
To simplify the description, we introduce the notation,

Integrals of Wilson loops in the large N limit
Fundamental integral and its behavior at large N The most fundamental integral in performing general unitary matrix integrals is ) where Wg(σ, N ) is the Weingarten function [29,30], , where χ λ (σ) is the character of S n corresponding to the partition λ and s λ,N is the Schur polynomial of λ.Since we are interested in the large N limit, we do not need the detail of this function but need only its asymptotic behavior in N → ∞, which is given by [31], where we have assumed σ ∈ S n is a product of cycles σ = c 1 • • • c k , |σ| denotes the minimal number of factors necessary to write σ as a product of transpositions, and Cat(m) = (2m)!m!(m+1)! is the Catalan number.From (4.2), we see that the leading contributions of the 1/N expansion of the integral (4.1) come from τ = σ and thus it can be written as Note that this is also obtained in [32][33][34].All the properties we show in this subsection follow from this expression.

Integral of a single Wilson loop
Let us consider a Wilson loop.Recalling that a Wilson loop is associated with an equivalent class of reduced cycles and a reduced cycle is consist of a power of a primitive cycle, the Wilson loop is written as Tr P C (U ) l for a primitive loop C in G and an integer l.We can then evaluate the integral of | Tr P C (U ) l | 2 at large N as To show it, we pick one of the unitary matrices included in P C (U ) and call it U .Then Tr P C (U ) can be written in general as where i = ±1 (i = 1, • • • , r) and A i is a product of U e (e ∈ E) except for U .Since C is a primitive cycle, all A i are never the same.In the following, we assume i = 1 for all of i = 1, • • • , r because, as we will see soon, the signatures are irrelevant in the following discussion.
If there is such a unitary matrix U that satisfies r = 1, the left-hand side of (4.4) is reduced to an integral over a single unitary matrix by absorbing all the other unitary matrices in P C (U ) to U via the left-right invariance of the Haar measure.Then (4.4) immediately follows from the formula given in [35,36], Tr where λ = (l m1 1 • • • l mr r ) is a partition, |λ| = r i=1 l i m i , and z λ is defined by (3.19).However, r is not always 1 in general and so we have to evaluate the integral more carefully.
Let us first show the case of l = 1, which includes all the essence.From (4.3), the left-hand side of (4.4) is written as where O(1/N ) includes terms with different permutations among i's and j's.For each σ, we obtain a product of traces of A i and A † i but the maximal number of the traces in the product is r from the construction.
The point is that the trace of A i and A † i gives Tr(A i A † i ) = N because A i 's are unitary matrices and this is the only chance to give an O(N ) contribution.Therefore the maximal power of N in the summation of σ ∈ S n is Tr( = N r and it occurs only once in the summation because C is now assumed to be a primitive cycle, and the other contributions are at most of the order of 1/N compared to this term.In fact, in the above expression, σ = 1 gives the contribution N r and other terms consist of a product of traces times a power of N less than r.Note that, if we take the signatures i into account, the permutation σ which gives N r is not σ = 1 in general, but the structure is the same: Only one element of S r gives the contribution N r and the other terms are O(1/N ) compared to it.This is the reason why the signatures i are irrelevant in this discussion.Therefore (4.7) becomes The remaining O(1/N ) terms include functions of other unitary matrices than U .However, they disappear in the large N limit.This is because the integral of such a function is of the order of unity in general and the number of terms appearing in this part is O(N 0 ), since we are considering sufficiently large N than r.For the same reason, we can ignore the O(1/N ) terms in (4.7).We can then conclude (4.4) for l = 1.We next consider the case of l > 1.In this case, the leading contribution (Tr( where ) and the indices of D are assumed to be written in mod l.Therefore, as the same reason of the case of l = 1, we can conclude (4.4) for l > 1.

Large N decomposition
We next consider two Wilson loops Tr(P C1 (U )) and Tr(P C2 (U )) corresponding to two (not necessarily primitive) cycles C 1 and C 2 , respectively.We can show that the integral of the product of the two Wilson loops are decomposed into the product of the integrals of the individual Wilson loops at large N : The proof is essentially the same as the proof of (4.4), or the proof of the large N decomposition of the expectation values of Wilson loops in QCD [37] (see also [38] and references therein).We assume that P C1 (U ) and P C2 (U ) share a unitary matrix U , otherwise (4.8) is trivial.We express Tr P C1 (U ) and Tr P C2 (U ) as as same as (4.5),where we have ignored the irrelevant signatures for simplicity.Then the left-hand side of (4.8) is written as From the same discussion as evaluating (4.7), the leading contribution is N r+s which comes from Tr(B j B † j ) and the other contributions disappear at large N .In this case, all of such terms are included in the summation of σ ∈ S r × S s ⊂ S r+s which gives |Tr P C1 (U )| This means that all the terms that mix Tr P C1 (U ) and Tr P C2 (U ) in the integral vanish at large N .Therefore we can conclude the decomposition (4.8) in the large N limit.

Integral of general Wilson loops
Combining (4.4) and (4.8), we see where the factor m! comes from the number of ways to choose m pairs of Tr P C (U ) l and Tr P C (U ) †l in Tr P C (U ) l 2m .
In general, we can associate a partition λ C = (l m1 1 l m2 2 • • • ) to a chiral primitive cycle C and consider the quantity, We can then evaluate the integral of a product of general Wilson loops at large N as where Π + is the set of the equivalence classes of the chiral primitive cycles in G, and z λ C is given by (3.19).
Comparing this result to the formula (4.6), we see that we can treat the Wilson loops corresponding to the chiral primitive cycles as if independent variables of the unitary integral on the graph in the large N limit.
4.2 Evaluation of Z G (q) at large N Let us then evaluate the matrix integral Z G (q) given in (3.12) at large N based on the properties we discussed so far.
Using the fact that Z G (q) counts all the Wilson loops with their length as weights and a reduced cycle is expressed as a power of a primitive cycle, we can rewrite Z G (q) as where we have used the large N decomposition in the last line.
Using the general relation, where we have used (4.9), p(n) is the number of the partition of n, and we have used the fact that ∞ i=1 (1 − q i ) is the generating function of p(n).Then we conclude that Z G (q) is written as an infinity product of (the square root of) the Ihara zeta function: This means that the tuned gKM model can be solved exactly in the large N limit.We see that the large N limit of (3.21) for |q| < 1 becomes (4.10) as expected.
The key to the derivation of the partition function (4.10) in the large N limit is that the integral of the unitary matrix in (3.12) is completely separated into the primitive cycles.The same is true for the expression (3.15) obtained by performing the U integral first.In fact, in [39], an equation satisfied by the eigenvalue distribution of the scalar fields of the original KM model has been derived through evaluating the HCIZ integral in the large N limit.Furthermore, an exact solution of this equation was obtained in [40] where the eigenvalue distribution of the scalar field becomes semi-circle, regardless of the dimension of the square lattice.Although our model cannot be compared directly with this exact solution of the original KM model because the parameters are adjusted so that the Ihara zeta function is realized, it is interesting that the scalar fields on the any dimensional square lattice behave in the same way as those on the one-dimensional square lattice.The one dimensional square lattice is nothing but the cycle graph discussed in Sec.3.4, and the partition function in large N limit becomes Z Cn (q) −→ 16) and (3.21).Combining this observation to (4.10), we see that the partition function of the gKM model on an arbitrary graph G can be rewritten as a product of the partition function of the model on the cycle graphs in the large N limit: Therefore, the gKM model in the large N limit on any graph G can be expected to behave similarly to the model on the cycle graph in general, which supports the results obtained with [40] in an extended way.
As another comment, it is suggestive that the expression (4.10) can be written as 3).Regarding it as an Ihara zeta function, the number of the reduced cycles of length k in the corresponding graph is counted by N k = d|k dN k/d and the free energy has the exact coupling expansion in q, in the large N limit of the gKM.The meaning of this expression will be discussed separately.

Finite N effect
We make a comment on Z G (q) on finite N .Looking at the results (3.18) and (4.10), one may think that (4.10) might hold even for finite N , but it is not the case in general.
As an example, we consider the graph showed in Fig. 3, which is the simplified version of the double triangle graph shown in Appendix A. Writing W k (U ) ≡ Tr W k U , we can easily show W 2i−1 (U ) = 0 (i ∈ N) for this graph and the expansion of the extended zeta function is expanded as After a straightforward calculation using (4.1) and (4.6), we obtain for N ≥ 3. Comparing it to the expansion for N ≥ 3, it is obvious that e dU e ζ G (q; U ) 2 at finite N .This is of course expected from the analysis made in Sect.4.1.

Conclusion and Discussion
In this paper, as an extension of the Ihara zeta function associated with graphs, we proposed a function that generates all non-collapsing Wilson loops with their lengths as weights.We also proposed the Kazakov-Migdal model on graphs (gKM model) and showed that, when the parameters of the model are appropriately tuned, the partition function is represented by the unitary matrix integral of the extended Ihara zeta function.Since the KM model has both Hermitian matrices and unitary matrices, there are two different ways to represent the same partition function, depending on which integral is performed first.Applying this duality to the tuned gKM model of a cyclic graph (single polygon), we exactly evaluated the unitary matrix integral dU exp ∞ n=1 q n n | Tr U n | 2 at finite N .We also discussed a possibility to analytically connect this integral to the region |q| > 1.We further showed that the partition function of the tuned gKM model corresponding to any graph can be computed exactly at large N and expressed as the infinite product of the square roots of the Ihara zeta functions.
The tuned gKM model is the minimal integral to count all the Wilson loops by unitary matrix integral and this is exactly the reason why the partition function of this model is expressed as a product of the square roots of the Ihara zeta function rather than the zeta function itself.The point is that an integral of a unitary matrix U vanishes unless the integrand contains equal numbers of U ij and U * ij .For this reason, the simple integration of a Wilson loop e∈E dU e Tr P C (U ) on the graph generally vanishes, except for such special Wilson loops that happen to contain equal numbers of U ij and U * ij .Therefore, in order that the U -integral of a Wilson loop always produces a non-zero value, we need to pair a Wilson loop with its Hermitian conjugate and evaluate the quantity e∈E dU e | Tr P C (U )| 2 .From the standpoint of the cycles on the graph, the Wilson loop and its Hermitian conjugate correspond to a reduced cycle C and its inverse cycle C −1 , respectively.Then, if we want to count up all the Wilson loops by unitary matrix integral, we have only to count up exactly half of the reduced cycles, which is nothing but the square root of the Ihara zeta function.This suggests that the gKM model on an arbitrary graph becomes effectively a direct product of the models on cycle graphs (one-dimensional lattices) in the large N limit.This is expected to be consistent with the result obtained in [40] that the exact solution of the original KM model exhibits universal behavior regardless of the lattice dimension.
As discussed in Sect.3.4, in evaluating the integral dU exp , we did not integrate the unitary matrix directly, but rather performed the Hermitian matrix integral (3.17) using the duality of the tuned gKM model.Interestingly enough, this Hermitian matrix integral is defined not only for |q| < 1 but also for |q| > 1 as in (3.21).On the other hand, the tuned gKM model is defined only for |q| < 1 since otherwise β becomes negative.It is also clear from the fact that the partition function at large N is represented by the Ihara zeta function, whose infinite product representation is defined only for q with sufficiently small |q|.Therefore, we would normally expect that only the result (3.21) in |q| < 1 would be meaningful.However, as we show explicitly in Appendix C, we can perform the unitary matrix integration directly for small N by using Cauchy's integral theorem both in the regions |q| < 1 and |q| > 1, and the results are identical with (3.21) and (3.23), respectively.This suggests that this unitary matrix integral can be analytically connected to the region of |q| > 1.This also opens an possibility for the tuned gKM model to be extended to the region |q| > 1.
It is also interesting to think of the tuned gKM model as an induced QCD as discussed in the original KM model.The original KM model does indeed add up all Wilson loops, but it includes also the contributions of backtrackings and the length of those Wilson loops does not match the power of q.By applying the parameter setting introduced in this paper, the relation between the length of the Wilson loop and the power of q is established, which is expected to make the connection to the gauge theory easier.Furthermore, although we mainly concentrate on Z G (q) in the partition function (3.12) in this paper, the partition function also includes the factor α − n V N 2

2
. In order to evaluate the large N behavior of the model, we would have to control the partition function by carefully tuning the value of the factor α. We may then see the large N behavior of QCD explicitly through the exact expression of the partition function obtained in this paper.
In [15], it is shown that HCIZ integral and KM model can be embedded into a supersymmetric model on the graph, where the matrix integral is evaluated as a result of the Duistermaat-Heckman localization formula (fixed point theorem).We expect that our (tuned) gKM model can also be embedded into a supersymmetric model and the related Ihara zeta function can be understood as a vev of the supersymmetric Wilson loops.On the other hand, the partition function of the supersymmetric gauge theory, including the quiver quantum mechanics (see e.g.[41,42]), gives various kinds of supersymmetric indices, which count the number of the gauge invariant operators (chiral rings) or BPS bound states.So it is interesting to find a relation between the BPS state counting and Ihara zeta function via the supersymmetric quiver gauge theories on the graph.Since the supersymmetric quiver quantum mechanics also counts the BPS bound states of the multi-centered black holes in the dual gravitational (closed string) picture [43], the Ihara zeta function, involving graph theory as well, has the potential to shed new light on gravity and string theory.

A Some examples of the Ihara zeta function A.1 Cycle graph
The cycle graph C n is shown in Fig. 2 (a).C n has n vertices and edges, namely n V = n E = n.The adjacency matrix is an n × n matrix and given by , and the degree matrix is given by D Cn = diag(2, 2, . . ., 2).Then the Ihara zeta function is given by an inverse of the determinant of the following matrix , which is nothing but the q-deformed Cartan matrix of the affine Lie algebra Ân−1 (the C n graph is the Dynkin diagram of Ân−1 ).Thus we get The square in this expression means that there are a primitive cycle of the length n and its inverse (anti-clockwise and clockwise cycle).So the square root of the Ihara zeta function picks up the primitive cycle of one direction only, which is the chiral primitive cycle.
The logarithm of the Ihara zeta function has the expansion, This means that N k = 2n iff k ≡ 0 (mod n) and then all the reduced cycles are expressed by the power of the one-directed primitive cycle of the length n.

A.2 Double triangle
Here we consider the graph depicted in Fig. 4.This graph contains two triangles, so we call this graph the double triangle (DT), which is also obtained by removing one edge from the complete graph K 4 .DT has four vertices and five edges (n V = 4 and n E = 5).There are infinitely many primitive cycles on DT, since there are infinite ways to go around the combination of the different cycles.So the original definition of the Ihara zeta function becomes an infinite product over the primitive cycles.However, Ihara's theorem says that the Ihara zeta function is expressed by a determinant made from the 4 × 4 adjacency matrix of DT , and the degree matrix D DT = diag (3,2,3,2).We obtain a rational function of q ζ DT (q) = 1 (1 − q 2 ) det (1 4 − qA DT + q 2 (D DT − 1 4 )) = 1 (1 − q 4 )(1 + q 2 − 2q 3 )(1 − q 2 − 2q 3 ) .

A.3 Tetrahedron
If we would like to consider the Wilson loops of the gauge theory defined on the surface of the polyhedra [17], it is useful to consider the polyhedra as graphs.For example, let us consider the tetrahedron.If we give connections between vertices by edges of the tetrahedron, we obtain a graph in the three-dimensional picture depicted in Fig. 5 (a).The tetrahedron graph is equivalent to the planner graph shown in Fig. 5 (b) and (c).This graph is also called the complete graph K 4 , which has four vertices and six edges (n V = 4 and n E = 6).The adjacency matrix of the complete graph K 4 is given by , and the degree matrix is D K4 = diag (3,3,3,3).Then we can obtain the Ihara zeta function, The logarithm of the Ihara zeta function for K 4 has the series expansion, log ζ K4 (q) = 8q 3 + 6q 4 + 16q 6 + 24q 7 + 21q 8 + 176q 9  3 + 120q 10 + 168q 11 + 330q 12 + • • • .
From this, we have The number of the reduced cycles increases more rapidly than the case of the DT.

B Formulas for symmetric polynomials
In this appendix, we give a proof of the formulas (3.20) and (3.22) 3 .
We first define the homogeneous symmetric function associated to n ∈ N by x ij , where x = {x i |i ∈ N}.It is easy to see that the generating function of the homogeneous symmetric functions is given by t n h n (x) .
Here we transform H(t) as is a sequence of the edges which satisfies t(e a ) = s(e a+1 ) (a = 1, • • • , −1), where is called the length of the path P , which is denoted as (P ).If two paths P = e 1 • • • e and P = e 1 • • • e satisfy t(e k ) = s(e 1 ), we can construct a new path of length + by connecting them as P P ≡ e 1 • • • e e 1 , • • • , e .A backtracking of P is a part of P which satisfies e −1 a+1 = e a .When the path P satisfies s(e 1 ) = t(e ), P is called a cycle of length .A cycle C = e 1 e 2 • • • e is called tailless when e −1 = e 1 , which is equivalent to that C 2 has no backtracking.A cycle C is called reduced when C has no backtracking nor tail.A reduced cycle C is called primitive when C satisfies C = B r (r ≥ 2) for any reduced cycle B.We call two cycles C = e 1 • • • e and C = e 1 • • • e equivalent when e a = e a+r for an integer r, namely the difference between two cycles is only the starting point.Using the equivalence of the cycles C ∼ C , we can define the equivalence class of the cycles, which is denoted by [C]; [C] = {e 1 e 2 • • • e , e 2 e 3 • • • e 1 , . . ., e e 1 • • • e −1 }.

Figure 1 :
Figure 1: Discretizations of a flat plane or torus T 2 can be obtained from graphs.The original KM model is defined on the square lattice graph (b).In this sence, the gKM model has possibilities for constructing various kind of discretized scalar field theories by choosing the graph structure.
the Barnes function, A and B are Hermitian matrices whose eigenvalues are (a 1 , • • • , a N ) and (b 1 , • • • , b N ), respectively, and ∆(a) and ∆(b) are the Vandermonde determinant with respect to A and B, respectively; 14) and(3.15)are difficult to carry out in general, but there are exceptions; it becomes a Gaussian integral when deg v = 2 for all the vertices, which is achieved when G is a cycle graph (a single polygon).Let us consider the case when G is a cycle graph C n (n-gon).(See Fig.2.)In this case, the number of the vertices and edges are the same, n V = n E = n, and we have the same number of unitary matrices U e and Hermitian matrices Φ v (e, v = 1, • • • , n).

Figure 2 :
Figure 2: A cycle graph C n has n vertices and edges.This is also called Ân−1 quiver chain.(b) is a simplest case of C 1 (1-gon).Z Cn reduces to Z C1 with the length n cycle.

Figure 3 :
Figure 3: An example of the graph.

Figure 4 :
Figure 4: This graph has four vertices and five edges.Two triangles are glued by the edge e 5 .

N 3 = 12 , N 4 = 8 ,
N 6 = 24, N 7 = 28, N 8 = 8, N 9 = 48, N 10 = 120, N 11 = 44, N 12 = 104, • • • ,and the others are zero.Note that N k are even numbers, since the equivalent class of the reduced cycle always contains the cycle of both directions.For example, N 3 /2 = 3 + 3 stands for the number of two chiral primitive cycles with length three on each triangle, N 4 /2 = 4 stands for the number of chiral primitive cycles around the outer rhombus with length four, N 6 /2 = 3 + 3 + 6 stands for two squares of length three chiral primitive cycles on each triangle, and a length six product of two chiral primitive cycles (triangles), and so on.

Figure 5 :
Figure 5: This graph has four vertices and six edges and forms a tetrahedron, where four triangles are gluing with each other.The three-dimensional tetrahedron graph is equivalent to a planner complete graph K 4 , since the graph has only the information on the connections between vertices.

p
p l (x) m l m l !l m l t lm l =where p n (x) is called the power sum symmetric function,p n (x) ≡ λ l (x) , for a partition λ = (λ 1 , • • • , λ r ) of N which satisfies λ 1 ≥ • • • ≥ λ r , |λ| ≡ r i=1 λ r = N .Comparing the coefficient of t n , we obtain the relation,