Entanglement of Free Fermions on Hadamard Graphs

Free Fermions on vertices of distance-regular graphs are considered. Bipartition are defined by taking as one part all vertices at a given distance from a reference vertex. The ground state is constructed by filling all states below a certain energy. Borrowing concepts from time and band limiting problems, algebraic Heun operators and Terwilliger algebras, it is shown how to obtain, quite generally, a block tridiagonal matrix that commutes with the entanglement Hamiltonian. The case of the Hadamard graphs is studied in details within that framework and the existence of the commuting matrix is shown to allow for an analytic diagonalization of the restricted two-point correlation matrix and hence for an explicit determination of the entanglement entropy.


Introduction
As a probe in particular of the correlations in quantum many body systems and field theories, the study of entanglement is of fundamental interest. We shall here be concerned with free Fermions on graphs, that is we shall explore systems of Fermions hopping between the vertices of distance-regular graphs with dynamics controlled by the adjacency matrix taken as the Hamiltonian. We shall develop general results and focus on their realization in Hadamard graphs. A key finding will be the identification of a block tridiagonal matrix that commutes with the entanglement Hamiltonian. As will be shown, this will allow for an analytic diagonalization of the chopped correlation matrix on Hadamard graphs. The commuting operator will be found by extending the algebraic Heun operator construction to the framework of the Terwilliger algebra arising in the context of association schemes. This will generalize previous entanglement studies of free Fermion chains [27,9,10,4,5], that exploited ideas borrowed from time and band limiting problems [31,21,32,13,15].
To talk about entanglement, one needs to take the system under consideration in a given state and to examine how a part of it is correlated to the rest. Often this is done in the ground state. In the case of Fermionic systems, this requires a certain filling of the Fermi sea, there is hence a limited set of energies that is thus selected. Determining a part of the systems in turn involves picking an ensemble of vertices.
Because of their simplicity, free chains have been fairly well studied. In these situations, owing to Wick's theorem, the reduced density matrix can be expressed in thermodynamic form in terms of a so-called entanglement Hamiltonian which is determined by the two-point correlation matrix restricted to the chosen part of the system (for reviews see [22,28]). Quite remarkably in some instances it has proven possible to find a tridiagonal matrix that commutes with the correlation matrix and thus shares with it joint eigenvectors. Either for numerical computation issues or asymptotic analyses, this fortuitous circumstance has proved important.
The existence of this commuting operator for free Fermionic chains could be suspected from the parallel with problems of time and band limiting in signal processing where one wishes to optimize the concentration in time of signals comprised of a limited frequency interval. (The restrictions in the ranges of two dual variables is at the root of the analogy.) The seminal work in the signal analysis context was done by Slepian and coworkers [31,21,32], who identified such a useful operator commuting with the limiting operator. It should be said that this feature has emerged in various other contexts among which random matrix theory or integrable models.
The bearing of all this on the examination of entanglement in free Fermionic chains was brought to light by Eisler and Peschel (among others see also [11]) who found for example in [10], the tridiagonal matrix that commutes with the truncated correlation matrix of a uniform finite chain and studied in detail how this commuting matrix compares with the entanglement Hamiltonian and how it may yield the spectrum of the latter.
The reasons behind the rather miraculous existence of the commuting operators have intrigued many and Grnbaum who has devoted much effort to understand this phenomenon has tied it to the field of bispectral problems that he has established [7]. This later led to the concept of algebraic Heun operator introduced in [15,16] (see also [24,25]) which was shown to naturally provide commuting operators. These methods were applied in [4,5] to various free Fermionic finite chains with non-uniform couplings and possessing a bispectral underpinning. More precisely, they were used to analyse chains described by a hopping matrix h tridiagonal in the characteristic (site) basis and admitting another "position" operator X tridiagonal in the eigenbasis of h. (Such a set of operators conforms to the definition a Leonard pair [34].) In situations where the limiting in energy and space involved sequences of consecutive values beginning with the initial one, the algebraic Heun operator -a bilinear expression in the two Leonard operators -was readily found in [4,5] to give the commuting matrix.
We here extend this study to systems of free Fermions on graphs and in particular on those of Hadamard. While the entanglement of such Fermionic configurations has been looked at in [19,20], our main interest will be to identify the commuting matrices in the case of distance-regular graphs. This will prove possible when the limitings will respect the distance structure of the graph. Because of the presence of degeneracies, it will call for an extension of the algebraic Heun construction to the framework of Terwilliger algebras [33] which have been introduced to study the properties of association schemes. Simply put these algebras are generated by the adjacency matrix which plays in our models the role of the Hamiltonian and the orthogonal projectors E * i on subspaces spanned by the characteristic vectors corresponding to the vertices at distance i from the reference one. (For an introduction to Terwilliger algebras see [12].) This will lead to commuting operators that are block-tridiagonal (instead of being simply tridiagonal in the case of the Fermionic chains described in the preceding paragraph). This should be put in parallel with the recent investigations of Grnbaun and collaborators to identify commuting operators stemming from matrix orthogonal polynomials [14].
The analysis will be carried out in details for Hadamard graphs and the usefulness of the commuting block-tridiagonal matrix will prove manifest in this case. The Hadamard graph of order n has diameter 4 and 4n vertices [1]; their set form a self-dual association scheme. The spectral properties of Hadamard graphs make them a prime candidates, amongst distance regular graphs, for a variety of applications; for example they are studied in the context of state transfer in quantum walks, see [3]. Given a reference vertex, the chopped correlation matrix will be expressed in terms of projectors made out of the primitive and dual idempotents of the scheme; it will therefore belong to the associated Terwilliger algebra. One usually want to obtain the eigenvalues of the correlation matrix to determine for instance the entanglement entropy. The essential merit of the commuting operator that shares common eigenvectors with the correlation matrix is that in contrast to the latter, it generally lends itself to a tractable numerical diagonalization. Here, quite remarkably, for Hadamard graphs, it will be seen to allow in fact for an analytic determination of the eigenvectors.
The paper will now unfold as follows. Section 2 will provide the description of free Fermions on graphs. The Hamiltonian will be diagonalized with the help of the primitive idempotents which will also be used to express the two-point correlation matrix in the ground state. Section 3 frames the discussion of entanglement. The system is split in two parts with one comprising the complete neighbourhoods that are consecutive to the reference vertex. That the entanglement Hamiltonian is determined by the restricted correlation matrix is reviewed and this central entity is seen to be given as the product of projectors. Attention will be focused on graphs that belong to association schemes that are P-and Q-polynomials. These will be surveyed in Section 4 where the Terwilliger algebras will be defined. Section 5 contains one of the main results and shows how a block-tridiagonal operator that commutes with the restricted correlation matrix is obtained from the algebraic Heun construction and the Terwilliger algebra. Section 6 shows in detail how the analysis applies in the case of Hadamard graphs. Free Fermions on these graphs will be shown to provide a system where the entanglement Hamiltonian can be diagonalized analytically thanks to the existence of the commuting operator. Exact expressions for the Von Neumann entropies will be obtained as a result. Concluding remarks will be found in the last section.

Free Fermions on graphs
Let G be a graph, V = {v 1 , . . . , v N } the set of its vertices and E the set of its edges. We define a vector space of dimension N with the orthonormal canonical basis {|1 , . . . , |N }, called position basis, where the vector |i is associated to the vertex v i . We define the adjacency matrix A as a N × N matrices with the components, for 1 ≤ i, j ≤ N , We restrict ourselves to non-oriented edges i.e.
In this case, the matrix A is symmetric. We introduce the fermionic operators c n and c † n (n = 1, 2, . . . , N ) satisfying We consider the following open quadratic free-Fermion Hamiltonian on the graph G It describes the hopping of Fermions on the underlying graph G. As usual, in free-Fermion models, to diagonalise the Hamiltonian H, one diagonalises the matrix A. As that matrix is symmetric, it can be diagonalised by an orthonormal basis |θ i k for k = 0, 1, . . . , δ and i = 1, 2 . . . f k , where θ k are the pairwise distinct eigenvalues of A and f k is the multiplicity of the eigenvalue θ k . The eigenfunctions φ j (θ i k ) are the components of the vector |θ i k in the position basis: We order the eigenvalues such that θ 0 < θ 1 < · · · < θ δ . We introduce E k , the primitive idempotent associated to the eigenspace with eigenvalue θ k : Let us remark that f k = rank(E k ). We define also the following projector on the first K + 1 eigenspaces: These projectors permit to write the adjacency matrix as follows Having diagonalised A, the Hamiltonian H can be written as follows where the transformed fermionic operatorsc i k and (c i k ) † are related to the first ones bỹ Then, the eigenvalues E of the Hamiltonian H are described by the set of integers k ∈ {0, 1, . . . f k } and are given by Let K ∈ {0, 1, . . . , δ} such that The number K can be modified by adding a matrix proportional to the identity at the adjacency matrix A. That corresponds to consider an additional constant magnetic field for the Hamiltonian H. We can also choose this magnetic field such that the inequalities in (12) are always strict which simplifies the latter discussion about the ground state and the entropy. The ground state of the Hamiltonian H is where the vacuum state |0 is annihilated by the fermionic operatorsc i k : The correlation matrix C in the ground state |Ψ 0 is the N ×N matrix with components, for 1 ≤ n, m ≤ N , By using the second relation in (10), we can rewrite the correlation matrix as follows By rewriting the ground state with its definition (13), by using the anticommutation relations of the fermionic operators and the properties (14), one gets Therefore the correlation matrix corresponds to the projector π 2 defined by (7): C = π 2 (K).

Entanglement entropy
In order to examine entanglement, we must first define a bipartition of our free-Fermionic graph. We denote by d(v, w) the distance between the two vertices v and w of the graph G.
From now on, we fix a vertex v of G and, without loss of generality, we take v = v 1 . We choose as subsystem 1 the vertex v 1 and all the vertices at distance at most from v 1 : Let us denote by N 1 the cardinal of the previous set i.e. the number of vertices in the subsystem 1. We shall find how it is intertwined with the rest of the graph in the ground state |Ψ 0 . To that end, we need the reduced density matrix where the subsystem 2 is the complement of the subsystem 1, from which one can compute for instance the von Neumann entropy The explicit computations of this entanglement entropy amounts to finding the eigenvalues of ρ 1 . It has been observed that this reduced density matrix ρ 1 is determined by the spatially "chopped" correlation matrix C , which is the following N 1 × N 1 submatrix of C: The argument which we take from [26] (see also [29]) goes as follows. Because the ground state of the Hamiltonian H is a Slater determinant, all correlations can be expressed in terms of the one-particle functions, i.e. in terms of the matrix elements of C. Restricting to observables associated to part 1, since A = tr(ρ 1 A), the factorization property will hold according to Wick's theorem if ρ 1 is of the form with the entanglement Hamiltonian H given by The hopping matrix h = |h mn | 1≤m,n≤N 1 is defined so that holds, and one finds through diagonalisation that We thus see that the 2 N 1 ×2 N 1 matrix ρ 1 is obtained from the N 1 ×N 1 matrix C or equivalently, the entanglement Hamiltonian H.
To describe this chopped correlation matrix C, we introduce the projector where d, the diameter of the graph G, is the maximum distance of any vertex from v 1 and I is the identity matrix. Let π 1 ( ) denote the projection onto the first neighbourhoods of v 1 : and n i be the number of vertices at distance i from v 1 . The operator π 1 ( ) projects on N 1 = i=0 n i sites. Finally, the chopped correlation matrix can be written as (see for instance [23,17]) To calculate the entanglement entropies one therefore has to compute the eigenvalues of C. As explained in [27], this is not easy to do numerically because the eigenvalues of that matrix are exponentially close to 0 and 1. This motivates the search for a matrix T such that and easier to diagonalise numerically than C. We denote the commutator bracket by [·, ·]; that is, [A, B] = AB − BA for matrices A, B. We shall hence find a T satisfying (29) by looking for a matrix commuting with both projectors: In [4,5], a such operator T , which is a tridiagonal matrix, for free-Fermion on a line has been constructed when A, defining the Hamiltonian H, belongs to a bispectral problem. In this case, the operator T turns out to be a Heun operator. Inspired by these results, we consider a special class of graphs whose associations schemes are P -and Q-polynomials (which we will provide background on in next section) which replaces the notion of bispectrality. Then, we show that in this case an operator T tridiagonal by block satisfying (30) exists. We may observe that the matrix D defined by would describe a dual entanglement situation where the vacuum state would be filled with excitations by the set {0, ..., }, and the subsystem would consist of the sites in the neighbourhood {0, . . . , K}. We recall that for n × n matrices M and N , the characteristic polynomials of M N and N M coincide, see Section 11.6 of [30]. Recalling that π i (j) is idempotent for i ∈ {1, 2}, we see that C and D both have the same eigenvalues as the following matrix: Since C and D have the same non-zero eigenvalues, the entanglement entropies will be the same in these two instances. Such duality generalizes the one studied in [2].
4 P -and Q-polynomial association schemes

Association schemes and distance regular graphs
A distance regular graph is a graph G of diameter d such that there exists constants p k ij for i, j, k ∈ {0, . . . , d} such that, for any pair of vertices x, y at distance k, there are p k ij vertices at distance i from x and distance j from y. We follow the standard text [1] and the dynamic survey [6] and defer to them for further background on distance regular graphs. For i = 0, . . . , d, the distance matrix, denoted A i , is given by The distance matrix A 1 corresponds to the adjacency matrix A used in the previous sections.
The following is the main theorem about the eigenspaces of distance regular graphs and can be found in [1]. In the following theorem, let N the number of vertices of G and J denotes the all ones matrix, I denotes the identity matrix and • denotes the Schur product, the entry-wise matrix multiplication.
More generally, if {A 0 , . . . , A d } are symmetric matrices such that Theorem 4.1 holds, then we say that they are the associate matrices of an (commutative) association scheme. For brevity, we sometimes refer to an commutative association scheme as a scheme. The constants p k ij are known as the intersection numbers of the scheme and the constants q k ij are known as the Krein parameters of the scheme.
The distance matrices of a distance regular graph form the basis of a (d + 1)-dimensional, semi-simple, commutative subalgebra M of the algebra of N × N matrices over C, also call the Bose-Mesner algebra of the graph. Since are both bases of M, there exists change of bases matrices between them. The eigenmatrices of the association scheme are d + 1 × d + 1 matrices P and Q such that Note that this implies that {P ij } d i=0 are the eigenvalues of A j . Since M is closed under multiplication and addition, there is a choice of the ordering of A 1 , . . . , A d such that A i is a polynomial in A 1 for all i. We say the scheme is P -polynomial if A i is a polynomial in A 1 of degree i, for each i = 0, . . . , d. It is not hard to see that, up to reordering the associate matrices, the condition of being P -polynomial is equivalent to requiring that p k ij = 0 whenever the sum of two of {i, j, k} is strictly smaller than the third element. The latter condition is also called metric and we may use P -polynomial and metric exchangeably. The class metric schemes are exactly those where A 1 is the adjacency matrix of a distance regular graph.
Similarly, we say a scheme is Q-polynomial or cometric if E i is a polynomial under Schur multiplication in E 1 of degree i, for each i = 0, . . . , d. Equivalent, a scheme is Q-polynomial if q k ij = 0 whenever the sum of two of {i, j, k} is strictly smaller than the other element.

Terwilliger algebra
We will now look at the Terwilliger algebra. We consider an association scheme A with associate matrices A 0 , . . . , A d . Throughout this section, let x ∈ V be a fixed vertex of A.
For i = 0, . . . , d, we will define the diagonal matrix E * i (x) as follows: We call E * i (x) the ith dual idempotent with respect of x. Similarly, for i = 0, . . . , d, we will consider diagonal matrices A * i (x) with entries as follows: We call A * i the ith dual distance matrix with respect to x. When the context is clear, we . Note, that if A 1 is the adjacency matrix of a distance regular graph, E * i is the diagonal characteristic matrix for the set of vertices at distance i from x, also known as the ith neighbourhood of x. We also have that Recalling the eigenmatrices of a scheme, we can see that The Terwilliger algebra T (x) is the subalgebra of Mat N ×N (C) generated by M and M * . The following lemma is found in [33].

Lemma. (Terwilliger 1992)
The above lemma suggest that the Terwilliger algebra is easier to work with when many of the intersection numbers or the Krein parameters of the scheme vanish. We will later look at the class of Hadamard graphs, which have many vanishing Krein parameters and intersection numbers.

Chopped correlation matrix for P -and Q-polynomials scheme
Let G be a distance regular graph of diameter d that is Q-polynomial. We recall that, in this case δ = d. We will consider the association scheme of G, with the idempotents in the Qpolynomial ordering. We choose as the matrix A defining the Hamiltonian (3) the adjacency matrix A 1 of this scheme. The eigenvalues of A denoted previously by θ i corresponds to the constants P i1 of the association scheme. As seen previously, the matrices A = A 1 and A * = A * 1 can be expand in the basis of idempotents and dual idempotents, respectively: However, due to the fact that we consider a graph in a P − and Q− polynomials scheme, A and A * satisfy supplementary relations. Indeed, for a P − and Q− polynomials scheme, one gets p j 1i = 0 and q j 1i = 0 for |i − j| > 1. Then, from Lemma 4.2 and from d i=0 E * i = 1, one deduces the following relation Similarly, by using d i=0 E i = 1, one obtains In the basis where E i is diagonal, A is diagonal and A * is tridiagonal by block. In contrary, in the basis where E * i is diagonal, A is tridiagonal by block and A * is diagonal by block. The choice of a graph G in a P − and Q− polynomials association scheme leads to interesting properties for the matrix A and A * which generalize the properties obtained from the fact to belong to a bispectral problem used in [4,5]. Let us emphasize that A and A * are not a tridiagonal pair [18] in general since there exists a subspace W of the space V spans by the vertices of the graph such that AW ⊆ W, A * W ⊆ W, W = 0 and W = V. However restricted to that space W, A and A * become a tridiagonal pair.
Recall that we are concerned with the eigenvalues of the chopped correlation matrix where we recall that π 1 ( ) = s=0 E * s and π 2 (K) = K k=0 E k . Note that matrix Π(j, ) is not in the Bose-Mesner algebra of G, but it is in the Terwilliger algebra of G.

Properties of particular chopped correlation matrix
We will state, without proof, some simple observations about Π(K, ) in the following lemma. We denote by 0 a block of 0s of the appropriate size. We will give the spectrum as a list of distinct eigenvalues θ with their multiplicities m θ in superscript in parentheses. For example, if the eigenvalues of M are a with multiplicity 1, 1 with multiplicity n and 0 with multiplicity 3n − 1, we write the spectrum of M as {a (1) , 1 (n) , 0 (3n−1) }. We denote by 0 m×n the m × n matrix of all zeroes and we omit the subscript when the order is clear from the context.    ( , 0 (N −1) }.

(Generalized) Heun operator
For µ, ν ∈ R, we define the Heun operator T as follows: where {·, ·} denotes the anti-commutator bracket; that is, {A, B} := AB + BA for matrices A, B. From relations (34)-(37) and the properties of the idempotents, we obtain two expressions for the Heun operator: We will give a choice of µ and ν such that the Heun operator commutes with the operator that we are considering.
Proof. By using (39), one gets With the value of ν given in the theorem the commutator vanishes. The second relation is proven similarly by using (40).

Formally self-dual association schemes
We say that an association scheme is formally self-dual if, for some ordering of the primitive idempotents, P = Q. It follows then that the intersection numbers and Krein parameters agree; more formally, for all 0 ≤ i, j, k ≤ d, p k ij = q k ij . The Hamming scheme and the Hadamard scheme are examples of formally self-dual schemes. Let A be a formally self-dual, metric and co-metric association scheme. We may assume that the primitive idempotents {E i } d i=0 , are ordered such that n i = f i . With this ordering, the P -and Q matrices of the scheme are equal. In particular Q i,1 = θ i .

Theorem.
If G is a graph whose adjacency matrix is the A 1 matrix of a formally self-dual, metric and co-metric association scheme, Π(K, ) is cospectral to Π( , K) for all K, = 0, . . . , d.
Proof. We have shown previously that D defined by (31) is cospectral with C = Π( , K). Now we show D is cospectral to Π(K, ).
For bipartite graph in a Q-polynomial and Q-bipartite scheme, one gets for all i = 0, . . . , d. By taking the forms (34) and (37) for A and A * , one can show that where we recall that θ i = P i,1 . Similarly, by taking the forms (35) and (36) for A and A * , one can show that where we recall that Q i,1 = θ i for formally self-dual association schemes. The Hadamard graph of order n (see below for the definition) is bipartite and formally self-dual with eigenvalues θ 0 = n, θ 1 = √ n, θ 2 = 0, θ 3 = − √ n, θ 4 = −n. Then, for τ = n and ρ = √ n, one gets θ 2 i − ρθ i θ i−1 + θ 2 i−1 − τ = 0 for i = 1, 2, 3, 4. One gets For Hamming graph, one gets θ i = −L + 2i (for i = 0, 1, . . . L). For τ = 4 and ρ = 2, one gets θ 2 i − ρθ i θ i−1 + θ 2 i−1 − τ = 0 and we recover the result [12]: 6 Hadamard graphs 6.1 Definition of the Hadamard graphs An example of a formally self-dual association scheme is the Bose-Mesner scheme of the Hadamard graphs. Hadamard graphs are exposited in Section 1.8 of [1]. We will give a full description here. The Hadamard graph of order n is a distance-regular, antipodal and bipartite graphs constructed as follows. Let H be an n × n Hadamard matrix. The graph H of H has two vertices c + and c − for each column of H and vertices r + and r − for each row of H. For 0 ≤ i, j ≤ n − 1, the graph H has edges (r + i , c + j ) and (r − i , c − j ) if H i,j = 1 and edges The resulting graph has diameter 4 and 4n vertices. For the ordering of the idempotents as above, the P and Q matrices are equal and are as follows: An equivalent definition for a distance-regular graph is a graph for which there exist integers  Figure 1 shows the Hadamard graph H 4 corresponding to H 4 , with its two color classes shown in white and gray. Up to isomorphism, H 4 is the unique Hadamard graph of order 4, on 16 vertices. We observe that H 4 is isomorphic to the hypercube graph Q 4 ; note that no other Hadamard graph is isomorphic to the hypercube on the same number of vertices, since the graphs will have different diameters.
Let H be Hadamard matrix of n. We may assume that the first row and first column of H all have entry 1. Let H be obtained from H by deleting the first column. We consider the Hadamard graph constructed from H and we will write the distance matrices as block matrices, partitioned by the distance partition with respect to vertex 0. The block sizes are 1, n, 2n − 2, n and 1 and the sets are: with the ordering of the vertices as given above.
Let R m be m × m back diagonal matrix given by R i,m+1−i = 1 for i = 1, . . . , m and all other entries are 0. Note that 0 n×m and J n×m denote the n × m block of 0s and 1s, respectively. We denote by 1 m a column vector of order m with all entries equal to 1. Recall that I m is the identity matrix of order m. Note that we will write our matrices as block matrices and the blocks are separated by lines. When the order of the submatrix is clear, we will omit the subscript of 0 n×m , J n×m , 1 m and I m . We let Note that M 1 + M 2 = J n,2n−2 and M 1 − M 2 = J n,n−1 −J n,n−1 . We observe that HH T = nI n − J n , so we may compute and With this setup, we may write the distance matrices more explicitly as block matrices: Let R m be the m × m back-diagonal matrix. Then For completeness, A 0 = I 4n .

Direct computations of the spectrum for chopped correlation matrices
We will find simpler expression for Π(j, ) for a Hadamard graph of order n using (50) and (51). In light of Lemma 5.1, we will focus on Π(j, ) for j, ∈ {1, 2, 3}.
6.1 Lemma. For any n, we consider Π(j, ) with respect to vertex 0. The following hold.
Proof. By using Theorem 5.4, the spectrum of Π(j, ) for j > is similar to the ones of Π( , j). We need to prove the above theorem only for Π(j, ) with j ≤ .

Spectrum of the chopped correlation matrix Π(2, 2)
Lemma 5.1 and Proposition 6.1 provide the spectrum for all the chopped correlation matrices expect for Π (2,2). We see that this chopped correlation matrix is quite complicated and a direct computation seems hopeless. We shall use the Heun operator introduced in Section 5.2 commuting with Π(2, 2): we find the eigenvectors of the Heun operator and let Π(2, 2) act on them to get the eigenvalues of Π(2, 2). As proved previously, the Heun operator (39) with ν = −Q 1,1 − Q 2,1 = − √ n and µ = −P 1,1 − P 2,1 = − √ n commutes with Π(2, 2) and reads as follows: where M 1 and M 2 are as in (50). Since T and Π(2, 2) commute, there exists a basis of R 4n which is a common eigenbasis for both matrices. Due to the block structure of both matrices, we need only diagonalise the first principal 3n − 1 × 3n − 1 submatrix of T . Proof. Since Π(2, 2) has a 3n − 1 × 3n − 1 block surrounded by zeroes, it has eigenvalues 0 with multiplicity (at least) n + 1. We will find the eigenvectors for the first 3n − 1 × 3n − 1 principal submatrix of Π(2, 2) by find an eigenbasis for the same principal submatrix of T given in (53). Since T and Π(2, 2) commute by Corollary 5.3, they share an eigenbasis. In general, not every eigenbasis of T need be an eigenbasis of Π(2, 2). we will find one that is. We will diagonalize the following matrix: Since each block of T has constant rows, we will find eigenvalues of the matrix: Consider an eigenvector v of B with eigenvalue λ. It is easy to see that the first entry of v is non-zero (else v = 0) so we may scale so that v = 1 a b T for some a, b. We obtain that Bv =   n 3 /2 + a n 2 + 2n 3 /2 n + 2n 1 /2 + an + 2bn 1 /2 (n − 1) Thus we have that We need only find the eigenvalues of B and obtain their eigenvector using (54). Observe that λ 1 = 2n 3 /2 + 2n, a(λ 1 ) = n − 1 /2 , b(λ 1 ) = 1 2n 1 /2 + 2 is one such solution. Let µ, ν be the two other eigenvalues of B. We have that Thus, µ, ν are the two solution of (x − µ)(x − ν) = x 2 − (µ + ν)x + µν, and thus We have thus the following three eigenvectors for T with eigenvalues λ 1 , µ and ν, respectively. Let f 0 , . . . , f n−2 be the elementary basis for R n−1 ; we mean that f i is the (n−1)-dimensional vector with 1 in the ith position and 0 elsewhere. We observe that is an eigenvector of T with eigenvector 0 for i = 0, . . . , n − 3. It is also easy to observe that {v i } n−3 i=0 is a linearly independent set since {f i − f n−2 } n−3 i=0 is linearly independent. Let e 0 , . . . , e n−2 be the elementary basis for R n−1 Consider the following: Observe that 2 √ n2 √ nM 1 M T 1 (e i − e n−1 ) = 2n(nI n − (n − 2)J n )(e i − e n−1 ) = 2n 2 (e i − e n−1 ).
For this, we can find vectors x which give eigenvectos of T . Let . Now we will find the eigenvalues of Π(2, 2). Let P be the first 3n − 1 × 3n − 1 principal submatrix of Π(2, 2).
We can compute that P u 0 = u 0 and the eigenvalues of P corresponding to u 1 and u 2 are From the results of Propositions 6.1 and 6.2, we can determine exactly the entanglement entropy for different configurations. In particular, the following exact asymptotics are found as a function of the size n of the Hadamard matrix: S 1,1 ∼ S 1,2 ∼ S 2,2 ∼ n→∞ 2 ln(2) − 3 4 ln(3) n .
These results can be qualitatively understood as follows. We should keep in mind that the adjacency matrix A 1 as a Hamiltonian is non-local. Now the Hadamard graph H n has 4n vertices. While a given site does not interact with all sites it is nevertheless coupled to n of them, that is to a number equal to the order of the graph. Let L = Im(Π 1 ( )) be one part of the bipartition and define its volume by V(L ) = Card(L ) (= N 1 ). The boundary of L is ∂L = {v i ∈ L : ∃ v j ∈ H n L with d(v 1 , v j ) = 1}. Figure 2 shows L 2 and its boundary on the Hadamard graph H 4 . The area S(L ) would be defined as S(L ) = Card(∂L ), see [8]. If l = 1, 2, 3, from the definition of the Hadamard graph, it is immediate to see that V(L ) and S(L ) are linear functions of n. The distinction between area and volume law in the thermodynamic limit is therefore pointless in this case.
Consider formula (58). The part L 3 of the graph contains all but one of its vertices. It is hence not a surprise that the distribution of quantum correlations does not depend on n in this situation. Looking at (59), the Fermi sea in this case extends to almost all modes. The correlation matrix C = Π 1 ( )Π 2 (K)Π 1 ( ) does then get closer and closer to 1 (like Π 2 (K)) as n becomes large. We thus realize that its eigenvalues will all approach 1 and that the entropy will go to zero as n grows. Finally, for = 1, 2, in view of the non-local character of the Hamiltonian, it is expected that the entropy will increase like the volume equivalently, like n, in the thermodynamic limit as is confirmed by (60).

Concluding Remarks
This paper has explored the entanglement of free Fermions on graphs of association schemes. It has shown how a block-tridiagonal operator that commutes with the entanglement Hamiltonian can be determined in certain instances. This relied on the extension of the algebraic Heun construct to the realm of Terwilliger algebras. The formalism has been applied to the study of Fermions on Hadamard graphs. The existence of the commuting operator has proved its use in this case as it allowed to diagonalize exactly the associated entanglement Hamiltonian -an achievement which is not so common.
Looking forward, the examination of the entanglement of free Fermions (or harmonic oscillators) on graphs of other associations schemes such as the Hamming and Johnson ones obviously comes to mind. With an eye to wider applications of algebraic Heun operators, we feel it would be most interesting to explore how such constructions could be further extended to tridiagonal pairs and algebras and how this could then also apply to entanglement studies. We plan on pursuing some of these questions.