MnLargeSymbols’164 MnLargeSymbols’171

UUITP-09/21

# Virasoro constraints revisited

Luca Cassia  Department of Physics and Astronomy, Uppsala University, Box 516,
SE-75120 Uppsala, Sweden<EMAIL_ADDRESS>,  Rebecca Lodin
Department of Physics and Astronomy, Uppsala University, Box 516, SE-75120
Uppsala, Sweden<EMAIL_ADDRESS>and  Maxim Zabzine  Department
of Physics and Astronomy, Uppsala University, Box 516, SE-75120 Uppsala,
Sweden<EMAIL_ADDRESS>
###### Abstract.

We revisit the Virasoro constraints and explore the relation to the Hirota
bilinear equations. We furthermore investigate and provide the solution to
non-homogeneous Virasoro constraints, namely those coming from matrix models
whose domain of integration has boundaries. In particular, we provide the
example of Hermitean matrices with positive eigenvalues in which case one can
find a solution by induction on the rank of the matrix model.

All authors are supported in part by the grant “Geometry and Physics” from the
Knut and Alice Wallenberg foundation.

###### Contents

  1. 1 Introduction
  2. 2 Solving homogeneous Virasoro constraints
    1. 2.1 Symmetry and universality of correlators
    2. 2.2 A digression into Hirota equations
  3. 3 Solving non-homogeneous Virasoro constraints
    1. 3.1 Gaussian constraints in the hypercube
    2. 3.2 Gaussian constraints in the orthant
  4. 4 Conclusion
  5. A Special functions
  6. B ABJ-like models
  7. C Cohomological description of the constraints

##  1\. Introduction

Matrix models have established their physical significance through their
appearance in the description of quantum gravity in two dimensions. They
appear when transforming the integral over all possible geometries and
topologies to its discrete analogue, in other words the summation over random
triangulations of surfaces of arbitrary genus (as reviewed in [1]). As such
they give a powerful tool to probe both simple models of 2d quantum gravity as
well as the theory of strings propagating in low dimensional spacetimes.

The simplest example of a matrix model is given by the Hermitean 1-matrix
model. This model can be shown to satisfy Virasoro constraints, which can be
derived as Ward identities corresponding to the fact that the matrix integral
is invariant under certain reparametrizations of the integration variables.
Upon introducing a set of auxiliary “time” variables, the Virasoro constraints
can be compactly written as linear differential equations for the generating
function [2]. Is is also often the case that this generating function of
infinitely-many times acts as a $\tau$-function of some integrable hierarchy,
thus satisfying an additional set of bilinear Hirota equations. It is then
clear that matrix models form a special class of solutions to both sets of
equations: Virasoro and Hirota. However, given a generic solution to the
Virasoro constraints it is not obvious (and in general not true) that this
also solves Hirota bilinear equations. One can conclude that the moduli space
of solutions of Virasoro constraints and that of Hirota relations intersect
along a subspace containing the matrix model solutions but the precise
relation between the two is still not completely understood (see [3] for a
review).

Inspired by this classical story of integrable models, we investigate such
phenomena for generalizations of the Hermitean matrix model (HMM). Namely, we
consider Virasoro constraints derived from HMMs with arbitrary polynomial
potentials and we study their most general space of solutions. In particular,
it is well known that in the case of a Gaussian potential the constraints can
be solved combinatorially through repeated application of cut-and-join
operators to some specific initial data [4]. In the case of higher polynomial
potentials this technique can be generalized but the solution one finds is no
longer unique, and one finds a larger space of possible initial data [5, 6].

A secondary type of generalization which is natural to consider consists in
deforming the domain of integration of the matrix model. In the HMM case one
usually integrates over all Hermitean matrices, or equivalently over all
possible values of the real eigenvalues (upon restriction of the integral to
diagonal matrices). It is then interesting to ask what kind of Ward identities
one can derive if the integral is restricted to a subspace of such matrices.
For simplicity one can consider subspaces invariant under the adjoint action
of $U(N)$ in such a way that the diagonalization procedure still applies and
the integral can again be rewritten as an integral over eigenvalues. One of
the simplest such deformations is the one in which we restrict each eigenvalue
to some finite interval on the real line. This defines an integral over an
$N$-dimensional box (hypercube) with boundary. This results in a boundary
contribution to the Virasoro constraints, rendering them non-homogeneous.
However, as we will show, the resulting constraints can still be solved.

The paper is organized as follows. In Section 2 we consider homogeneous
Virasoro constraints, where we firstly review the Hermitean matrix model and
the Virasoro constraints it satisfies. The most generic form of their solution
is derived for potentials of low polynomial degrees. Secondly, we digress into
the subject of classical integrability and integrable hierarchies of equations
and investigate the conditions under which our solutions to the Virasoro
constraints also satisfy Hirota bilinear equations. In Section 3 we move on to
discuss non-homogeneous Virasoro constraints, where we provide the solution
when the eigenvalue model is allowed to have a domain of integration with non-
empty boundary. We discuss the specific examples in which the domain of
integration is taken to be an hypercube or an orthant in $\mathbb{R}^{N}$. In
Section 4 we summarize our results and discuss further directions. In the
Appendices we provide definitions for some special functions, we discuss the
case of Virasoro constraints for ABJ-like matrix models and finally we provide
a cohomological description of the constraints in terms of Lie algebra
cohomology.

##  2\. Solving homogeneous Virasoro constraints

We begin by reviewing the classical Virasoro constraints satisfied by the
Hermitean 1-matrix model, as inspired by [4] in the case of a Gaussian
potential and [6] in the case of a more general polynomial potential. We will
here extend the results of [6] by providing additional details to the
solutions of the Virasoro constraints, in particular by specifying
normalisations (i.e. empty correlators) and also providing an explicit example
for the case of cubic potential.

We consider the $\beta$-deformed Hermitean matrix model given by the
generating function

$\tau_{N}(\mathsf{a};t)=\frac{1}{N!}\int_{\Gamma^{N}}\,\prod_{i=1}^{N}\operatorname{d}x_{i}\prod_{1\leq
i<j\leq
N}|x_{i}-x_{j}|^{2\beta}\operatorname{e}^{-\sum_{i=1}^{N}V(x_{i})+\sum_{s=1}^{\infty}t_{s}\sum_{i=1}^{N}x_{i}^{s}}$
(2.1)

where $N\in\mathbb{Z}_{\geq 0}$ labels the rank of the matrices and
$\beta\in\mathbb{C}$ is a generic deformation parameter [7]. We take the
potential $V(x)$ to be a polynomial function of the form

$V(x)=\sum_{k=1}^{\mathsf{p}}\frac{\mathsf{a}_{k}x^{k}}{k}$ (2.2)

with coefficients or couplings $\mathsf{a}_{k}\in\mathbb{C}$. The contour
$\Gamma$ is chosen among all such that the integral is convergent as a
functions of the parameters $\mathsf{a}_{k}$. It is also important that all
eigenvalues $x_{i}$ have the same integration contour $\Gamma$ so that the
Weyl symmetry of the model is preserved, i.e. symmetry under permutations of
the integration variables. This type of polynomial potential can equivalently
be obtained as shift in the first $\mathsf{p}$ times

$t_{s}\mapsto t_{s}-\mathsf{a}_{s}/s,\quad\quad\quad
s=1,\dots,\mathsf{p}\leavevmode\nobreak\ .$ (2.3)

In formula (2.1) we use the condensed notations
$\mathsf{a}=\\{\mathsf{a}_{1},\dots,\mathsf{a}_{\mathsf{p}}\\}$ to denote
dependence on the polynomial coefficients and $t=\\{t_{1},t_{2},\dots\\}$ to
denote dependence on an infinite set of formal time variables. Here it should
be noted that the infinite set of times $\\{t_{s}\\}$ are formal variables, in
the sense that the power series expression of the generating function in these
variables is not required to be convergent, only the coefficients in this
expansion, i.e. the correlators, are required to converge. On the other hand,
these correlators are well defined analytical functions of the coupling
constants $\\{\mathsf{a}_{k}\\}$.

The model in (2.1) satisfies the (homogeneous) Virasoro constraints

$\left(\sum_{k=1}^{\mathsf{p}}\mathsf{a}_{k}\frac{\partial}{\partial
t_{k+n}}+\delta_{n,-1}\mathsf{a}_{1}N-L_{n}\right)\tau_{N}(\mathsf{a};t)=0,\qquad
n\geq-1$ (2.4)

with the Virasoro generators $L_{n}$ being given by

$\displaystyle L_{n>0}$ $\displaystyle=2\beta N\frac{\partial}{\partial
t_{n}}+\beta\sum_{\mu+\nu=n}\frac{\partial^{2}}{\partial t_{\mu}\partial
t_{\nu}}+(1-\beta)(n+1)\frac{\partial}{\partial
t_{n}}+\sum_{s>0}st_{s}\frac{\partial}{\partial t_{s+n}}$ (2.5) $\displaystyle
L_{0}$ $\displaystyle=\beta
N^{2}+(1-\beta)N+\sum_{s>0}st_{s}\frac{\partial}{\partial t_{s}}$
$\displaystyle L_{-1}$
$\displaystyle=Nt_{1}+\sum_{s>0}st_{s}\frac{\partial}{\partial
t_{s-1}}\leavevmode\nobreak\ ,$

In the literature it is often common to introduce an additional time $t_{0}$
and to let the generating function be proportional to
$\operatorname{e}^{Nt_{0}}$ so that derivatives w.r.t. $t_{0}$ act as
multiplications by factors of $N$. From this point of view we can further
write the generating function as a linear combination of functions at
different values of $N$ without modifying the Virasoro constrains. However,
because the $t_{0}$-derivative operator commutes with all Virasoro generators,
we can restrict ourselves to a given eigenspace with fixed eigenvalue $N$
without any loss of generality. We also remark that from the point of view of
the matrix model integral in (2.1) the parameter $N$ is required to be an
integer, however once we consider the Virasoro constraints in (2.4) this
parameter does not need to satisfy any special condition and in fact it can be
taken to be any complex number. As we will see in the end of this section, the
solution to the constraints can always be analytically continued in $N$.

Because of (2.3) we also have the additional constraints

$\left(\frac{\partial}{\partial
t_{k}}+k\frac{\partial}{\partial\mathsf{a}_{k}}\right)\tau_{N}(\mathsf{a};t)=0\,,\qquad
k=1,\dots,\mathsf{p}$ (2.6)

whose solution can be trivially written as

$\tau_{N}(\mathsf{a};t)\equiv\tau_{N}(\\{t_{s}-\mathsf{a}_{s}/s\\})$ (2.7)

We also remark that equations (2.4) and (2.6) are invariant under the action
of the symmetry

$\displaystyle\beta\to\frac{1}{\beta},\quad\quad\quad N\to-\beta
N,\quad\quad\quad
t_{s}\to-\frac{1}{\beta}t_{s},\quad\quad\quad\mathsf{a}_{k}\to-\frac{1}{\beta}\mathsf{a}_{k}\leavevmode\nobreak\
.$ (2.8)

where we regard $\beta$ and $N$ as arbitrary complex parameters. Therefore we
can expect the solution to satisfy the relation

$\tau_{N}^{\beta}(\mathsf{a};t)=\tau_{-\beta
N}^{1/\beta}\left(-\frac{1}{\beta}\mathsf{a};-\frac{1}{\beta}t\right)$ (2.9)

for an appropriate choice of initial data on both sides. In the context of AGT
correspondence [8, 9] the parameter $\beta$ is the ratio of Nekrasov
parameters $\epsilon_{1}$, $\epsilon_{2}$ and the symmetry
$\beta\leftrightarrow 1/\beta$ corresponds to the exchange
$\epsilon_{1}\leftrightarrow\epsilon_{2}$.

Starting from the Virasoro constraints in (2.4) one can recursively solve this
in order to obtain correlators $c_{\lambda}(\mathsf{a})$ defined by the
expansion

$\tau_{N}(\mathsf{a};t)=\sum_{\lambda}\frac{1}{|\operatorname{Aut}(\lambda)|}c_{\lambda}(\mathsf{a})\prod_{s\in\lambda}t_{s}\leavevmode\nobreak\
,$ (2.10)

as was shown in [6], where the correlators $c_{\lambda}(\mathsf{a})$ are only
determined up to an overall normalization factor corresponding to the empty
correlator $c_{\emptyset}(\mathsf{a})$. Here the sum runs over all integer
partitions $\lambda$, with $\operatorname{Aut}(\lambda)$ being the
automorphism group of the partition, i.e. the group of permutations that leave
the partition $\lambda$ invariant.

The case $\mathsf{p}=0$ is the most degenerate in the sense that the only
solution is the constant generating function $\tau_{N}(t)\equiv\tau_{N}(0)$
and all correlators are identically zero. This is also reflected by the fact
that the corresponding matrix model correlators are expressed by integrals
which do not converge, due to the absence of an appropriate potential.

In the case of $\mathsf{p}=1,2$, it has been shown in [4, 6] that the Virasoro
constraints can be re-summed in order to express the generating function
through cut-and-join operators. We review here the derivation of such
solutions. For $\mathsf{p}=1$ one can multiply each Virasoro constraint in
(2.4) by a weight factor $(n+1)t_{n+1}$ and re-sum them all together from
$n=0$ to $n=\infty$ to obtain a single equation

$\mathsf{a}_{1}\,D\tau_{N}(\mathsf{a}_{1};t)=W_{-1}\tau_{N}(\mathsf{a}_{1};t)\leavevmode\nobreak\
.$ (2.11)

where $D$ is the dilatation operator written as

$D=\sum_{n=1}^{\infty}nt_{n}\frac{\partial}{\partial t_{n}}$ (2.12)

and $W_{-1}$ is the cut-and-join operator given by the series

$W_{-1}=\sum_{n=0}^{\infty}(n+1)t_{n+1}L_{n}\leavevmode\nobreak\ .$ (2.13)

Observe that the dilatation operator $D$ acts diagonally on monomials with
eigenvalues given by the total degree of the monomial 111We define a notion of
$t$-degree such that $\deg(t_{s})=s$ and
$\deg(t_{a}t_{b})=\deg(t_{a})+\deg(t_{b})$. Then $D$ can equivalently be
regarded as the corresponding grading operator.. Using the fact that $W_{-1}$
is homogeneous of degree 1 w.r.t. $D$, we can solve (2.11) as follows. First
we notice that $D$ annihilates the constant term in $\tau_{N}$. With some
simple manipulations we can write

$\left(\mathsf{a}_{1}\,D-W_{-1}\right)\left[\tau_{N}(\mathsf{a}_{1};t)-\tau_{N}(\mathsf{a}_{1};0)\right]=W_{-1}\tau_{N}(\mathsf{a}_{1};0).$
(2.14)

Having removed the constant term from the generating function in the l.h.s.,
we now have that $D$ acts with non-zero eigenvalues and it becomes invertible

$\left(1-\frac{1}{\mathsf{a}_{1}}D^{-1}W_{-1}\right)\left[\tau_{N}(\mathsf{a}_{1};t)-\tau_{N}(\mathsf{a}_{1};0)\right]=\frac{1}{\mathsf{a}_{1}}D^{-1}W_{-1}\tau_{N}(\mathsf{a}_{1};0)\leavevmode\nobreak\
.$ (2.15)

Finally we use the identity $\frac{1}{1-z}=\sum_{s=0}^{\infty}z^{s}$ to invert
the operator on the l.h.s. and we get

$\displaystyle\tau_{N}(\mathsf{a}_{1};t)$
$\displaystyle=\tau_{N}(\mathsf{a}_{1};0)+\sum_{s=0}^{\infty}\left(\frac{D^{-1}W_{-1}}{\mathsf{a}_{1}}\right)^{s+1}\tau_{N}(\mathsf{a}_{1};0)$
(2.16)
$\displaystyle=\exp\left(\frac{W_{-1}}{\mathsf{a}_{1}}\right)\tau_{N}(\mathsf{a}_{1};0)\leavevmode\nobreak\
.$

where in the last step we used that $\deg(W_{-1})=1$ so that all factors of
$D^{-1}$ can be substituted with actual numbers that combine to give a factor
$1/(s+1)!$ in the sum. Notice that all the time dependence of the generating
function comes from the exponential of the cut-and-join operator, while
$\tau_{N}(\mathsf{a}_{1};0)\equiv c_{\emptyset}(\mathsf{a}_{1})$ can be
regarded as initial data.

For $\mathsf{p}=2$ we can perform a similar computation to get the cut-and-
join expression of the solution. In this case we re-sum the Virasoro
constraints in (2.4) with weight $(n+2)t_{n+2}$ from $n=-1$ to $n=\infty$ and
obtain

$\mathsf{a}_{2}\,D\,\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2};t)=\left(W_{-2}-\mathsf{a}_{1}\,L_{-1}\right)\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2};t)$
(2.17)

where

$W_{-2}=\sum_{n=-1}^{\infty}(n+2)t_{n+2}L_{n}\leavevmode\nobreak\
,\quad\quad\quad\deg(W_{-2})=2\leavevmode\nobreak\ .$ (2.18)

In this case the cut-and-join operator in the r.h.s. is not homogeneous (in
fact $\deg(L_{-1})=1$) however, using that $[W_{-2},L_{-1}]=0$, we can still
get a compact formula for the solution,

$\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2};t)=\exp\left(\frac{W_{-2}}{2\mathsf{a}_{2}}-\frac{\mathsf{a}_{1}}{\mathsf{a}_{2}}L_{-1}\right)c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})\leavevmode\nobreak\
.$ (2.19)

For $\mathsf{a}_{1}=0$, $\mathsf{a}_{2}=1$ and $\beta=1$ this matches exactly
with the original solution of [4].

For $\mathsf{p}=3$ and higher we do encounter a problem. Namely, in order to
get a dilatation operator in the l.h.s. of the equation we need to multiply by
the weight $(n+\mathsf{p})t_{n+\mathsf{p}}$ and then sum from $n+\mathsf{p}=1$
to infinity, which means $n\geq 1-\mathsf{p}$. However we immediately see that
for $n<-1$ there is no corresponding Virasoro constraint to sum. This clearly
first happens when $\mathsf{p}=3$ in which case we would need to sum starting
from the $L_{-2}$ constraint which is not defined. The best we can do in this
case is to start the sum from the lowest constarint, i.e. $n=-1$. The
resulting equation is

$\mathsf{a}_{\mathsf{p}}\hat{D}\tau_{N}(\mathsf{a};t)=\hat{W}\tau_{N}(\mathsf{a};t)$
(2.20)

with

$\hat{D}=\sum_{n=-1}^{\infty}(n+\mathsf{p})t_{n+\mathsf{p}}\frac{\partial}{\partial
t_{n+\mathsf{p}}}$ (2.21)

and

$\hat{W}=\sum_{n=-1}^{\infty}(n+\mathsf{p})t_{n+\mathsf{p}}L_{n}-\sum_{k=1}^{\mathsf{p}-1}\mathsf{a}_{k}\sum_{n=1-\mathsf{p}}^{\infty}(n+\mathsf{p})t_{n+\mathsf{p}}\frac{\partial}{\partial
t_{n+k}}-\mathsf{a}_{1}(\mathsf{p}-1)t_{\mathsf{p}-1}N\leavevmode\nobreak\ .$
(2.22)

Because $\mathsf{p}\geq 3$ the sum in (2.21) is missing some of the initial
terms, which means that the operator $\hat{D}$ has a larger kernel and it is
no longer a good grading operator. In particular, the operator $\hat{D}$ is
blind to all times $t_{s}$ with $s<\mathsf{p}-1$.

Let us consider the example $\mathsf{p}=3$. The kernel of $\hat{D}$ in this
case contains the constant term as well as all powers of $t_{1}$. As a
consequence of this fact, all correlators of the form $c_{\\{1,1,\dots,1\\}}$
cannot be fixed by the recursion and they have to be regarded as additional
initial data for the equation. For convenience we package all such missing
data in a generating function of such correlators,

$\displaystyle\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t_{1})$
$\displaystyle:=\left.\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)\right|_{t_{s>1}=0}$
(2.23) $\displaystyle=\sum_{k=0}^{\infty}c_{\\{\underbrace{\scriptstyle
1,\dots,1}_{k}\\}}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})\frac{t_{1}^{k}}{k!}$

in other words we drop all higher times from the generating function and we
only keep the series in $t_{1}$. We can now repeat the derivation of the
previous cases and find the partial solution

$\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)=\sum_{s=0}^{\infty}\left(\frac{\hat{D}^{-1}\hat{W}}{\mathsf{a}_{3}}\right)^{s}\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t_{1})$
(2.24)

Unfortunately, we can no longer give a closed formula for the series. First we
observe that while $\hat{D}$ is still invertible on the image of $\hat{W}$,
the cut-and-join operator $\hat{W}$ is a sum of non-commuting terms of
different degrees with respect to $\hat{D}$. This implies that even if we
substitute each factor of $\hat{D}^{-1}$ with its numerical value all such
terms no longer combine to give a nice combinatorial factor. We are not aware
of any nice resummation property of this series.

As we have shown in the previous paragraphs a solution to the Virasoro
constraints (2.4) can always be found combinatorially provided a certain
amount of initial data is provided. One should regard this data as boundary
conditions to fix in order to solve the PDEs obtained from the constraints.
The choice of such initial data however seems to be completely arbitrary. For
some specific choice one recovers the generating function of the matrix model
in (2.1), but for more general choices the solution the we constructed might
not have a matrix integral representation. Furthermore, there are cases in
which the solution does admit a matrix model representation but not of the
type of (2.1), for instance the ABJ-like model which we review in Appendix B.

The last step is to impose the additional equations (2.6). This will restrict
further the possible choices of initial data (as a function of the coupling
constants).

####  2.0.1. Degree $\mathsf{p}=1$

The Virasoro constraints in (2.4) imply an infinite set of recursion relations
between the coefficients of the power series of the generating function in
(2.10). The additional constraint in (2.6) gives instead an infinite set of
relations between correlators and their $\mathsf{a}_{1}$-derivatives. Of all
such relations we just need to consider the first non-trivial one, i.e. the
one we get by setting all times to zero because all other relations will
follow by our recursive solution of the Virasoro constraints. Hence we are led
to consider the ODE

$\left.\frac{\partial}{\partial
t_{1}}\tau_{N}(\mathsf{a}_{1};t)\right|_{t=0}=-\left.\frac{\partial}{\partial\mathsf{a}_{1}}\tau_{N}(\mathsf{a}_{1};t)\right|_{t=0}$
(2.25)

where setting the times to zero corresponds to restricting to the initial data
only. Explicitly we can then write this equation as

$c_{\\{1\\}}(\mathsf{a}_{1})=-\frac{\partial}{\partial\mathsf{a}_{1}}c_{\emptyset}(\mathsf{a}_{1})\leavevmode\nobreak\
.$ (2.26)

But now we can use that the Virasoro constraints already tell us how to
compute the correlator $c_{\\{1\\}}$,

$c_{\\{1\\}}(\mathsf{a}_{1})=\frac{N(\beta(N-1)+1)}{\mathsf{a}_{1}}\,c_{\emptyset}(\mathsf{a}_{1})\leavevmode\nobreak\
.$ (2.27)

Equating the r.h.s. of (2.26) to that of (2.27) we get a differential equation
for $c_{\emptyset}$ as a function of $\mathsf{a}_{1}$, which we can solve as

$c_{\emptyset}(\mathsf{a}_{1})=k_{N,\beta}\cdot\mathsf{a}_{1}^{-N(\beta(N-1)+1)}$
(2.28)

where $k_{N,\beta}$ is an integration constant which might depend on $N$ and
$\beta$ but not on the coupling. This fixes the initial data essentially
uniquely as a function of the coupling constant $\mathsf{a}_{1}$ and therefore
also the full generating function.

####  2.0.2. Degree $\mathsf{p}=2$

This case follows similarly. The additional constraint (2.6) now gives rise to
2 independent conditions, one for each coupling constant:

$c_{\\{1\\}}(\mathsf{a}_{1},\mathsf{a}_{2})=-\frac{\partial}{\partial\mathsf{a}_{1}}c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})\leavevmode\nobreak\
,$ (2.29)
$c_{\\{2\\}}(\mathsf{a}_{1},\mathsf{a}_{2})=-2\frac{\partial}{\partial\mathsf{a}_{2}}c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})\leavevmode\nobreak\
.$ (2.30)

Substituting the correlators $c_{\\{1\\}}$, $c_{\\{2\\}}$ in the l.h.s. with
the expressions we get by expanding (2.19),

$c_{\\{1\\}}(\mathsf{a}_{1},\mathsf{a}_{2})=-\frac{\mathsf{a}_{1}N}{\mathsf{a}_{2}}\,c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})\leavevmode\nobreak\
,$ (2.31)
$c_{\\{2\\}}(\mathsf{a}_{1},\mathsf{a}_{2})=\frac{N(\mathsf{a}_{2}(\beta(N-1)+1)+\mathsf{a}_{1}^{2})}{\mathsf{a}_{2}^{2}}\,c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})\leavevmode\nobreak\
,$ (2.32)

we end up with 2 independent PDEs for the function
$c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})$. The solution can be computed
exactly and it gives

$c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2})=k_{N,\beta}\cdot\mathsf{a}_{2}^{-\frac{1}{2}N(\beta(N-1)+1)}\exp\left(\frac{N\mathsf{a}_{1}^{2}}{2\mathsf{a}_{2}}\right)$
(2.33)

so that we can fix completely the dependence on the couplings up to some
integration constant $k_{N,\beta}$.

####  2.0.3. Degree $\mathsf{p}=3$

In this case we have an unbalance between the number of coupling constants and
that of non-trivial PDEs coming from the constraints. In particular, we
observe that we can write 3 independent differential equations coming from the
additional constraint (2.6), however only two of them are meaningful once we
substitute the solution of the Virasoro constraints. This is because the
Virasoro recursion does not allow us to solve for all correlators. In other
words we do not just have one function $c_{\emptyset}$ to determine but an
infinite number of functions, i.e. all the coefficients of (2.23).

The additional constraint can be split into the three equations

$-\frac{\partial}{\partial\mathsf{a}_{1}}c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})=c_{\\{1\\}}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})\leavevmode\nobreak\
,$ (2.34)
$-2\frac{\partial}{\partial\mathsf{a}_{2}}c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})=\left(-\frac{\mathsf{a}_{1}N}{\mathsf{a}_{3}}+\frac{\mathsf{a}_{2}}{\mathsf{a}_{3}}\frac{\partial}{\partial\mathsf{a}_{1}}\right)c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})\leavevmode\nobreak\
,$ (2.35)
$-3\frac{\partial}{\partial\mathsf{a}_{3}}c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})=\left(\frac{(1-\beta(N-1))N}{\mathsf{a}_{3}}-\frac{\mathsf{a}_{2}^{2}}{\mathsf{a}_{3}^{2}}\frac{\partial}{\partial\mathsf{a}_{1}}+\frac{\mathsf{a}_{1}\mathsf{a}_{2}N}{\mathsf{a}_{3}^{2}}+\frac{\mathsf{a}_{1}}{\mathsf{a}_{3}}\frac{\partial}{\partial\mathsf{a}_{1}}\right)c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})\leavevmode\nobreak\
,$ (2.36)

where in the first equation we are not able to substitute the explicit
expression of $c_{\\{1\\}}$ as it cannot be determined using Virasoro. In
fact, the correlator $c_{\\{1\\}}$ is itself part of the initial data. The
other two equations however do contain non-trivial information and can be
solved w.r.t. the couplings $\mathsf{a}_{2}$ and $\mathsf{a}_{3}$:

$c_{\emptyset}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3})=\exp\left(-\frac{\mathsf{a}_{2}\left(\mathsf{a}_{2}^{2}-6\mathsf{a}_{1}\mathsf{a}_{3}\right)N}{12\mathsf{a}_{3}^{2}}\right)\mathsf{a}_{3}^{-\frac{1}{3}(1-\beta(N-1))N}g\left(\frac{\mathsf{a}_{2}^{2}-4\mathsf{a}_{1}\mathsf{a}_{3}}{2\mathsf{a}_{3}^{4/3}}\right)\leavevmode\nobreak\
.$ (2.37)

We observe that the solution is no longer unique up to normalization, and in
fact it depends upon an arbitrary choice of function $g(z)$. Fixing this
function explicitly then fixes completely the generating function in (2.23)
using the equation in (2.34). To summarize, for $\mathsf{p}=3$ we cannot give
a unique solution to the system of constraints however we managed to reduce
all of the indeterminacy of the solution to a finite amount of information
carried by the univariate function $g$. Moreover, $g$ might depend also on $N$
and $\beta$ but we do not write this dependence explicitly.

If we want to further restrict the solution then we need more conditions to
impose on the function $g(z)$. For instance, if we insist that our solution is
of matrix model type, then we need to impose other constraints that one might
naturally obtain from matrix integrals. Namely, one can impose non-trivial
relations between correlators that express the fact that traces of finite size
matrices are not all independent. In fact, given a matrix $A$ of size $N$,
only the first $N$ traces $\operatorname{Tr}(A^{s})$, $s=1,\dots,N$ are
linearly independent while all higher traces can be expressed as polynomial
combinations of those222From the group theoretical point of view this
corresponds to the fact that the algebra of $U(N)$ has $N$ linearly
independent Casimir functions.. In the language of the eigenvalue model, this
means that all power-sum variables $p_{s}=\sum_{i}x_{i}^{s}$ of degree $s$
higher than the rank $N$ can be written in terms of $p_{1},\dots,p_{N}$ and
therefore similar relations must follow for the correlators. Let us consider
the simple case of rank $N=1$ as an example to illustrate this point. For
1-dimensional matrices we have that the trace commutes with taking powers,
therefore

$p_{s}=p_{1}^{s}\leavevmode\nobreak\ ,$ (2.38)

and therefore two correlators must be equal if their partitions have the same
size. It is straightforward then to obtain the corresponding Ward identity on
the generating function $\tau_{1}(\mathsf{a};t)$,

$\frac{\partial}{\partial
t_{s}}\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)=\left(\frac{\partial}{\partial
t_{1}}\right)^{s}\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)\leavevmode\nobreak\
.$ (2.39)

The only independent equation is that for $s=2$ which combined with the
additional constraint (2.6) becomes

$-2\frac{\partial}{\partial\mathsf{a}_{2}}\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)=\left(-\frac{\partial}{\partial\mathsf{a}_{1}}\right)^{2}\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};t)\leavevmode\nobreak\
.$ (2.40)

Applying to the solution in (2.37) we finally find the equation

$zg(z)-8g^{\prime\prime}(z)=0,\quad\quad\quad
z=\frac{\mathsf{a}_{2}^{2}-4\mathsf{a}_{1}\mathsf{a}_{3}}{2\mathsf{a}_{3}^{4/3}}$
(2.41)

which (up to rescaling of $z$) is the Airy ODE (see (A.3)), so that

$\left.g(z)\right|_{N=1}=k_{A}\operatorname{Ai}(z/2)+k_{B}\operatorname{Bi}(z/2)$
(2.42)

for arbitrary coefficients $k_{A},k_{B}$. Similar considerations apply to rank
$N>1$.

####  2.0.4. Comments about degree $\mathsf{p}>3$

The situation in degree $\mathsf{p}>3$ is similar to that of $\mathsf{p}=3$,
in the sense that the solution is not unique and the initial data spans an
infinite-dimensional space generated by all correlators whose partitions only
contain integers that are smaller than $\mathsf{p}-2$. These are exactly the
correlators that couple to monomials in times which are in the kernel of the
partial dilatation operator $\hat{D}$ in (2.21). Of the $\mathsf{p}$
additional contraints only $2$ can be used to obtain non-trivial equations for
$c_{\emptyset}(\mathsf{a})$. These equations allow to reduce the number of
functional variables of $c_{\emptyset}(\mathsf{a})$ from $\mathsf{p}$ to
$\mathsf{p}-2$.

###  2.1. Symmetry and universality of correlators

Now that we found explicit solutions to the constraints we can go back to the
symmetry equation in (2.9) and check explicitly that this is satisfied for our
generating functions. We notice first that this is indeed the case when the
solution is unique ($\mathsf{p}=1,2$) provided we rescale both sides of the
equation (2.9) by the corresponding empty correlator, which is the only piece
of initial data required. In the cases for which the solution is not unique
($\mathsf{p}\geq 3$) then more assumptions on the initial data are required.

Assuming now that the initial data for the solution is chosen to be compatible
with the symmetry $\beta\to 1/\beta$, one can show that there is a special
value of $\beta$ such that the correlators simplify and the solution becomes
homogeneous in $N$, $\mathsf{a}_{k}$ and $t_{s}$. Namely, for $\beta=-1/N$ the
generating function satisfies the property

$\frac{\tau_{N}(\mathsf{a};t)}{\tau_{N}(\mathsf{a};0)}=\frac{\tau_{N/\kappa}(\kappa\,\mathsf{a};\kappa\,t)}{\tau_{N/\kappa}(\kappa\,\mathsf{a};0)}\quad\quad\quad\text{for
}\beta=-\frac{1}{N}$ (2.43)

for an arbitrary scalar $\kappa\neq 0$. This implies that $(N,\mathsf{a},t)$
are not independent parameters of the generating function. Using this property
we can scale away the rank by choosing $\kappa=N$, hence we obtain the
universal relation

${\tau_{N}(\mathsf{a};t)}=\frac{\tau_{1}(N\mathsf{a};Nt)}{\tau_{1}(N\mathsf{a};0)}{\tau_{N}(\mathsf{a};0)}\quad\quad\quad\text{for
}\beta=-\frac{1}{N}$ (2.44)

which expresses the rank $N$ generating function as a scalar multiple of the
rank $1$ (normalized) generating function. As an immediate consequence of this
fact we have that all correlators in rank $N$ become universal in the sense
that they only depend on the size of the corresponding partition and not on
the actual shape (exactly as in the rank 1 case). Moreover, these correlators
do depend on the rank $N$ as simple functions, namely their dependence only
enters through a power $N^{\ell(\lambda)}$, where $\ell(\lambda)$ is the
length of the partition.

###  2.2. A digression into Hirota equations

In this section we analyze whether our solutions of Virasoro constraints are
actual $\tau$-functions, by which we mean generating functions in times which
satisfy a bilinear equation analogous to the Plücker identities for a finite
dimensional Grassmannian. This integrability condition is called Hirota
equation and, following work of Sato, it can be regarded as the infinite
dimensional analogue of Plücker relations for the Sato Grassmannian (for a
review see references [10, 11, 12]).

Hirota bilinear equations are written as the integral identity

$\oint{\mathrm{d}z}\exp\left(2\sum_{s=1}^{\infty}z^{s}v_{s}\right)\exp\left(-\sum_{s=1}^{\infty}\frac{z^{-s}}{s}\frac{\partial}{\partial
v_{s}}\right)\tau(u+v)\tau(u-v)=0$ (2.45)

where $u_{s}$, $v_{s}$ are two independent sets of times, $z$ is an auxiliary
complex variable and $\tau$ is a generic generating function. The integration
in $z$ is used to project out the coefficient of $z^{-1}$ in the Laurent
series expansion of the integrand. If this residue is zero then we say that
$\tau$ is a $\tau$-function.

Equation (2.45) is a rather compact way to express Hirota relation but not
very explicit in terms of what it implies for the correlators. In order to get
a better understanding of its meaning we can expand the generating function
$\tau$ as a power series in the times $v_{s}$ as

$\tau(u\pm v)=\sum_{d=0}^{\infty}[\tau(u\pm
v)]_{(d)}=\sum_{d=0}^{\infty}\left[\sum_{\lambda\vdash
d}\frac{1}{|\mathrm{Aut}(\lambda)|}\tau_{\lambda}(u)\prod_{l\in\lambda}\left(\pm
v_{l}\right)\right]\leavevmode\nobreak\ ,$ (2.46)

where $[\tau(u\pm v)]_{(d)}$ is used to denote the part of $\tau(u\pm v)$
which is of degree $d$ in the times $v_{s}$ and $\lambda\vdash d$ denotes that
$\lambda$ is an integer partition of $d\in\mathbb{Z}_{\geq 0}$. Here we also
use the notation
$\tau_{\lambda}(u)=\prod_{a\in\lambda}\frac{\partial}{\partial u_{a}}\tau(u)$.
Similarly, we expand the exponentials using Cauchy’s identity

$\exp\left(\sum_{s=1}^{\infty}\frac{z^{s}p_{s}}{s}\right)=\sum_{m=0}^{\infty}z^{m}\mathsf{Schur}_{\\{m\\}}(p_{s})\leavevmode\nobreak\
,$ (2.47)

where $\mathsf{Schur}_{\\{m\\}}(p_{s})$ is the symmetric Schur polynomial of
order $m$. We can therefore rewrite (2.45) as

$\displaystyle\frac{1}{2\pi\operatorname{i}}\oint{\mathrm{d}z}\sum_{k=0}^{\infty}\sum_{l=0}^{\infty}z^{k-l}\mathsf{Schur}_{\\{k\\}}\left(p_{s}=2sv_{s}\right)\mathsf{Schur}_{\\{l\\}}\left(p_{s}=-\frac{\partial}{\partial
v_{s}}\right)\tau(u+v)\tau(u-v)=$ (2.48)
$\displaystyle=\sum_{k=0}^{\infty}\mathsf{Schur}_{\\{k\\}}\left(p_{s}=2sv_{s}\right)\mathsf{Schur}_{\\{k+1\\}}\left(p_{s}=-\frac{\partial}{\partial
v_{s}}\right)\tau(u+v)\tau(u-v)=$
$\displaystyle=\sum_{d=0}^{\infty}\left[\sum_{k=0}^{d}\mathsf{Schur}_{\\{k\\}}\left(p_{s}=2sv_{s}\right)\mathsf{Schur}_{\\{k+1\\}}\left(p_{s}=-\frac{\partial}{\partial
v_{s}}\right)\sum_{j=0}^{d+1}[\tau(u+v)]_{(j)}[\tau(u-v)]_{(d+1-j)}\right]=0\leavevmode\nobreak\
.$

The expression inside the square brackets in the last line represents the
projection of (2.45) onto its degree $d$ part w.r.t. the times $v_{s}$. Since
each degree $d$ must be independent from the others, we get an infinite set of
relations. Moreover, by explicit computation one can check that for degree
$0\leq d<3$ this equation is trivially satisfied, while the first non-trivial
condition on the $\tau_{\lambda}(u)$ appears at degree $d=3$.

In order to get an equation for the correlators we need to set the times
$u_{s}$ to zero so that

$\tau_{\lambda}(0)\equiv c_{\lambda}\leavevmode\nobreak\ .$ (2.49)

For instance, in degree $d=3$ we get the identity

$\begin{split}3c_{\\{1,1\\}}^{2}&+3c_{\emptyset}c_{\\{2,2\\}}-4c_{\emptyset}c_{\\{3,1\\}}-4c_{\\{1\\}}c_{\\{1,1,1\\}}+c_{\emptyset}c_{\\{1,1,1,1\\}}-3c_{\\{2\\}}^{2}+4c_{\\{1\\}}c_{\\{3\\}}=0\leavevmode\nobreak\
,\end{split}$ (2.50)

and in degree $d=4$

$3c_{\\{1,1\\}}c_{\\{2,1\\}}+2c_{\emptyset}c_{\\{3,2\\}}-3c_{\\{1\\}}c_{\\{2,1,1\\}}+c_{\emptyset}c_{\\{2,1,1,1\\}}+3c_{\\{1\\}}c_{\\{4\\}}-3c_{\emptyset}c_{\\{4,1\\}}-c_{\\{2\\}}c_{\\{1,1,1\\}}-2c_{\\{2\\}}c_{\\{3\\}}=0\leavevmode\nobreak\
.$ (2.51)

In degree 5 and higher there are more than one independent bilinear equation
(one for each partition of size $d$).

Matrix models (at $\beta=1$) are known to provide solutions to Hirota
equations and the reason is that the generating function can be written as a
determinant of a $N\times N$ symmetric matrix (known as the Hankel matrix).
This property follows from Andréief’s integration formula [13] and it is known
to fail when a $\beta$-deformation is turned on. Restricting to the case
$\beta=1$, Andréief’s integration formula allows us to write the determinantal
identity

$\tau_{N}(\mathsf{a};t)=\underset{1\leq i,j\leq
N}{\det}\left[\left(\frac{\partial}{\partial
t_{1}}\right)^{i+j-2}\tau_{1}(\mathsf{a};t)\right]\leavevmode\nobreak\ .$
(2.52)

This way we are able to express the rank $N$ matrix integral as the
determinant of a Hankel matrix constructed out of simpler rank 1 integrals,

$\tau_{1}(\mathsf{a};t)=\int_{\Gamma}\operatorname{d}x\,\operatorname{e}^{-V(x)+\sum_{s=1}^{\infty}t_{s}x^{s}}\leavevmode\nobreak\
.$ (2.53)

Now, assuming the form of the potential $V(x)$ given in (2.2), one can trade
the derivatives with respect to $t_{1}$ with derivatives in the coupling
constant $\mathsf{a}_{1}$ using the additional constraint in (2.6). Evaluating
the generating function at zero times then gives a similar determinant
expression for the empty correlator of the matrix model (i.e. the partition
function)

$\tau_{N}(\mathsf{a};0)=\underset{1\leq i,j\leq
N}{\det}\left[\left(-\frac{\partial}{\partial\mathsf{a}_{1}}\right)^{i+j-2}\tau_{1}(\mathsf{a};0)\right]$
(2.54)

For arbitrary values of the $\beta$-parameter there are no known
generalizations of Hirota bilinear equations. However, for $\beta=1$ the
matrix model generating functions do satisfy Hirota equations by construction
and formula (2.52) gives a combinatorial way to compute all correlators in
arbitrary finite rank. In the following we try to match the Virasoro solutions
that we found in the previous section against explicit matrix model
expressions. Generically speaking, the matching boils down to the choice of
initial data for the recursion. Initial data compatible with matrix model
integrals will automatically lead to Virasoro solutions which are also
solutions to Hirota equations. More general choices of initial data, when
plugged into Hirota equations give rise to additional non-trivial differential
equations that one needs to impose. For the sake of concreteness we explicitly
consider the cases $\mathsf{p}=1,2,3$ and analyze them one by one.

####  2.2.1. $\mathsf{p}=1$

Because the solution in (2.28) allows us to rewrite all correlators as scalar
multiples of the empty one, once we plug them into the bilinear equation
(2.48) there is an overall factor of $c_{\emptyset}^{2}$ which factors out of
the equation. This means that Hirota equations are blind to the choice of
initial data for the $\mathsf{p}=1$ Virasoro solution. Since the overall
normalization of the correlators is no longer important we can rescale
$c_{\emptyset}$ to match with that of a matrix model with linear potential.
Hence the matrix integral solution is essentially the unique one. For
arbitrary $\beta$ this solution does not satisfy Hirota however when $\beta=1$
it does automatically. This can be checked using the rank $N=1$ integral

$\tau_{1}(\mathsf{a}_{1};0)=\int_{0}^{\infty}\operatorname{e}^{-\mathsf{a}_{1}x}\operatorname{d}x=\frac{1}{\mathsf{a}_{1}},\quad\quad\quad\Re(\mathsf{a}_{1})>0$
(2.55)

together with the determinantal formula (2.54) which becomes

$\displaystyle\tau_{N}(\mathsf{a}_{1};0)$ $\displaystyle=\underset{1\leq
i,j\leq N}{\det}\left[\Gamma(i+j-1)\mathsf{a}_{1}^{1-i-j}\right]$ (2.56)
$\displaystyle=\operatorname{BarnesG}(N+1)^{2}\cdot\mathsf{a}_{1}^{-N^{2}}\leavevmode\nobreak\
.$

This reproduces exactly the solution (2.28) at $\beta=1$ once we appropriately
choose the value of the integration constant
$k_{N,1}=\operatorname{BarnesG}(N+1)^{2}$.

####  2.2.2. $\mathsf{p}=2$

A similar discussion holds for the case $\mathsf{p}=2$. The solution is also
essentially unique and can be rescaled to match the matrix model one. At
$\beta\neq 1$ Hirota does not hold but for $\beta=1$ we can use Andréief’s
integration formula to express the matrix model generating function as a
determinant. An explicit computation gives

$\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2};0)=\int_{-\infty}^{\infty}\operatorname{e}^{-\mathsf{a}_{1}x-\mathsf{a}_{2}\frac{x^{2}}{2}}\operatorname{d}x=\sqrt{\frac{2\pi}{\mathsf{a}_{2}}}\operatorname{e}^{\frac{\mathsf{a}_{1}^{2}}{2\mathsf{a}_{2}}},\quad\quad\quad\Re(\mathsf{a}_{2})>0$
(2.57)

so that

$\displaystyle\tau_{N}(\mathsf{a}_{1},\mathsf{a}_{2};0)$
$\displaystyle=\underset{1\leq i,j\leq
N}{\det}\left[\left(-\frac{\partial}{\partial
a_{1}}\right)^{i+j-2}\sqrt{\frac{2\pi}{\mathsf{a}_{2}}}\operatorname{e}^{\frac{\mathsf{a}_{1}^{2}}{2\mathsf{a}_{2}}}\right]$
(2.58)
$\displaystyle=(2\pi)^{\frac{N}{2}}\operatorname{BarnesG}(N+1)\cdot\mathsf{a}_{2}^{-\frac{N^{2}}{2}}\operatorname{e}^{\frac{N\mathsf{a}_{1}^{2}}{2\mathsf{a}_{2}}}$

This reproduces the solution (2.33) at $\beta=1$ for an appropriate choice of
integration constant
$k_{N,1}=(2\pi)^{\frac{N}{2}}\operatorname{BarnesG}(N+1)$.

####  2.2.3. $\mathsf{p}=3$

When the potential is of degree higher than 2 the correlators are not all
proportional to the empty one and the actual Virasoro solution depends on more
initial data. For $\mathsf{p}=3$ all correlators can be rewritten as
polynomials in the totally anti-symmetric correlators contained in the
generating function (2.39). Using the additional constraint in (2.6) we can
rewrite each one of them using $\mathsf{a}$-derivatives acting on the empty
correlator as

$c_{\\{\underbrace{\scriptstyle
1,\dots,1}_{k}\\}}(\mathsf{a})=\left(-\frac{\partial}{\partial\mathsf{a}_{1}}\right)^{k}c_{\emptyset}(\mathsf{a})\leavevmode\nobreak\
.$ (2.59)

If we then plug this Virasoro solution inside of Hirota equations we get non-
linear differential equations for $c_{\emptyset}(\mathsf{a})$. In degree $d=3$
for instance, equation (2.50) becomes

$\displaystyle\mathsf{a}_{3}^{2}\left(c_{\emptyset}\partial_{\mathsf{a}_{1}}^{4}c_{\emptyset}-4\left(\partial_{\mathsf{a}_{1}}^{3}c_{\emptyset}\right)\left(\partial_{\mathsf{a}_{1}}c_{\emptyset}\right)+3\left(\partial_{\mathsf{a}_{1}}^{2}c_{\emptyset}\right)^{2}\right)+(\mathsf{a}_{2}^{2}-4\mathsf{a}_{1}\mathsf{a}_{3})$
$\displaystyle\left(\left(\partial_{\mathsf{a}_{1}}c_{\emptyset}\right)^{2}-c_{\emptyset}\partial_{\mathsf{a}_{1}}^{2}c_{\emptyset}\right)+$
(2.60)
$\displaystyle-2\mathsf{a}_{3}c_{\emptyset}\partial_{\mathsf{a}_{1}}c_{\emptyset}+\mathsf{a}_{2}Nc_{\emptyset}^{2}=0$

which is clearly a non-trivial constraint on the function
$c_{\emptyset}(\mathsf{a})$ and it means that Hirota equations are no longer
automatically satisfied not even in the $\beta=1$ case. Using the explicit
expression of $c_{\emptyset}$ that we found in (2.37) we then get an equation
for the unknown function $g$, namely

$2\left(zg^{\prime}(z)^{2}-8g^{\prime\prime\prime}(z)g^{\prime}(z)+6g^{\prime\prime}(z)^{2}\right)+g(z)\left(g^{\prime}(z)-2zg^{\prime\prime}(z)+4g^{\prime\prime\prime\prime}(z)\right)=0$
(2.61)

Again, this is a non-trivial differential constraint on the initial data
$g(z)$ and therefore we have that some of our Virasoro solutions do not
satisfy Hirota333Here we only considered the $d=3$ Hirota equation, however
one should check for all equations in all degrees $d>3$. even for $\beta=1$.

Clearly, there is one obvious set of solutions which is provided by matrix
models at $\beta=1$, however in this case it is not so simple to write down a
closed formula for the function $g(z)$ by matching with the integral
expression of the generating function (especially at higher rank $N>1$). As an
example we consider the simple case of $N=1$. The cubic potential $V(x)$ can
always be rewritten as a depressed cubic via a change of variable
$x=(\mathsf{a}_{3}^{-1/3}x^{\prime}-\mathsf{a}_{2}/2\mathsf{a}_{3})$, so that
we can write the rank-1 matrix integral as

$\displaystyle\tau_{1}(\mathsf{a}_{1},\mathsf{a}_{2},\mathsf{a}_{3};0)$
$\displaystyle=\int_{\Gamma}\exp\left(-\frac{1}{3}\mathsf{a}_{3}x^{3}-\frac{1}{2}\mathsf{a}_{2}x^{2}-\mathsf{a}_{1}x\right)\operatorname{d}x$
(2.62)
$\displaystyle=\mathsf{a}_{3}^{-1/3}\exp\left(-\frac{\mathsf{a}_{2}(\mathsf{a}_{2}^{2}-6\mathsf{a}_{1}\mathsf{a}_{3})}{12\mathsf{a}_{3}^{2}}\right)\int_{\Gamma^{\prime}}\exp\left(-\frac{x^{\prime
3}}{3}+\left(\frac{\mathsf{a}_{2}^{2}-4\mathsf{a}_{1}\mathsf{a}_{3}}{4\mathsf{a}_{3}^{4/3}}\right)x^{\prime}\right)\operatorname{d}x^{\prime}\leavevmode\nobreak\
.$

By comparison with (2.37) we have

$\left.g(z)\right|_{N=1}=\int_{\Gamma^{\prime}}\exp\left(-\frac{x^{\prime
3}}{3}+\frac{z}{2}x^{\prime}\right)\operatorname{d}x^{\prime}\leavevmode\nobreak\
,$ (2.63)

which, for an appropriate choice of contour $\Gamma^{\prime}$, is a linear
combination of Airy functions. It is then straightforward to show that (2.61)
and all higher degree Hirota equations follow from the Airy differential
equation (2.41). We should mention however, that it is not obvious to us that
the converse is true as well, i.e. that Hirota equations of the form (2.61)
imply the Airy ODE on $g(z)$, which would imply that solution to Hirota plus
Virasoro at $\mathsf{p}=3$ are essentially uniquely given by a 1-matrix
integrals (at least for rank 1).

##  3\. Solving non-homogeneous Virasoro constraints

Now we wish to generalize the above matrix model computations and allow for
the domain of integration of the eigenvalues to have a boundary. As a concrete
example we focus on the case in which each eigenvalue $x_{i}$ belongs to the
finite interval $[a,b]$ on the real line444Note here that $a$ is the limit of
integration which is not to be confused with the couplings $\mathsf{a}_{k}$
appearing in the potential $V(x)$ in (2.2).. The generating function of the
matrix model is then given by

$\tau_{N}(\mathsf{a};t)=\frac{1}{N!}\int_{[a,b]^{N}}\prod_{i=1}^{N}\operatorname{d}x_{i}\prod_{1\leq
i<j\leq
N}|x_{i}-x_{j}|^{2\beta}\operatorname{e}^{-\sum_{i=1}^{N}V(x_{i})+\sum_{s=1}^{\infty}t_{s}\sum_{i=1}^{N}x_{i}^{s}}$
(3.1)

where the compact domain $[a,b]^{N}$ is essentially an hypercube inside of
$\mathbb{R}^{N}$. Furthermore, it is important that all integrals have the
same domain of integration $[a,b]$ because we want to preserve the symmetry of
the integral under permutation of the eigenvalues (this is the action of the
Weyl group of $U(N)$ on the Cartan subalgebra). The choice of domain becomes
part of the parameters of the generating function $\tau_{N}(\mathsf{a};t)$,
however we do not explicitly write this dependence in order not to clutter too
much the notation.

The matrix model in (3.1) no longer satisfies standard Virasoro constraints
because none of the infinitesimal deformations generated by the vector fields

$\xi_{n}=\sum_{i=1}^{N}x_{i}^{n+1}\frac{\partial}{\partial
x_{i}},\quad\quad\quad n\geq-1$ (3.2)

preserve the domain of integration $[a,b]^{N}$. This means that the integral
is no longer invariant under such deformations. However, we can compute
exactly what is the variation of the integral and write down a modified
version of the Virasoro constraints. As it turns out, these equations are
modified by a boundary term which can be rewritten as a matrix integral of
rank $N-1$. This boundary term can be directly associated to a matrix model
living on the $(N-1)$-dimensional faces of the hypercube and it turns the
Virasoro constraints into non-homogeneous equations. More concretely, we can
write the equations

$\left(\sum_{k=1}^{\mathsf{p}}\mathsf{a}_{k}\frac{\partial}{\partial
t_{k+n}}+\mathsf{a}_{1}N\delta_{n,-1}-L_{n}\right)\tau_{N}(\mathsf{a};t)=B_{n}(t),\qquad
n\geq-1$ (3.3)

where the l.h.s. is the same as that of the standard Virasoro constraints,
while the r.h.s. of the equation is the boundary term mentioned above. As
before, the generating function satisfies the additional constraints

$\left(\frac{\partial}{\partial
t_{k}}+k\frac{\partial}{\partial\mathsf{a}_{k}}\right)\tau_{N}(\mathsf{a};t)=0\,,\qquad
k=1,\dots\mathsf{p}$ (3.4)

due to the presence of a non-trivial polynomial potential.

We now compute each such boundary term $B_{n}(t)$ explicitly. Let us call
$\omega(t)\in\Omega^{N}(\mathbb{R}^{N})\otimes\mathbb{C}[[t_{1},t_{2},\dots]]$
the top degree differential form in the integrand of (3.1), so that
$\tau_{N}(t)=\int\omega(t)$. Virasoro constraints are obtained by integrating
the variation of $\omega(t)$ with respect to the infinitesimal diffeomorphism
generated by the vector field $\xi_{n}$. We can therefore write this variation
as a Lie derivative along the vector which becomes

$B_{n}(t)=-\int_{\mathcal{D}}\mathcal{L}_{\xi_{n}}\omega(t)=-\int_{\mathcal{D}}\operatorname{d}\iota_{\xi_{n}}\omega(t)=-\int_{\partial\mathcal{D}}\iota_{\xi_{n}}\omega(t)$
(3.5)

where $\mathcal{D}\subset\mathbb{R}^{N}$ is a generic domain of integration
and we used that $\omega(t)$ is closed555The form $\operatorname{d}\omega$
might have some legs in the $t$-directions however these are projected out
when contracted with $\xi_{n}$. because of dimensional reasons. Substituting
the explicit expressions for $\omega(t)$, $\xi_{n}$ and $\mathcal{D}$ we have

$\iota_{\xi_{n}}\omega(t)=\frac{1}{N!}\sum_{i=1}^{N}(-1)^{i+1}x_{i}^{n+1}\prod_{j\neq
i}\operatorname{d}x_{i}\prod_{1\leq j<k\leq
N}|x_{j}-x_{k}|^{2\beta}\operatorname{e}^{-\sum_{j=1}^{N}V(x_{j})+\sum_{s=1}^{\infty}t_{s}\sum_{j=1}^{N}x_{j}^{s}}$
(3.6)

and integrating over $\partial[a,b]^{N}$ we get

$B_{n}(t)=-\left\llangle\prod_{j=1}^{N-1}|z-x_{j}|^{2\beta}\right\rrangle_{N-1}z^{n+1}\left.\exp\left(-V(z)+\sum_{s=1}^{\infty}t_{s}z^{s}\right)\right|_{z=a}^{z=b}$
(3.7)

where we introduced the notation

$\left.f(z)\right|_{z=a}^{z=b}=f(b)-f(a)$ (3.8)

and

$\left\llangle
f(x)\right\rrangle_{N}:=\frac{1}{N!}\int_{[a,b]^{N}}\prod_{i=1}^{N}\operatorname{d}x_{i}\,f(x)\prod_{1\leq
i<j\leq
N}|x_{i}-x_{j}|^{2\beta}\exp\left({-\sum_{i=1}^{N}V(x_{i})+\sum_{s=1}^{\infty}t_{s}\sum_{i=1}^{N}x_{i}^{s}}\right)$
(3.9)

so that $\llangle\,\rrangle_{N}$ is the time-dependent expectation value in a
matrix model of rank $N$.

Before discussing the solution of these non-homogeneous Virasoro constraints
we remark here that the choice of contour $[a,b]^{N}$ has a direct
interpretation in terms of matrices and not just eigenvalues. Namely, one can
define the function

$\sigma:H_{N}\to\mathbb{R}^{N}/S_{N}$ (3.10)

which assigns to each Hermitean matrix in $H_{N}$ its spectrum which can be
represented as a $N$-tuple of real numbers up to permutations (hence the
quotient by the symmetric group $S_{N}$). This map is surjective and moreover
it is clearly invariant under the adjoint action of $U(N)$ on $H_{N}$. This
means that the preimage of any subset of the codomain is a disjoint union of
adjoint orbits of $U(N)$. The last step consists of identifying the
appropriate subset of the codomain. We choose the compact subspace
$[a,b]^{N}/S_{N}\hookrightarrow\mathbb{R}^{N}/S_{N}$ so that this choice
matches the domain in (3.1) (up to a numerical factor equal to the order of
the Weyl group). The preimage under $\sigma$ of this region defines a subset
of $H_{N}$ invariant under the $U(N)$ action. This is the domain of
integration of the actual matrix model.

###  3.1. Gaussian constraints in the hypercube

For the sake of concreteness, we now consider the particular case of a
Gaussian potential with $\mathsf{p}=2$ and coupling constants
$\mathsf{a}_{1}=0$ and $\mathsf{a}_{2}=1$, in other words
$V_{\mathrm{G}}(x)=\frac{x^{2}}{2}$. Since the coupling constants are fixed,
we avoid writing them as variables of the generating function and simply
denote it as $\tau_{N}(t)$. Other choices of polynomial potentials are not
conceptually more difficult to treat except for a slightly more involved
discussion of the initial data. Here we chose a Gaussian potential so that we
can demonstrate the main ideas for solving non-homogeneous constraints without
having to discuss the possible non-uniqueness of the solution.

The non-homogeneous Gaussian Virasoro constraints are given by

$\left(\frac{\partial}{\partial
t_{n+2}}-L_{n}\right)\tau_{N}(t)=B_{n}(t)\leavevmode\nobreak\ ,\quad\quad\quad
n\geq-1\leavevmode\nobreak\ ,$ (3.11)

with the boundary term $B_{n}(t)$ as in (3.7). Re-summing the constraints with
weight $(n+2)t_{n+2}$ from $n=-1$ to $n=\infty$ we find the equation

$\left(D-W_{-2}\right)\tau_{N}(t)=B(t)$ (3.12)

where $D$ and $W_{-2}$ are given in (2.12) and (2.18) respectively, while

$B(t):=\sum_{n=-1}^{\infty}(n+2)t_{n+2}B_{n}(t)\leavevmode\nobreak\ .$ (3.13)

This can be solved inductively on the degree as

$\tau_{N}^{(d)}(t)=\frac{1}{d}\left(W_{-2}\tau_{N}^{(d-2)}(t)+B^{(d)}(t)\right)$
(3.14)

where $d$ is the degree with respect to the dilatation operator $D$ and we
assume the expansion

$B(t)=\sum_{d=1}^{\infty}B^{(d)}(t)\leavevmode\nobreak\
,\quad\quad\quad[D,B^{(d)}(t)]=d\cdot B^{(d)}(t)\leavevmode\nobreak\ .$ (3.15)

Observe that the sum over $d$ starts from one and not from zero because from
(3.7) it follows that $B^{(0)}(t)=B(0)=0$. Then we have

$\begin{split}\tau_{N}(t)=&\tau_{N}(0)+(D-W_{-2})^{-1}(W_{-2}\tau_{N}(0)+B(t))\\\
=&\sum_{s=0}^{\infty}(D^{-1}W_{-2})^{s}\tau_{N}(0)+\sum_{s=0}^{\infty}(D^{-1}W_{-2})^{s}D^{-1}B(t)\\\
=&\exp\left(W_{-2}/2\right)\tau_{N}(0)+\sum_{s=0}^{\infty}(D^{-1}W_{-2})^{s}D^{-1}B(t)\\\
=&\exp(W_{-2}/2)\tau_{N}(0)+\sum_{s=0}^{\infty}\sum_{d=1}^{\infty}\frac{W_{-2}^{s}B^{(d)}(t)}{2^{s}d\left(\frac{d}{2}+1\right)_{s}}\leavevmode\nobreak\
,\end{split}$ (3.16)

with $(x)_{n}$ being the Pochhammer symbol defined in (A.2). The result above
generalizes the solution to the homogeneous Virasoro constraints given in
(2.19) for $\mathsf{a}_{k}=\delta_{k,2}$. Observe that the non-homogeneous
solution naturally factorizes into the sum of a homogeneous solution and a
boundary term depending on $B(t)$. See Appendix C for a cohomological
interpretation of this fact.

In the following we provide explicit solutions to the non-homogeneous
constraints in terms of correlators for the cases of low rank $N=1,2$ as
examples. By doing so we also elaborate on the fact that the rank $N$ solution
can be quite explicitly constructed from knowledge of the one in rank $N-1$,
so that there are two types of recursions at play. The first is the usual
recursion on correlators for fixed rank and the second is the recursion in the
rank itself.

####  3.1.1. Examples

Starting with the case of $N=1$, the boundary term in (3.7) becomes

$B_{n}(t)=a^{n+1}\exp\left(-\frac{a^{2}}{2}+\sum_{s=1}^{\infty}t_{s}a^{s}\right)-b^{n+1}\exp\left(-\frac{b^{2}}{2}+\sum_{s=1}^{\infty}t_{s}b^{s}\right)$
(3.17)

so that the degree $d$ term in the expansion (3.15) can be written as

$B^{(d)}(t)=\sum_{n=-1}^{d-2}(n+2)t_{n+2}\sum_{\lambda\vdash(d-n-2)}\frac{1}{|\textrm{Aut}(\lambda)|}\left(a^{d-1}\operatorname{e}^{-\frac{a^{2}}{2}}-b^{d-1}\operatorname{e}^{-\frac{b^{2}}{2}}\right)\prod_{\ell\in\lambda}t_{\ell}\leavevmode\nobreak\
.$ (3.18)

Plugging in (3.16) we can explicitly compute the first few
correlators666Notice that in rank $N=1$ the correlators only depend on the
size of the partitions and not on the actual shape.

$\begin{split}c^{N=1}_{\emptyset}&=\tau_{1}(0)\\\
c^{N=1}_{\\{1\\}}&=\operatorname{e}^{-\frac{a^{2}}{2}}-\operatorname{e}^{-\frac{b^{2}}{2}}\\\
c^{N=1}_{\\{2\\}}&=\tau_{1}(0)+a\operatorname{e}^{-\frac{a^{2}}{2}}-b\operatorname{e}^{-\frac{b^{2}}{2}}\\\
c^{N=1}_{\\{3\\}}&=\left(a^{2}+2\right)\operatorname{e}^{-\frac{a^{2}}{2}}-\left(b^{2}+2\right)\operatorname{e}^{-\frac{b^{2}}{2}}\\\
c^{N=1}_{\\{4\\}}&=3\tau_{1}(0)+a\left(a^{2}+3\right)\operatorname{e}^{-\frac{a^{2}}{2}}-b\left(b^{2}+3\right)\operatorname{e}^{-\frac{b^{2}}{2}}\\\
c^{N=1}_{\\{5\\}}&=\left(a^{4}+4a^{2}+8\right)\operatorname{e}^{-\frac{a^{2}}{2}}-\left(b^{4}+4b^{2}+8\right)\operatorname{e}^{-\frac{b^{2}}{2}}\leavevmode\nobreak\
.\end{split}$ (3.19)

The empty correlator in this case can be computed explicitly by evaluating the
integral and we get

$\tau_{1}(0)=\int_{a}^{b}\exp(-x^{2}/2)\operatorname{d}x=\sqrt{\frac{\pi}{2}}\left(\operatorname{erf}\left(\frac{b}{\sqrt{2}}\right)-\operatorname{erf}\left(\frac{a}{\sqrt{2}}\right)\right)\leavevmode\nobreak\
,$ (3.20)

then we can write a closed formula for all correlators in rank 1 as

$c^{N=1}_{\lambda}=F_{|\lambda|}$ (3.21)

with

$\displaystyle F_{s}:=$
$\displaystyle\frac{(1+(-1)^{s})}{2}(s-1)!!\,\tau_{1}(0)+$ (3.22)
$\displaystyle+\left(\sum_{k=0}^{\left\lfloor\frac{s-1}{2}\right\rfloor}\frac{(s-1)!!}{(s-1-2k)!!}a^{s-1-2k}\right)\operatorname{e}^{-\frac{a^{2}}{2}}-\left(\sum_{k=0}^{\left\lfloor\frac{s-1}{2}\right\rfloor}\frac{(s-1)!!}{(s-1-2k)!!}b^{s-1-2k}\right)\operatorname{e}^{-\frac{b^{2}}{2}}\leavevmode\nobreak\
.$

Moving on to the case of rank $N=2$ the boundary term takes the form

$B_{n}(t)=\left\llangle|a-x_{1}|^{2\beta}\right\rrangle_{1}a^{n+1}\exp\left(-\frac{a^{2}}{2}+\sum_{s=1}^{\infty}t_{s}a^{s}\right)-\left\llangle|b-x_{1}|^{2\beta}\right\rrangle_{1}b^{n+1}\exp\left(-\frac{b^{2}}{2}+\sum_{s=1}^{\infty}t_{s}b^{s}\right)\leavevmode\nobreak\
.$ (3.23)

For arbitrary values of $\beta$ the expectation values that appear are not
polynomial in the integration variables of the matrix model therefore this
boundary term is a complicated function for which we do not have a closed
formula. For integer $\beta$ however they become expectation values of
polynomial functions and thus they can be computed combinatorially. For
instance, let us consider the simple case $\beta=1$. Then we can write the
boundary term of degree $d$ as

$\displaystyle B^{(d)}$
$\displaystyle(t)=\sum_{n=-1}^{d-2}(n+2)t_{n+2}\sum_{\begin{subarray}{c}\lambda,\rho\\\
|\lambda|+|\rho|=d-(n+2)\end{subarray}}\frac{1}{|\operatorname{Aut}(\lambda)|}\frac{1}{|\operatorname{Aut}(\rho)|}\prod_{\ell\in\lambda}t_{\ell}\prod_{r\in\rho}t_{r}\times$
(3.24)
$\displaystyle\times\left[\left(F_{2+|\lambda|}-2aF_{1+|\lambda|}+a^{2}F_{|\lambda|}\right)a^{n+1+|\rho|}\operatorname{e}^{-\frac{a^{2}}{2}}-\left(F_{2+|\lambda|}-2bF_{1+|\lambda|}+b^{2}F_{|\lambda|}\right)b^{n+1+|\rho|}\operatorname{e}^{-\frac{b^{2}}{2}}\right]\leavevmode\nobreak\
.$

We can now use (3.16) to explicitly find the first few correlators as

$\displaystyle c^{N=2}_{\emptyset}$ $\displaystyle=\tau_{2}(0)$ (3.25)
$\displaystyle c^{N=2}_{\\{1\\}}$
$\displaystyle=\operatorname{e}^{-\frac{a^{2}}{2}}\left(a^{2}F_{0}-2aF_{1}+F_{2}\right)-\operatorname{e}^{-\frac{b^{2}}{2}}\left(b^{2}F_{0}-2bF_{1}+F_{2}\right)$
$\displaystyle c^{N=2}_{\\{1,1\\}}$
$\displaystyle=2\tau_{2}(0)+\operatorname{e}^{-\frac{a^{2}}{2}}\left(a^{3}F_{0}-a^{2}F_{1}-aF_{2}+F_{3}\right)+\operatorname{e}^{-\frac{b^{2}}{2}}\left(b^{3}F_{0}-b^{2}F_{1}-bF_{2}+F_{3}\right)$
$\displaystyle c^{N=2}_{\\{2\\}}$
$\displaystyle=4\tau_{2}(0)+a\operatorname{e}^{-\frac{a^{2}}{2}}\left(F_{0}a^{2}-2F_{1}a+F_{2}\right)-b\operatorname{e}^{-\frac{b^{2}}{2}}\left(F_{0}b^{2}-2F_{1}b+F_{2}\right)$
$\displaystyle\dots&$

Similarly, one can use the data from rank 1 and 2 to construct the rank 3
solution and higher.

###  3.2. Gaussian constraints in the orthant

We now consider the special case of a domain $\mathcal{D}=[0,\infty)^{N}$,
i.e. the orthant obtained as the limit in which $a\to 0$ and $b\to\infty$. The
generating function on the orthant is given by

$\tau_{N}(t)=\frac{1}{N!}\int_{[0,\infty)^{N}}\prod_{i=1}^{N}\operatorname{d}x_{i}\prod_{1\leq
k<l\leq
N}|x_{k}-x_{l}|^{2\beta}\operatorname{e}^{-\sum_{i=1}^{N}V_{\mathrm{G}}(x_{i})+\sum_{s=1}^{\infty}t_{s}\sum_{i=1}^{N}x_{i}^{s}}$
(3.26)

where we again have a Gaussian potential $V_{\mathrm{G}}(x)=\frac{x^{2}}{2}$.
Since one of the boundaries has been pushed to infinity we are left with only
one boundary contribution per integration variable, namely the contribution
coming from the boundary at zero. By explicitly computing formula (3.7) we
immediately see that all such terms vanish if $n\geq 0$, while for $n=-1$ we
have the only non-trivial boundary term

$B_{n}(t)=\delta_{n,-1}\left\llangle\prod_{i=1}^{N-1}x_{i}^{2\beta}\right\rrangle_{N-1}\leavevmode\nobreak\
.$ (3.27)

Therefore the Virasoro constraints are almost homogeneous with the exception
of the $n=-1$ constraint, also known as the string equation [14]. This
suggests that we can use the following trick. We compute the solution to the
homogeneous $n\geq 0$ constraints independently of the string equation and
then solve the string equation in a second step. In order to perform the first
step, we resum all homogeneous constraints with weight $(n+2)t_{n+2}$ to get

$\hat{D}\tau_{N}(t)=\hat{W}_{-2}\tau_{N}$ (3.28)

with

$\hat{D}=\sum_{s=2}^{\infty}st_{s}\partial_{s}\quad\quad\quad\text{and}\quad\quad\quad\hat{W}_{-2}=\sum_{n=0}^{\infty}(n+2)t_{n+2}L_{n}\leavevmode\nobreak\
.$ (3.29)

Observe that $\hat{D}$ is the same partial dilatation operator that we
encountered in the cubic case $\mathsf{p}=3$ in (2.21). This means that the
solution of (3.28) is not unique and it must depend on a choice of initial
data given by generating function of anti-symmetric correlators, i.e.
$\tau_{N}(t_{1})$ as defined in (2.23). The solution of the recursion then can
be written as

$\tau_{N}(t)=\sum_{s=0}^{\infty}(\hat{D}^{-1}\hat{W}_{-2})^{s}\,\tau_{N}(t_{1})$
(3.30)

This fixes all correlators except for those of the form $c_{\\{1,1,\dots\\}}$
which can only be fixed in a second step by analyzing the string equation.
Equivalently, using (3.30) we can fix completely the dependence of the
generating function upon all times $t_{s}$ with $s>1$. The dependence on the
time $t_{1}$ is entirely contained in the function $\tau_{N}(t_{1})$ and can
only be fixed via the string equation. Setting all the higher times to zero we
can simplify the string equation as

$\left(\partial_{1}-Nt_{1}\right)\tau_{N}(t_{1})=B_{-1}(t_{1})\leavevmode\nobreak\
.$ (3.31)

We define the operator on the l.h.s. as $M_{N}:=\partial_{1}-Nt_{1}$. Since
$M_{N}$ is invertible on the complement of the homogeneous solution (i.e. the
kernel of $M_{N}$), the solution becomes

$\tau_{N}(t_{1})=\tau_{N}^{\mathrm{h}}(t_{1})+M_{N}^{-1}B_{-1}(t_{1})$ (3.32)

with $\tau_{N}^{\mathrm{h}}=\operatorname{e}^{Nt_{1}^{2}/2}\tau_{N}(0)$ being
a solution to the homogeneous constraints. Furthermore, we can study the
boundary term $B_{-1}(t_{1})$ and observe that for integer $\beta$ the
boundary term is the expectation value of a polynomial function and therefore
it can be rewritten as the action of a differential operator on the generating
function in rank one less. We call this operator $S_{N-1}(t_{1},\partial_{1})$
and define it by the equation

$B_{-1}(t_{1})=:S_{N-1}\,\tau_{N-1}(t_{1})\leavevmode\nobreak\ .$ (3.33)

Finally, we can write the solution to the non-homogeneous string equation as

$\tau_{N}(t_{1})=\operatorname{e}^{\frac{Nt_{1}^{2}}{2}}\tau_{N}(0)+M_{N}^{-1}S_{N-1}\tau_{N-1}(t_{1})\leavevmode\nobreak\
,$ (3.34)

whose computation can now be done inductively w.r.t. the rank $N$. We remark
here that the operator $S_{N}$ is well-defined only for integer $\beta$ when
the argument of the expectation value in (3.27) is polynomial. For arbitrary
$\beta$ one should compute the non-homogeneous contribution as in (3.16).

For $\beta=1$ one can compute the operator $S_{N}$ for the first few values of
$N$ by making use of the Virasoro constraints. For instance,

$S_{0}(t_{1},\partial_{1})=1$ (3.35)
$S_{1}(t_{1},\partial_{1})=\partial_{1}^{2}$ (3.36)
$S_{2}(t_{1},\partial_{1})=\frac{1}{4}\left(\partial_{1}^{4}-12\partial_{1}^{2}-2t_{1}\partial_{1}^{3}+11t_{1}\partial_{1}+t_{1}^{2}\partial_{1}^{2}+24\right)$
(3.37)
$\begin{split}S_{3}(t_{1},\partial_{1})=&\frac{1}{36}(104t_{1}^{3}\partial_{1}+4t_{1}^{4}\partial_{1}^{2}-264t_{1}^{2}\partial_{1}^{2}-12t_{1}^{3}\partial_{1}^{3}+13t_{1}^{2}\partial_{1}^{4}+576t_{1}^{2}\\\
&-1260t_{1}\partial_{1}+639\partial_{1}^{2}+213t_{1}\partial_{1}^{3}-6t_{1}\partial_{1}^{5}-54\partial_{1}^{4}+\partial_{1}^{6}-540)\end{split}$
(3.38)

The claim is thus that using the above procedure the model can be solved
completely in any rank, up to initial data given by the collection of all
empty correlators $\tau_{N}({0})$.

####  3.2.1. Examples

Let us consider first the case of rank $N=1$. The boundary terms in the
Virasoro constraints on the orthant are

$B_{n}(t)=\delta_{n,-1}\quad\quad\quad B(t)=t_{1}$ (3.39)

Plugging in (3.16) we get

$\tau_{1}(t)=\exp\left(W_{-2}/2\right)\tau_{1}(0)+\sum_{s=0}^{\infty}\frac{W_{-2}^{s}}{2^{s}(3/2)_{s}}t_{1}\leavevmode\nobreak\
.$ (3.40)

The string equation is

$(\partial_{1}-t_{1})\tau_{1}(t_{1})=1$ (3.41)

whose solution is

$\tau_{1}(t_{1})=\operatorname{e}^{\frac{t_{1}^{2}}{2}}\tau_{1}(0)+\sqrt{\frac{\pi}{2}}\operatorname{e}^{\frac{t_{1}^{2}}{2}}\operatorname{erf}\left(\frac{t_{1}}{\sqrt{2}}\right)$
(3.42)

which is exactly what we get from (3.40) upon setting all higher times to zero
(so that $W_{-2}=t_{1}^{2}$). Here $\tau(0)$ is an integration constant which
we can fix by computing the following integral explicitly

$\tau_{1}(0)=\int_{0}^{\infty}\exp(-x^{2}/2)\operatorname{d}x=\sqrt{\frac{\pi}{2}}\leavevmode\nobreak\
,$ (3.43)

which we can use to give a closed formula for each correlator as

$c_{\lambda}^{N=1}=2^{\frac{|\lambda|-1}{2}}\Gamma\left(\frac{|\lambda|+1}{2}\right)\leavevmode\nobreak\
.$ (3.44)

This is consistent with (3.21) upon taking the limit $a\rightarrow 0$ and
$b\rightarrow\infty$.

For rank $N=2$ we have

$B_{-1}(t)=\left\llangle
x^{2\beta}\right\rrangle_{1}=\int_{0}^{\infty}x^{2\beta}\exp\left(-\frac{x^{2}}{2}+\sum_{s=1}^{\infty}t_{s}x^{s}\right)\operatorname{d}x$
(3.45)

and

$B^{(d)}(t)=\sum_{|\lambda|=d-1}\frac{1}{|\mathrm{Aut}(\lambda)|}2^{\frac{|\lambda|-1+2\beta}{2}}\Gamma\left(\frac{|\lambda|+1+2\beta}{2}\right)t_{1}\prod_{a\in\lambda}t_{a}\leavevmode\nobreak\
.$ (3.46)

In this case we cannot get a closed formula for the correlators, however we
can still get explicit expressions for individual correlators by using (3.16).
The first few correlators are

$\displaystyle c^{N=2}_{\emptyset}$ $\displaystyle=\tau_{2}(0)$ (3.47)
$\displaystyle c^{N=2}_{\\{1\\}}$
$\displaystyle=2^{\beta-\frac{1}{2}}\Gamma\left(\beta+\frac{1}{2}\right)$
$\displaystyle c^{N=2}_{\\{1,1\\}}$
$\displaystyle=2\tau_{2}(0)+2^{\beta}\Gamma(\beta+1)$ $\displaystyle
c^{N=2}_{\\{2\\}}$ $\displaystyle=2(\beta+1)\tau_{2}(0)$ $\displaystyle
c^{N=2}_{\\{1,1,1\\}}$
$\displaystyle=2^{\beta-\frac{1}{2}}(2\beta+5)\Gamma\left(\beta+\frac{1}{2}\right)$
$\displaystyle c^{N=2}_{\\{2,1\\}}$
$\displaystyle=\frac{2^{\beta+\frac{3}{2}}\Gamma\left(\beta+\frac{5}{2}\right)}{2\beta+1}$
$\displaystyle c^{N=2}_{\\{3\\}}$
$\displaystyle=2^{\beta+\frac{1}{2}}(\beta+1)\Gamma\left(\beta+\frac{1}{2}\right)\leavevmode\nobreak\
.$

In general, one can compute the rank $N$ solution using knowledge of the rank
$N-1$ correlators.

##  4\. Conclusion

In the present article we studied various versions of Virasoro constraints
obtained as Ward identities for Hermitean 1-matrix models with polynomial
potential. First we analyzed the moduli space of solutions for the standard
Virasoro constraints and compared the solution on a generic point in this
space to that given by the Hermitean matrix model itself. In this regard we
showed that there must exist points in the moduli space which do not satisfy
Hirota bilinear relations and we investigated the minimal set of assumptions
to define the subspace of Hirota solutions.

An interesting question that arises then is whether there are solutions to
Hirota that do not come from solving Virasoro constraints. One possible set of
candidates for such $\tau$-functions are those obtained from $q,t$-models at
$t=q$ (or $\beta=1$) for which the matrix integral is substituted by a Jackson
integral. These model satisfy $q$-Virasoro constraints instead of usual ones,
however they do still admit a determinant representation via Andréief’s
identity and therefore they satisfy the Hirota equation. These and other
discrete generalizations of the Hermitean matrix models certainly do merit
further examination.

Next, we studied non-homogeneous Virasoro constraint obtained from certain
matrix models whose integration domains has non-empty boundary. We showed that
a solution can be algebraically constructed in a similar way as in the
homogeneous case and we further argued that such solutions can be obtained by
induction on the rank of the model.

It would be interesting to study more general cases of non-homogeneous
Virasoro constraints. In particular, one could consider domains of integration
defined by inequalities of the type

$p_{s}=\sum_{i=1}^{N}x_{i}^{s}\leq C$ (4.1)

for some fixed constant $C\in\mathbb{R}$ (which generalizes the role of our
$a$ and $b$ parameters). This leads to a well defined matrix model because the
power-sums are Casimirs of the Lie algebra of $U(N)$ hence they define
conjugation invariant functions on the space of Hermitean matrices. Moreover,
each Virasoro vector field $\xi_{n}$ is proportional to the gradient of the
power-sum function $p_{n+2}$ which means that every subspace defined by (4.1)
is maximally deformed by one of the Virasoro generators. Correspondingly, at
least one constraints will be non-homogenous and the boundary term should have
some geometric interpretation in terms of the theory of coadjoint orbits for
$U(N)$. We leave the study of these models for future investigations.

##  Appendix A Special functions

In this section we will summarize the special functions used throughout. The
first such function is the Gamma function $\Gamma(z)$ defined as the analytic
continuation of the integral

$\Gamma(z)=\int_{0}^{\infty}t^{z-1}\operatorname{e}^{-t}\operatorname{d}t\,,\qquad\Re(z)>0$
(A.1)

from which we can define the Pochhammer symbol $(z)_{n}$,

$(z)_{n}=\frac{\Gamma(z+n)}{\Gamma(z)}\leavevmode\nobreak\ .$ (A.2)

We also need the Airy functions which are defined as being the linearly
independent solutions to the differential equation

$f^{\prime\prime}(z)-zf(z)=0$ (A.3)

i.e.

$f(z)=k_{A}\,\operatorname{Ai}(z)+k_{B}\,\operatorname{Bi}(z)\leavevmode\nobreak\
.$ (A.4)

with $k_{A},k_{B}$ arbitrary complex coefficients. The Airy functions can be
given in terms of the contour integrals

$\operatorname{Ai}(z)=\frac{1}{2\pi\operatorname{i}}\int_{\infty\operatorname{e}^{-\frac{\pi\operatorname{i}}{3}}}^{\infty\operatorname{e}^{\frac{\pi\operatorname{i}}{3}}}\exp\left(\frac{t^{3}}{3}-tz\right)\,\operatorname{d}t$
(A.5)

and

$\operatorname{Bi}(z)=\frac{1}{2\pi}\int_{-\infty}^{\infty\operatorname{e}^{-\frac{\pi\operatorname{i}}{3}}}\exp\left(\frac{t^{3}}{3}-tz\right)\,\operatorname{d}t\,+\,\frac{1}{2\pi}\int_{-\infty}^{\infty\operatorname{e}^{\frac{\pi\operatorname{i}}{3}}}\exp\left(\frac{t^{3}}{3}-tz\right)\,\operatorname{d}t\leavevmode\nobreak\
.$ (A.6)

Furthermore, we make use of the error function $\operatorname{erf}(z)$ which
we define by the series

$\operatorname{erf}(z)=\frac{2}{\sqrt{\pi}}\sum_{n=0}^{\infty}\frac{(-1)^{n}z^{2n+1}}{n!(2n+1)}\leavevmode\nobreak\
,$ (A.7)

and the Barnes G-function

$\operatorname{BarnesG}(z+1)=(2\pi)^{z/2}\exp\left(-{\frac{z+z^{2}(1+\gamma)}{2}}\right)\,\prod_{k=1}^{\infty}\left\\{\left(1+{\frac{z}{k}}\right)^{k}\exp\left({\frac{z^{2}}{2k}}-z\right)\right\\}\leavevmode\nobreak\
,$ (A.8)

for $\gamma$ the Euler-Mascheroni constant. For positive integers
$n\in\mathbb{Z}$ the Barnes G-function evaluates to

$\operatorname{BarnesG}(n+1)=\prod_{k=1}^{n-1}k!\leavevmode\nobreak\ .$ (A.9)

##  Appendix B ABJ-like models

We consider here a matrix model of non-Hermitean type which does satisfy the
same Virasoro (homogeneous) constraints as those considered in the bulk of
this article. Let us define a $\beta$-deformed ABJ-like matrix model by the
matrix integral generating function

$\displaystyle\tau_{N,M}(t)=\frac{1}{N!}\frac{1}{M!}\int_{\Gamma_{x}^{N}}\,\prod_{i=1}^{N}\operatorname{d}x_{i}\int_{\Gamma_{y}^{M}}\,\prod_{a=1}^{M}\operatorname{d}y_{a}$
$\displaystyle\frac{\prod_{1\leq i<j\leq N}|x_{i}-x_{j}|^{2\beta}\prod_{1\leq
a<b\leq
M}|y_{a}-y_{b}|^{2/\beta}}{\prod_{i=1}^{N}\prod_{a=1}^{M}|x_{i}-y_{a}|^{2}}\times$
(B.1)
$\displaystyle\times\exp\left(\sum_{s=1}^{\infty}t_{s}\left(\sum_{i=1}^{N}x_{i}^{s}-\frac{1}{\beta}\sum_{a=1}^{M}y_{a}^{s}\right)\right)$

where we observe that there are two sets of eigenvalues variables
$\\{x_{i}\\}$ and $\\{y_{a}\\}$, each integrated over the corresponding
contour $\Gamma_{x}$ or $\Gamma_{y}$, respectively. Moreover, the familiar
Vandermonde determinant has been replaced by a ($\beta$-deformed) Cauchy
determinant in the first line. This model can also be interpreted as the
matrix model corresponding to an integral over supermatrices in the algebra of
the supergroup $U(N|M)$ [15] which plays an important role in computing the
localized partition function of the ABJ theory [16, 17] and also for the
analytically continued Chern-Simons partition function on lens space [18].

The expression in (B.1) is to be regarded as a formal integral due to the
absence of a potential which would make the integral convergent, however we
remark that a polynomial potential can always be introduced by shifting the
time variables as in (2.3). The Virasoro constraints can then be derived by
acting with the Lie derivative along the vectors

$\xi_{n}=\sum_{i=1}^{N}x_{i}^{n+1}\frac{\partial}{\partial
x_{i}}+\sum_{a=1}^{M}y_{a}^{n+1}\frac{\partial}{\partial y_{a}}$ (B.2)

and they can be rewritten as

$\left(\sum_{k=1}^{\mathsf{p}}\mathsf{a}_{k}\frac{\partial}{\partial
t_{k+n}}+\delta_{n,-1}\mathsf{a}_{1}\left(N-\frac{1}{\beta}M\right)-L^{\mathrm{ABJ}}_{n}\right)\tau_{N,M}(\mathsf{a};t)=0\leavevmode\nobreak\
,\quad\quad\quad n\geq-1$ (B.3)

with operators $L_{n}^{\mathrm{ABJ}}$ given by

$\displaystyle L^{\mathrm{ABJ}}_{n>0}$
$\displaystyle=2\beta\left(N-\frac{1}{\beta}M\right)\frac{\partial}{\partial
t_{n}}+\beta\sum_{\mu+\nu=n}\frac{\partial^{2}}{\partial t_{\mu}\partial
t_{\nu}}+\left(1-\beta\right)(n+1)\frac{\partial}{\partial
t_{n}}+\sum_{s>0}st_{s}\frac{\partial}{\partial t_{s+n}}$ (B.4) $\displaystyle
L^{\mathrm{ABJ}}_{0}$
$\displaystyle=\beta\left(N-\frac{1}{\beta}M\right)^{2}+\left(1-\beta\right)\left(N-\frac{1}{\beta}M\right)+\sum_{s>0}st_{s}\frac{\partial}{\partial
t_{s}}$ $\displaystyle L^{\mathrm{ABJ}}_{-1}$
$\displaystyle=\left(N-\frac{1}{\beta}M\right)t_{1}+\sum_{s>0}st_{s}\frac{\partial}{\partial
t_{s-1}}\leavevmode\nobreak\ .$

It is then clear that if we define the effective rank

$N_{\mathrm{eff}}:=N-\frac{1}{\beta}M$ (B.5)

the constraints in (B.3) became equivalent to those in (2.4) for the HMM.
Therefore, for an appropriate choice of initial data, the solutions to the two
sets of equations will be equivalent. In other words,

$\tau_{N,M}\left(\mathsf{a};t\right)\cong\tau_{N_{\mathrm{eff}}}\left(\mathsf{a};t\right)$
(B.6)

This shows that the same Virasoro constraints do have different solutions from
the HMM one, in this case not of matrix model type but supermatrix model.

Finally, we observe that the ABJ generating function in (B.1) and the
constraints in (B.3) are both invariant under the symmetry transformation

$\beta\to\frac{1}{\beta},\quad\quad\quad N_{\mathrm{eff}}\to-\beta
N_{\mathrm{eff}},\quad\quad\quad t_{s}\to-\frac{1}{\beta}t_{s},\quad\quad\quad
x_{i}\leftrightarrow y_{a},\quad\quad\quad N\leftrightarrow M.$ (B.7)

Notice that this is a generalization of the symmetry in (2.8) which can easily
be recovered for $M=0$. Moreover, there is a special value of $\beta$ for
which the solution becomes homogeneous w.r.t. a rescaling of effective rank
and times, namely for

$\beta=-\frac{1-M}{N}$ (B.8)

we have that the generating function satisfies

$\frac{\tau_{N,M}(\mathsf{a};t)}{\tau_{N,M}(\mathsf{a};0)}=\frac{\tau_{N/\kappa,M}(\kappa\mathsf{a};\kappa\,t)}{\tau_{N/\kappa,M}(\kappa\mathsf{a};0)}=\frac{\tau_{N,\kappa(M-1)+1}(\kappa\mathsf{a};\kappa\,t)}{\tau_{N,\kappa(M-1)+1}(\kappa\mathsf{a};0)}\leavevmode\nobreak\
,\quad\quad\quad\kappa\in\mathbb{C}^{\times}$ (B.9)

which then implies the universal relation

$\tau_{N,M}(\mathsf{a};t)=\frac{\tau_{1,0}\left(\frac{N}{1-M}\mathsf{a};\frac{N}{1-M}t\right)}{\tau_{1,0}\left(\frac{N}{1-M}\mathsf{a};0\right)}\tau_{N,M}(\mathsf{a};0)\leavevmode\nobreak\
.$ (B.10)

##  Appendix C Cohomological description of the constraints

In this section we use Lie algebra cohomology (with coefficients in a module)
to describe the relation between homogeneous and non-homogeneous Virasoro
constraints. In particular we restrict to the Gaussian case, i.e.
$\mathsf{p}=2$ and $\mathsf{a}_{1}=0$, $\mathsf{a}_{2}=1$.

We consider the parabolic subalgebra $\mathfrak{p}\subset\mathrm{Vir}$
generated by $\\{U_{-1},U_{0},U_{1},\dots\\}$, with
$U_{n}:=L_{n}-\partial_{n+2}$. Let $V_{N}=\mathbb{C}[[t_{1},t_{2},\dots]]$ be
the ring of formal series in times together with a $\mathfrak{p}$-module
structure given by the action of the Virasoro generators as differential
operators in times. Observe that all $V_{N}$ for $N\in\mathbb{Z}$ are
isomorphic to each other as vector spaces but not as $\mathfrak{p}$-modules.
They differ by the action of the Virasoro generators $U_{n}$ which explicitly
depends on the integer $N$.

We define cochains
$C^{n}(\mathfrak{p},V_{N})\cong\mathrm{Hom}(\Lambda^{n}\mathfrak{p},V_{N})$
and a differential $\delta:C^{n}(\mathfrak{p},V_{N})\to
C^{n+1}(\mathfrak{p},V_{N})$,

$\displaystyle(\delta f)(X_{1},\dots,X_{n+1})=$ $\displaystyle\sum_{1\leq
s<t\leq
n+1}(-1)^{s+t-1}f([X_{s},X_{t}],X_{1},\dots,\hat{X}_{s},\dots,\hat{X}_{t},\dots,X_{n+1})$
(C.1) $\displaystyle+\sum_{s=1}^{n+1}(-1)^{s}X_{s}\cdot
f(X_{1},\dots,\hat{X}_{s},\dots,X_{n+1})$

with $f\in C^{n}(\mathfrak{p},V_{N})$ and $X_{s}\in\mathfrak{p}$.

A 1-cochain $\alpha$ is a 1-cocycle if it is $\delta$-closed, i.e.

$\delta\alpha(X_{1},X_{2})=-\alpha([X_{1},X_{2}])-X_{1}\cdot\alpha(X_{2})+X_{2}\cdot\alpha(X_{1})=0$
(C.2)

and it is exact if it is the $\delta$-image of a 0-cochain $f\in
C^{0}(\mathfrak{g},V_{N})\cong V_{N}$,

$\alpha(X)=\delta f(X)=-X\cdot f$ (C.3)

If we define a 1-cochain $\alpha\in C^{1}(\mathfrak{g},V_{N})$ by
$\alpha(U_{n}):=B_{n}$, with $B_{n}\in V_{N}$ defined as in (3.7), then we can
interpret equation (3.3) as saying that $\tau_{N}\in V_{N}$ is a 0-cochain and
that $\alpha$ is its differential,

$\delta\tau_{N}=\alpha$ (C.4)

so that $\alpha$ is trivial in cohomology and $\tau_{N}$ is a trivialization.
In the Lie algebra cohomology language one says that $\alpha$ is an inner
derivation of $\mathfrak{p}$. Moreover, two different trivializations differ
by a 0-cocycle $Z_{N}\in H^{0}(\mathfrak{p},V_{N})\subset V_{N}$, i.e. a
solution to the homogeneous Virasoro constraints

$\delta Z_{N}(U_{n})=-U_{n}\cdot Z_{N}=0,\quad\quad\quad
n\geq-1\leavevmode\nobreak\ .$ (C.5)

Therefore we have the familiar statement that non-homogeneous solutions form
an affine space modeled over the space of homogeneous solutions
$H^{0}(\mathfrak{p},V_{N})$.

Because $\alpha$ is exact then it must also be closed and its cocycle
condition can be written as

$\delta\alpha(U_{m},U_{n})=-\alpha([U_{m},U_{n}])-U_{m}\cdot\alpha(U_{n})+U_{n}\cdot\alpha(U_{m})=0\leavevmode\nobreak\
,$ (C.6)

which we can rewrite using the explicit expression of the commutator in
$\mathfrak{p}$ to get

$(m-n)B_{m+n}=-U_{m}\cdot B_{n}+U_{n}\cdot B_{m}\leavevmode\nobreak\ .$ (C.7)

This means that (C.7) is a necessary (but not sufficient) condition on the
non-homogeneous terms $B_{n}$ in order for the constraints to have a solution.

## References

  * [1] P. H. Ginsparg and G. W. Moore, “Lectures on 2-D gravity and 2-D string theory,” in Theoretical Advanced Study Institute (TASI 92): From Black Holes and Strings to Particles, pp. 277–469.  10, 1993.  arXiv:hep-th/9304011.
  * [2] A. Mironov and A. Morozov, “On the origin of Virasoro constraints in matrix models: Lagrangian approach,” Phys. Lett. B252 (1990) 47–52.
  * [3] A. Mironov, A. Morozov, and Z. Zakirova, “Discrete Painlevé equation, Miwa variables and string equation in 5d matrix models,” JHEP 19 (2020) 227, arXiv:1908.01278 [hep-th].
  * [4] A. Morozov and S. Shakirov, “Generation of Matrix Models by W-operators,” JHEP 04 (2009) 064, arXiv:0902.2627 [hep-th].
  * [5] C. Cordova, B. Heidenreich, A. Popolitov, and S. Shakirov, “Orbifolds and Exact Solutions of Strongly-Coupled Matrix Models,” Commun. Math. Phys. 361 no. 3, (2018) 1235–1274, arXiv:1611.03142 [hep-th].
  * [6] L. Cassia, R. Lodin, and M. Zabzine, “On matrix models and their $q$-deformations,” JHEP 10 (2020) 126, arXiv:2007.10354 [hep-th].
  * [7] A. Morozov, “Challenges of beta-deformation,” Theor. Math. Phys. 173 (2012) 1417–1437, arXiv:1201.4595 [hep-th].
  * [8] L. F. Alday, D. Gaiotto, and Y. Tachikawa, “Liouville Correlation Functions from Four-dimensional Gauge Theories,” Lett. Math. Phys. 91 (2010) 167–197, arXiv:0906.3219 [hep-th].
  * [9] N. Wyllard, “A(N-1) conformal Toda field theory correlation functions from conformal N = 2 SU(N) quiver gauge theories,” JHEP 11 (2009) 002, arXiv:0907.2189 [hep-th].
  * [10] V. G. Kac, A. K. Raina, and N. Rozhkovskaya, Bombay lectures on highest weight representations of infinite dimensional Lie algebras.  Advanced Series in Mathematical Physics. World Scientific, Singapore, 2013\.  http://cds.cern.ch/record/1601455.
  * [11] T. Miwa, E. Date, and M. Jimbo, Solitons: Differential Equations, Symmetries and Infinite Dimensional Algebras, vol. 135 of Cambridge Tracts in Mathematics.  4, 2020.
  * [12] A. Alexandrov and A. Zabrodin, “Free fermions and tau-functions,” Journal of Geometry and Physics 67 (May, 2013) 37–80. http://dx.doi.org/10.1016/j.geomphys.2013.01.007.
  * [13] C. Andréief, “Note sur une relation entre les intégrales définies des produits des fonctions,” Mém. Soc. Sci. Phys. Nat. Bordeaux (3) 2 (1886) 1–14. https://kms.univer.kharkov.ua/pdf/1879_1883/1883_1/1883_1_2_Andreev.pdf.
  * [14] R. Dijkgraaf, H. L. Verlinde, and E. P. Verlinde, “Loop equations and Virasoro constraints in nonperturbative 2-D quantum gravity,” Nucl. Phys. B348 (1991) 435–456. [,435(1990)].
  * [15] N. Drukker and D. Trancanelli, “A Supermatrix model for N=6 super Chern-Simons-matter theory,” JHEP 02 (2010) 058, arXiv:0912.3006 [hep-th].
  * [16] O. Aharony, O. Bergman, D. L. Jafferis, and J. Maldacena, “N=6 superconformal Chern-Simons-matter theories, M2-branes and their gravity duals,” JHEP 10 (2008) 091, arXiv:0806.1218 [hep-th].
  * [17] O. Aharony, O. Bergman, and D. L. Jafferis, “Fractional M2-branes,” JHEP 11 (2008) 043, arXiv:0807.4924 [hep-th].
  * [18] M. Marino and P. Putrov, “Exact Results in ABJM Theory from Topological Strings,” JHEP 06 (2010) 011, arXiv:0912.3074 [hep-th].