Abstract
We attempt to systematically develop the matrix-model/quantum-theory correspondence by working out explicitly various non-trivial examples.
The path integral and scalar field theory
Path integral
We consider a single free particle moving in one dimension. The solution $|{\psi}_s(t)>$ of the Schrodinger equation is then of the form
\begin{eqnarray}
|{\psi}_s(t)>&=&\int dx <x|{\psi}_s(t)>|x>\nonumber\\
&=&\int dx \int dx_0 G(x,t;x_0,t_0)<x_0|{\psi}_s(t_0)>.
\end{eqnarray}
The Green function $G(x,t;x_0,t_0)$ is the transition amplitude from the point $x_0$ at time $t_0$ to the point $x$ at time $t$ which is the most basic object in the quantum theory. It is defined by
\begin{eqnarray}
G(x,t;x_0,t_0)&=&<x,t|x_0,t_0>\nonumber\\
&=&<x|e^{-\frac{i}{\hbar}H(t-t_0)} |x_0>\nonumber\\
&=&\int {\cal D}p{\cal D}x~e^{\frac{i}{\hbar}\int_{t_0}^{t} ds(p\dot{x}-H(p,x))}\nonumber\\
&=&{\cal N}\int {\cal D}x~e^{\frac{i}{\hbar}S[x]}.
\end{eqnarray}
This fundamental result holds for any Hamiltonian of the form $H=p^2/2m+ V(x)$ and thus $S[x]$ is the corresponding action given in terms of the Lagrangian $L(x,\dot{x})$ by the usual formula $S[x]=\int dt~L(x,\dot{x})=\int dt~(m\dot{x}^2/2-V(x))$.
The generalization of the above result to matrix elements of operators is given by (with $T$ being the time-ordering operator)
\begin{eqnarray}
<x,t|T(X(t_1)...X(t_n))|x_0,t_0>
&=&{\cal N}\int {\cal D}x~x(t_1)...x(t_n)~e^{\frac{i}{\hbar}S[x]}.\label{me}
\end{eqnarray}
In the limit $t_0\longrightarrow -\infty$ and $t\longrightarrow \infty$ only the ground state contributes. Thus, we obtain in this limit the matrix elements $<0|T(X(t_1)...X(t_n))|0>$ given by
\begin{eqnarray}
<0|T(X(t_1)...X(t_n))|0>
&=&\frac{\int {\cal D}x~x(t_1)...x(t_n)~e^{\frac{i}{\hbar}S[x]}}{\int {\cal D}x~e^{\frac{i}{\hbar}S[x]}}.
\end{eqnarray}
We introduce the path integral $Z[J]$ in the presence of a source $J(t)$ by
\begin{eqnarray}
Z[J]
&=&\int {\cal D}x~e^{\frac{i}{\hbar}S[x]+\frac{i}{\hbar}\int dt J(t)x(t)}.
\end{eqnarray}
This path integral is the generating functional of all the matrix elements $<0|T(X(t_1)...X(t_n))|0>$ where $X(t)$ is the coordinate operator. Indeed, we have
\begin{eqnarray}
<0|T(X(t_1)...X(t_n))|0>
&=&\frac{1}{Z[0]}\bigg(\frac{\hbar}{i}\bigg)^n\frac{{\delta}^nZ[J]}{{\delta}J(t_1)...\delta J(t_n)}|_{J=0}.
\end{eqnarray}
From this discussion $Z[J]$ is seen to be the vacuum-to-vacuum amplitude in the presence of the source $J(t)$.
Scalar phi-four theory
Next, we consider a phi-four field theory in a $d-$dimensional Minkowski spacetime. The action is given explicitly by
\begin{eqnarray}
S[\phi]=\int d^dx\bigg[\frac{1}{2}{\partial}_{\mu}\phi{\partial}^{\mu}\phi -\frac{1}{2}m^2{\phi}^2-\frac{\lambda}{4}{\phi}^4\bigg].\label{scalar}
\end{eqnarray}
This is the only renormalizable interacting scalar field theory in $d=4$ dimensions.
By analogy with the case of the point particle the path integral of this scalar field theory is formally given by the formula
\begin{eqnarray}
Z[J]
&=&\int {\cal D}\phi~e^{\frac{i}{\hbar}S[\phi]+\frac{i}{\hbar}\int d^dx J(x)\phi(x)}.
\end{eqnarray}
This path integral is the generating functional of all the matrix elements $<0|T(\Phi(x_1)...\Phi(x_n))|0>$ (also called $n-$point functions).
Indeed, we have
\begin{eqnarray}
<0|T(\Phi(x_1)...\Phi(x_n))|0>
&=&\frac{1}{Z[0]}\bigg(\frac{\hbar}{i}\bigg)^n\frac{{\delta}^nZ[J]}{{\delta}J(x_1)...\delta J(x_n)}|_{J=0}\nonumber\\
&=&\frac{\int {\cal D}\phi~\phi(x_1)...\phi(x_n)~e^{\frac{i}{\hbar}S[\phi]}}{\int {\cal D}\phi~e^{\frac{i}{\hbar}S[\phi]}}.
\end{eqnarray}
Euclidean rotation
Euclidean spacetime is obtained from Minkowski spacetime via the so-called Wick rotation. This is also called the imaginary time formulation which is obtained by the substitutions $t\longrightarrow -i\tau$, $x^0=ct\longrightarrow-ix^d=-ic\tau$, ${\partial}_0\longrightarrow i{\partial}_d$. Hence ${\partial}_{\mu}\phi{\partial}^{\mu}\phi\longrightarrow -({\partial}_{\mu}\phi)^2$ and $iS\longrightarrow -S_E$ where
\begin{eqnarray}
S_E[\phi]=\int d^dx\bigg[\frac{1}{2}({\partial}_{\mu}\phi)^{2} +\frac{1}{2}m^2{\phi}^2+\frac{\lambda}{4}{\phi}^4\bigg].
\end{eqnarray}
The Euclidean path integral and the Euclidean $n-$point functions become
\begin{eqnarray}
Z_E[J]=\int {\cal D}\phi ~e^{-\frac{1}{\hbar}S_E[\phi]+\frac{1}{\hbar}\int d^d x J(x)\phi(x)}.
\end{eqnarray}
\begin{eqnarray}
<0|T({\Phi}(x_1)...{\Phi}(x_n))|0>_E&=&\frac{\hbar^n}{Z[0]}\frac{{\delta}^nZ_E[J]}{{\delta}J(x_1)...{\delta}J(x_n)}|_{J=0}\nonumber\\
&=&\frac{\int {\cal D}\phi~{\phi}(x_1)...{\phi}(x_n) ~e^{-\frac{1}{\hbar}S_E[\phi]}}{\int {\cal D}\phi ~e^{-\frac{1}{\hbar}S_E[\phi]}}.
\end{eqnarray}
Thermal field theory
In the case of the point particle we have derived the formula (we set $\hbar=1$)
\begin{eqnarray}
<x|e^{-i H(t-t_0)} |x_0>={\cal N}\int {\cal D}x~e^{i S[x]}.
\end{eqnarray}
After Euclidean rotation we obtain
\begin{eqnarray}
<x|e^{-H(\tau-\tau_0)} |x_0>={\cal N}\int {\cal D}x~e^{-S_E[x]}.
\end{eqnarray}
The thermodynamical partition function is given by (we set $k_B=1$, i.e. $\beta=1/T$)
\begin{eqnarray}
Z&=&Tre^{-\beta H}\nonumber\\
&=&\int dx <x|e^{-\beta H} |x>\nonumber\\
&=&{\cal N}\int_{\rm periodic} {\cal D}x~e^{-S_E[x]}.\label{per}
\end{eqnarray}
In other words, we have
\begin{eqnarray}
\beta=\tau-\tau_0.
\end{eqnarray}
And
\begin{eqnarray}
x\equiv x_0\iff x(\tau_0+\beta)= x(\tau_0).
\end{eqnarray}
Hence, temperature is introduced into quantum field theory by restricting the Euclidean (or imaginary) time to a finite interval $[0,\beta]$ and imposing periodic boundary conditions on the coordinate degrees freedom $x(\tau)$. Indeed, the integration measure in (\ref{per}) is over all closed paths with a time period equal the inverse temperature while the Euclidean action $S_E$ is given explicitly by
\begin{eqnarray}
S_E=\int_0^{\beta} d\tau L_E(x,\dot{x}).
\end{eqnarray}
Vector scalar phi-four theory
The vector scalar phi-four theory is a generalization of the scalar phi-four theory (\ref{scalar}) to an $O(N)-$symmetric theory given by the Euclidean action
\begin{eqnarray}
S[\phi]=-\int d^dx \bigg(\frac{1}{2}(\partial_{\mu}\phi^i)^2+\frac{1}{2}m^2\phi^i\phi^i+\frac{\lambda}{4}(\phi^i\phi^i)^2\bigg).\label{contaction}
\end{eqnarray}
Lattice regularization
Field-theoretic calculations are made more explicit and more rigorous by working on an Euclidean lattice spacetime. In fact, the Euclidean lattice provides a concrete non-perturbative definition of the theory.
Thus, we will employ lattice regularization in which $x=an$, $\int d^dx =a^d\sum_n$, $\phi^i(x)=\phi_n^i$ and $\partial_{\mu}\phi^i=(\phi_{n+\hat{\mu}}^i-\phi_n^i)/a$.
The lattice action (\ref{contaction}) becomes then
\begin{eqnarray}
S[\phi]&=&\sum_n \bigg(a^{d-2}\sum_{\mu}\phi_n^i\phi_{n+\hat{\mu}}^i-\frac{a^{d-2}}{2}(m^2a^2+2d)\phi_n^i\phi_n^i-\frac{a^d\lambda}{4}(\phi_n^i\phi_n^i)^2\bigg)\nonumber\\
&=&\sum_n\bigg(2\kappa \sum_{\mu}\Phi_n^i\Phi_{n+\hat{\mu}}^i-{\Phi}_n^i\Phi_n^i-g(\Phi_n^i\Phi_n^i-1)^2 \bigg).
\end{eqnarray}
The mass parameter $m^2$ is replaced by the so-called hopping parameter $\kappa$ and the coupling constant $\lambda$ is replaced by the coupling constant $g$ where
\begin{eqnarray}
m^2a^2=\frac{1-2g}{\kappa}-2d~,~\frac{\lambda}{a^{d-4}}=\frac{g}{\kappa^2}.
\end{eqnarray}
The fields $\phi_n^i$ and $\Phi_n^i$ are related by
\begin{eqnarray}
\phi_n^i=\sqrt{\frac{2\kappa}{a^{d-2}}}\Phi_n^i.
\end{eqnarray}
The partition function is given by
\begin{eqnarray}
Z&=&\int \prod_{n,i} d\Phi_n^{i}~e^{S[\phi]}\nonumber\\
&=&\int d\mu(\Phi)~e^{2\kappa\sum_n\sum_{\mu}\Phi_n^i\Phi_{n+\hat{\mu}}^i}.
\end{eqnarray}
The measure $d\mu(\phi)$ is given explicitly by
\begin{eqnarray}
d\mu(\Phi)&=&\prod_{n,i} d\Phi_n^i ~e^{-\sum_n\big({\Phi}_n^i\Phi_n^i+g(\Phi_n^i\Phi_n^i-1)^2 \big)}\nonumber\\
&=&\prod_n \bigg(d^N\vec{\Phi}_n ~e^{-\vec{\Phi}_n^2-g(\vec\Phi_n^2-1)^2 }\bigg)\nonumber\\
&\equiv &\prod_n d\mu(\Phi_n).
\end{eqnarray}
This is a generalized Ising model. Indeed, the limit $g\longrightarrow \infty $ of the $O(1)$ model is precisely the Ising model in $d$ dimensions. The limit $g\longrightarrow \infty $ of the $O(3)$ model corresponds to the Heisenberg model in $d$ dimensions. The $O(N)$ models on the lattice are thus intimately related to spin models.
Matrix scalar field theory
Matrix phi-four theory is a very well known model which depends on a single hermitian matrix $M$ given by the action
\begin{eqnarray}
V&=&B{ Tr} M^2+C{Tr} M^4\nonumber\\
&=&\frac{N}{g}\bigg(-{ Tr} M^2+\frac{1}{4} { Tr} M^4\bigg).\label{quarticM}
\end{eqnarray}
The model depends actually on a single coupling $g$ such that (by scaling the field $M$ as $M\longrightarrow \alpha M$)
\begin{eqnarray}
B\alpha^2=-\frac{N}{g}~,~C\alpha^4=\frac{N}{4g}\Rightarrow B^2=\frac{4NC}{g}.
\end{eqnarray}
The partition function (path integral) is given by
\begin{eqnarray}
Z=\int d M ~e^{-V}.
\end{eqnarray}
We can now diagonalize the scalar matrix $M$ as
\begin{eqnarray}
M=U\Lambda U^{-1}.
\end{eqnarray}
The measure is therefore given by
\begin{eqnarray}
d M= d\Lambda dU \Delta^2(\Lambda).
\end{eqnarray}
The $dU$ is the usual Haar measure over the group $SU(N)$ which is normalized such that $\int dU=1$ whereas the Jacobian $\Delta^2(\Lambda)$ is precisely the so-called Vandermonde determinant given explicitly by
\begin{eqnarray}
\Delta^2(\Lambda)= \prod_{i>j}(\lambda_i-\lambda_j)^2.
\end{eqnarray}
The partition function becomes
\begin{eqnarray}
Z=\int d \Lambda~\Delta^2(\Lambda) ~\exp\big(-{ Tr}\big(B{\Lambda}^2+C{\Lambda}^4\big)\big).
\end{eqnarray}
We are therefore dealing with an effective potential given by
\begin{eqnarray}
V_{\rm eff}=B\sum_{i=1}\lambda_i^2+C\sum_{i=1}\lambda_i^4-\frac{1}{2}\sum_{i\neq j}\ln (\lambda_i-\lambda_j)^2.
\end{eqnarray}
In the large $N$ limit we can apply the saddle point method and thus the dominant configuration of eigenvalues $\Lambda=(\lambda_1,...,\lambda_N)$ must be a solution of the equation of motion
\begin{eqnarray}
\frac{dV_{\rm eff}}{d\lambda_i}=0.
\end{eqnarray}
There are two stable phases in this model given in terms of an eigenvalue distribution $\rho(\lambda)$ as follows \cite{Shimamune:1981qf}:
-
Disordered phase (one-cut) for $g\geq g_c$
This is characterized by the eigenvalues distribution of the matrix $M$ given by
\begin{eqnarray}
\rho(\lambda)&=&\frac{1}{N\pi}(2C\lambda^2+B+C\delta^2)\sqrt{\delta^2-\lambda^2}\nonumber\\
&=&\frac{1}{g\pi}\bigg(\frac{1}{2}\lambda^2-1+r^2\bigg)\sqrt{4r^2-\lambda^2}.\label{pred1}
\end{eqnarray}
This is a single cut solution with the cut defined by
\begin{eqnarray}
-2r\leq \lambda\leq 2r.
\end{eqnarray}
\begin{eqnarray}
r=\frac{1}{2}\delta.
\end{eqnarray}
\begin{eqnarray}
\delta^2&=&\frac{1}{3C}(-B+\sqrt{B^2+12 NC})=\frac{1}{3}(1+\sqrt{1+3g}).
\end{eqnarray}
-
Non-uniform ordered phase (two-cut) for $g\le g_c$
This is characterized by the eigenvalues distribution of the matrix $M$ given by
\begin{eqnarray}
\rho(\lambda)&=&\frac{2C|\lambda|}{N\pi}\sqrt{(\lambda^2-\delta_1^2)(\delta_2^2-\lambda^2)}\nonumber\\
&=&\frac{|\lambda|}{2g\pi}\sqrt{(\lambda^2-r_{-}^2)(r_{+}^2-\lambda^2)}.\label{pred2}
\end{eqnarray}
Here there are two cuts defined by
\begin{eqnarray}
r_{-}\leq |\lambda|\leq r_{+}.
\end{eqnarray}
\begin{eqnarray}
r_{-}=\delta_1~,~r_{+}=\delta_2.
\end{eqnarray}
\begin{eqnarray}
r_{\mp}^2&=&\frac{1}{2C}(-B\mp 2\sqrt{NC})\nonumber\\
&=&2(1\mp \sqrt{g}).
\end{eqnarray}
\begin{eqnarray}
g_c=1\leftrightarrow B_c^2=4NC \leftrightarrow B_c=-2\sqrt{NC}.
\end{eqnarray}
The scaled parameters are determined in the large $N$ limit to be given by
\begin{eqnarray}
\tilde{B}=\frac{B}{N^{3/2}}~,~\tilde{C}=\frac{C}{N^2}.
\end{eqnarray}
The critical line in the plane $(x,y)=(\tilde{C},-\tilde{B})$ lies at $-\tilde{B}_c=2\sqrt{\tilde{C}}$. This line for a fixed value of $\tilde{C}$ separates the disordered phase with $\langle M\rangle=0$ for $-\tilde{B}\lt -\tilde{B}_c$ from the non-uniform ordered phase with $\langle M\rangle\ne 0$ for $-\tilde{B}\gt -\tilde{B}_c$.
For $C=0$ (perturbative regime) the eigenvalues distribution should be given by the celebrated Wigner semi-circle law, viz
\begin{eqnarray}
\rho(\lambda)&=&\frac{\tilde{B}}{\pi}\sqrt{\delta^2-\lambda^2}~,~\delta^2=\frac{2}{\tilde{B}}.
\end{eqnarray}
We will use the Metropolis algorithm to study numerically the physics of this model. Under the change $\lambda_i\longrightarrow \lambda_i+h$ of the eigenvalue $\lambda_i$ the above effective potential changes as $V_{\rm eff}\longrightarrow V_{\rm eff}+\Delta V_{i,h}$ where
\begin{eqnarray}
\Delta V_{i,h}=B\Delta S_2+C\Delta S_4+\Delta S_{\rm Vand}.
\end{eqnarray}
The variations $\Delta S_2$, $\Delta S_4$ and $\Delta S_{\rm Vand}$ are given explicitly by
\begin{eqnarray}
\Delta S_2=h^2+2h\lambda_i.
\end{eqnarray}
\begin{eqnarray}
\Delta S_4=6h^2\lambda_i^2+4h\lambda_i^3+4h^3\lambda_i+h^4.
\end{eqnarray}
\begin{eqnarray}
\Delta S_{\rm Vand}=-2\sum_{j \ne i}\ln|1+\frac{h}{\lambda_i-\lambda_j}|.
\end{eqnarray}
The metropolis accept/reject step is simply given by the transition probability
\begin{eqnarray}
W(\lambda_i\longrightarrow \lambda_i+h)={\rm min}(1,\exp(-\Delta V_{i,h})).
\end{eqnarray}
As a test for our simulations we measure an exact result given by the Schwinger-Dyson identity
\begin{eqnarray}
\langle 2BTr M^2+4CTrM^4\rangle=N^2.
\end{eqnarray}
We also measure the energy given by the average value of the action $\langle V\rangle$, the specific heat $C_v= \langle(V-\langle V\rangle)^2\rangle$, the magnetization $m=|Tr M|$ and the suceptibility $\chi=\langle(m-\langle m\rangle)^2\rangle$.
The M-(atrix) BFSS theory (or lattice string theory)
The so-called M-(atrix) theory is a matrix quantum mechanics given by
the BFSS matrix model \cite{Banks:1996vh} which admits supergravity in
$11$ dimensions as a low energy limit and thus can provide the UV
completion of the $11-$dimensional M-theory which unifies all five
superstring theories in $10$ dimensions and provides their
non-perturbative formulation.
The degrees of freedom of the BFSS matrix model are: a hermitian $N\times N$ time-dependent bosonic matrix (representing the gauge field), $d=9$ hermitian $N\times N$ time-dependent bosonic matrices $\Phi_i$ (representing the coordinates) and $16$ hermitian $N\times N$ time-dependent fermionic matrices $\psi_{\alpha}$ (in other words $\psi$ is a Majorana-Weyl fermion in $d+1=10$ dimensions implementing supersymmetry). The BFSS action is given explicitly by
\begin{eqnarray}
S=\frac{1}{g^2}\int_0^{\beta}dt{\rm Tr}\bigg[\frac{1}{2}(D_t\Phi_i)^2-\frac{1}{4}[\Phi_i,\Phi_j]^2+\frac{1}{2}\psi_{\alpha}D_t\psi_{\alpha}-\frac{1}{2}\psi_{\alpha}(\gamma_i)_{\alpha\beta}[\Phi_i,\psi_{\beta}]\bigg].\label{BFSS}
\end{eqnarray}
This action arises from the reduction of $10-$dimensional supersymmetric Yang-Mills gauge theory to one dimension and it describes the dynamics of $N$ coincident D0-branes in type IIA string theory. In some precise sense, this action can thus be thought of as a non-perturbative definition of M-theory hence its alternative name M-(atrix) theory.
The energy of the bosonic truncation of the above BFSS matrix model (\ref{BFSS}) is defined by (where $Z(\beta)$ is the partition function at temperature $T=1/\beta$)
\begin{eqnarray}
E=-\frac{1}{Z(\beta)}\frac{Z(\beta^{'})-Z(\beta)}{\Delta\beta}~,~\Delta\beta=\beta^{'}-\beta\longrightarrow 0.\label{formu}
\end{eqnarray}
We compute immediately \cite{Hotta:1998en}
\begin{eqnarray}
\frac{E}{N^2}=\frac{3T}{N^2}\langle {\rm commu}\rangle~,~{\rm commu}=-\frac{1}{4g^2}\int_{0}^{\beta}dt{Tr}[{\Phi}_i^{},{\Phi}_j^{}]^2.
\end{eqnarray}
The proof goes as follows. We relate the partition functions $Z(\beta^{'})$ and $Z(\beta)$ by the following scalings
\begin{eqnarray}
\frac{t^{\prime}}{t}=\frac{\beta^{'}}{\beta}~,~\frac{A^{\prime}}{A}=\frac{\beta}{\beta^{'}}~,~\frac{\Phi^{\prime}}{\Phi}=\sqrt{\frac{\beta^{'}}{\beta}}.
\end{eqnarray}
This guarantees that the kinetic term is fully invariant. We assume that the measures over $\Phi_i(t)$ and $A(t)$ are also invariant under these scalings. Then from the non-invariance of the Yang-Mills term we can derive the above formula of the energy.
Other important observers (which act as order parameters) are the Polyakov line $\langle |P|\rangle$ and the extent of the space $R^2$ defined below.
At high temperatures the bosonic part of the BFSS quantum mechanics reduces to the bosonic part of the IKKT model \cite{Kawahara:2007ib}. The leading behavior of the various observables of interest at high temperatures can be obtained in terms of the corresponding expectation values in the IKKT model. We get then
\begin{eqnarray}R^2=\sqrt{T}\chi_1.
\end{eqnarray}
\begin{eqnarray}
\langle |P|\rangle=1-\frac{1}{2d}{T}^{-3/2}\chi_1.
\end{eqnarray}
\begin{eqnarray}
\frac{E}{N^2}=\frac{3}{4}{T}\chi_2~,~\chi_2=(d-1)(1-\frac{1}{N^2}).
\end{eqnarray}
The coefficient $\chi_1$ for various $d$ and $N$ can be read off from table $1$ of \cite{Kawahara:2007ib} whereas the coefficient $\chi_2$ was determined exactly from the Schwinger-Dyson equation.
This behavior can be used to calibrate our Monte Carlo simulations at high temperatures.
The degrees of freedom of the BFSS matrix model are: a hermitian $N\times N$ time-dependent bosonic matrix (representing the gauge field), $d=9$ hermitian $N\times N$ time-dependent bosonic matrices $\Phi_i$ (representing the coordinates) and $16$ hermitian $N\times N$ time-dependent fermionic matrices $\psi_{\alpha}$ (in other words $\psi$ is a Majorana-Weyl fermion in $d+1=10$ dimensions implementing supersymmetry). The BFSS action is given explicitly by
\begin{eqnarray}
S=\frac{1}{g^2}\int_0^{\beta}dt{\rm Tr}\bigg[\frac{1}{2}(D_t\Phi_i)^2-\frac{1}{4}[\Phi_i,\Phi_j]^2+\frac{1}{2}\psi_{\alpha}D_t\psi_{\alpha}-\frac{1}{2}\psi_{\alpha}(\gamma_i)_{\alpha\beta}[\Phi_i,\psi_{\beta}]\bigg].\label{BFSS}
\end{eqnarray}
This action arises from the reduction of $10-$dimensional supersymmetric Yang-Mills gauge theory to one dimension and it describes the dynamics of $N$ coincident D0-branes in type IIA string theory. In some precise sense, this action can thus be thought of as a non-perturbative definition of M-theory hence its alternative name M-(atrix) theory.
The energy of the bosonic truncation of the above BFSS matrix model (\ref{BFSS}) is defined by (where $Z(\beta)$ is the partition function at temperature $T=1/\beta$)
\begin{eqnarray}
E=-\frac{1}{Z(\beta)}\frac{Z(\beta^{'})-Z(\beta)}{\Delta\beta}~,~\Delta\beta=\beta^{'}-\beta\longrightarrow 0.\label{formu}
\end{eqnarray}
We compute immediately \cite{Hotta:1998en}
\begin{eqnarray}
\frac{E}{N^2}=\frac{3T}{N^2}\langle {\rm commu}\rangle~,~{\rm commu}=-\frac{1}{4g^2}\int_{0}^{\beta}dt{Tr}[{\Phi}_i^{},{\Phi}_j^{}]^2.
\end{eqnarray}
The proof goes as follows. We relate the partition functions $Z(\beta^{'})$ and $Z(\beta)$ by the following scalings
\begin{eqnarray}
\frac{t^{\prime}}{t}=\frac{\beta^{'}}{\beta}~,~\frac{A^{\prime}}{A}=\frac{\beta}{\beta^{'}}~,~\frac{\Phi^{\prime}}{\Phi}=\sqrt{\frac{\beta^{'}}{\beta}}.
\end{eqnarray}
This guarantees that the kinetic term is fully invariant. We assume that the measures over $\Phi_i(t)$ and $A(t)$ are also invariant under these scalings. Then from the non-invariance of the Yang-Mills term we can derive the above formula of the energy.
Other important observers (which act as order parameters) are the Polyakov line $\langle |P|\rangle$ and the extent of the space $R^2$ defined below.
At high temperatures the bosonic part of the BFSS quantum mechanics reduces to the bosonic part of the IKKT model \cite{Kawahara:2007ib}. The leading behavior of the various observables of interest at high temperatures can be obtained in terms of the corresponding expectation values in the IKKT model. We get then
\begin{eqnarray}R^2=\sqrt{T}\chi_1.
\end{eqnarray}
\begin{eqnarray}
\langle |P|\rangle=1-\frac{1}{2d}{T}^{-3/2}\chi_1.
\end{eqnarray}
\begin{eqnarray}
\frac{E}{N^2}=\frac{3}{4}{T}\chi_2~,~\chi_2=(d-1)(1-\frac{1}{N^2}).
\end{eqnarray}
The coefficient $\chi_1$ for various $d$ and $N$ can be read off from table $1$ of \cite{Kawahara:2007ib} whereas the coefficient $\chi_2$ was determined exactly from the Schwinger-Dyson equation.
This behavior can be used to calibrate our Monte Carlo simulations at high temperatures.
The matrix harmonic oscillator
For quenched fermions, low temperatures and large number of dimensions $d\longrightarrow\infty$ the BFSS matrix model is equivalent to a matrix harmonic oscillator problem given by the following simple matrix scalar field theory
\begin{eqnarray}
S[\Phi]=\frac{1}{g^2}\int_0^{\beta} dt Tr\bigg[\frac{1}{2}(\partial_t\Phi_i)^2+\frac{1}{2}m^2(\Phi_i)^2\bigg].\label{ac1}
\end{eqnarray}
The mass $m$ is given by
\begin{eqnarray}
m=d^{1/3}.
\end{eqnarray}
This matrix harmonic oscillator is a thermal field theory in one dimension and thus $t$ must be an imaginary time, viz $\tilde{t}=-it$ is the real time. The fields are periodic with period $\beta=1/T$ where $T$ is the Hawking temperature as
\begin{eqnarray}
\Phi_i(t+\beta)=\Phi_i(\beta).
\end{eqnarray}
Of course, we will study the system in the t'Hooft limit given by
\begin{eqnarray}
\lambda=g^2N.
\end{eqnarray}
There seems to be therefore two independent coupling constants $\lambda$ and $T$. However, $g^2$ can always be rescaled away. By dimensional analysis we find that $\Phi_i$ behaves as inverse length and $\lambda$ as inverse length cubed whereas $T$ behaves as inverse length and as a consequence the dimensionless coupling constant must be given by
\begin{eqnarray}
\tilde{\lambda}=\frac{\lambda}{T^3}.
\end{eqnarray}
We will choose $\lambda=1$, i.e. $g^2=1/N$.
First, we need to put the action on a lattice. We define $t=(n-1)a$, $n=1,...,\Lambda$ with $\beta=\Lambda a$. The periodicity condition becomes $\Phi_i(n+\Lambda)=\Phi_i(n)$. The time derivative is then given by the difference
\begin{eqnarray}
\partial_t\Phi_i(t)=\frac{\Phi_i(n+1)-\Phi_i(n)}{a}~,~a\longrightarrow 0.
\end{eqnarray}
The lattice action is then given by
\begin{eqnarray}
S=\frac{N}{a}\sum_{n=1}^{\Lambda}{Tr}\bigg[\Phi_i^2(n)-\Phi_i(n+1)\Phi_i(n)+\frac{m^2a^2}{2}\Phi_i^2(n)\bigg].
\end{eqnarray}
The exact variation of this action (due to the variation of the matrix $\Phi_i(n)$) is trivial given by
\begin{eqnarray}
\Delta S_i(n)= \frac{N}{a}{Tr}\bigg[2(1+\frac{m^2a^2}{2})\Phi_i(n)\delta \Phi_i(n)-(\Phi_i(n+1)+\Phi_i(n-1))\delta\Phi_i(n)+(1+\frac{m^2a^2}{2})\delta\Phi_i^2(n)\bigg].
\end{eqnarray}
We choose the update
\begin{eqnarray}
(\delta\Phi_i(n))_{pq}=h\delta_{pa}\delta_{qb}+h^*\delta_{qa}\delta_{pb}.
\end{eqnarray}
The variation is then given explicitly by
\begin{eqnarray}
\Delta S_i(n)^{ab}=\frac{N}{a}\bigg[4(1+\frac{m^2a^2}{2}){\rm Re}(h^*\phi_i(n))_{ab}-2{\rm Re}(h^*\Phi_i(n+1)+h^*\Phi_i(n-1))_{ab}+(1+\frac{m^2a^2}{2})(2hh^*+(h^2+h^{*2})\delta_{ab})\bigg].
\end{eqnarray}
The metropolis accept/reject step is simply given by the transition probability
\begin{eqnarray}
W((\Phi_i(n))_{ab}\longrightarrow (\Phi_i(n))_{ab}+h)={\rm min}(1,\exp(-\Delta S_{i}(n)^{ab})).
\end{eqnarray}
A one-dimensional gauge field
Interaction of the scalar fields $\Phi_i(t)$ with a one-dimensional U(N) gauge field $A(t)$ is implemented in the usual way through minimal coupling, viz
\begin{eqnarray}
\partial_t\Phi_i\longrightarrow D_t\Phi_i=\partial_t\Phi_i-i[A,\Phi_i]. \end{eqnarray}
Obviously, the scalar fields $\Phi_i(t)$ transform in the adjoint representation of the gauge group. The action (\ref{ac1}) becomes then
\begin{eqnarray}
S=\frac{1}{g^2}\int_0^{\beta}dt{Tr}\bigg[\frac{1}{2}(D_t\Phi_i)^2+\frac{1}{2}m^2(\Phi_i)^2\bigg].
\end{eqnarray}
The above action is invariant under the $U(N)$ gauge transformations
\begin{eqnarray}
\psi_{\alpha}\longrightarrow U\psi_{\alpha}U^{\dagger}~,~\Phi_i\longrightarrow U\Phi_iU^{\dagger}~,~A\longrightarrow UAU^{\dagger}-\frac{1}{i}U\partial_tU^{\dagger}.
\end{eqnarray}
We will now gauge-fix this symmetry non-perturbatively \cite{Filev:2015hia}.
The gauge field $A$ becomes the link variable
\begin{eqnarray}
U_{n,n+1}={\cal P}\exp(-i\int_{na}^{(n+1)a}dt A(t))=1-iaA(n)+O(a^2)~,~a\longrightarrow 0.\label{link}
\end{eqnarray}
The covariant derivative becomes then
\begin{eqnarray}
\frac{U_{n,n+1}\Phi_i(n+1)U_{n+1,n}-\Phi_i(n)}{a}=D_t\Phi_i~,~a\longrightarrow 0.
\end{eqnarray}
The action becomes given by
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^2(n)-\frac{1}{a}U_{n,n+1}\Phi_i(n+1)U_{n+1,n}\Phi_i(n)+\frac{1}{2}m^2a\Phi_i(n)^2\bigg].
\end{eqnarray}
We have a local $U(N)$ symmetry at each lattice site. We use this symmetry to rotate the link variables as follows
\begin{eqnarray}
&&\Phi^{\prime}_i(1)=\Phi_i(1)\nonumber\\
&&\Phi^{\prime}_i(2)=U_{1,2}\Phi_i(2)U_{1,2}^{\dagger}\nonumber\\
&&\Phi^{\prime}_i(3)=U_{1,2}U_{2,3}\Phi_i(3)(U_{1,2}U_{2,3})^{\dagger}\nonumber\\
&&..\nonumber\\
&&..\nonumber\\
&&\Phi^{\prime}_i(\Lambda)=U_{1,2}...U_{\Lambda-1,\Lambda}\Phi_i(\Lambda)(U_{1,2}...U_{\Lambda-1,\Lambda})^{\dagger}.
\end{eqnarray}
The action becomes then
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^{\prime 2}(n)+\frac{1}{2}a^2m^2\Phi_i^{\prime 2}(n)\bigg]-\frac{1}{g^2a}\bigg[\sum_{n=1}^{\Lambda-1}{Tr}\Phi_i^{\prime}(n)\Phi_i^{\prime}(n+1)+\Phi_i^{\prime}(\Lambda+1)W^{\dagger}\Phi_i^{\prime}(\Lambda)W\bigg].\nonumber\\
\end{eqnarray}
The matrix $W$ is defined by
\begin{eqnarray}
W=U_{1,2}...U_{\Lambda-1,\Lambda}U_{\Lambda,\Lambda+1}.
\end{eqnarray}
We can diagonalize the matrix $W$ in the usual way as $W=UDU^{\dagger}$ where $D={\rm diag}(\exp(i\theta_1),...,\exp(i\theta_N))$.
The action becomes (with $\tilde{\Phi}_i(n)=U^{\dagger}\Phi_i^{\prime}(n)U$)
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{Tr}\bigg[\frac{1}{a}\tilde{\Phi}_i^{ 2}(n)+\frac{1}{2}a^2m^2\tilde{\Phi}_i^{ 2}(n)\bigg]-\frac{1}{g^2a}\bigg[\sum_{n=1}^{\Lambda-1}{ Tr}\tilde{\Phi}_i^{}(n)\tilde{\Phi}_i^{}(n+1)+\tilde{\Phi}_i^{}(\Lambda+1)D^{\dagger}\tilde{\Phi}_i^{}(\Lambda)D\bigg].\nonumber\\
\end{eqnarray}
Alternatively, this action can also be rewritten in terms of $D_{\Lambda}=D^{1/\Lambda}$ and $\bar{\Phi}_i(n)=h_n^{\dagger}\tilde{\Phi}_i(n)h_n$, where $h_n=D_{\Lambda}^n$, as
\begin{eqnarray}
S&=&\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\bar{\Phi}_i^{ 2}(n)+\frac{1}{2}a^2m^2\bar{\Phi}_i^{ 2}(n)\bigg]\nonumber\\
&-&\frac{1}{g^2a}\sum_{n=1}^{\Lambda}{Tr}\bar{\Phi}_i^{}(n)D_{\Lambda}\bar{\Phi}_i^{}(n+1)D_{\Lambda}^{\dagger}.
\end{eqnarray}
The measure $\prod_{i,n}dX_i(n)$ is invariant under all the above unitary transformations. Finally, we reduce the measure over the transporter fields as follows
\begin{eqnarray}
\prod_{n=1}^{\Lambda}{\cal D}U_{n,n+1}&=&\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}U_{1,2}\nonumber\\
&\sim &\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}W\nonumber\\
&\sim &\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}U{\cal D}D\Delta^2(D)\nonumber\\
&\sim&\prod_{a=1}^Nd\theta_a.\prod_{a>b}|e^{i\theta_a}-e^{i\theta_b}|^2\nonumber\\
&\sim&\prod_{a=1}^Nd\theta_a.\prod_{a>b}\sin^2\frac{\theta_a-\theta_b}{2}\nonumber\\
&\sim &\prod_{a=1}^Nd\theta_a.\exp(-S_{\rm FP}).
\end{eqnarray}
The Faddeev-Popov gauge-fixing action is then given explicitly by
\begin{eqnarray}
S_{\rm FP}=-\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^2(n)-\frac{1}{a}U_{n,n+1}\Phi_i(n+1)U_{n+1,n}\Phi_i(n)+\frac{1}{2}m^2a\Phi_i(n)^2\bigg].\label{BFSS_appr}
\end{eqnarray}
The adopted gauge-fixing condition is effectively the so-called static gauge given by
\begin{eqnarray}
U_{n,n+1}=D_{\Lambda}={\rm diag}(\exp(i\frac{\theta_1}{\Lambda}),...,\exp(i\frac{\theta_N}{\Lambda}))~,~\forall n.
\end{eqnarray}
In terms of the gauge field $A$ this condition reads (by means of equation (\ref{link}))
\begin{eqnarray}
A(t)=-\frac{1}{\beta}{\rm diag}(\theta_1,...,\theta_N).
\end{eqnarray}
In summary, we are interested in the total action
\begin{eqnarray}
S_{\rm total}&=&N\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}{\Phi}_i^{ 2}(n)+\frac{1}{2}am^2{\Phi}_i^{ 2}(n)\bigg]\nonumber\\
&-&\frac{N}{a}\sum_{n=1}^{\Lambda}{Tr}{\Phi}_i^{}(n)D_{\Lambda}{\Phi}_i^{}(n+1)D_{\Lambda}^{\dagger}-\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
Before we discuss the variations of this action we present the most important observables. The Polyakov line is the order parameter of the Hagedorn transition in string theory and the deconfinment transition in gauge theory which is associated with the spontaneous breakdown of the $U(1)$ symmetry $A(t)\longrightarrow A(t)+C.{\bf 1}$ (see the extensive list of references in \cite{Kawahara:2007fn}). The Polyakov line is defined in terms of the holonomy matrix $U$ by the relation
\begin{eqnarray}
P=\frac{1}{N}Tr U~,~U={\cal P}\exp(-i\int_0^{\beta} dt A(t)).
\end{eqnarray}
We compute \cite{Kawahara:2007fn}
\begin{eqnarray}
U&=&U_{1,2}U_{2,3}...U_{\Lambda-1,\Lambda}U_{\Lambda,1}=D_{\Lambda}^{\Lambda}={\rm diag}(\exp(i\theta_1),...,\exp(i\theta_N).
\end{eqnarray}
Hence
\begin{eqnarray}
P=\frac{1}{N}\sum_ae^{i\theta_a}.
\end{eqnarray}
Another important observable is the radius (or extent of space or more precisely the extent of the eigenvalue distribution) is defined by
\begin{eqnarray}
R^2=\frac{a}{\Lambda N^2}\langle {\rm radius}\rangle~,~{\rm radius}=\frac{N}{a}\sum_{n=1}^{\Lambda}{Tr}{\Phi}_i^{2}(n).
\end{eqnarray}
The energy in the present model (\ref{BFSS_appr}) is given effectively by the extent of space. Indeed, we compute (by using the formula (\ref{formu})) for the model (\ref{BFSS_appr}) the energy
\begin{eqnarray}
\frac{E}{N^2}=\frac{a^2Tm^2}{N^{2}}\langle{\rm radius}\rangle=m^2R^2.
\end{eqnarray}
We also measure the eigenvalues distribution of the holonomy matrix $U$ and the eigenvalues distribution of the bosonic matrices $\Phi_i(n)$.
In our present model ($d$ gauged matrix harmonic oscillators) the variation of the bosonic matrices $\Phi_i(n)$ is straightforward (because of the non-zero mass term and the absence of the Yang-Mills term the so-called problem of flat directions does not pose itself here). The variation of the angles $\theta_a$ requires however a careful treatment due to the center of mass $\theta_{\rm cm}=\sum_{a=1}^N\theta_a/N$ which does not appear in the action and behaves as a random walker.
We compute as before the variation under the change (for example for $a\neq b$)
\begin{eqnarray}
(\Phi_i(n))_{ab}\longrightarrow (\Phi_i(n))_{ab}+h
\end{eqnarray}
to find the result
\begin{eqnarray}
\Delta S_{{\rm total}~i}(n)^{ab}=\frac{N}{a}\bigg[4(1+\frac{m^2a^2}{2}){\rm Re}(h^*\phi_i(n))_{ab}-2{\rm Re}(f^*\Phi_i(n+1)+g^*\Phi_i(n-1))_{ab}+(1+\frac{m^2a^2}{2})(2hh^*+(h^2+h^{*2})\delta_{ab})\bigg]
\end{eqnarray}
where $f^*$ and $g^*$ are given by
\begin{eqnarray}
f^*=h^*\exp(\frac{i}{\Lambda}(\theta_a-\theta_b))~,~g^*=h^*\exp(-\frac{i}{\Lambda}(\theta_a-\theta_b))
\end{eqnarray}
But there is also the variation under the change of the angle $\theta_c$ of the holonomy matrix, viz
\begin{eqnarray}
\theta_c\longrightarrow\theta_c^{\prime}=\theta_c+\alpha.\label{theta_var}
\end{eqnarray}
The relevant action here is given by
\begin{eqnarray}
S(\theta)&=&-\frac{N}{a}\sum_{n=1}^{\Lambda}\sum_{a,b}e^{-\frac{i}{\Lambda}(\theta_a-\theta_b)}({\Phi}_i(n))_{ab}({\Phi}_i(n+1))_{ba}\nonumber\\
&-&\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
The variation of this action is given explicitly by
\begin{eqnarray}
\Delta S_c(\theta)&=&-\frac{2N}{a}{\rm Re}\bigg(\sum_{n=1}^{\Lambda}\sum_{a\ne c}e^{\frac{i}{\Lambda}(\theta_a-\theta_c)}(e^{-i\frac{\alpha}{\Lambda}}-1)({\Phi}_i(n+1))_{ac}({\Phi}_i(n))_{ca}\bigg)\nonumber\\
&-&\sum_{a\ne c}\ln\sin^2\frac{\theta_a-\theta_c-\alpha}{2}+\sum_{a\ne c}\ln\sin^2\frac{\theta_a-\theta_c}{2}.
\end{eqnarray}
An extremely important remark is now in order. The action $S(\theta)$ does not depend on the center of mass $\theta_{\rm cm}=\sum_{a=1}^N\theta_a/N$. Indeed, the functional integration over $\theta_a$ is (almost) identical to the functional integration over $\tilde{\theta}_a=\theta_a-\theta_{\rm cm}$ which satisfies $\sum_{a}\tilde{\theta}_a=0$.
Furthermore, by means of the $U(1)$ gauge transformation $A(t)\longrightarrow A(t) +C .{\bf 1}$ (a large gauge transformation with a non-zero winding number) we can choose the static gauge $A(t)=-{\rm diag}(\theta_1,...,\theta_N)/\beta$ in such a way that (see for example \cite{Anagnostopoulos:2007fw})
\begin{eqnarray}
{\rm max}(\tilde{\theta}_a)-{\rm min}(\tilde{\theta}_a)\lt 2\pi.
\end{eqnarray}
We have then explicitly
\begin{eqnarray}
\int d\theta F(\Delta\theta)=\int d\tilde{\theta}\delta(\sum_a\tilde{\theta}_a)F(\Delta\tilde{\theta})\int d\theta_{\rm cm}.
\end{eqnarray}
Since $\theta=\tilde{\theta}+\theta_{\rm cm}$ and $-\pi\lt \theta\le +\pi$ while $ {\rm min}(\tilde{\theta}_a)\le\tilde{\theta}\le {\rm max}(\tilde{\theta}_a)$ we conclude immediately that the center of mass $\theta_{\rm cm}$ must be in the interval $]-\pi-{\rm max}(\tilde{\theta}_a),\pi-{\rm min}(\tilde{\theta}_a)]$. Hence the above integral becomes (with $\mu={\rm max}(\tilde{\theta}_a)-{\rm min}(\tilde{\theta}_a)$)
\begin{eqnarray}
\int d\theta F(\Delta\theta)=\int d\tilde{\theta}\delta(\sum_a\tilde{\theta}_a)F(\Delta\tilde{\theta})(2\pi+\mu).
\end{eqnarray}
Clearly, this is true as long as $\mu\le 2\pi$ while for $\mu\gt 2\pi$ the additional Boltzmann weight is identically zero. We have then the extra Boltzmann weight \cite{code_Hanada}
\begin{eqnarray}
&&w(\mu)=2\pi+\mu~,~\mu\le 2\pi\nonumber\\
&&w(\mu)=0~,~\mu\gt 2\pi.
\end{eqnarray}
In other words, we can replace the functional integration over $\theta_a$ with the functional integration over $\tilde{\theta}_a$ with an additional Boltzmann weight $w(\mu)$.
Thus, the variation (\ref{theta_var}) should be thought of as $\tilde{\theta}_c\longrightarrow \theta_c^{\prime}=\tilde{\theta}_c+\alpha$ where $\alpha$ is uniformly distributed in the range $]-\pi-{\rm max}(\tilde{\theta}_a),\pi-{\rm min}(\tilde{\theta}_a)]$.
\bibitem{Shimamune:1981qf}
Y.~Shimamune,
``On the Phase Structure of Large $N$ Matrix Models and Gauge Models,''
Phys.\ Lett.\ {\bf 108B}, 407 (1982).
doi:10.1016/0370-2693(82)91223-0
%%CITATION = doi:10.1016/0370-2693(82)91223-0;%%
%47 citations counted in INSPIRE as of 09 Dec 2019
%\cite{Banks:1996vh}
\bibitem{Banks:1996vh}
T.~Banks, W.~Fischler, S.~H.~Shenker and L.~Susskind,
``M theory as a matrix model: A Conjecture,''
Phys.\ Rev.\ D {\bf 55}, 5112 (1997) doi:10.1103/PhysRevD.55.5112
[hep-th/9610043].
%%CITATION = doi:10.1103/PhysRevD.55.5112;%% %2782 citations counted in INSPIRE as of 02 Oct 2019
%\cite{Kawahara:2007fn}
\bibitem{Kawahara:2007fn}
N.~Kawahara, J.~Nishimura and S.~Takeuchi,
``Phase structure of matrix quantum mechanics at finite temperature,''
JHEP {\bf 0710}, 097 (2007)
doi:10.1088/1126-6708/2007/10/097
[arXiv:0706.3517 [hep-th]].c
%%CITATION = doi:10.1088/1126-6708/2007/10/097;%%
%37 citations counted in INSPIRE as of 26 Apr 2017
%\cite{Filev:2015hia}
\bibitem{Filev:2015hia}
V.~G.~Filev and D.~O'Connor,
``The BFSS model on the lattice,''
JHEP {\bf 1605} (2016) 167
doi:10.1007/JHEP05(2016)167
[arXiv:1506.01366 [hep-th]].c
%%CITATION = doi:10.1007/JHEP05(2016)167;%%
%14 citations counted in INSPIRE as of 25 Dec 2016
%\cite{Anagnostopoulos:2007fw}\bibitem{Anagnostopoulos:2007fw}
K.~N.~Anagnostopoulos, M.~Hanada, J.~Nishimura and S.~Takeuchi,
``Monte Carlo studies of supersymmetric matrix quantum mechanics with sixteen supercharges at finite temperature,''
Phys.\ Rev.\ Lett.\ {\bf 100}, 021601 (2008)
doi:10.1103/PhysRevLett.100.021601
[arXiv:0707.4454 [hep-th]].
%%CITATION = doi:10.1103/PhysRevLett.100.021601;%%
%171 citations counted in INSPIRE as of 20 Dec 2019
%\cite{code_Hanada}
\bibitem{code_Hanada}
Masanori Hanada, ``Users’ Manual of Monte Carlo Code for the Matrix Model of M-theory'', https://sites.google.com/site/hanadamasanori/home/mmmm.
%\cite{Kawahara:2007ib}\bibitem{Kawahara:2007ib}
N.~Kawahara, J.~Nishimura and S.~Takeuchi,
``High temperature expansion in supersymmetric matrix quantum mechanics,''
JHEP {\bf 0712}, 103 (2007)
doi:10.1088/1126-6708/2007/12/103
[arXiv:0710.2188 [hep-th]].
%%CITATION = doi:10.1088/1126-6708/2007/12/103;%%
%41 citations counted in INSPIRE as of 28 Dec 2019
\partial_t\Phi_i\longrightarrow D_t\Phi_i=\partial_t\Phi_i-i[A,\Phi_i]. \end{eqnarray}
Obviously, the scalar fields $\Phi_i(t)$ transform in the adjoint representation of the gauge group. The action (\ref{ac1}) becomes then
\begin{eqnarray}
S=\frac{1}{g^2}\int_0^{\beta}dt{Tr}\bigg[\frac{1}{2}(D_t\Phi_i)^2+\frac{1}{2}m^2(\Phi_i)^2\bigg].
\end{eqnarray}
The above action is invariant under the $U(N)$ gauge transformations
\begin{eqnarray}
\psi_{\alpha}\longrightarrow U\psi_{\alpha}U^{\dagger}~,~\Phi_i\longrightarrow U\Phi_iU^{\dagger}~,~A\longrightarrow UAU^{\dagger}-\frac{1}{i}U\partial_tU^{\dagger}.
\end{eqnarray}
We will now gauge-fix this symmetry non-perturbatively \cite{Filev:2015hia}.
The gauge field $A$ becomes the link variable
\begin{eqnarray}
U_{n,n+1}={\cal P}\exp(-i\int_{na}^{(n+1)a}dt A(t))=1-iaA(n)+O(a^2)~,~a\longrightarrow 0.\label{link}
\end{eqnarray}
The covariant derivative becomes then
\begin{eqnarray}
\frac{U_{n,n+1}\Phi_i(n+1)U_{n+1,n}-\Phi_i(n)}{a}=D_t\Phi_i~,~a\longrightarrow 0.
\end{eqnarray}
The action becomes given by
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^2(n)-\frac{1}{a}U_{n,n+1}\Phi_i(n+1)U_{n+1,n}\Phi_i(n)+\frac{1}{2}m^2a\Phi_i(n)^2\bigg].
\end{eqnarray}
We have a local $U(N)$ symmetry at each lattice site. We use this symmetry to rotate the link variables as follows
\begin{eqnarray}
&&\Phi^{\prime}_i(1)=\Phi_i(1)\nonumber\\
&&\Phi^{\prime}_i(2)=U_{1,2}\Phi_i(2)U_{1,2}^{\dagger}\nonumber\\
&&\Phi^{\prime}_i(3)=U_{1,2}U_{2,3}\Phi_i(3)(U_{1,2}U_{2,3})^{\dagger}\nonumber\\
&&..\nonumber\\
&&..\nonumber\\
&&\Phi^{\prime}_i(\Lambda)=U_{1,2}...U_{\Lambda-1,\Lambda}\Phi_i(\Lambda)(U_{1,2}...U_{\Lambda-1,\Lambda})^{\dagger}.
\end{eqnarray}
The action becomes then
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^{\prime 2}(n)+\frac{1}{2}a^2m^2\Phi_i^{\prime 2}(n)\bigg]-\frac{1}{g^2a}\bigg[\sum_{n=1}^{\Lambda-1}{Tr}\Phi_i^{\prime}(n)\Phi_i^{\prime}(n+1)+\Phi_i^{\prime}(\Lambda+1)W^{\dagger}\Phi_i^{\prime}(\Lambda)W\bigg].\nonumber\\
\end{eqnarray}
The matrix $W$ is defined by
\begin{eqnarray}
W=U_{1,2}...U_{\Lambda-1,\Lambda}U_{\Lambda,\Lambda+1}.
\end{eqnarray}
We can diagonalize the matrix $W$ in the usual way as $W=UDU^{\dagger}$ where $D={\rm diag}(\exp(i\theta_1),...,\exp(i\theta_N))$.
The action becomes (with $\tilde{\Phi}_i(n)=U^{\dagger}\Phi_i^{\prime}(n)U$)
\begin{eqnarray}
S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{Tr}\bigg[\frac{1}{a}\tilde{\Phi}_i^{ 2}(n)+\frac{1}{2}a^2m^2\tilde{\Phi}_i^{ 2}(n)\bigg]-\frac{1}{g^2a}\bigg[\sum_{n=1}^{\Lambda-1}{ Tr}\tilde{\Phi}_i^{}(n)\tilde{\Phi}_i^{}(n+1)+\tilde{\Phi}_i^{}(\Lambda+1)D^{\dagger}\tilde{\Phi}_i^{}(\Lambda)D\bigg].\nonumber\\
\end{eqnarray}
Alternatively, this action can also be rewritten in terms of $D_{\Lambda}=D^{1/\Lambda}$ and $\bar{\Phi}_i(n)=h_n^{\dagger}\tilde{\Phi}_i(n)h_n$, where $h_n=D_{\Lambda}^n$, as
\begin{eqnarray}
S&=&\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\bar{\Phi}_i^{ 2}(n)+\frac{1}{2}a^2m^2\bar{\Phi}_i^{ 2}(n)\bigg]\nonumber\\
&-&\frac{1}{g^2a}\sum_{n=1}^{\Lambda}{Tr}\bar{\Phi}_i^{}(n)D_{\Lambda}\bar{\Phi}_i^{}(n+1)D_{\Lambda}^{\dagger}.
\end{eqnarray}
The measure $\prod_{i,n}dX_i(n)$ is invariant under all the above unitary transformations. Finally, we reduce the measure over the transporter fields as follows
\begin{eqnarray}
\prod_{n=1}^{\Lambda}{\cal D}U_{n,n+1}&=&\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}U_{1,2}\nonumber\\
&\sim &\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}W\nonumber\\
&\sim &\prod_{n=2}^{\Lambda}{\cal D}U_{n,n+1}{\cal D}U{\cal D}D\Delta^2(D)\nonumber\\
&\sim&\prod_{a=1}^Nd\theta_a.\prod_{a>b}|e^{i\theta_a}-e^{i\theta_b}|^2\nonumber\\
&\sim&\prod_{a=1}^Nd\theta_a.\prod_{a>b}\sin^2\frac{\theta_a-\theta_b}{2}\nonumber\\
&\sim &\prod_{a=1}^Nd\theta_a.\exp(-S_{\rm FP}).
\end{eqnarray}
The Faddeev-Popov gauge-fixing action is then given explicitly by
\begin{eqnarray}
S_{\rm FP}=-\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
Observables and the variation of the holonomy angles
The initial gauge-invariant lattice action is
\begin{eqnarray}S=\frac{1}{g^2}\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}\Phi_i^2(n)-\frac{1}{a}U_{n,n+1}\Phi_i(n+1)U_{n+1,n}\Phi_i(n)+\frac{1}{2}m^2a\Phi_i(n)^2\bigg].\label{BFSS_appr}
\end{eqnarray}
The adopted gauge-fixing condition is effectively the so-called static gauge given by
\begin{eqnarray}
U_{n,n+1}=D_{\Lambda}={\rm diag}(\exp(i\frac{\theta_1}{\Lambda}),...,\exp(i\frac{\theta_N}{\Lambda}))~,~\forall n.
\end{eqnarray}
In terms of the gauge field $A$ this condition reads (by means of equation (\ref{link}))
\begin{eqnarray}
A(t)=-\frac{1}{\beta}{\rm diag}(\theta_1,...,\theta_N).
\end{eqnarray}
In summary, we are interested in the total action
\begin{eqnarray}
S_{\rm total}&=&N\sum_{n=1}^{\Lambda}{ Tr}\bigg[\frac{1}{a}{\Phi}_i^{ 2}(n)+\frac{1}{2}am^2{\Phi}_i^{ 2}(n)\bigg]\nonumber\\
&-&\frac{N}{a}\sum_{n=1}^{\Lambda}{Tr}{\Phi}_i^{}(n)D_{\Lambda}{\Phi}_i^{}(n+1)D_{\Lambda}^{\dagger}-\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
Before we discuss the variations of this action we present the most important observables. The Polyakov line is the order parameter of the Hagedorn transition in string theory and the deconfinment transition in gauge theory which is associated with the spontaneous breakdown of the $U(1)$ symmetry $A(t)\longrightarrow A(t)+C.{\bf 1}$ (see the extensive list of references in \cite{Kawahara:2007fn}). The Polyakov line is defined in terms of the holonomy matrix $U$ by the relation
\begin{eqnarray}
P=\frac{1}{N}Tr U~,~U={\cal P}\exp(-i\int_0^{\beta} dt A(t)).
\end{eqnarray}
We compute \cite{Kawahara:2007fn}
\begin{eqnarray}
U&=&U_{1,2}U_{2,3}...U_{\Lambda-1,\Lambda}U_{\Lambda,1}=D_{\Lambda}^{\Lambda}={\rm diag}(\exp(i\theta_1),...,\exp(i\theta_N).
\end{eqnarray}
Hence
\begin{eqnarray}
P=\frac{1}{N}\sum_ae^{i\theta_a}.
\end{eqnarray}
Another important observable is the radius (or extent of space or more precisely the extent of the eigenvalue distribution) is defined by
\begin{eqnarray}
R^2=\frac{a}{\Lambda N^2}\langle {\rm radius}\rangle~,~{\rm radius}=\frac{N}{a}\sum_{n=1}^{\Lambda}{Tr}{\Phi}_i^{2}(n).
\end{eqnarray}
The energy in the present model (\ref{BFSS_appr}) is given effectively by the extent of space. Indeed, we compute (by using the formula (\ref{formu})) for the model (\ref{BFSS_appr}) the energy
\begin{eqnarray}
\frac{E}{N^2}=\frac{a^2Tm^2}{N^{2}}\langle{\rm radius}\rangle=m^2R^2.
\end{eqnarray}
We also measure the eigenvalues distribution of the holonomy matrix $U$ and the eigenvalues distribution of the bosonic matrices $\Phi_i(n)$.
In our present model ($d$ gauged matrix harmonic oscillators) the variation of the bosonic matrices $\Phi_i(n)$ is straightforward (because of the non-zero mass term and the absence of the Yang-Mills term the so-called problem of flat directions does not pose itself here). The variation of the angles $\theta_a$ requires however a careful treatment due to the center of mass $\theta_{\rm cm}=\sum_{a=1}^N\theta_a/N$ which does not appear in the action and behaves as a random walker.
We compute as before the variation under the change (for example for $a\neq b$)
\begin{eqnarray}
(\Phi_i(n))_{ab}\longrightarrow (\Phi_i(n))_{ab}+h
\end{eqnarray}
to find the result
\begin{eqnarray}
\Delta S_{{\rm total}~i}(n)^{ab}=\frac{N}{a}\bigg[4(1+\frac{m^2a^2}{2}){\rm Re}(h^*\phi_i(n))_{ab}-2{\rm Re}(f^*\Phi_i(n+1)+g^*\Phi_i(n-1))_{ab}+(1+\frac{m^2a^2}{2})(2hh^*+(h^2+h^{*2})\delta_{ab})\bigg]
\end{eqnarray}
where $f^*$ and $g^*$ are given by
\begin{eqnarray}
f^*=h^*\exp(\frac{i}{\Lambda}(\theta_a-\theta_b))~,~g^*=h^*\exp(-\frac{i}{\Lambda}(\theta_a-\theta_b))
\end{eqnarray}
But there is also the variation under the change of the angle $\theta_c$ of the holonomy matrix, viz
\begin{eqnarray}
\theta_c\longrightarrow\theta_c^{\prime}=\theta_c+\alpha.\label{theta_var}
\end{eqnarray}
The relevant action here is given by
\begin{eqnarray}
S(\theta)&=&-\frac{N}{a}\sum_{n=1}^{\Lambda}\sum_{a,b}e^{-\frac{i}{\Lambda}(\theta_a-\theta_b)}({\Phi}_i(n))_{ab}({\Phi}_i(n+1))_{ba}\nonumber\\
&-&\frac{1}{2}\sum_{a\ne b}\ln\sin^2\frac{\theta_a-\theta_b}{2}.
\end{eqnarray}
The variation of this action is given explicitly by
\begin{eqnarray}
\Delta S_c(\theta)&=&-\frac{2N}{a}{\rm Re}\bigg(\sum_{n=1}^{\Lambda}\sum_{a\ne c}e^{\frac{i}{\Lambda}(\theta_a-\theta_c)}(e^{-i\frac{\alpha}{\Lambda}}-1)({\Phi}_i(n+1))_{ac}({\Phi}_i(n))_{ca}\bigg)\nonumber\\
&-&\sum_{a\ne c}\ln\sin^2\frac{\theta_a-\theta_c-\alpha}{2}+\sum_{a\ne c}\ln\sin^2\frac{\theta_a-\theta_c}{2}.
\end{eqnarray}
An extremely important remark is now in order. The action $S(\theta)$ does not depend on the center of mass $\theta_{\rm cm}=\sum_{a=1}^N\theta_a/N$. Indeed, the functional integration over $\theta_a$ is (almost) identical to the functional integration over $\tilde{\theta}_a=\theta_a-\theta_{\rm cm}$ which satisfies $\sum_{a}\tilde{\theta}_a=0$.
Furthermore, by means of the $U(1)$ gauge transformation $A(t)\longrightarrow A(t) +C .{\bf 1}$ (a large gauge transformation with a non-zero winding number) we can choose the static gauge $A(t)=-{\rm diag}(\theta_1,...,\theta_N)/\beta$ in such a way that (see for example \cite{Anagnostopoulos:2007fw})
\begin{eqnarray}
{\rm max}(\tilde{\theta}_a)-{\rm min}(\tilde{\theta}_a)\lt 2\pi.
\end{eqnarray}
We have then explicitly
\begin{eqnarray}
\int d\theta F(\Delta\theta)=\int d\tilde{\theta}\delta(\sum_a\tilde{\theta}_a)F(\Delta\tilde{\theta})\int d\theta_{\rm cm}.
\end{eqnarray}
Since $\theta=\tilde{\theta}+\theta_{\rm cm}$ and $-\pi\lt \theta\le +\pi$ while $ {\rm min}(\tilde{\theta}_a)\le\tilde{\theta}\le {\rm max}(\tilde{\theta}_a)$ we conclude immediately that the center of mass $\theta_{\rm cm}$ must be in the interval $]-\pi-{\rm max}(\tilde{\theta}_a),\pi-{\rm min}(\tilde{\theta}_a)]$. Hence the above integral becomes (with $\mu={\rm max}(\tilde{\theta}_a)-{\rm min}(\tilde{\theta}_a)$)
\begin{eqnarray}
\int d\theta F(\Delta\theta)=\int d\tilde{\theta}\delta(\sum_a\tilde{\theta}_a)F(\Delta\tilde{\theta})(2\pi+\mu).
\end{eqnarray}
Clearly, this is true as long as $\mu\le 2\pi$ while for $\mu\gt 2\pi$ the additional Boltzmann weight is identically zero. We have then the extra Boltzmann weight \cite{code_Hanada}
\begin{eqnarray}
&&w(\mu)=2\pi+\mu~,~\mu\le 2\pi\nonumber\\
&&w(\mu)=0~,~\mu\gt 2\pi.
\end{eqnarray}
In other words, we can replace the functional integration over $\theta_a$ with the functional integration over $\tilde{\theta}_a$ with an additional Boltzmann weight $w(\mu)$.
Thus, the variation (\ref{theta_var}) should be thought of as $\tilde{\theta}_c\longrightarrow \theta_c^{\prime}=\tilde{\theta}_c+\alpha$ where $\alpha$ is uniformly distributed in the range $]-\pi-{\rm max}(\tilde{\theta}_a),\pi-{\rm min}(\tilde{\theta}_a)]$.
Phase structure
https://badisydri.blogspot.com/2020/01/phases-of-gauge-theory-in-lower.htmlReferences
%\cite{Shimamune:1981qf}\bibitem{Shimamune:1981qf}
Y.~Shimamune,
``On the Phase Structure of Large $N$ Matrix Models and Gauge Models,''
Phys.\ Lett.\ {\bf 108B}, 407 (1982).
doi:10.1016/0370-2693(82)91223-0
%%CITATION = doi:10.1016/0370-2693(82)91223-0;%%
%47 citations counted in INSPIRE as of 09 Dec 2019
%\cite{Banks:1996vh}
\bibitem{Banks:1996vh}
T.~Banks, W.~Fischler, S.~H.~Shenker and L.~Susskind,
``M theory as a matrix model: A Conjecture,''
Phys.\ Rev.\ D {\bf 55}, 5112 (1997) doi:10.1103/PhysRevD.55.5112
[hep-th/9610043].
%%CITATION = doi:10.1103/PhysRevD.55.5112;%% %2782 citations counted in INSPIRE as of 02 Oct 2019
%\cite{Kawahara:2007fn}
\bibitem{Kawahara:2007fn}
N.~Kawahara, J.~Nishimura and S.~Takeuchi,
``Phase structure of matrix quantum mechanics at finite temperature,''
JHEP {\bf 0710}, 097 (2007)
doi:10.1088/1126-6708/2007/10/097
[arXiv:0706.3517 [hep-th]].c
%%CITATION = doi:10.1088/1126-6708/2007/10/097;%%
%37 citations counted in INSPIRE as of 26 Apr 2017
%\cite{Filev:2015hia}
\bibitem{Filev:2015hia}
V.~G.~Filev and D.~O'Connor,
``The BFSS model on the lattice,''
JHEP {\bf 1605} (2016) 167
doi:10.1007/JHEP05(2016)167
[arXiv:1506.01366 [hep-th]].c
%%CITATION = doi:10.1007/JHEP05(2016)167;%%
%14 citations counted in INSPIRE as of 25 Dec 2016
%\cite{Anagnostopoulos:2007fw}\bibitem{Anagnostopoulos:2007fw}
K.~N.~Anagnostopoulos, M.~Hanada, J.~Nishimura and S.~Takeuchi,
``Monte Carlo studies of supersymmetric matrix quantum mechanics with sixteen supercharges at finite temperature,''
Phys.\ Rev.\ Lett.\ {\bf 100}, 021601 (2008)
doi:10.1103/PhysRevLett.100.021601
[arXiv:0707.4454 [hep-th]].
%%CITATION = doi:10.1103/PhysRevLett.100.021601;%%
%171 citations counted in INSPIRE as of 20 Dec 2019
%\cite{code_Hanada}
\bibitem{code_Hanada}
Masanori Hanada, ``Users’ Manual of Monte Carlo Code for the Matrix Model of M-theory'', https://sites.google.com/site/hanadamasanori/home/mmmm.
%\cite{Kawahara:2007ib}\bibitem{Kawahara:2007ib}
N.~Kawahara, J.~Nishimura and S.~Takeuchi,
``High temperature expansion in supersymmetric matrix quantum mechanics,''
JHEP {\bf 0712}, 103 (2007)
doi:10.1088/1126-6708/2007/12/103
[arXiv:0710.2188 [hep-th]].
%%CITATION = doi:10.1088/1126-6708/2007/12/103;%%
%41 citations counted in INSPIRE as of 28 Dec 2019
ReplyDelete*how we introduce phi-bar- in action (70)?
*between the action (70) and (73) did we use phi or phi-bar- at final ?
and thanks Mr.
phi-bar is only change of notation (a rescaling); after it is done it is irrelevant. so I think it is phi-bar used in variation; Check again though because many new equations have been added since your question.
Delete