# One-matrix differential reformulation of two-matrix models

Joren Brunekreef<EMAIL_ADDRESS>Luca Lionni<EMAIL_ADDRESS>IMAPP,
Radboud University, Nijmegen, The Netherlands.  Johannes Thürigen
<EMAIL_ADDRESS>Mathematisches Institut der Westfälischen
Wilhelms-Universität Münster
Einsteinstr. 62, D-48419 Münster, EU  Institut für Mathematik/Institut für
Physik der Humboldt-Universität zu Berlin
Unter den Linden 6, D-10099 Berlin, EU

###### Abstract

Differential reformulations of field theories are often used for explicit
computations. We derive a one-matrix differential formulation of two-matrix
models, with the help of which it is possible to diagonalize the one- and two-
matrix models using a formula by Itzykson and Zuber that allows diagonalizing
differential operators with respect to matrix elements of Hermitian matrices.
We detail the equivalence between the expressions obtained by diagonalizing
the partition function in differential or integral formulation, which is not
manifest at first glance. For one-matrix models, this requires transforming
certain derivatives to variables. In the case of two-matrix models, the same
computation leads to a new determinant formulation of the partition function,
and we discuss potential applications to new orthogonal polynomials methods.

## I. Introduction

Two-matrix models have been widely studied in the literature since their
introduction by Itzykson and Zuber [1], see for instance [2, 3, 4, 5, 6, 7].
Their Feynman expansions count ribbon graphs with two kinds of vertices. They
match, for different choices of potentials, the partition functions of some
statistical models such as the Ising model [8], colored triangulations [9],
hard particles [10] etc. on random surfaces (thereby relating to 2D quantum
gravity).

In quantum field theory, differential formulations such as

$\begin{split}\frac{\mathcal{Z}}{\mathcal{Z}_{0}}&=\frac{1}{\mathcal{Z}_{0}}\int\mathrm{D}\phi\,\textrm{e}^{-\int\mathrm{d}^{d}x\,\mathrm{d}^{d}y\,\frac{1}{2}\phi(x)K(x,y)\phi(y)-\int\mathrm{d}^{d}x\,V\left(\phi(x)\right)}\\\
&=\left[\textrm{e}^{-\int\mathrm{d}^{d}w\,V\left(\frac{\partial}{\partial{J(w)}}\right)}\textrm{e}^{\int\mathrm{d}^{d}x\,\mathrm{d}^{d}y\,J(x)K^{-1}(x,y)J(y)}\right]_{J=0}\,,\end{split}$
(1)

are widely used (Zee coins it the “central identity of quantum field theory”
in his famous book [11]). Here the $\phi$ are fields on a $D$-dimensional
spacetime, $\int\mathrm{D}\phi$ is a functional integration, the interaction
potential $V$ is a polynomial in the fields, $K$ is the propagator, and
$\mathcal{Z}_{0}$ is a normalization. Such formulations allow expressing
Gaussian expectation values and Wick’s formula, the key ingredient of
perturbative expansions, as well-defined algebraic expressions that do not
suffer the potential divergences of the combinatorially equivalent integrals
(see for instance [12, 13, 14, 15]). They moreover provide a practical tool
for computations of Feynman diagrams.

In particular, the differential formulation provides well-defined expressions
for the Gaussian expectation values of _two-matrix models_ , while the
integrals are not well-defined and only understood formally through Wick’s
theorem.

Remarkably, as we will show in the paper, two-matrix models can be expressed
in terms of a single matrix in a differential formulation: the two matrices
necessary to control the interaction between the two kinds of vertices (or
equivalently impose bipartiteness of the ribbon graphs) are no longer needed
in the differential formulations. Since the necessity of using two matrices
complicates the resolution of two-matrix models, the existence of a one-matrix
differential formulation raises the question whether this allows for a simpler
resolution.

With the help of the one-matrix differential reformulation, it is possible to
diagonalize the two-matrix models directly in the differential formulation.
The steps leading to the expression of the partition function as a Slater
determinant and thereby to the resolution by biorthogonal polynomials involve
elegant and interesting equations whose equivalence is not always manifest. We
detail their relations and take the opportunity to review and clarify some
aspects of the differential formulation that can be found in the matrix-model
literature, as well as the formula needed to diagonalize the differential
operators with respect to matrix elements for Hermitian matrices, implicitly
shown by Itzykson and Zuber in [1] and stated by Zuber in [16]. These
computations do not seem to indicate any simplification in the resolution of
two-matrix models in differential formulation as initially hoped: the role
played by the two matrices is now played by the matrix and the derivatives
with respect to the matrix, leading to the same resolution by biorthogonal
polynomials.

On the other hand, the analysis of the equivalence between the different
formulas for one-matrix models does provide new insights on two-matrix models.
Indeed, the proof of equivalence between two of the diagonalized differential
formulations of the one-matrix models involves transforming certain
derivatives to variables in the Slater determinant formulation. The
computation is not trivial and is detailed in the last section. The same
computation can be applied with little modifications to two-matrix models, and
leads to a new determinant formulation of their partition functions. Building
on this, we comment on a potential new resolution using orthogonal polynomials
instead of biorthogonal polynomials.

In order to show that the different expressions of the partition functions
indeed coincide, we verify that all the prefactors and normalizations match.
These constants are often omitted in the literature as they are sometimes
tedious to compute and do not contribute to the critical behaviors or the
combinatorial interpretations. An advantage of the computations in
differential formulations is that they most often do not need an overall
normalization, as seen for instance in (1).

Differential formulations are often used in the literature with little or
without justifications or references, and in addition to the developments
presented here, the present paper is also intended to provide a text of
reference by giving an exhaustive review in the context of random matrix
models.

For the benefit of the reader, we first present a diagrammatic overview of all
the different representations of the one- and two-matrix models that we
consider in this article, with lines connecting two representations when we
have a direct computation for translating between them. The expressions are
shown in a schematic way, omitting normalization terms, traces, factors of
$N$, and the like. The full form of each expression can be found in the main
text of the article. Function arguments and integration variables are omitted
where possible. If confusion may arise, the arguments of the potential
functions $V$ (for one-matrix models) and $V_{1},V_{2}$ (for two-matrix
models) are indicated by superscripts. The action is understood to be
$S=N\mathrm{Tr}(\frac{1}{2}M^{2}-V(M))$ for one-matrix models, and
$S=N\mathrm{Tr}(AB-V_{1}(A)-V_{2}(B))$ for two-matrix models. Elements of the
matrices appearing inside the determinants are labeled $i,j$. In the external
expressions (261) and (262), the $r,s,p$ are related to the $i$ of (240) and
(203) as $i=(p-1)r+s$ for a polynomial potential
$V_{1}(x)=\sum_{k}^{p}\frac{\alpha_{k}}{k}x^{k}$ of degree $p$.

The structure is the same in both diagrams. In the bottom layer, we have all
the integral representations of the matrix model under consideration. In the
top layer, we show the corresponding expressions in differential form.
$\Delta$ denotes the Vandermonde determinant. The front layers gather the
expressions obtained after diagonalizing in integral form, while the back
layers gather the expressions obtained after diagonalizing in differential
form.

Several lines connecting parts within a diagram share common features,
indicated by a common greek letter. The lines $\alpha$ symbolize the direct
equivalence between integral and differential formulations. Lines $\beta$
symbolize the equivalence between the one- and two-matrix differential
formulations. Lines $\delta$ involve a diagonalization procedure in integral
formulation, leading to the appearance of Vandermonde determinants. Similarly,
lines $\eta$ refer to a diagonalization in the differential formulation. When
moving from left to right along lines marked $\phi$, the common strategy is to
expand the Vandermonde determinants, thereby moving the integration or
differentiation inside the matrix elements of the determinant. The appearance
of $\kappa$ indicates that a heat kernel method was used to transform between
an integral and differential expression. The steps marked $\sigma$ absorb the
$1/\Delta$ prefactor by extracting the symmetric part of the determinant.
Finally, for vertical lines labeled $\lambda$, it is shown how to transform
some of the derivatives to variables and vice-versa, while $\tilde{\lambda}$
corresponds to the equivalent computation in integral formulation. These are
the computations leading to our newly found determinantal expressions for two-
matrix models, that could potentially lead to a new resolution in terms of
orthogonal polynomials.

(82) Figure 1: Diagrammatic overview of the various expressions for one-matrix
models. (185) Figure 2: Diagrammatic overview of the various expressions for
two-matrix models.

## II. Brief presentation of the models

We start with a reminder on the two models presenting all relevant details, in
particular their perturbative expansions and diagonalizations as well as the
reformulation of two-matrix models in terms of biorthogonal polynomials.

### II. 1. One-matrix models

#### (a) The model

The partition function of the Hermitian one-matrix model with potential $V$ is

$Z_{V}=\int_{\mathbb{H}_{N}}\frac{\mathrm{d}M}{a_{N}}\textrm{e}^{-N\mathrm{Tr}(\frac{1}{2}M^{2}-V(M))},$
(186)

where $\mathbb{H}_{N}$ is the set of Hermitian matrices of size $N\times N$,
the potential is $V(x)=\sum_{k\geq 1}\frac{\lambda_{k}}{k}x^{k}$, the measure
is $\mathrm{d}M=\prod_{i=1}^{N}\mathrm{d}M_{ii}\prod_{i<j}\mathrm{d}\Re
M_{ij}\mathrm{d}\Im M_{ij}$, and the normalization
$a_{N}={2^{-N/2}}(\pi/N)^{-{N^{2}}/2}$ is such that $Z_{V=0}=1$ (by explicit
computation of the Gaussian integral). The free energy of the model is
$F_{V}=\log Z_{V}$.

#### (b) Diagrammatic expansion

By expanding the exponential of the potential $V$, exchanging the summation
and the integration, and using Wick’s theorem, one obtains a Feynman expansion
of the one-matrix model over ribbon graphs, which can be seen as graphs with
an additional cyclic ordering of the edges around each vertex. Given a ribbon
graph, a cyclic sequence of edges such that every edge in the sequence also
precedes the edge that follows in the ordering around a vertex of the graph is
called a _face_ (for more details see [17, 18, 19, 20, 21]). In the expansion
of the partition function, the vertices of valency $k$ (that is, with $k$
edges attached) are counted with a weight $N\frac{\lambda_{k}}{k}$ and the
edges with a weight $1/N$, as they correspond to propagators

$\int_{\mathbb{H}_{N}}\frac{\mathrm{d}M}{a_{N}}\textrm{e}^{-\frac{N}{2}\mathrm{Tr}(M^{2})}M_{ij}M_{kl}=\frac{1}{N}\delta_{i,l}\delta_{j,k}.$
(187)

To each face corresponds a trace of these Kronecker symbols, resulting in an
additional factor $N$ per face. The Feynman expansion of the free energy
involves only connected graphs, which is not true for the partition function.
Denoting by $V_{k}(\Gamma)$ the number of vertices of valency $k$ of a ribbon
graph $\Gamma$, the weight (or amplitude) of a graph in the perturbative
expansion is therefore

$\mathcal{A}(\Gamma)=\frac{N^{\chi(\Gamma)}}{|\mathrm{Aut}(\Gamma)|}\prod_{k\geq
1}\lambda_{k}^{V_{k}(\Gamma)},$ (188)

where $\chi(\Gamma)$ is the Euler characteristic

$\chi=V-E+F=2K-2g,$ (189)

in which $V=\sum_{k}V_{k}$ is the total number of vertices, $E$ is the total
number of edges, $F$ the total number of faces, and $K$ the number of
connected components of the ribbon graph ($K=1$ for the graph expansion of the
free energy), and $g$, the genus of the ribbon graph, is a non-negative
integer. The combinatorial factor $|\mathrm{Aut}(\Gamma)|$ is the number of
automorphisms of the ribbon graph $\Gamma$ (see [17, 18, 19]).

#### (c) Clarification regarding matrix integrals

All matrix integrals in the paper are understood as formal integrals: They are
perturbative expansions labelled by Feynman diagrams, that is, formal
generating functions for ribbon graphs (see [20, 22, 19] and references
therein). Every equality between partition functions is understood as equality
between their perturbative expansions, term by term, regardless of the
analytic properties of these series. The questions of whether the matrix
integrals are well-defined, whether the series are summable when the coupling
constants are in some domain and in some certain limits (e.g. large $N$), and
whether an equality between these different quantities holds non-
perturbatively are not addressed in this paper.

#### (d) Diagonalization

One of the standard ways of solving this model is to first diagonalize the
Hermitian matrices as $M=U^{\dagger}XU$, where $U\in\textrm{U}(N)$ and
$X=\textrm{diag}(\mathbf{x})$,
$\mathbf{x}=(x_{1},\ldots,x_{N})\in\mathbb{R}^{N}$. Due to the unitary
invariance, the integration over unitary matrices factors out the volume of
the unitary group, and the partition function is expressed as an integral over
the $N$ real eigenvalues,

$Z_{V}=\int_{\mathbb{R}^{N}}\frac{\mathrm{d}\mathbf{x}}{b_{N}}\,\Delta^{2}(\mathbf{x})\,\textrm{e}^{-N(\frac{1}{2}\mathbf{x}^{2}-V(\mathbf{x}))},$
(190)

where the measure is $\mathrm{d}\mathbf{x}=\prod_{i=1}^{N}\mathrm{d}x_{i}$,
the normalization is $b_{N}={(2\pi)^{N/2}\prod_{j=1}^{N}j!}\,{N^{-N^{2}/2}}$
(this can be computed from Mehta’s integral Eq. (3.3.10) in [23]),
$\Delta(\mathbf{x})=\prod_{i<j}(x_{j}-x_{i})$ is the Vandermonde determinant,
$\mathbf{x}^{2}=\sum_{i=1}^{N}x_{i}^{2}$, and we use the notation

$V(\mathbf{x})\equiv\sum_{i=1}^{N}V(x_{i}).$ (191)

To solve the model one may apply for example the saddle point method, loop or
Schwinger-Dyson equations, or the method of orthogonal polynomials [19]. Here
we focus on the steps leading to the resolution by orthogonal polynomials
relying on the diagonalization.

### II. 2. Two-matrix models

#### (a) Partition function

The models we are interested in have the following form:

$Z_{V_{1},V_{2}}=\int_{\mathbb{H}_{N}\times\mathbb{H}_{N}}\frac{\mathrm{d}A\,\mathrm{d}B}{c_{N}}\textrm{e}^{-N\mathrm{Tr}(AB-
V_{1}(A)-V_{2}(B))}$ (192)

where $V_{1}(x)=\sum_{k\geq 1}\frac{\alpha_{k}}{k}x^{k}$ and
$V_{2}(x)=\sum_{k\geq 1}\frac{\beta_{k}}{k}x^{k}$, and the normalization
$c_{N}$ is chosen such that for $V_{1}=V_{2}=0$, $Z_{0,0}=1$. The Gaussian
part is formal and should be understood as follows. For two real variables, we
formally have for $\alpha\in\mathbb{R}$:

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\alpha
Nxy}=\frac{2\pi}{\imath\alpha N},$ (193)

and more generally for $n,m\geq 0$:

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\alpha
Nxy}\,x^{n}y^{m}=\frac{\delta_{n,m}}{(-N)^{n}}\frac{\partial^{n}}{\partial\alpha^{n}}\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\alpha
Nxy}=\frac{\delta_{n,m}}{(-N)^{n}}\frac{\partial^{n}}{\partial\alpha^{n}}\frac{2\pi}{\imath\alpha
N},$ (194)

so that

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\alpha
Nxy}\,x^{n}y^{m}=\delta_{n,m}\,\frac{2\pi}{\imath}\,\left(\frac{1}{\alpha
N}\right)^{n+1}n!.$ (195)

This leads to the normalization:

$c_{N}=\int_{\mathbb{H}_{N}\times\mathbb{H}_{N}}\mathrm{d}A\,\mathrm{d}B\,\textrm{e}^{-N\mathrm{Tr}(AB)}=2^{N}\left(\frac{\pi}{\imath
N}\right)^{N^{2}},$ (196)

and to the propagator

$\int_{\mathbb{H}_{N}\times\mathbb{H}_{N}}\frac{\mathrm{d}A\,\mathrm{d}B}{c_{N}}\,\textrm{e}^{-N\mathrm{Tr}(AB)}\,A_{ab}B_{cd}=\frac{1}{N}\delta_{a,d}\delta_{b,c}.$
(197)

The integrals in Eq. (193), Eq. (195), Eq. (196) and Eq. (197) are not well-
defined. On the other hand, the values we set formally here are the ones
needed in order to recover the correct combinatorial expansions with the
correct overall normalizations of the partition functions.111The formal
propagator Eq. (195) can be found in [[9], Eq. (4.2)], but we have added the
factor ${2\pi}/{\imath}$ to agree with the value of Eq. (192) when
$V_{1}(x)=V_{2}(x)=-\epsilon x^{2}/2$, in the limit where $\epsilon\rightarrow
0$ see e.g. [24], App. A. In the differential formulation, the equivalent
expressions are well-defined and automatically normalized, as will be detailed
in the rest of the text.

#### (b) Diagrammatic expansion

The Feynman graph expansion of a two-matrix model involves a summation over
ribbon graphs that have two kinds of vertices, respectively associated to the
traces of the matrices $A$ and $B$, so that the edges only link vertices of
different kinds (the ribbon graphs are bipartite 222As a convention, if there
are quadratic terms in $V_{1}$ and $V_{2}$, we interpret them as interactions
corresponding to bi-valent vertices in the graphs.). Denoting by
$V^{A}_{k}(\Gamma)$ and $V^{B}_{k}(\Gamma)$ the number of vertices of valency
$k$ of each type of vertices of a ribbon graph $\Gamma$, the weight of a
Feynman graph is given by

$\mathcal{A}(\Gamma)=\frac{N^{\chi(\Gamma)}}{|\mathrm{Aut}(\Gamma)|}\prod_{k\geq
1}\alpha_{k}^{V^{A}_{k}(\Gamma)}\beta_{k}^{V^{B}_{k}(\Gamma)},$ (198)

the graphs contributing to the free-energy being connected.

#### (c) Diagonalization

The first step in the resolution is to diagonalize $A$ and then use the
Harish-Chandra–Itzykson–Zuber formula [25, 1] to diagonalize $B$:

$\int_{\textrm{U}(N)}\mathrm{d}U\,\textrm{e}^{\alpha\mathrm{Tr}(AUBU^{\dagger})}={\alpha^{-\frac{N(N-1)}{2}}}\prod_{k=1}^{N-1}k!\
\frac{\det\\{\textrm{e}^{\alpha\,a_{i}b_{j}}\\}_{1\leq i,j\leq
N}}{\Delta(\mathbf{a})\Delta(\mathbf{b})},$ (199)

where $\mathrm{d}U$ is the normalized Haar measure on the unitary group
$\textrm{U}(N)$ and $\alpha$ is a complex coefficient. As shown in [2], for
$\alpha=N$ this leads to

$Z_{V_{1},V_{2}}=\int_{\mathbb{R}^{2N}}\frac{\mathrm{d}\mathbf{a}\,\mathrm{d}\mathbf{b}}{N!\,d_{N}}\,\Delta(\mathbf{a})\Delta(\mathbf{b})\,\textrm{e}^{-N(\mathbf{a}\cdot\mathbf{b}-V_{1}(\mathbf{a})-V_{2}(\mathbf{b}))},$
(200)

where $d_{N}=(2\pi/\imath)^{N}\prod_{j=1}^{N-1}j!\,N^{-N(N+1)/2}$.

#### (d) Determinant form

Let us detail the steps leading to the usual resolution using biorthogonal
polynomials [2]. The Leibniz determinant formula

$\det(\\{f_{i,j}\\}_{1\leq i,j\leq N})=\sum_{\sigma\in
S_{N}}(-1)^{\sigma}\prod_{i=1}^{N}\,f_{i,\sigma(i)}$ (201)

where $(-1)^{\sigma}$ is the sign of the permutation $\sigma$, is then used to
develop the Vandermonde determinants
$\Delta(\mathbf{x})=\det\\{x_{i}^{j-1}\\}_{1\leq i,j\leq N}$, so that Eq.
(200) reads:

$\displaystyle Z_{V_{1},V_{2}}$
$\displaystyle=\frac{1}{N!\,d_{N}}\sum_{\sigma,\sigma^{\prime}\in
S_{N}}(-1)^{{\sigma+\sigma^{\prime}}}\prod_{i=1}^{N}\int_{\mathbb{R}}\mathrm{d}x\int_{\mathbb{R}}\mathrm{d}y\,x^{\sigma(i)-1}y^{\sigma^{\prime}(i)-1}\,\textrm{e}^{-N(xy-
V_{1}(x)-V_{2}(y))}$ $\displaystyle=\frac{1}{d_{N}}\sum_{\sigma\in
S_{N}}(-1)^{{\sigma}}\prod_{i=1}^{N}\int\mathrm{d}x\,\mathrm{d}y\,x^{\sigma(i)-1}y^{i-1}\,\textrm{e}^{-N(xy-
V_{1}(x)-V_{2}(y))},$ (202)

leading to the following expression of the partition function as a Hankel
determinant (again by Eq. 201):

$Z_{V_{1},V_{2}}=\frac{1}{d_{N}}\det\left\\{\int\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-Nxy}x^{i}\textrm{e}^{NV_{1}(x)}y^{j}\textrm{e}^{NV_{2}(y)}\right\\}_{0\leq
i,j\leq N-1}.$ (203)

#### (e) Biorthogonal polynomials

Using the properties of the determinant, one can subtract from each row or
column a linear combination of rows or columns of lower indices without
affecting the result, so that one can replace $x^{i}$ and $y^{j}$ by any
polynomials $P_{i}(x)$, $Q_{j}(y)$ of degrees $i$ and $j$ respectively with
unit leading coefficients (such polynomials are said to be monic):

$Z_{V_{1},V_{2}}=\frac{1}{d_{N}}\det\left(\left\\{\langle
P_{i}|Q_{j}\rangle\right\\}_{0\leq i,j\leq N-1}\right),\qquad\langle
P_{i}|Q_{j}\rangle=\int\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-Nxy}P_{i}(x)\textrm{e}^{NV_{1}(x)}Q_{j}(y)\textrm{e}^{NV_{2}(y)}.$
(204)

In particular, the monic polynomials $P_{i},Q_{j}$ can be chosen to be
orthogonal for the formal symmetric bilinear form on the right-hand side of
Eq. (204), in which case the determinant Eq. (203) is the product of the
diagonal terms:

$\langle P_{i}|Q_{j}\rangle=h_{i}\delta_{i,j},\qquad\qquad
Z_{V_{1},V_{2}}=\frac{1}{d_{N}}\prod_{j=1}^{N}h_{i}.$ (205)

Two sequences of polynomials satisfying this orthogonality relation are said
to be biorthogonal. They are determined recursively for specific choices of
the potentials $V_{1}$ and $V_{2}$ [2, 3, 9, 10, 6].

## III. Differential reformulation of matrix integrals

In this section we show how to reformulate two-matrix models in terms of a
single matrix and its derivatives. We first review the usual differential
formulations of one- and two-matrix models. Then we present the new one-matrix
differential formulation of the latter. We finally show how to diagonalize the
models in this formulation.

### III. 1. The differential formulation

It is known (see e.g. [14]) that a Hermitian one-matrix model admits the
following reformulation:

$Z_{V}=\int_{\mathbb{H}_{N}}\frac{\mathrm{d}M}{a_{N}}\textrm{e}^{-N\mathrm{Tr}(\frac{1}{2}M^{2}-V(M))}=\left[\textrm{e}^{+\frac{1}{2N}\mathrm{Tr}(\frac{\partial}{\partial
M})^{2}}\textrm{e}^{N\mathrm{Tr}\,V(M)}\right]_{M=0},$ (206)

where $\mathrm{Tr}\bigl{(}\frac{\partial}{\partial
M}\bigr{)}^{2}=\sum_{a,b=1}^{N}\frac{\partial}{\partial
M_{ab}}\frac{\partial}{\partial M_{ba}}$. The role of the matrices and
derivatives can be exchanged:

$Z_{V}=\left[\textrm{e}^{N\mathrm{Tr}V(\frac{\partial}{\partial
M})}\textrm{e}^{\frac{1}{2N}\mathrm{Tr}\,(M^{2})}\right]_{M=0},$ (207)

which is the formulation analogous to Eq. (1) for a one-matrix model. Note
that the differential expressions are clearly well-normalized, as for $V=0$,
the only contributing term in the series-expansion of the exponential is 1.

As mentioned in Sec. II. 1. (c), the sign “$=$” means that the equality holds
at the perturbative level: the Feynman expansions coincide term by term. While
the perturbative expansion of the matrix model is obtained by expanding the
exponential of $V$ in series and exchanging the summation and the integration
(Sec. II. 1. (b)), in the differential formulation, the perturbative expansion
is obtained by exchanging the summation and the evaluation $M=0$.333In the
case where both sides of Eq. (206) are well-defined non-perturbatively, one
may wonder whether the two coincide.

More precisely, recalling that $V(x)=\sum_{k\geq
1}\frac{\lambda_{k}}{k}x^{k}$, the perturbative expansion of the right-hand
side of Eq. (206) reads:

$\begin{split}\left[\textrm{e}^{\frac{1}{2N}\mathrm{Tr}(\frac{\partial}{\partial
M})^{2}}\textrm{e}^{N\mathrm{Tr}\,V(M)}\right]_{M=0}&=\sum_{\\{n_{k}\geq
0\\}}\prod_{k}\frac{(N\lambda_{k}/k)^{n_{k}}}{n_{k}!}\left[\textrm{e}^{\frac{1}{2N}\mathrm{Tr}(\frac{\partial}{\partial
M})^{2}}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}}\right]_{M=0}\\\
&=\sum_{\\{n_{k}\geq
0\\}}\prod_{k}\frac{(N\lambda_{k}/k)^{n_{k}}}{n_{k}!}\sum_{i\geq
0}\frac{1}{i!}\left[\left(\frac{1}{2N}\mathrm{Tr}(\frac{\partial}{\partial
M})^{2}\right)^{i}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}}\right]_{M=0}\\\
&=\sum_{\\{n_{k}\\}}\prod_{k}\frac{(N\lambda_{k}/k)^{n_{k}}}{n_{k}!}\frac{1}{({\ell(\mathbf{n})}/2)!}\Biggl{[}\left(\frac{1}{2N}\mathrm{Tr}(\frac{\partial}{\partial
M})^{2}\right)^{\frac{\ell(\mathbf{n})}{2}}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}}\Biggr{]},\end{split}$
(208)

where $\ell(\mathbf{n})=\sum_{k}kn_{k}$. Developing the term originating from
the propagator:

$\left(\mathrm{Tr}\frac{\partial}{\partial M}\frac{\partial}{\partial
M}\right)^{\frac{\ell(\mathbf{n})}{2}}=\sum_{p=1}^{{\ell(\mathbf{n})}/2}\sum_{i_{p},j_{p}=1}^{N}\prod_{p=1}^{{\ell(\mathbf{n})}/2}\frac{\partial}{\partial
M_{i_{p}j_{p}}}\frac{\partial}{\partial M_{j_{p}i_{p}}},$ (209)

we see that the expression between brackets in the last line of Eq. (208) is a
sum over all possible ways to pair the derivatives $\frac{\partial}{\partial
M_{ij}}$ and the matrix elements $M_{kl}$ originating from the interaction
potential, with:

$\frac{1}{2N}\mathrm{Tr}\frac{\partial}{\partial
M_{i_{p}j_{p}}}\frac{\partial}{\partial
M_{j_{p}i_{p}}}M_{ab}M_{cd}=\frac{1}{2N}(\delta_{i_{p},a}\delta_{j_{p},b}\delta_{j_{p},c}\delta_{i_{p},d}+\delta_{i_{p},c}\delta_{j_{p},d}\delta_{j_{p},a}\delta_{i_{p},b}).$
(210)

Since each $i_{p}$ and $j_{p}$ appear only in one term of the sort, the sums
over $i_{p},j_{p}$ can be carried out with the following simplifications:

$\sum_{i_{p},j_{p}}\frac{1}{2N}\mathrm{Tr}\frac{\partial}{\partial
M_{i_{p}j_{p}}}\frac{\partial}{\partial
M_{j_{p}i_{p}}}M_{ab}M_{cd}=\frac{1}{N}\delta_{a,d}\delta_{b,c},$ (211)

where we recognize the propagator Eq. (187). The combinatorial factor
$\frac{1}{({\ell(\mathbf{n})}/2)!}$ compensates for the number of permutations
of the traces of pairs of derivatives (the index $p$), as all permutations
give the same result. In other words, the term in square brackets in the last
line of Eq. (208) coincides with a Gaussian expectation through Wick’s
theorem:

$\frac{1}{({\ell(\mathbf{n})}/2)!}\left(\frac{1}{2N}\mathrm{Tr}\frac{\partial}{\partial
M}\frac{\partial}{\partial
M}\right)^{{\ell(\mathbf{n})}/2}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}}=\int_{\mathbb{H}_{N}}\frac{\mathrm{d}M}{a_{N}}\textrm{e}^{-\frac{N}{2}\mathrm{Tr}(M^{2})}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}},$
(212)

which is expressed as a sum over ribbon graphs with $n_{k}$ vertices of
valency $k$ and ${\ell(\mathbf{n})}/2=\frac{1}{2}{\sum_{k}kn_{k}}$ edges. We
thus recover the usual perturbative expansion of Eq. (186), which proves Eq.
(206) at the perturbative level. Eq. (207) is shown the same way, with:

$\frac{1}{({\ell(\mathbf{n})}/2)!}\prod_{k}\left(\mathrm{Tr}\left(\frac{\partial}{\partial
M}\right)^{k}\right)^{n_{k}}\left(\frac{1}{2N}\mathrm{Tr}(M^{2})\right)^{\frac{\ell(\mathbf{n})}{2}}=\int_{\mathbb{H}_{N}}\frac{\mathrm{d}M}{a_{N}}\textrm{e}^{-\frac{N}{2}\mathrm{Tr}(M^{2})}\prod_{k}\mathrm{Tr}(M^{k})^{n_{k}}.$
(213)

### III. 2. Differential formulation of two-matrix models

For a two-matrix model the standard differential formulation is simply

$Z_{V_{1},V_{2}}=\int_{\mathbb{H}_{N}^{2}}\frac{\mathrm{d}A\,\mathrm{d}B}{c_{N}}\,\textrm{e}^{-N\mathrm{Tr}(AB-
V_{1}(A)-V_{2}(B))}=\left[\textrm{e}^{\frac{1}{N}\mathrm{Tr}\,(\frac{\partial}{\partial
A}\frac{\partial}{\partial
B})}\textrm{e}^{N\mathrm{Tr}\left(V_{1}(A)+V_{2}(B)\right)}\right]_{A,B=0}.$
(214)

Again, the differential expression is clearly well-normalized. Eq. (214) is
shown perturbatively using the fact that now the propagator corresponds to:

$\sum_{i_{p},j_{p}}\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A_{i_{p}j_{p}}}\frac{\partial}{\partial
B_{j_{p}i_{p}}}A_{ab}B_{cd}=\frac{1}{N}\delta_{a,d}\delta_{b,c},$ (215)

so that through Wick’s theorem, denoting by
$\ell(\mathbf{a}+\mathbf{b})=\sum_{k}k(a_{k}+b_{k})$, we have

$\frac{1}{\frac{\ell(\mathbf{a}+\mathbf{b})}{2}!}\left(\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}\right)^{\frac{\ell(\mathbf{a}+\mathbf{b})}{2}}\prod_{k}\mathrm{Tr}(A^{k})^{a_{k}}\mathrm{Tr}(B^{k})^{b_{k}}=\int\frac{\mathrm{d}A\mathrm{d}B}{c_{N}}\textrm{e}^{-N\mathrm{Tr}(AB)}\prod_{k}\mathrm{Tr}(A^{k})^{a_{k}}\mathrm{Tr}(B^{k})^{b_{k}},$
(216)

the other steps being the same as for one-matrix models. This last quantity is
expressed as a sum over bipartite ribbon graphs with $a_{k}$ and $b_{k}$
vertices of valency $k$ associated to the matrix $A$ respectively $B$ and
$\sum_{k}ka_{k}=\sum_{k}kb_{k}$ edges.

While the right-hand side of Eq. (216) is formal and understood via Wick’s
theorem and the formal propagator Eq. (197), the left-hand side is well-
defined.

In the differential formulation however, it is no longer necessary to use two
matrices: this is important in the integral formulation to impose
bipartiteness of the Feynman ribbon graphs. This bipartiteness can also be
implemented in a differential formulation using the fact that derivatives only
act on matrices. The resulting differential formulation thereby only involves
a single matrix:

$\boxed{Z_{V_{1},V_{2}}=\left[\textrm{e}^{N\mathrm{Tr}\,V_{1}(\frac{1}{\sqrt{N}}\frac{\partial}{\partial
M})}\textrm{e}^{N\mathrm{Tr}\,V_{2}(\frac{1}{\sqrt{N}}M)}\right]_{M=0}.}$
(217)

This can be proven by showing that this differential formula generates the
same ribbon graphs, together with the same combinatorial weights Eq. (198). It
can also be proven in the differential formulation starting from the right-
hand side of Eq. (214):

$Z_{V_{1},V_{2}}=\left[\textrm{e}^{\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}}\,\textrm{e}^{N\mathrm{Tr}\,V_{1}(A)}\,\textrm{e}^{N\mathrm{Tr}\,V_{2}(B)}\right]_{A=B=0}=\biggl{[}\left[\textrm{e}^{\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}}\,\textrm{e}^{N\mathrm{Tr}\,V_{1}(A)}\right]_{A=0}\textrm{e}^{N\mathrm{Tr}\,V_{2}(B)}\biggr{]}_{B=0},$
(218)

where, recalling that $\ell(\mathbf{a})=\sum_{k}ka_{k}$:

$\displaystyle\left[\textrm{e}^{\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}}\,\textrm{e}^{N\mathrm{Tr}\,V_{1}(A)}\right]_{A=0}$
$\displaystyle=\sum_{\\{a_{k}\geq 0\\}}\sum_{i\geq
0}\prod_{k}\frac{(N\alpha_{k}/k)^{a_{k}}}{a_{k}!}\frac{1}{i!}\left[\Bigl{(}\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}\Bigr{)}^{i}\prod_{k}(\mathrm{Tr}A^{k})^{a_{k}}\right]_{A=0}$ (219)
$\displaystyle=\sum_{\\{a_{k}\geq
0\\}}\prod_{k}\frac{(N\alpha_{k}/k)^{a_{k}}}{a_{k}!}\frac{1}{\ell(\mathbf{a})!}\Bigl{(}\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}\Bigr{)}^{\ell(\mathbf{a})}\prod_{k}(\mathrm{Tr}A^{k})^{a_{k}}.$

Assuming that we have proven that:

$\Bigl{(}\mathrm{Tr}\frac{\partial}{\partial A}\frac{\partial}{\partial
B}\Bigr{)}^{\ell(\mathbf{a})}\prod_{k}\Bigl{(}\mathrm{Tr}\left(A^{k}\right)\Bigr{)}^{a_{k}}=\ell(\mathbf{a})!\prod_{k}\left(\mathrm{Tr}\left(\frac{\partial}{\partial
B}\right)^{k}\right)^{a_{k}},$ (220)

we obtain

$\left[\textrm{e}^{\frac{1}{N}\mathrm{Tr}\frac{\partial}{\partial
A}\frac{\partial}{\partial
B}}\,\textrm{e}^{N\mathrm{Tr}\,V_{1}(A)}\right]_{A=0}=\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\partial}{\partial
B})},$ (221)

and thereby

$Z_{V_{1},V_{2}}=\left[\textrm{e}^{N\mathrm{Tr}\,V_{1}(\frac{1}{N}\frac{\partial}{\partial
B})}\textrm{e}^{N\mathrm{Tr}\,V_{2}(B)}\right]_{B=0}.$ (222)

The sought formula Eq. (217) is obtained by change of variable
$M=\sqrt{N}B$.444Note that in this differential formulation, changes of
variables are made without “Jacobian”, since both sides are clearly
normalized.

Let us now prove Eq. (220). We write

$\left(\mathrm{Tr}\frac{\partial}{\partial A}\frac{\partial}{\partial
B}\right)^{\ell(\mathbf{a})}=\sum_{p=1}^{\ell(\mathbf{a})}\sum_{i_{p},j_{p}=1}^{N}\prod_{p=1}^{\ell(\mathbf{a})}\,\frac{\partial}{\partial
A_{i_{p}j_{p}}}\frac{\partial}{\partial B_{j_{p}i_{p}}}.$

For every value taken by $\\{i_{p},j_{p}\\}_{p}$, the quantity
$\prod_{p=1}^{\ell(\mathbf{a})}\frac{\partial}{\partial
A_{i_{p}j_{p}}}\frac{\partial}{\partial
B_{j_{p}i_{p}}}\prod_{k}\left(\mathrm{Tr}\left(A^{k}\right)\right)^{a_{k}}$ is
the sum over all possible ways of assigning all the derivatives
$\frac{\partial}{\partial B_{j_{p}i_{p}}}\frac{\partial}{\partial
A_{i_{p}j_{p}}}$ to the matrix elements $A_{kl}$ in the product of traces.
Every time a term of the form $\frac{\partial}{\partial
B_{j_{p}i_{p}}}\frac{\partial}{\partial A_{i_{p}j_{p}}}$ acts on a matrix
element $A_{kl}$, the latter is replaced in the product of traces by
$\frac{\partial}{\partial B_{lk}}\delta_{i_{p},k}\delta_{j_{p},l}.$ All the
indices $i_{p},j_{p}$ belong to a single Kronecker symbol, so that the sum
$\sum_{p=1}^{\ell(\mathbf{a})}\sum_{i_{p},j_{p}=1}^{N}$ simply eliminates all
Kronecker symbols. We are left with

$\Bigl{(}\mathrm{Tr}\frac{\partial}{\partial A}\frac{\partial}{\partial
B}\Bigr{)}^{\ell(\mathbf{a})}\prod_{k}\mathrm{Tr}(A^{k})^{a_{k}}=\sum_{\begin{subarray}{c}{\textrm{possible}}\\\
{\textrm{
assignments}}\end{subarray}}\prod_{k}\left(\mathrm{Tr}\Bigl{(}\bigl{(}\frac{\partial}{\partial
B}\bigr{)}^{T}\Bigr{)}^{k}\right)^{a_{k}}=\ell(\mathbf{a})!\,\prod_{k}\left(\mathrm{Tr}\bigl{(}\frac{\partial}{\partial
B}\bigr{)}^{k}\right)^{a_{k}},$ (223)

where the last equality holds because there are $(\sum_{k}ka_{k})!$ ways of
assigning the derivatives $\frac{\partial}{\partial
B_{j_{p}i_{p}}}\frac{\partial}{\partial A_{i_{p}j_{p}}}$ to the matrix
elements $A_{kl}$ in the product of traces, and the result for each assignment
is independent of the specific assignment chosen. This concludes the proof.

### III. 3. Diagonalization in the differential formulation

It is known [1, 16] that for Hermitian matrices $A$, the differential operator
$D_{k}(\frac{\partial}{\partial A})=\mathrm{Tr}(\frac{\partial}{\partial
A})^{k}$ acting on $U(N)$-invariant functions555One may of course wonder if
the developments of the paper still hold for matrix models involving real
symmetric matrices. The diagonalization would require a formula analogous to
Eq. (224), but for the derivatives of a real symmetric matrix-valued
$O(N)$-invariant function. Such a formula is not known, because the equivalent
of the HCIZ formula Eq. (199) with the integration over $U(N)$ matrices
replaced with an integration over $O(N)$ matrices is known for _skew-_
symmetric matrices $A,B$, whereas the case of symmetric matrices is much more
involved [16].  diagonalizes as

$D_{k}\Bigl{(}\frac{\partial}{\partial
a_{i}}\Bigr{)}=\frac{1}{\Delta(\mathbf{a})}\sum_{i=1}^{N}\frac{\partial^{k}}{\partial
a_{i}^{k}}\Delta(\mathbf{a}).$ (224)

This implies in particular that the differential operator in Eq. (217) has the
diagonal form

$\mathrm{Tr}\,V_{1}\left(\frac{1}{\sqrt{N}}\frac{\partial}{\partial
M}\right)=\sum_{k\geq
1}\frac{\alpha_{k}}{k\sqrt{N}^{k}}D_{k}\Bigl{(}\frac{\partial}{\partial
M}\Bigr{)}=\frac{1}{\Delta({\mathbf{m}})}\sum_{i=1}^{N}V_{1}\left(\frac{1}{\sqrt{N}}\frac{\partial}{\partial
m_{i}}\right)\Delta({\mathbf{m}}),$ (225)

since it acts on $\exp(V_{2}(M))$, which depends on traces of $M$ and is thus
$U(N)$-invariant.

Since this is only briefly mentioned in the final remarks of [16], let us make
the argument more precise here. One defines the operator $\widehat{D}_{k}$ for
some positive integer $k$ as the operator which acts on the exponential trace
of the product of two Hermitian matrices $A$ and $B$ as

$\widehat{D}_{k}\textrm{e}^{\alpha\mathrm{Tr}AB}\equiv\alpha^{k}\mathrm{Tr}(B^{k})\textrm{e}^{\alpha\mathrm{Tr}AB}$
(226)

and $\alpha$ is a complex coefficient. From this definition it follows that
this operator has the explicit representation as a functional of matrix
derivatives

$D_{k}\Bigl{(}\frac{\partial}{\partial
A}\Bigr{)}=\mathrm{Tr}\Bigl{(}\frac{\partial}{\partial
A}\Bigr{)}^{k}=\sum_{i_{1},...,i_{k}}\frac{\partial}{\partial
A_{i_{1}i_{2}}}\frac{\partial}{\partial
A_{i_{2}i_{3}}}\dots\frac{\partial}{\partial A_{i_{k}i_{1}}}\,$ (227)

on functions of Hermitian matrices $f(A)$ whose function space is spanned by
the basis $\\{\textrm{e}^{\alpha\mathrm{Tr}AB}\\}_{B\in\mathbb{H}_{N}}$. For
$\alpha=\imath$, any Hermitian matrix-valued function belongs to that space by
Fourier transformation on Hermitian matrices (see e.g. [26])

$f(A)=\frac{1}{2^{N}\pi^{N^{2}}}\int_{\mathbb{H}_{N}}\mathrm{d}B\,\tilde{f}(B)\,\textrm{e}^{\imath\mathrm{Tr}AB}\,.$
(228)

Consider now the representation of the operator $\hat{D}_{k}$ on the subspace
of $\textrm{U}(N)$-invariant functions $f(A)=f(UAU^{-1})$, or equivalently

$f(A)=\int_{\textrm{U}(N)}\mathrm{d}U\,f(UAU^{-1})=\frac{1}{2^{N}\pi^{N^{2}}}\int_{\mathbb{H}_{N}}\mathrm{d}B\tilde{f}(B)\int_{\textrm{U}(N)}\mathrm{d}U\,\textrm{e}^{\alpha\mathrm{Tr}(UAU^{-1}B)}\,.$
(229)

The domain of such invariant functions is $N$-dimensional, and the functions
can be represented on diagonal matrices. The defining equation Eq. (226) on
the invariant subspace thus has the form

$\widehat{D}_{k}\int_{\textrm{U}(N)}\mathrm{d}U\textrm{e}^{\alpha\mathrm{Tr}(UAU^{-1}B)}\equiv\alpha^{k}\mathrm{Tr}(B^{k})\int\mathrm{d}U\textrm{e}^{\alpha\mathrm{Tr}(UAU^{-1}B)}$
(230)

Using the Itzykson-Zuber integral formula Eq. (199), it follows directly from
a Laplace expansion of the determinant,

$\displaystyle\sum_{l}\Bigl{(}\frac{\partial}{\partial
a_{l}}\Bigr{)}^{k}\det(\textrm{e}^{\alpha a_{m}b_{n}})_{m,n}$
$\displaystyle=\sum_{l}\Bigl{(}\frac{\partial}{\partial
a_{l}}\Bigr{)}^{k}\sum_{j}(-1)^{l+j}\textrm{e}^{\alpha
a_{l}b_{j}}\det(\textrm{e}^{\alpha a_{m}b_{n}})_{\begin{subarray}{c}m\neq
l,\\\ n\neq j\end{subarray}}$ (231) $\displaystyle=\sum_{j}(\alpha
b_{j})^{k}\sum_{l}(-1)^{j+l}\textrm{e}^{\alpha
a_{l}b_{j}}\det(\textrm{e}^{\alpha a_{m}b_{n}})_{\begin{subarray}{c}m\neq
l,\\\ n\neq
j\end{subarray}}=\alpha^{k}\mathrm{Tr}(B^{k})\det(\textrm{e}^{\alpha
a_{m}b_{n}})_{m,n}\,,$

that the functional of eigenvalue derivatives $D_{k}(\frac{\partial}{\partial
a_{i}})$ Eq. (224) is the representation of $\widehat{D}_{k}$ on
$\textrm{U}(N)$-invariant functions in diagonalized variables, that is on
eigenvalues.

Another equivalent way to prove this formula is to see the function $f$ as a
function of the matrix elements $\\{A_{i,j}\\}$ in the canonical basis and
consider the eigenvalues $\\{a_{1},\ldots,a_{N}\\}$ of $A\in\mathbb{H}_{N}$,
together with $N(N-1)$ variables that determine uniquely $U\in\textrm{U}(N)$
such that $A=U\textrm{Diag}(\\{a_{i}\\})U^{\dagger}$, as new variables after a
change of variables (see the detailed discussion in Chapter 3 of [23]).
Differentiating $\sum_{l,k}U^{\dagger}_{il}A_{lk}U_{kj}$ with respect to the
$a_{p}$, we find $\sum_{l,k}U^{\dagger}_{il}\frac{\partial A_{lk}}{\partial
a_{p}}U_{kj}=\delta_{ij}\delta_{jp},$ which is inverted as $\frac{\partial
A_{ij}}{\partial a_{p}}=U_{ip}U^{\dagger}_{pj}.$ From this, the derivative of
any function of $A$, $\frac{\partial}{\partial a_{p}}f(\\{A_{ij}\\})$ is well-
defined, and we may verify explicitly that both the right-hand-side $D^{(1)}$
of Eq. (227) acting on a unitary invariant function expressed as Eq. (228),
and the right-hand-side $D^{(2)}$ of Eq. (224) acting on Eq. (229) expressed
in eigenvalue-variables using Eq. (199) give the same expression:

$D^{(1)}f(A)=D^{(2)}f(A)=\frac{\imath^{k}}{2^{N}\pi^{N^{2}}}\int_{\mathbb{H}_{N}}\mathrm{d}B\tilde{f}(B)\mathrm{Tr}(B^{k})\textrm{e}^{\imath\mathrm{Tr}(AB)}.$
(232)

This requires using the fact that the Fourier transform $\tilde{f}(B)$ is also
unitary invariant.

We use this formula to diagonalize the differential formulations. For a one-
matrix model, applying Eq. (224) to Eq. (206), we obtain

$Z_{V}=\left[\textrm{e}^{\frac{1}{2N}\mathrm{Tr}\frac{\partial^{2}}{\partial
M^{2}}}\textrm{e}^{NV(M)}\right]_{M=0}=\left[\frac{1}{\Delta(\mathbf{x})}\textrm{e}^{\frac{1}{2N}(\frac{\partial}{\partial\mathbf{x}})^{2}}\Delta(\mathbf{x})\,\textrm{e}^{NV(\mathbf{x})}\right]_{\mathbf{x}=0}$
(233)

where
$(\frac{\partial}{\partial\mathbf{x}})^{2}=\sum_{i}\frac{\partial^{2}}{\partial
x_{i}^{2}}$. The formula is clearly well-normalized, since
$(\frac{\partial}{\partial\mathbf{x}})^{2}\Delta(\mathbf{x})=0$. Note that the
differential reformulation of the usual integral over eigenvalues Eq. (190)
reads

$Z_{V}=\int_{\mathbb{R}^{N}}\frac{\mathrm{d}\mathbf{x}}{b_{N}}\textrm{e}^{-\frac{N}{2}\mathbf{x}^{2}}\Delta^{2}(\mathbf{x})\,\textrm{e}^{NV(\mathbf{x})}=e_{N}\left[\textrm{e}^{\frac{1}{2N}(\frac{\partial}{\partial\mathbf{x}})^{2}}\Delta^{2}(\mathbf{x})\,\textrm{e}^{NV(\mathbf{x})}\right]_{\mathbf{x}=0},$
(234)

where the normalization $e_{N}=N^{N(N-1)/2}/\prod_{j=1}^{N}j!$ is computed in
Appendix A.

For the diagonalization of two-matrix models it is crucial that we have a
differential formulation using only one matrix such as Eq. (217). It is not
possible to simultaneously diagonalize differential operators involving two
matrices such as $\mathrm{Tr}\,(\frac{\partial}{\partial
A}\frac{\partial}{\partial B})$ in Eq. (214), at least not by the arguments
used here. However, the differential formulation in terms of a single matrix
allows us to diagonalize two-matrix models.

Applying Eq. (225) to Eq. (217) yields

$\boxed{Z_{V_{1},V_{2}}=\left[\frac{1}{\Delta(\mathbf{x})}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\partial}{\partial\mathbf{x}})}\Delta(\mathbf{x})\,\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}\mathbf{x})}\right]_{\mathbf{x}=0}}\,.$
(235)

For comparison, the differential formulation of the usual integral over
eigenvalues Eq. (200) is

$Z_{V_{1},V_{2}}=e_{N}\left[\textrm{e}^{\frac{1}{N}\sum_{i=1}^{N}\frac{\partial}{\partial
a_{i}}\frac{\partial}{\partial
b_{i}}}\,\Delta(\mathbf{a})\Delta(\mathbf{b})\textrm{e}^{NV_{1}(\mathbf{a})}\,\textrm{e}^{NV_{2}(\mathbf{b})}\right]_{\mathbf{a}=\mathbf{b}=0}$
(236)

where the normalization is the same as for the one-matrix model (see Appendix
A).

The equivalence between Eq. (233) and Eq. (234) on one hand, and between Eq.
(235) and Eq. (236) on the other hand is not manifest in this form, however it
becomes clear when formulated in terms of Slater determinants, as detailed
later in Sec. IV. 1. (c) for the two-matrix case and in Sec. IV. 3. (a) for
the one-matrix model.

## IV. Expressions as Slater determinants and applications

### IV. 1. Equivalent formulations of a two-matrix model

In this first subsection, we detail the steps leading to the Slater
determinant formulation of the partition functions of two-matrix models, in
differential formulation. The same holds for the one-matrix models, as the
special case $V_{1}(x)=\frac{x^{2}}{2}$ (cf. Sec. IV. 3. (a)).

#### (a) Determinant form

Starting from Eq. (235):

$Z_{V_{1},V_{2}}=\left[\frac{1}{\Delta(\mathbf{x})}\left(\prod_{i=1}^{N}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\partial}{\partial
x_{i}})}\right)\Delta(\mathbf{x})\left(\prod_{i=1}^{N}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x_{i})}\right)\right]_{\mathbf{x}=0},$
(237)

we incorporate the products on the left and right of the Vandermonde
determinant into a Slater determinant:

$Z_{V_{1},V_{2}}=\left[\frac{1}{\Delta(\mathbf{x})}\det\left\\{\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\mathrm{d}}{\mathrm{d}x_{i}})}x_{i}^{j-1}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x_{i})}\right\\}_{1\leq
i,j\leq N}\right]_{\mathbf{x}=0}.$ (238)

We can extract the symmetric part of the determinant by using the following
identity ([27], Thm. 1.2.4 p. 24):

$\left[\frac{1}{\Delta(\mathbf{x})}\det_{i,j}\bigl{(}f_{j}(x_{i})\bigr{)}\right]_{\forall
i,\,x_{i}=x}=f_{N}\,\det_{i,j}\left(f_{j}^{(i-1)}(x)\right),$ (239)

with $f_{N}=1/\prod_{j=1}^{N-1}j!$.666The factor $(-1)^{\frac{N(N-1)}{2}}$ in
the reference is due to a different convention in the definition of the
Vandermonde determinant. This leads to the following Slater determinant:

$Z_{V_{1},V_{2}}=f_{N}\det\left\\{\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\mathrm{d}}{\mathrm{d}x})}x^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1}.$ (240)

This is again well-normalized, since for $V_{1}=V_{2}=0$ the determinant is
$\prod_{i=1}^{N-1}i!$.

#### (b) Biorthogonal polynomials

As in Sec. II. 2. (e), $\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}$ and $x^{j}$
can be replaced in Eq. (240) by any monic polynomials
$P_{i}(\frac{\mathrm{d}}{\mathrm{d}x})$, $Q_{j}(x)$ respectively of degrees
$i$ and $j$:

$Z_{V_{1},V_{2}}=f_{N}\det\left\\{\left[P_{i}\Bigl{(}\frac{\mathrm{d}}{\mathrm{d}x}\Bigr{)}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{d}{dx})}Q_{j}(x)\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1},$ (241)

the differential formulation of the bilinear form of Eq. (204) being:

$\langle
f|g\rangle=c_{1}\left[f\biggl{(}\frac{1}{\sqrt{N}}\frac{\mathrm{d}}{\mathrm{d}x}\biggr{)}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\mathrm{d}}{\mathrm{d}x})}g\biggl{(}\frac{x}{\sqrt{N}}\biggr{)}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x)}\right]_{x=0}.$
(242)

From this, one can apply the method of biorthogonal polynomials.

#### (c) Equivalence between the diagonalized formulations for two-matrix
models

Expanding the Vandermonde determinants in Eq. (236) as in Sec. II. 2. (d), we
obtain:

$Z_{V_{1},V_{2}}=e_{N}N!\det\left\\{\left[\textrm{e}^{\frac{1}{N}\frac{\partial}{\partial
a}\frac{\partial}{\partial
b}}a^{i}\textrm{e}^{NV_{1}(a)}b^{j}\textrm{e}^{NV_{2}(b)}\right]_{a=b=0}\right\\}_{0\leq
i,j\leq N-1}.$ (243)

This is seen to be equivalent to Eq. (240), since, from Eq. (221) which also
holds when $A,B$ are real variables (and seeing $N$ as a real coefficient):

$\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{d}{dx})}x^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}x)}\right]_{x=0}=N^{i/2}\left[\textrm{e}^{\frac{1}{\sqrt{N}}\frac{\partial}{\partial
a}\frac{\partial}{\partial
b}}a^{i}\textrm{e}^{NV_{1}(a)}b^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}b)}\right]_{a=b=0},$
(244)

which in turn by the change of variable $b=\sqrt{N}b^{\prime}$ is equal to
$N^{\frac{i+j}{2}}[\textrm{e}^{\frac{1}{N}\frac{\partial}{\partial
a}\frac{\partial}{\partial
b}}a^{i}\textrm{e}^{NV_{1}(a)}b^{j}\textrm{e}^{NV_{2}(b)}]_{a=b=0}.$ The
normalizations do agree since $f_{N}=e_{N}N!/N^{N(N-1)/2}$.

Eq. (243) is in turn seen to be equivalent to the Hankel determinant Eq.
(203), as from Eq. (214) for $N=1$:

$\left[\textrm{e}^{\frac{1}{N}\frac{\partial}{\partial
a}\frac{\partial}{\partial
b}}a^{i}\textrm{e}^{NV_{1}(a)}b^{j}\textrm{e}^{NV_{2}(b)}\right]_{a=b=0}=\int_{\mathbb{R}^{2}}\frac{\mathrm{d}x\,\mathrm{d}y}{c_{1}}\,\textrm{e}^{-Nxy}x^{i}\textrm{e}^{NV_{1}(x)}y^{j}\textrm{e}^{NV_{2}(y)},$
(245)

where we recall that $c_{1}=\frac{2\pi}{iN}$, so that
$e_{N}N!/(c_{1})^{N}=N^{\frac{N(N-1)}{2}}f_{N}/(c_{1})^{N}=1/d_{N}$.

#### (d) Expansion over Schur functions

Performing steps similar to Sec. II. 2. (d) but in the differential
formulation and the opposite direction, we may re-express Eq. (240) as:

$Z_{V_{1},V_{2}}=\frac{f_{N}}{N!}\left[{\Delta\biggl{(}\frac{\partial}{\partial\mathbf{x}}\biggr{)}}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\partial}{\partial\mathbf{x}})}\Delta(\mathbf{x})\,\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}\mathbf{x})}\right]_{\mathbf{x}=0}.$
(246)

As explained in ([28], Eq (2.2.4)), the expression Eq. (246) is easily seen —
by its action on Schur functions — to be a scalar product on symmetric
functions, with

$Z_{V_{1},V_{2}}=f_{N}\left\langle\exp\left(NV_{1}\left(\frac{1}{\sqrt{N}}\,\mathchoice{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\displaystyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\textstyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\scriptstyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\scriptscriptstyle\bullet$}}}}}\,\right)\right),\exp\left(NV_{2}\left(\frac{1}{\sqrt{N}}\,\mathchoice{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\displaystyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\textstyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\scriptstyle\bullet$}}}}}{\mathbin{\vbox{\hbox{\scalebox{0.5}{$\scriptscriptstyle\bullet$}}}}}\,\right)\right)\right\rangle.$
(247)

Using the orthogonality of Schur functions for this scalar product and the
Cauchy-Littlewood formula, one obtains directly an explicit double-series
expansion over Schur functions (Eq. (2.2.13) in [28]), which guarantees the
fact that $Z_{V_{1},V_{2}}$ is a tau function of the Toda hierarchy.

### IV. 2. Integration over one variable in the determinant form of two-matrix
models

In order to verify the equivalence — in differential formulation — between Eq.
(233) and Eq. (234), it is required to go to Slater determinant form, and then
transform certain derivatives to variables. We do this for general $V_{1}$, as
the computation is the same in differential formulation for general $V_{1}$ or
for $V_{1}(x)=x^{2}/2$ (corresponding to the one-matrix models). We also
detail the computation in the integral formulation for general $V_{1}$, which
is naturally more involved than for $V_{1}(x)=x^{2}/2$, in which case it is a
simple Gaussian integration. This leads to a new determinant formulation of
the partition functions of two-matrix models, and we comment on its potential
use for resolutions involving orthogonal polynomials instead of biorthogonal
polynomials.

#### (a) General computation

We start from Eq. (240):

$Z=f_{N}\det_{i,j}\left\\{\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\mathrm{d}}{\mathrm{d}x})}x^{j}\textrm{e}^{NV_{2}(\frac{x}{\sqrt{N}})}\right]_{x=0}\right\\}=f_{N}\det_{i,j}\left\\{\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}\right\\},$
(248)

where $0\leq i,j\leq N-1$, by change of variable.777The factors
$\sqrt{N}^{j-i}$ from the change of variable compensate when factorized out of
the determinant In this subsection, we assume $V_{1}$ to be a polynomial:

$V_{1}(x)=\sum_{k=1}^{p}\frac{\alpha_{k}}{k}x^{k},\qquad p\geq 1.$ (249)

We write, as in Eq. (221) but for real variables:

$\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}=\left[\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}y^{i}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}.$
(250)

We first notice that:

$\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}=(N^{1-p}\alpha_{p})^{r}S_{(p-1)r}(y)\,\textrm{e}^{NV_{1}(y)},$
(251)

where $S_{(p-1)r}(y)$ is a monic polynomial of degree $(p-1)r$.

Let us consider the row $i=(p-1)r+s$, where $s\in\\{0,1,\ldots,p-1\\}$. We can
replace in the determinant the right-hand side of Eq. (250) by
$\Bigl{[}\Bigl{[}\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}y^{s}S_{(p-1)r}(y)\textrm{e}^{NV_{1}(\frac{1}{N}y)}\Bigr{]}_{y=0}x^{j}\textrm{e}^{NV_{2}(x)}\Bigr{]}_{x=0}$
and therefore by

$Z=f_{N}\det\left\\{(N^{1-p}\alpha_{p})^{-r(i)}\left[\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}y^{s(i)}\frac{\mathrm{d}^{r(i)}}{\mathrm{d}y^{r(i)}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1}.$ (252)

We now show that:

$\left[\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}y^{s}\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}=\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}x^{j}\textrm{e}^{NV_{2}(x)}\Bigr{]}\right]_{x=0},$
(253)

that is, just as for Eq. (221) or Eq. (250), the result of the inner bracket
$[...]_{y=0}$ on the left hand side is obtained by “replacing” the $y$’s by
$\frac{\mathrm{d}}{\mathrm{d}x}$’s, but also the
$\frac{\mathrm{d}}{\mathrm{d}y}$’s by $y$’s, while reversing the left-to-right
ordering.

Developing first the right hand side, we see that:

$\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\left[x^{j}\textrm{e}^{NV_{2}(x)}\right]\right]_{x=0}=\sum_{\\{n_{k}\geq
0\\}_{k}}\prod_{k=1}^{p}\frac{(N^{1-k}\alpha_{k}/k)^{n_{k}}}{n_{k}!}\left[\frac{\mathrm{d}^{\ell(\mathbf{n})}}{\mathrm{d}x^{\ell(\mathbf{n})}}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\left[x^{j}\textrm{e}^{NV_{2}(x)}\right]\right]_{x=0},$
(254)

where $\ell(\mathbf{n})=\sum_{k=1}^{p}kn_{k}$. The bracket on the right hand
side of this equation is non-vanishing only if
$\ell(\mathbf{n})=\sum_{k=1}^{p}kn_{k}\geq r$ and $r$ of the derivatives act
on $x^{r}$, in which case:

$\left[\frac{\mathrm{d}^{\ell(\mathbf{n})}}{\mathrm{d}x^{\ell(\mathbf{n})}}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\left[x^{j}\textrm{e}^{NV_{2}(x)}\right]\right]_{x=0}=r!\binom{\ell(\mathbf{n})}{r}\left[\frac{\mathrm{d}^{\ell(\mathbf{n})-r+s}}{\mathrm{d}x^{\ell(\mathbf{n})-r+s}}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}.$
(255)

The right hand side of Eq. (253) is therefore equal to:

$\sum_{\bigl{\\{}n_{k}\geq
0\;\mid\;\ell(\mathbf{n})\geq\,r\bigr{\\}}}\prod_{k=1}^{p}\frac{(N^{1-k}\alpha_{k}/k)^{n_{k}}}{n_{k}!}\frac{\ell(\mathbf{n})!}{(\ell(\mathbf{n})-r)!}\left[\frac{\mathrm{d}^{\ell(\mathbf{n})-r+s}}{\mathrm{d}x^{\ell(\mathbf{n})-r+s}}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}.$
(256)

We now focus on the left hand side of Eq. (253):

$\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}y^{s}\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]=\sum_{q\geq
0}\frac{1}{q!}\left[\frac{\mathrm{d}^{q}}{\mathrm{d}y^{q}}y^{s}\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]\frac{\mathrm{d}^{q}}{\mathrm{d}x^{q}}=\sum_{q\geq
s}\frac{s!}{q!}\binom{q}{s}\left[\frac{\mathrm{d}^{q+r-s}}{\mathrm{d}y^{q+r-s}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]\frac{\mathrm{d}^{q}}{\mathrm{d}x^{q}},$
(257)

in which the three brackets are evaluated in $y=0$. We develop:

$\displaystyle\left[\frac{\mathrm{d}^{q+r-s}}{\mathrm{d}y^{q+r-s}}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}$
$\displaystyle=\sum_{\\{n_{k}\geq
0\\}_{k}}\prod_{k=1}^{p}\frac{(N^{1-k}\alpha_{k}/k)^{n_{k}}}{n_{k}!}\left[\frac{\mathrm{d}^{q+r-s}}{\mathrm{d}y^{q+r-s}}y^{\ell(\mathbf{n})}\right]_{y=0}$
(258) $\displaystyle=\sum_{\begin{subarray}{c}{\\{n_{k}\geq 0\\}\textrm{
s.t.}}\\\\[1.50694pt]
{\ell(\mathbf{n})=q+r-s}\end{subarray}}\prod_{k=1}^{p}\frac{(N^{1-k}\alpha_{k}/k)^{n_{k}}}{n_{k}!}\ell(\mathbf{n})!.$
(259)

Replacing this in Eq. (257) and then Eq. (257) in the left hand side of Eq.
(253) leads to the same expansion Eq. (256). We have shown so far that:
$Z_{V_{1},V_{2}}=f_{N}\det\left\\{M_{i,j}\right\\}_{0\leq i,j\leq N-1},$ (260)
where: $\displaystyle M_{(p-1)r+s\;,\;j}$
$\displaystyle=({N^{p-1}}/{\alpha_{p}})^{r}\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}x^{j}\textrm{e}^{NV_{2}(x)}\Bigr{]}\right]_{x=0}$
(261)
$\displaystyle=({N^{p-1}}/{\alpha_{p}})^{r}\int_{\mathbb{R}^{2}}\frac{\mathrm{d}x\,\mathrm{d}y}{2\pi/\imath}\,\textrm{e}^{-xy}\textrm{e}^{NV_{1}(\frac{1}{N}y)}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}x^{j}\textrm{e}^{NV_{2}(x)}\Bigr{]}.$
(262)

The equivalence between Eq. (261) and Eq. (262) is the usual one-matrix
differential formulation of two-matrix models of Sec. III. 2 but for $A,B$ two
real variables, and seeing
$x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}[x^{j}\textrm{e}^{NV_{2}(x)}]$ as
a function of $x$.

#### (b) Computation in the integral formulation

Starting from Eq. (203), we change variable as $y^{\prime}=Ny$ and include the
normalization by $2\pi/\imath$:

$Z_{V_{1},V_{2}}=f_{N}\det\left\\{\int_{\mathbb{R}^{2}}\frac{\mathrm{d}x\,\mathrm{d}y}{2\pi/\imath}\,\textrm{e}^{-xy}x^{i}\textrm{e}^{NV_{1}(\frac{1}{N}x)}y^{j}\textrm{e}^{NV_{2}(y)}\right\\}_{0\leq
i,j\leq N-1},$ (263)

where we have used the fact that
$f_{N}=(2\pi/\imath)^{N}d_{N}^{-1}N^{-\frac{N(N+1)}{2}}$. The steps up to
equation Eq. (252) are precisely the same in integral form, so we do not
detail them again. While a perturbative proof of Eq. (253) is certainly
possible, it does not seem to us that the proof above can directly be
translated in integral form (in particular the steps Eq. (257) to Eq. (259)),
so we propose a different proof. In integral form, Eq. (252) reads

$Z_{V_{1},V_{2}}=f_{N}\det\left\\{(N^{1-p}\alpha_{p})^{-r(i)}\int_{\mathbb{R}^{2}}\frac{\mathrm{d}x\,\mathrm{d}y}{2\pi/\imath}\,\textrm{e}^{-xy}y^{s(i)}\left[\frac{\mathrm{d}^{r(i)}}{\mathrm{d}y^{r(i)}}\textrm{e}^{NV_{1}(y)}\right]x^{j}\textrm{e}^{NV_{2}(x)}\right\\}_{0\leq
i,j\leq N-1}.$ (264)

By performing $r$ integrations by part:

$\int\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-xy}y^{s}\left[\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{NV_{1}(y)}\right]x^{j}\textrm{e}^{NV_{2}(x)}=(-1)^{r}\int\mathrm{d}x\,\mathrm{d}y\,\left[\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\textrm{e}^{-xy}y^{s}\right]\textrm{e}^{NV_{1}(y)}x^{j}\textrm{e}^{NV_{2}(x)},$
(265)

where the boundary terms are assumed to vanish every time due to the
exponential terms. We use the following identity:

$\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\bigl{[}(ay)^{s}\textrm{e}^{axy}\bigr{]}=\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\bigl{[}(ax)^{r}\textrm{e}^{axy}\bigr{]}.$
(266)

Indeed, using Leibniz formula:

$\displaystyle\frac{\mathrm{d}^{r}}{\mathrm{d}y^{r}}\left[(ay)^{s}\textrm{e}^{axy}\right]=\sum_{k=0}^{r}\binom{r}{k}a^{s}\frac{\mathrm{d}^{k}}{\mathrm{d}y^{k}}\left[y^{s}\right]\frac{\mathrm{d}^{r-k}}{\mathrm{d}y^{r-k}}\left[\textrm{e}^{axy}\right]=\textrm{e}^{axy}\sum_{k=0}^{\min(r,s)}\frac{1}{k!}\frac{s!}{(s-k)!}\frac{r!}{(r-k)!}\frac{(ay)^{s}(ax)^{r}}{(axy)^{k}}.$
(267)

The same expression is obtained developing
$\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\left[(ax)^{r}\textrm{e}^{axy}\right]$.
Using this identity, Eq. (265) is equal to:

$(-1)^{s}\int\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{NV_{1}(y)}\left[\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\textrm{e}^{-xy}x^{r}\right]x^{j}\textrm{e}^{NV_{2}(x)}=\int\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-xy}y^{s}\textrm{e}^{NV_{1}(y)}\left[\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}x^{j}\textrm{e}^{NV_{2}(x)}\right],$
(268)

where the equality is obtained by doing again $s$ integrations by part. We
thus recover Eq. (262).

### IV. 3. Potential application to new orthogonal polynomial method

#### (a) Equivalence between the diagonalized differential formulations for
one-matrix models

Performing the steps of Sec. IV. 1. (a), but starting from the differential
formulation Eq. (233) of one-matrix models leads to:

$Z_{V}=f_{N}\,\det\left\\{\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}x^{j}\textrm{e}^{NV(x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1}.$ (269)

This is also the differential formulation of the Slater determinant expression
Eq. (240) for a two-matrix model with $V_{1}(x)=\frac{x^{2}}{2}$.

On the other hand, performing steps similar to Sec. II. 2. (d) to express Eq.
(234) in determinant form, we obtain:

$Z_{V}=N!e_{N}\,\det\left\\{\left[\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}x^{i+j}\textrm{e}^{NV(x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1}.$ (270)

From the simpler formulation Eq. (270), one can use _orthogonal_ polynomials
satisfying

$\left[\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}P_{i}(x)P_{j}(x)\textrm{e}^{NV(x)}\right]_{x=0}=\delta_{ij}p_{i}$
(271)

instead of _biorthogonal_ polynomials Eq. (204) satisfying

$\left[P_{i}\Bigl{(}\frac{\mathrm{d}}{\mathrm{d}x}\Bigr{)}\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}Q_{j}(x)\textrm{e}^{NV(x)}\right]_{x=0}=\delta_{ij}p_{i},$
(272)

a priori leading to a simpler resolution.

The computations of Sec. IV. 2 in the differential formulation, in the case
where $V_{1}(x)=\frac{x^{2}}{2}$ ($p=2$, $\alpha_{2}=1$), prove the
equivalence between Eq. (269) and Eq. (270). Indeed for $p=2$, $r(i)=i$ so
that the coefficients $({N^{p-1}}/{\alpha_{p}})^{r(i)}=N^{i}$ in the
determinant just provide an overall factor $N^{N(N-1)/2}$, and
$f_{N}N^{\frac{N(N-1)}{2}}=N!e_{N}$.

In integral form, the proof is simpler, as the Gaussian integration for
$V_{1}(x)=\frac{x^{2}}{2}$ can be carried out explicitly after changing
variables for $x^{\prime}=Nx/\imath$.888There is an important subtlety in
changing variables from a real variable to a pure complex variable: in
addition to the Jacobian of the change of variables, an additional factor $-1$
has to be taken into account. This is detailed in Appendix B, see Eq. (306).
The problem does not appear for the one-matrix model in the computations in
differential formulation of Sec. IV. 2, or using the heat kernel Eq. (275).
Equivalently, one may recover directly the Slater determinant form of the
partition function of one-matrix models in the integral formulation, that is

$Z_{V_{2}}=\frac{N!}{b_{N}}\det\left\\{\int\mathrm{d}y\,\textrm{e}^{-N\frac{y^{2}}{2}}\textrm{e}^{NV_{2}(y)}y^{i+j}\right\\}_{0\leq
i,j\leq N-1},$ (273)

from the differential formulation of the Slater determinant Eq. (240) (after
the change of variables Eq. (248)),

$Z_{{x^{2}}/2,V_{2}}=f_{N}\det\left\\{\left[\frac{\mathrm{d}^{i}}{\mathrm{d}x^{i}}\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}x^{j}\textrm{e}^{NV_{2}(x)}\right]_{x=0}\right\\}_{0\leq
i,j\leq N-1},$ (274)

by using the _heat kernel_ formulation:

$\textrm{e}^{t\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}F(x)=\int\mathrm{d}y\,K_{t}(x-y)F(y),\qquad
K_{t}(x-y)=\frac{1}{\sqrt{4\pi t}}e^{-\frac{(x-y)^{2}}{4t}},$ (275)

for $t=\frac{1}{2N}$ and $F(x)=x^{j}\textrm{e}^{NV_{2}(x)}$. We get:

$Z_{{x^{2}}/2,V_{2}}=f_{N}\det\left\\{\sqrt{\frac{N}{2\pi}}\int\mathrm{d}y\left[\frac{\partial^{i}}{\partial
x^{i}}\textrm{e}^{-\frac{N(x-y)^{2}}{2}}\right]_{x=0}y^{j}\textrm{e}^{NV_{2}(y)}\right\\}_{0\leq
i,j\leq N-1},$ (276)

where the bracket evaluating $x$ at 0 has been moved inside the integral, as
only the derivative of the kernel depends on $x$. One can verify that

$\left[\frac{\partial^{i}}{\partial
x^{i}}\textrm{e}^{-\frac{N(x-y)^{2}}{2}}\right]_{x=0}=N^{i}\tilde{Q}_{i}(y)\textrm{e}^{-\frac{N}{2}y^{2}},$
(277)

where $\tilde{Q}_{i}(y)$ is a monic polynomial of degree $i$, so we do recover
Eq. (273) as
$f_{N}(\frac{N}{2\pi})^{\frac{N}{2}}N^{\frac{N(N-1)}{2}}=N!/b_{N}$.

One can also use the heat kernel to find the integral representation of Eq.
(233) and the one-matrix version of Eq. (238), again obtained by setting
$V_{1}(x)=\frac{x^{2}}{2}$. Applying the $N$-dimensional version of Eq. (275)
to Eq. (233) for each $x_{i}$ we find

$Z_{V}=\left(\frac{N}{2\pi}\right)^{\frac{N}{2}}\left[\frac{1}{\Delta(\mathbf{x})}\int\mathrm{d}\mathbf{y}\,\textrm{e}^{-\frac{N(\mathbf{x}-\mathbf{y})^{2}}{2}}\Delta(\mathbf{y})\,e^{NV(\mathbf{y})}\right]_{\mathbf{x}=0}.$
(278)

Similarly, in the one-matrix version of Eq. (238) we can use the heat kernel
to find

$Z_{x^{2}/2,V_{2}}=\left[\frac{1}{\Delta(\mathbf{x})}\det\left\\{\frac{1}{\sqrt{2\pi}}\int\mathrm{d}y\,\textrm{e}^{-\frac{(x_{i}-y)^{2}}{2}}y^{j-1}\textrm{e}^{NV_{2}\left(\frac{1}{\sqrt{N}}y\right)}\right\\}_{1\leq
i,j\leq N}\right]_{\mathbf{x}=0}.$ (279)

#### (b) Application to new orthogonal polynomial methods

As we have seen at the beginning of the present section, the computations of
Sec. IV. 2 in the case where $V_{1}(x)=\frac{x^{2}}{2}$, prove the equivalence
in the differential formulation between Eq. (269) and Eq. (270), and the
latter expression allows using orthogonal polynomials Eq. (271) instead of
biorthogonal polynomials Eq. (272). One may naturally wonder if for potentials
$V_{1}$ of degree higher than two, our new expressions Eq. (261), Eq. (262)
could also allow the use of orthogonal polynomials instead of biorthogonal
polynomials for two-matrix models.

To this aim, we reorganize the lines of the matrix $M$ in Eq. (261) according
to their remainder modulo $p-1$, that is:

$Z_{V_{1},V_{2}}=h_{N}\det\left\\{\tilde{M}_{i,j}\right\\}_{0\leq i,j\leq
N-1},\quad\mathrm{where}\quad\tilde{M}=\left(\begin{array}[]{c}T_{0}\\\
\hline\cr\vdots\\\ \hline\cr
T_{p-1}\end{array}\right)\quad\mathrm{and}\quad(T_{s})_{r,j}=M_{(p-1)r+s\;,\;j}.$
(280)

where $h_{N}=f_{N}\,\mathrm{sgn}(N,p)$, $\mathrm{sgn}(N,p)$ being the parity
of the permutation of rows, and we recall (Eq. (261)) that in differential
formulation,

$M_{(p-1)r+s\;,\;j}=({N^{p-1}}/{\alpha_{p}})^{r}\left[\textrm{e}^{NV_{1}(\frac{\mathrm{d}}{\mathrm{d}x})}x^{r}\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}x^{j}\textrm{e}^{NV_{2}(x)}\Bigr{]}\right]_{x=0}.$
(281)

If $s_{0}$ and $r_{0}$ are respectively the remainder and the quotient of the
Euclidean division of $N-1$ by $p-1$, then for $0\leq s\leq s_{0}$, the matrix
$T_{s}$ has $r_{0}+1$ lines, while for $s_{0}<s<p-1$, $T_{s}$ has $r_{0}$
lines.

By reorganizing the columns of $\tilde{M}$, and the lines of $T_{s}$
independently for each $s$, we may replace
$({N^{p-1}}/{\alpha_{p}})^{-r}(T_{s})_{r,j}=({N^{p-1}}/{\alpha_{p}})^{-r}M_{(p-1)r+s\;,\;j}$
in the determinant by

$\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}P_{r}^{(s)}(x)\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}Q_{j}(x)\textrm{e}^{NV_{2}(x)}\Bigr{]}\right]_{x=0}=\int\frac{\mathrm{d}x\,\mathrm{d}y}{2\pi/\imath}\,\textrm{e}^{-xy}\textrm{e}^{NV_{1}(\frac{1}{N}y)}P_{r}^{(s)}(x)\frac{\mathrm{d}^{s}}{\mathrm{d}x^{s}}\Bigl{[}Q_{j}(x)\textrm{e}^{NV_{2}(x)}\Bigr{]},$
(282)

where for $0\leq s<p-1$, the $P_{r}^{(s)}(x)$ are monic polynomials of degree
$r$, and $Q_{j}$ is a monic polynomial of degree $j$.

One may for instance choose $P_{r}^{(s)}=Q_{r}$ for all $0\leq s<p-1$. We use
the following notation for the determinant obtained this way:

$Z_{V_{1},V_{2}}=h_{N}\det\left\\{W_{i,j}\right\\}_{0\leq i,j\leq
N-1},\quad\mathrm{where}\quad W=\left(\begin{array}[]{c}J_{0}\\\
\hline\cr\vdots\\\ \hline\cr J_{p-1}\end{array}\right).$ (283)

For $s=0$, the elements $(J_{0})_{r,j}$ of the resulting matrix $J_{0}$ define
the symmetric bilinear form

$\langle Q_{r}\mid\mid
Q_{j}\rangle=\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}\textrm{e}^{NV_{2}(x)}Q_{r}(x)Q_{j}(x)\right]_{x=0}=\frac{\imath}{2\pi}\int\mathrm{d}x\,Q_{r}(x)Q_{j}(x)\textrm{e}^{NV_{2}(x)}\int\mathrm{d}y\,\textrm{e}^{-xy}\textrm{e}^{NV_{1}(\frac{1}{N}y)},$
(284)

and requiring $\langle Q_{r}\mid\mid Q_{j}\rangle=h_{r}\delta_{r,j}$ defines a
family of orthogonal polynomials for the following weight:

$\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}\textrm{e}^{NV_{2}(x)}\quad\leftrightarrow\quad\textrm{e}^{NV_{2}(x)}\int\mathrm{d}y\,\textrm{e}^{-xy}\textrm{e}^{NV_{1}(\frac{1}{N}y)}.$
(285)

In addition to the usual three-terms recurrence
$Q_{n+1}(x)=(x-\beta_{n})Q_{n}-\frac{h_{n}}{h_{n-1}}Q_{n-1}$ satisfied by any
family of orthogonal polynomials, to obtain families of recurrence relations
on the coefficients $\beta_{n}$ and $h_{n}$, one may for instance use the fact
that

$\displaystyle\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}V_{1}^{\prime}\Bigl{(}\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x}\Bigr{)}F(x)\right]_{x=0}$
$\displaystyle=\left[\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}V_{1}^{\prime}\Bigl{(}\frac{1}{N}y\Bigr{)}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}F(x)\right]_{x=0}$
(286)
$\displaystyle=\left[\left[\textrm{e}^{\frac{\mathrm{d}}{\mathrm{d}x}\frac{\mathrm{d}}{\mathrm{d}y}}\frac{\mathrm{d}}{\mathrm{d}y}\textrm{e}^{NV_{1}(\frac{1}{N}y)}\right]_{y=0}F(x)\right]_{x=0},$
(287)

so that

$\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}V_{1}^{\prime}\Bigl{(}\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x}\Bigr{)}F(x)\right]_{x=0}=\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}xF(x)\right]_{x=0},$
(288)

applied to $F(x)=Q_{r}(x)Q_{j}(x)\textrm{e}^{NV_{2}(x)}$ for any $r,j$. For
instance for the case $V_{1}(x)=x^{3}/3$, for which the weight for the
orthogonal polynomials includes the Airy function as a factor, this equation
leads to the family of equations

$\frac{2}{N}\langle V_{2}^{\prime}Q_{r}^{\prime}\mid\mid
Q_{j}\rangle+\frac{2}{N}\langle V_{2}^{\prime}Q_{r}\mid\mid
Q_{j}^{\prime}\rangle+\frac{1}{N^{2}}\Bigl{(}\langle
Q_{r}^{\prime\prime}\mid\mid Q_{j}\rangle+\langle Q_{r}\mid\mid
Q_{j}^{\prime\prime}\rangle+2\langle Q_{r}^{\prime}\mid\mid
Q_{j}^{\prime}\rangle\Bigr{)}+\langle R\,Q_{r}\mid\mid Q_{j}\rangle=0,$ (289)

with $R(x)=\frac{1}{N}V_{2}^{\prime\prime}(x)+(V_{2}^{\prime}(x))^{2}-x$.

The rest of the determinant Eq. (283) for $s>0$, obtained from Eq. (282), also
involves sums of terms of the form $\langle g\,Q_{r}\mid\mid
Q_{j}^{(k)}\rangle$ for some polynomials $g(x)$. Computing the determinant
therefore requires knowledge about the derivatives of the families of
orthogonal polynomials with weight Eq. (285). It is not clear at this point
whether the family of orthogonal polynomials can indeed be constructed and
whether the rest of the determinant can be computed this way.

On the other hand, everything seems to simplify drastically at large $N$: for
the choice $P_{r}^{(s)}=Q_{r}$ for all $0\leq s<p-1$, the leading
contributions in $N$ to the matrix elements Eq. (282) for the part of the
determinant corresponding to $s>0$ are obtained when all $s$ derivatives act
on $e^{NV_{2}(x)}$, raising a factor $N^{s}(V_{2}^{\prime}(x))^{s}$, so that
the elements of the determinant simplify to

$(J_{s})_{r,j}\sim_{N\rightarrow\infty}\frac{N^{(p-1)r+s}}{\alpha_{p}^{r}}\left\langle(V_{2}^{\prime})^{s}Q_{r}\mid\mid
Q_{j}\right\rangle.$ (290)

In the same way, the leading contribution to the family of relations Eq. (288)
seems to be given by

$\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}e^{NV_{2}(x)}V_{1}^{\prime}\bigl{(}V_{2}^{\prime}(x)\bigr{)}Q_{r}(x)Q_{j}(x)\right]_{x=0}=\left[\textrm{e}^{NV_{1}(\frac{1}{N}\frac{\mathrm{d}}{\mathrm{d}x})}xF(x)\right]_{x=0},\quad\Leftrightarrow\quad\langle
UQ_{r}\mid\mid Q_{j}\rangle=0,$ (291)

with $U(x)=V_{1}^{\prime}\bigl{(}V_{2}^{\prime}(x)\bigr{)}-x$. However, this
already goes beyond the original scope of the paper, and we leave this for
future work.

## Acknowledgements

The work of JB is supported by a Projectruimte grant of the Foundation for
Fundamental Research on Matter (FOM, now defunct), financially supported by
the Netherlands Organisation for Scientific Research (NWO). The work of JT is
funded by the German Research Foundation (DFG) primarily under the project
number 418838388 and further under Germany’s Excellence Strategy EXC
2044–390685587, Mathematics Münster: Dynamics–Geometry–Structure. LL
acknowledges support of the START-UP 2018 programme with project number
740.018.017, which is financed by the Dutch Research Council (NWO). LL thanks
R. Gurau and V. Rivasseau for mentioning references using the differential
formulation, and V. Bonzom for pointing out the work of Orlov.

##  Appendix A Normalizations of the differential formulations

#### Normalization of the differential formulation of diagonalized one-matrix
models

The normalization $e_{N}$ in Eq. (234) is computed developing the Vandermonde
determinants as in Sec. II. 2. (d):

$\begin{split}\frac{1}{e_{N}}&=\left[\textrm{e}^{\frac{1}{2N}(\frac{\partial}{\partial\mathbf{x}})^{2}}\Delta^{2}(\mathbf{x})\right]_{\mathbf{x}=0}\\\
&=N!\det_{i,j}\left\\{\left[\textrm{e}^{\frac{1}{2N}\frac{\mathrm{d}^{2}}{\mathrm{d}x^{2}}}x^{i+j}\right]_{x=0}\right\\}\\\
&=N!\det_{i,j}\Biggl{\\{}\sum_{n\geq
0}\frac{(2N)^{-n}}{n!}\left[\frac{\mathrm{d}^{2n}}{\mathrm{d}x^{2n}}x^{i+j}\right]_{x=0}\Biggr{\\}},\end{split}$
(292)

where the indices in the determinants range from 0 to $N-1$. The bracket in
the rightmost expression is non-vanishing only if $2n=i+j$, in which case
$\left[\frac{\mathrm{d}^{2n}}{\mathrm{d}x^{2n}}x^{i+j}\right]_{x=0}=(i+j)!$.
Therefore, $1/{e_{N}}=N!\det(R)(2N)^{-N(N-1)/2},$ where the $N\times N$ matrix
$R$ is such that $R_{ij}=0$ for $i+j$ odd, and $R_{ij}=(i+j)!/[(i+j)/2]!$ for
$i+j$ even. This determinant can be computed to be
$\det(R)=2^{N(N-1)/2}\prod_{j=1}^{N-1}j!$, so that
$e_{N}=N^{N(N-1)/2}/\prod_{j=1}^{N}j!.$

#### Normalization of the differential formulation of diagonalized two-matrix
models

The normalization is the same as for the one-matrix model, since:

$\displaystyle\left[\textrm{e}^{\frac{1}{N}\sum_{i=1}^{N}\frac{\partial}{\partial
a_{i}}\frac{\partial}{\partial
b_{i}}}\,\Delta(\mathbf{a})\Delta(\mathbf{b})\right]_{\mathbf{a}=\mathbf{b}=0}$
$\displaystyle=N!\det_{i,j}\left\\{\left[\textrm{e}^{\frac{1}{N}\sum_{i=1}^{N}\frac{\partial}{\partial
a}\frac{\partial}{\partial b}}\,a^{i}b^{j}\right]_{a=b=0}\right\\}$ (293)
$\displaystyle=N!\det_{i,j}\Biggl{\\{}\sum_{n\geq
0}\frac{(N)^{-n}}{n!}\left[\frac{d^{n}}{da^{n}}a^{i}\right]_{a=0}\left[\frac{d^{n}}{db^{n}}b^{j}\right]_{b=0}\Biggr{\\}}$
(294)
$\displaystyle=N!\det_{i,j}\Biggl{\\{}\delta_{i,j}i!N^{-i}\Biggr{\\}}=\frac{1}{e_{N}}.$
(295)

##  Appendix B Alternative representations using generalized heat kernels

In Sec. IV. 3. (a), we have detailed how, for one-matrix models, the Slater
determinant form Eq. (273) in the integral formulation can be recovered
directly from that in differential formulation Eq. (274) by using the heat
kernel formula Eq. (275).

In this appendix, we detail for completeness the analogous computation for
two-matrix models, that is, the integral analogues of expressions Eq. (235),
Eq. (238), and Eq. (240) using generalized heat kernels. Together these form
the back lower edge of the diagram shown in Fig. 2.

Our starting point is the differential representation Eq. (235):

$Z_{V_{1},V_{2}}=\left[\frac{1}{\Delta(\mathbf{x})}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\frac{\partial}{\partial\mathbf{x}})}\Delta(\mathbf{x})\,\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}\mathbf{x})}\right]_{\mathbf{x}=0}\,.$
(296)

In Eq. (275), we used the heat kernel $K_{t}(x)$ to transform a differential
representation of the one- model to an integral form. For two-matrix models,
we can use the following analogous identity (see e.g. [29]):

$\textrm{e}^{W(\frac{d}{d\mathbf{x}})}F(\mathbf{x})=\int_{\mathbb{R}^{N}}d\mathbf{y}\,K_{W}(\mathbf{x}-\mathbf{y})F(\mathbf{y}),$
(297)

where $W$ is a polynomial, and $K_{W}$ is a generalized heat kernel with the
following Fourier representation:

$K_{W}(\mathbf{x}-\mathbf{y})=\frac{1}{(2\pi)^{N}}\int_{\mathbb{R}^{N}}d\mathbf{p}\,\textrm{e}^{\imath\mathbf{p}\cdot(\mathbf{x}-\mathbf{y})+W(\imath\mathbf{p})}.$
(298)

If we set $W(\mathbf{x})=NV_{1}\left(\frac{1}{\sqrt{N}}\mathbf{x}\right)$ and
$F(\mathbf{x})=\Delta(\mathbf{x})e^{NV_{2}\left(\frac{1}{\sqrt{N}}\mathbf{x}\right)}$,
we can substitute Eq. (297) in Eq. (296) to find

$Z_{V_{1},V_{2}}=\frac{1}{(2\pi)^{N}}\left[\frac{1}{\Delta(\mathbf{x})}\int_{\mathbb{R}^{N}}\mathrm{d}\mathbf{y}\int_{\mathbb{R}^{N}}\mathrm{d}\mathbf{p}\,\textrm{e}^{\imath\mathbf{p}\cdot(\mathbf{x}-\mathbf{y})+NV_{1}(\frac{1}{\sqrt{N}}\imath\mathbf{p})}\Delta(\mathbf{y})\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}\mathbf{y})}\right]_{\mathbf{x}=0}.$
(299)

This is the integral analogue of Eq. (235). Subsequently, we can rewrite this
expression as a Slater determinant by expanding the Vandermonde determinant,
just like we did when rewriting Eq. (235) to Eq. (238), thus obtaining:

$Z_{V_{1},V_{2}}=\frac{1}{(2\pi)^{N}}\left[\frac{1}{\Delta(\mathbf{x})}\det\left\\{\int_{\mathbb{R}}\mathrm{d}y\int_{\mathbb{R}}\mathrm{d}p\,\textrm{e}^{\imath
p(x_{i}-y)+NV_{1}(\frac{1}{\sqrt{N}}\imath
p)}y^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}y)}\right\\}_{0\leq i,j\leq
N-1}\right]_{\mathbf{x}=0}.$ (300)

Finally, we can use identity Eq. (239) to obtain

$Z_{V_{1},V_{2}}=\frac{f_{N}}{(2\pi)^{N}}\det\left\\{\int_{\mathbb{R}}\mathrm{d}y\int_{\mathbb{R}}\mathrm{d}p\,\left[\frac{\partial^{i}}{\partial
x^{i}}\textrm{e}^{\imath
p(x-y)}\right]_{x=0}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\imath
p)}y^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}y)}\right\\}_{0\leq i,j\leq
N-1},$ (301)

where the bracket evaluating $x$ at 0 could now be moved inside the
determinant and the integral. This equation mirrors Eq. (240) using the
generalized heat kernel. Carrying out the differentiation and evaluating $x$
at 0:

$Z_{V_{1},V_{2}}=\frac{f_{N}}{(2\pi)^{N}}\det\left\\{\int_{\mathbb{R}}\mathrm{d}y\int_{\mathbb{R}}\mathrm{d}p\,(\imath
p)^{i}\textrm{e}^{-\imath py}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\imath
p)}y^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}y)}\right\\}_{0\leq i,j\leq
N-1}.$ (302)

To verify that the equation is well-normalized we may use Eq. (195) to carry
out the integration of the matrix elements. There is a subtlety however: Eq.
(193) and Eq. (195) have to be modified for $\alpha\in\imath\mathbb{R}$. For
$n=m=0$ and $\alpha\in\imath\mathbb{R}$, Eq. (193) should be replaced with:

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\imath\alpha
Nxy}=\frac{2\pi}{\alpha N},$ (303)

and for positive $n$ or $m$:

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\imath\alpha
Nxy}\,x^{n}y^{m}=\frac{\delta_{n,m}}{(-\imath
N)^{n}}\frac{\partial^{n}}{\partial\alpha^{n}}\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\imath\alpha
Nxy}=\frac{\delta_{n,m}}{(-\imath
N)^{n}}\frac{\partial^{n}}{\partial\alpha^{n}}\frac{2\pi}{\alpha N},$ (304)

so Eq. (195) has to be modified for:

$\int_{\mathbb{R}^{2}}\mathrm{d}x\,\mathrm{d}y\,\textrm{e}^{-\imath\alpha
Nxy}\,x^{n}y^{m}=\delta_{n,m}\,2\pi\imath\,\left(\frac{1}{\imath\alpha
N}\right)^{n+1}n!.$ (305)

These expressions allow showing that

$\int_{\mathbb{R}^{2}}\mathrm{d}y\,\mathrm{d}p\,\textrm{e}^{-\imath py}(\imath
p)^{i}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}\imath
p)}y^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}y)}=\imath\int_{\mathbb{R}^{2}}\mathrm{d}y\,\mathrm{d}x\,\textrm{e}^{-xy}x^{i}\textrm{e}^{NV_{1}(\frac{1}{\sqrt{N}}x)}y^{j}\textrm{e}^{NV_{2}(\frac{1}{\sqrt{N}}y)},$
(306)

in the sense that their perturbative expansions match, where Eq. (305) has to
be used to expand the left hand side and Eq. (195) for the right hand side,
which explains the factor $\imath$, while $\imath^{-1}$ would naively be
expected as a Jacobian by directly changing variables for $x=\imath p$ in the
integral. Finally, changing variables for $\sqrt{N}x^{\prime}=x$ and
$\sqrt{N}y^{\prime}=y$, we recover the integral Slater determinant form Eq.
(203), since $\frac{f_{N}}{(2\pi)^{N}}N^{N(N+1)/2}\imath^{N}=\frac{1}{d_{N}}$.

## References

  * [1] C. Itzykson and J.B. Zuber, _The planar approximation. II_ , _J. Math. Phys._ 21 (1980) 411.
  * [2] M.L. Mehta, _A method of integration over matrix variables_ , _Comm. Math. Phys._ 79 (1981) 327.
  * [3] B. Eynard, _Eigenvalue distribution of large random matrices, from one matrix to several coupled matrices_ , _Nuclear Physics_ 506 (1997) 633.
  * [4] B. Eynard, _Large-N expansion of the 2-matrix model_ , _JHEP_ 2003 (2003) 051.
  * [5] M. Bertola and B. Eynard, _Mixed correlation functions of the two-matrix model_ , _J. Phys. A: Math. Gen._ 36 (2003) 7733.
  * [6] M. Bertola, _Biorthogonal polynomials for two-matrix models with semiclassical potentials_ , _Journal of Approximation Theory_ 144 (2007) 162.
  * [7] M. Duits, A.B.J. Kuijlaars and M.Y. Mo, _The Hermitian Two Matrix Model with an Even Quartic Potential_ , _Memoirs of the American Mathematical Society_ 217 (2012) 1.
  * [8] V. Kazakov, _Ising model on a dynamical planar random lattice: Exact solution_ , _Physics Letters A_ 119 (1986) 140 .
  * [9] P. Di Francesco, B. Eynard and E. Guitter, _Coloring random triangulations_ , _Nuclear Physics B_ 516 (1998) 543 .
  * [10] J. Bouttier, P.D. Francesco and E. Guitter, _Critical and tricritical hard objects on bicolourable random lattices: exact solutions_ , _Journal of Physics A: Mathematical and General_ 35 (2002) 3821.
  * [11] A. Zee, _Quantum Field Theory in a Nutshell_ , Nutshell handbook, Princeton University Press, Princeton, NJ (2003).
  * [12] M. Salmhofer, _Renormalization: An introduction_ , Theoretical and Mathematical Physics, Springer-Verlag Berlin Heidelberg (2, 1999).
  * [13] R. Bauerschmidt, D.C. Brydges and G. Slade, _Introduction to a renormalisation group method_ , vol. 2242, Springer (2019), 10.1007/978-981-32-9593-3, [1907.05474].
  * [14] R. Gurau and T. Krajewski, _Analyticity results for the cumulants in a random matrix model_ , _Annales de l’Institut Henri Poincaré D: Combinatorics, Physics and their Interactions_ 2 (2015) 169.
  * [15] R. Gurau, V. Rivasseau and A. Sfondrini, _Renormalization: an advanced overview_ , 1401.5003.
  * [16] J.-B. Zuber, _The large-N limit of matrix integrals over the orthogonal group_ , _J. Phys. A: Math. Theor._ 41 (2008) 382001.
  * [17] A. Zvonkin, _Matrix integrals and map enumeration: An accessible introduction_ , _Mathematical and Computer Modelling_ 26 (1997) 281 .
  * [18] J. Bouttier, _Enumeration of maps_ , in _The Oxford Handbook of Random Matrix Theory_ , G. Akemann, J. Baik and P.D. Francesco, eds., (Oxford), pp. 534–560, Oxford University Press (2015).
  * [19] B. Eynard, T. Kimura and S. Ribault, _Random matrices_ , 1510.04430.
  * [20] B. Eynard, _Counting Surfaces_ , vol. 70 of _Progress in Mathematical Physics_ , Springer (2016), 10.1007/978-3-7643-8797-6.
  * [21] J. Thürigen, _Renormalization in combinatorially non-local field theories: the Hopf algebra of 2-graphs_ , 2102.12453.
  * [22] N. Orantin, _Chain of matrices, loop equations, and topological recursion_ , in _The Oxford Handbook of Random Matrix Theory_ , G. Akemann, J. Baik and P.D. Francesco, eds., (Oxford), pp. 329–355, Oxford University Press (2015).
  * [23] M.L. Mehta, _Random Matrices_ , Elsevier Science (2004).
  * [24] B. Eynard and A. Prats Ferrer, _2-Matrix versus Complex Matrix Model, Integrals over the Unitary Group as Triangular Integrals_ , _Comm. Math. Phys._ 277 (2008) 861.
  * [25] Harish-Chandra, _Differential operators on a semisimple lie algebra_ , _American Journal of Mathematics_ 79 (1957) 87.
  * [26] Zhang, J., Kieburg, M. and Forrester, P.J., _Harmonic analysis for rank-1 randomised Horn problems_ , _Letters in Mathematical Physics_ 111 (2021)  [1911.11316].
  * [27] L.K. Hua, _Harmonic Analysis of Functions of Several Complex Variables in the Classical Domains_ , Translations of Mathematical Monographs, American Mathematical Society (1963).
  * [28] A.Y. Orlov, _Tau Functions and Matrix Integrals_ , math-ph/0210012.
  * [29] D. Robinson, _Basic semigroup theory_ , _Proceedings of the Centre for Mathematics and its Applications_ (1996) 1 .