Color Confinement and Random Matrices
— A random walk down group manifold toward Casimir scaling —

Georg Bergnera, Vaibhav Gautamb,c, and Masanori Hanadab

a University of Jena, Institute for Theoretical Physics
Max-Wien-Platz 1, D-07743 Jena, Germany
b School of Mathematical Sciences, Queen Mary University of London
Mile End Road, London, E1 4NS, United Kingdom
c Department of Mathematics, University of Surrey
Guildford, Surrey, GU2 7XH, United Kingdom

Abstract

We explain the microscopic origin of linear confinement potential with the
Casimir scaling in generic confining gauge theories. In the low-temperature
regime of confining gauge theories such as QCD, Polyakov lines are slowly
varying Haar random modulo exponentially small corrections with respect to the
inverse temperature, as shown by one of the authors (M. H.) and Watanabe. With
exact Haar randomness, computation of the two-point correlator of Polyakov
loops reduces to the problem of random walk on group manifold. Linear
confinement potential with approximate Casimir scaling except at short
distances follows naturally from slowly varying Haar randomness. With
exponentially small corrections to Haar randomness, string breaking and loss
of Casimir scaling at long distance follow. Hence we obtain the Casimir
scaling which is only approximate and holds only at intermediate distance,
which is precisely needed to explain the results of lattice simulations. For
$(1+1)$-dimensional theories, there is a simplification that admits the
Casimir scaling at short distances as well.

###### Contents

  1. 1 Introduction
  2. 2 Microscopic origin of slowly varying Haar randomness
    1. 2.1 $\mathrm{S}_{N}$ invariance and Bose-Einstein condensation
    2. 2.2 $\mathrm{SU}(N)$ invariance and confinement
    3. 2.3 Polyakov lines are slowly varying Haar random
  3. 3 Casimir scaling from random walk
    1. 3.1 Vanilla random walk
    2. 3.2 Random walk with gradually changing velocity
    3. 3.3 Simplification for $(1+1)$-dimensional theories
  4. 4 String breaking and loss of Casimir scaling
  5. 5 Discussion

##  1 Introduction

In the low-temperature regime of confining gauge theory, the Polyakov line
must be slowly varying Haar random up to small corrections and gauge
transformations [1, 2, 3].111 In this paper, we consider bare Polyakov loop. A
multiplicative renormalization factor (see e.g., Ref. [4]) does not affect our
argument on string tension.  We show how linear confinement potential with
Casimir scaling [5, 6, 7, 8, 9]222 In this paper, we consider the Casimir
scaling of string tension in the confined phase. The Casimir scaling was also
studied for the one-point function of the Polyakov loop in the deconfined
phase [4, 10, 11].  follows from this fact. The property of a Haar random
distribution is explained below in more detail, repeating some of the
arguments in Ref. [1]. It is not based on a strong coupling limit in lattice
gauge theory but rather a property that emerges for the Polyakov line at a
certain range of energy scales.

For concreteness, we consider SU($N$) gauge theory, although the same argument
applies to any gauge group and any matter content. We quantize this theory via
Euclidean path integral at temperature $T=\beta^{-1}$. At each spatial point
$\vec{x}$, we can define Polyakov line $P_{\vec{x}}$ as

$\displaystyle P_{\vec{x}}=\textrm{Path\
ordering}\left[e^{i\int_{0}^{\beta}dtA_{t}(t,\vec{x})}\right]\,.$ (1)

We choose a reference point $\vec{x}_{0}$ and define the Wilson line
connecting these two points,

$\displaystyle U(\vec{x}_{0},{\vec{x}})=\textrm{Path\
ordering}\left[e^{i\int_{C}d\vec{x}\vec{A}}\right]\,.$ (2)

where $C$ is a contour connecting $\vec{x}_{0}$ and $\vec{x}$ at fixed
Euclidean time $t=0$. Then, $P^{\prime}_{\vec{x}}$ defined by

$\displaystyle P^{\prime}_{\vec{x}}\equiv U(\vec{x}_{0},\vec{x})\cdot
P_{\vec{x}}\cdot U(\vec{x}_{0},\vec{x})^{-1}$ (3)

(see Fig. 1) must be close to slowly-varying Haar random. Namely,
statistically, $P^{\prime}_{\vec{x}}$ is Haar random at each $\vec{x}$, and it
is slowly varying as a function of $\vec{x}$ in each typical field
configuration that has a large weight in the path integral. In Sec. 2, we
review the previous work and explain why the Polyakov line has such a
property. In Sec. 2.3, we estimate the rate of change (i.e., the meaning of
‘slow’) based on a lattice regularization. Note that slowly-varying Haar
randomness does not apply to the abelianized regime where $\mathrm{SU}(N)$
breaks to $\mathrm{U}(1)^{N-1}$ because the correction to slowly-varying Haar
randomness is too large. Ref. [12] confirmed that the Casimir scaling does not
hold in such a theory.

Figure 1: Illustration of $P^{\prime}_{\vec{x}}\equiv
U(\vec{x}_{0},\vec{x})\cdot P_{\vec{x}}\cdot(U(\vec{x}_{0},\vec{x}))^{-1}$
defined in (3) in terms of the path on compact space.

Along a straight line $\vec{x}=\vec{x}_{0}+L\hat{u}$ with some unit vector
$\hat{u}$ and length $L>0$, the change of $P_{\vec{x}}$ is a smooth analog of
random walk on the group manifold.333An earlier application of a random walk
on group manifold to confinement can be found in Refs. [13, 14, 15, 16]. These
references introduced a random walk based on phenomenological observation from
lattice simulations. We provide a more complete picture based on the
theoretical analysis of gauge theory.  For each typical field configuration in
path integral, we can take smooth functions $v_{\alpha}(L)$
($\alpha=1,2,\cdots,N^{2}-1$) and write the change of Polyakov line as

$\displaystyle
P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}=\textrm{Path\
ordering}\left[e^{i\sum_{\alpha}\int_{0}^{L}dL^{\prime}v_{\alpha}(L^{\prime})T_{\alpha}}\right]\,.$
(4)

We can do it for any representation r:

$\displaystyle
R^{\rm(r)}\left(P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}\right)=\textrm{Path\
ordering}\left[e^{i\sum_{\alpha}\int_{0}^{L}dL^{\prime}v_{\alpha}(L^{\prime})T^{\rm(r)}_{\alpha}}\right]\,.$
(5)

Character $\chi_{\rm r}$ associated with irreducible representation r is
defined by taking the trace of the representation matrix. Specifically,

$\displaystyle\chi_{\rm r}(P_{\vec{x}})=\mathrm{Tr}_{\rm
r}\left(R^{\rm(r)}(P_{\vec{x}})\right)$ (6)

is the Polyakov loop in representation r at point $\vec{x}$. Because of

$\displaystyle\chi_{\rm r}(P_{\vec{x}})=\chi_{\rm r}(P^{\prime}_{\vec{x}})$
(7)

and

$\displaystyle\chi_{\rm r^{\prime}}(P_{\vec{x}+L\hat{u}})=\chi_{\rm
r^{\prime}}(P^{\prime}_{\vec{x}+L\hat{u}})=\chi_{\rm
r^{\prime}}\left(P_{\vec{x}}^{\prime}\cdot
P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}\right)\,,$ (8)

the two-point function of Polyakov loops can be written as

$\displaystyle\left\langle\left(\chi_{\rm
r}(P_{\vec{x}})\right)^{\ast}\cdot\chi_{\rm
r^{\prime}}(P_{\vec{x}+L\hat{u}})\right\rangle=\left\langle\left(\chi_{\rm
r}(P^{\prime}_{\vec{x}})\right)^{\ast}\cdot\chi_{\rm
r^{\prime}}\left(P_{\vec{x}}^{\prime}\cdot
P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}\right)\right\rangle$ (9)

So far, $\langle\ \rangle$ meant the path integral. From here on, we interpret
$\langle\ \rangle$ as the average over slowly-varying Haar random
configurations $\\{P^{\prime}_{\vec{x}}\\}$, neglecting corrections to slowly-
varying Haar randomness. (As we will see, this approximation fails if $L$ is
too large.)

We consider random walk from $P^{\prime}_{\vec{x}}$ to
$P^{\prime}_{\vec{x}+L\hat{u}}$. The starting point $P^{\prime}_{\vec{x}}$ is
Haar random, and it is not correlated to the displacement
$P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}$. Therefore, we can treat
$P^{\prime}_{\vec{x}}$ and
$P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}$ as independent variables
and integrate out $P^{\prime}_{\vec{x}}$ by using the orthogonality condition
of representations under the Haar-random average,

$\displaystyle\frac{1}{{\rm Vol}(\mathrm{SU}(N))}\int dP(R^{({\rm
r})}_{ij}(P))^{\ast}R^{({\rm r^{\prime}})}_{kl}(P)=d_{\rm r}^{-1}\delta_{\rm
rr^{\prime}}\delta_{il}\delta_{kj}\,$ (10)

to obtain

$\displaystyle\langle(\chi_{\rm r}(P_{\vec{x}}))^{\ast}\chi_{\rm
r^{\prime}}(P_{\vec{x}+L\hat{u}})\rangle$ $\displaystyle=d_{\rm
r}^{-1}\delta_{rr^{\prime}}\left\langle\chi_{\rm
r}(P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}})\right\rangle\,.$ (11)

Combined with (4) and (5), the computation of the two-point function reduces
to a sort of random walk on group manifold. The difference from the vanilla
random walk problem is that velocity changes slowly.

In the rest of the paper, we will show that confinement with Casimir scaling
follows naturally from such a random walk on a group manifold. Our argument is
based only on the slowly varying Haar randomness, and the choice of gauge
group or matter content is irrelevant. However, corrections to Haar randomness
depend on these specifications of the theory, and they are relevant for long-
distance physics where the Casimir scaling is lost.

Our convention is as follows. We normalize generators of the Lie algebra
$T_{\alpha}$ ($\alpha=1,2,\cdots,N^{2}-1$) in the fundamental representation
as

$\displaystyle\mathrm{Tr}(T_{\alpha}T_{\beta})=2\delta_{\alpha\beta}\,.$ (12)

For SU(2) and SU(3), we can use Pauli matrices and Gell-Mann matrices,
respectively. For irreducible representation r, we define the quadratic
Casimir invariant $C_{\rm r}$ by

$\displaystyle\sum_{\alpha}\left(T^{\rm(r)}_{\alpha}\right)^{2}=4C_{\rm
r}\textbf{1}\,,$ (13)

where 1 is the identity. For spin-$j$ representation of SU(2), $C_{{\rm
spin}\mathchar 45\relax j}=j(j+1)$. For SU(3), $C_{\rm fund}=\frac{4}{3}$,
$C_{\rm adj}=3$, $C_{\rm 2\mathchar 45\relax sym}=\frac{10}{3}$, $C_{\rm
3\mathchar 45\relax sym}=6$ for fundamental, adjoint, rank-2 symmetric, and
rank-3 symmetric representations. Casimir scaling means that string tension
associated with representation r is proportional to $C_{\rm r}$. Specifically,
the Polyakov loop correlator is

$\displaystyle\left\langle\chi_{\rm r}(P_{\vec{x}})\cdot\chi_{\rm
r}(P_{\vec{x}+L\hat{u}})\right\rangle\sim e^{-\beta C_{\rm r}\sigma_{\rm
string}^{(0)}L}\,,$ (14)

with $\sigma_{\rm string}^{(0)}$ independent of $r$ in the regime where
Casimir scaling holds. This scaling does not hold at very long distances where
strings can break due to pair creation of particles such as quarks or gluons.

This paper is organized as follows. Sec. 2 explains how slowly varying Haar
randomness of the Polyakov loop emerges from a generic mechanism of
confinement. Sec. 3, the main part of this paper, derives Casimir scaling from
a random walk. The simplest, ‘vanilla’ random walk is studied in Sec. 3.1. It
is solvable, and Casimir scaling can be analytically derived. A more
sophisticated variant of the random walk takes into account the slowly-varying
condition and is studied in Sec. 3.2. We find a neat correspondence with a
lattice regularization and confirm Casimir scaling and consistency with the
microscopic mechanism of confinement based on the similarity to Bose-Einstein
condensation [1, 2, 3]. Note, however, that we only model the random walk and
details that can affect the short-distance behavior might depend on the
specific model. We expect that different confining theories can show different
short-distance behaviors because of the dependence on such details. Larger
distance behavior, on the other hand, emerges quite universally from the
different models of the random walk. In Sec. 4, it is discussed that
corrections to exact Haar randomness lead to a breaking of Casimir scaling.
Sec. 5 is devoted to discussion.

##  2 Microscopic origin of slowly varying Haar randomness

In this section, we review the origin of the slowly varying Haar randomness of
the Polyakov line [1]. As a first step, we rewrite the partition function in
terms of an extended Hilbert space ($\mathcal{H}_{\text{ext}}$), including
also non-gauge invariant states, and a projection to the gauge-invariant part
($\mathcal{H}_{\text{inv}}$) by integration over all gauge variations (gauge
orbit). For a gauge theory with gauge group $G$, the canonical partition
function can be written as

$\displaystyle
Z(T)=\mathrm{Tr}_{\mathcal{H}_{\text{inv}}}e^{-\hat{H}/T}=\frac{1}{\text{vol}G}\int_{G}dg\mathrm{Tr}_{\mathcal{H}_{\text{ext}}}\left(\hat{g}e^{-\hat{H}/T}\right)\;.$
(15)

where the integral is taken over the Haar measure. A simple but crucial fact
is that $\hat{g}$ in the expression above corresponds to the Polyakov line
[3].

###  2.1 $\mathrm{S}_{N}$ invariance and Bose-Einstein condensation

To understand the origin of slowly varying Haar randomness, it is instructive
to review the connection between the phenomenon of Bose-Einstein Condensation
(BEC) and confinement in gauge theories [1, 3]. In this case the gauge group
is represented by the permutations $\mathrm{S}_{N}$ of $N$ particles, the
group integral becomes a sum over all possible group elements.

Consider a system of $N$ indistinguishable bosons with a harmonic potential.
The Hamiltonian is given by

$\displaystyle\hat{H}=\frac{1}{2}\sum_{i=1}^{N}\left(\hat{p}^{2}_{xi}+\hat{p}^{2}_{yi}+\hat{p}^{2}_{zi}+\hat{x}^{2}_{i}+\hat{y}^{2}_{i}+\hat{z}^{2}_{i}\right)\,.$
(16)

This system has a gauged $\mathrm{S}_{N}$ permutation group symmetry. Namely,
states connected by $\mathrm{S}_{N}$ transformation (exchange of labels
$1,2,\cdots,N$) must be identified. A generic state spanning the extended
Hilbert space $\mathcal{H}_{\text{ext}}$ is given in the Fock basis as

$\displaystyle|\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}\rangle=\prod_{i=x,y,z}\frac{a^{{\dagger}n_{1i}}_{1i}}{\sqrt{n_{1i}!}}\frac{a^{{\dagger}n_{2i}}_{2i}}{\sqrt{n_{2i}!}}\ldots\frac{a^{{\dagger}n_{Ni}}_{Ni}}{\sqrt{n_{Ni}!}}|\textbf{0}\rangle$
(17)

where $\vec{n}_{i}=(n_{ix},n_{iy},n_{iz})$ denotes the occupation numbers for
each boson. The Fock vacuum is described as
$|\textbf{0}\rangle=|0,0,\ldots,0\rangle$.

For the system of $N$ bosons under consideration here with a gauged $S_{N}$
symmetry, the partition function takes the form,

$\displaystyle Z(T)$ $\displaystyle=$
$\displaystyle\frac{1}{N!}\sum_{\sigma\in
S_{N}}\mathrm{Tr}_{\mathcal{H}_{\text{ext}}}\left(\hat{\sigma}e^{-\hat{H}/T}\right)$
(18) $\displaystyle=$ $\displaystyle\frac{1}{N!}\sum_{\sigma\in
S_{N}}\sum_{\vec{n}_{N},\vec{n}_{N},\ldots,\vec{n}_{N}}\langle\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}|\hat{\sigma}e^{-\hat{H}/T}|\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}\rangle$
$\displaystyle=$
$\displaystyle\frac{1}{N!}\sum_{\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}}e^{-(E_{\vec{n}_{1}}+\ldots+E_{\vec{n}_{N}})/T}\sum_{\sigma\in
S_{N}}\langle\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}|\vec{n}_{\sigma(1)},\vec{n}_{\sigma(2)},\ldots,\vec{n}_{\sigma(N)}\rangle$
$\displaystyle=$
$\displaystyle\frac{1}{N!}\sum_{\vec{n_{1}},\vec{n_{2}},\ldots,\vec{n_{N}}}e^{-(E_{\vec{n}_{1}}+\ldots+E_{\vec{n}_{N}})/T}\sum_{\sigma\in
S_{N}}\delta_{\vec{n}_{1}\vec{n}_{\sigma(1)}}\delta_{\vec{n}_{2}\vec{n}_{\sigma(2)}}\ldots\delta_{\vec{n}_{N}\vec{n}_{\sigma(N)}}$
$\displaystyle=$
$\displaystyle\frac{1}{N!}\sum_{\vec{n_{1}},\vec{n_{2}},\ldots,\vec{n_{N}}}e^{-(E_{\vec{n}_{1}}+\ldots+E_{\vec{n}_{N}})/T}V_{|\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}\rangle}$

Here,
$V_{|\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}\rangle}\equiv\sum_{\sigma\in
S_{N}}\delta_{\vec{n}_{1}\vec{n}_{\sigma(1)}}\delta_{\vec{n}_{2}\vec{n}_{\sigma(2)}}\ldots\delta_{\vec{n}_{N}\vec{n}_{\sigma(N)}}$
is the volume of the stabilizer subgroup of the state
$|\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}\rangle$ i.e. all the elements of
$S_{N}$ that leave the state invariant. If all the particles have different
occupation numbers ($\vec{n}_{i}\neq\vec{n}_{j}$ for all $i\neq j$), then the
stabilizer of such a state is just the identity element $\sigma=\mathbf{1}$.
The contribution of such a state to the partition function thus remains
suppressed by a factor of $1/N!$. However, if all the particles have the same
occupation number ($\vec{n}_{1}=\ldots=\vec{n}_{N}=\vec{n}$), then the
stabilizer of such a state is the entire permutation group $S_{N}$ with
$V_{|\vec{n},\ldots,\vec{n}\rangle}=N!$. The suppression factor in the
partition function thus experiences an enhancement and becomes equal to 1 and
such a state is then favoured. This is exactly what happens in the ground
state $|0\rangle$ ($\vec{n}_{1}=\ldots=\vec{n}_{N}=\vec{0}$). More generally,
if the energy is distributed to a part of particles and the rest remains in
the one-particle ground state (say, $\vec{n}_{1},\cdots,\vec{n}_{M}$ becomes
nonzero while $n_{M+1}=\cdots\vec{N}=\vec{0}$), a large enhancement factor
$(N-M)!$ is obtained. In general, the optimal value of $M$ between $0$ and $N$
is realized so that the combination of such an enhancement factor and the
entropy of the excited sector is maximized. This is the BEC of $N-M$
particles. This mechanism works with finite interactions as well, and unbroken
$\mathrm{S}_{N-M}$ symmetry in the extended Hilbert space characterizes BEC.
All ‘Polyakov line’ $\sigma\in\mathrm{S}_{N}$ contribute equally, which means
the Haar-randomness.

###  2.2 $\mathrm{SU}(N)$ invariance and confinement

Coming back from BEC to our original question about gauge theories including
$\mathrm{SU}(N)$ pure Yang-Mills theory, a similar enhancement can be
observed. In general, an $\mathrm{SU}(M)$ subgroup can be deconfined, while
the rest of the degrees of freedom remain confined [17, 18, 19, 20] as
illustrated in Fig. 2. Roughly speaking, such partial deconfinement is
characterized by an unbroken $\mathrm{SU}(N-M)$ symmetry. Since this is
considered in the extended Hilbert space, there is no conflict with gauge
invariance; simply, all embeddings of SU($M$) into SU($N$) are equivalent. The
completely confined phase we are discussing in this paper is characterized by
unbroken $\mathrm{SU}(N)$ symmetry in the extended Hilbert space. Let us
elaborate on this point for the case of Yang-Mills theory and see how the
slowly varying Haar randomness of the Polyakov line emerges. The differences
from the example of BEC are that the gauge group is bigger (SU($N$) instead of
SN) and that we need to deal with spatial coordinates because we deal with QFT
rather than quantum mechanics.

Figure 2:  Partial deconfinement in SU($N$) QCD with $N_{f}$ fundamental
quarks. [Left] Gauge field. $M\times M$ sub-matrix is deconfined. [Right]
Quarks. $M$ components are deconfined. Any embeddings of SU($M$) into SU($N$)
are equivalent because of gauge symmetry. This figure is taken from Ref. [20].

For concreteness, we consider pure Yang-Mills theory and use a lattice
regularization with lattice spacing $a$ and unitary link variables
$U_{\vec{n},\mu}$ which is related to the Hermitian gauge field $A_{\mu}$ as
$U_{\vec{n},\mu}=e^{iag_{\rm YM}A_{\vec{n},\mu}}$. Here, $\vec{n}$ is an
integer-valued vector that specifies a spatial site on the lattice, and $\mu$
is a spatial direction $x,y,\ldots$. In the extended-Hilbert space picture,
wave functions are defined on
$\prod_{\mu,\vec{n}}[\mathrm{SU}(N)]_{\mu,\vec{n}}$, where
$[\mathrm{SU}(N)]_{\mu,\vec{n}}$ corresponds to each link. Specifically, on
each link, we consider a linear combination of coordinate eigenstates
$|U_{\vec{n},\mu}\rangle$. The confined vacuum in or sufficiently close to the
continuum limit is a wave packet localized around $U_{\vec{n},\mu}=\mathbf{1}$
up to gauge transformation. (As we will see shortly, a large enhancement
factor is associated with such a state, which justifies its dominance in the
partition function.) Wave packets connected by a (spacial) gauge
transformation $\Omega\in SU(N)$ are equivalent. Therefore, we should consider
wave packets localized around
$U_{\vec{n},\mu}=\Omega_{\vec{n}}^{-1}\Omega_{\vec{n}+\hat{\mu}}$.

Corresponding to
$\langle\vec{n}_{1},\vec{n}_{2},\ldots,\vec{n}_{N}|\vec{n}_{\sigma(1)},\vec{n}_{\sigma(2)},\ldots,\vec{n}_{\sigma(N)}\rangle$,
we have a factor $\langle\Phi|\hat{g}|\Phi\rangle$ in the partition function,
where we can take $|\Phi\rangle$ to be a wave packet. By definition, this
factor measures overlap between $|\Phi\rangle$ and $\hat{g}|\Phi\rangle$; see
Fig. 3. Because the Polyakov lines correspond to $\hat{g}$, such Polyakov
lines that do not move the wave packet too much can contribute to the
enhancement factor. Let us consider the action of the Polyakov loop
corresponding to ‘global’ $SU(N)$ transformation,
$P_{\vec{n}}\equiv\Omega_{\vec{n}}^{-1}V\Omega_{\vec{n}}$. It is easy to see
that the wave packet under consideration does not move regardless of the
choice of $V$, because

$\displaystyle
P_{\vec{n}}^{-1}(\Omega_{\vec{n}}^{-1}\Omega_{\vec{n}+\hat{\mu}})P_{\vec{n}+\hat{\mu}}=\Omega_{\vec{n}}^{-1}V^{-1}\Omega_{\vec{n}}(\Omega_{\vec{n}}^{-1}\Omega_{\vec{n}+\hat{\mu}})\Omega_{\vec{n}+\hat{\mu}}^{-1}V\Omega_{\vec{n}+\hat{\mu}}=\Omega_{\vec{n}}^{-1}\Omega_{\vec{n}+\hat{\mu}}\,.$
(19)

Thus, the vacuum is invariant under this ‘global’ $SU(N)$ transformation,
leading to an enhancement factor as in the BEC case. The enhancement factor
scales with $N$ as $e^{N^{2}}$. For partially-deconfined states shown in Fig.
2,

$\displaystyle V=\left(\begin{array}[]{cc}\textbf{1}_{M}&\textbf{0}\\\
\textbf{0}&\tilde{V}\end{array}\right)$ (22)

with any $\tilde{V}\in\textrm{SU}(N-M)$ could leave the states invariant. This
leads to an enhancement factor $\sim e^{(N-M)^{2}}$.

Let us elaborate on the argument above and show that the enhancement factor
can increase exponentially with volume. (Roughly speaking, we just promote a
‘global’ transformation to ‘slowly varying’ transformation.) Let us start with
the case of the confined vacuum ($M=0$). We consider a ‘local’ $SU(N)$
transformation with $V$ now being dependent on the site $\vec{n}$ with the
Polyakov line taking the form
$P_{\vec{n}}\equiv\Omega_{\vec{n}}^{-1}V_{\vec{n}}\Omega_{\vec{n}}$:

$\displaystyle
P_{\vec{n}}^{-1}U_{\vec{n},\mu}P_{\vec{n}+\hat{\mu}}=\Omega_{\vec{n}}^{-1}V_{\vec{n}}^{-1}\Omega_{\vec{n}}(\Omega_{\vec{n}}^{-1}\Omega_{\vec{n}+\hat{\mu}})\Omega_{\vec{n}+\hat{\mu}}^{-1}V_{\vec{n}+\hat{\mu}}\Omega_{\vec{n}+\hat{\mu}}=\Omega_{\vec{n}}^{-1}V_{\vec{n}}^{-1}V_{\vec{n}+\hat{\mu}}\Omega_{\vec{n}+\hat{\mu}}$
(23)

If $V_{\vec{n}}$ is slowly varying, then
$V_{\vec{n}}^{-1}V_{\vec{n}+\hat{\mu}}$ is close to $\mathbf{1}$, and we can
still have a significant enhancement factor as before. For the confined
vacuum, we can choose $V_{\vec{n}}$ arbitrarily as long as it is slowly
varying, indeed leading to a large enhancement factor that increases with
spatial volume. Intuitively, we can split the volume into multiple pieces that
are large enough so that the correlation between Polyakov lines in different
pieces is negligible, and then, we can associate $\sim e^{N^{2}}$ to each
piece, which leads to the scaling of the enhancement factor
$e^{N^{2}\times\mathrm{volume}}$.

Figure 3:  For a wave packet $|\Phi\rangle$, such $g$ that does not move
$|\Phi\rangle$ too much can lead to an enhancement factor and give a large
contribution to the partition function.

###  2.3 Polyakov lines are slowly varying Haar random

So far, we have seen that a large enhancement factor is indeed associated with
the confined vacuum. We have also seen that Polyakov lines contributing to the
enhancement factor are written as
$P_{\vec{n}}\equiv\Omega_{\vec{n}}^{-1}V_{\vec{n}}\Omega_{\vec{n}}$, where
$V_{\vec{n}}$ can be arbitrary as long as it is slowly varying. Typical
Polyakov lines dominating the partition function (equivalently, typical
Polyakov lines obtained from dominant configurations in the Euclidean path
integral) are the ones that contribute to the enhancement factor. Therefore,
typical Polyakov lines are slowly varying SU($N$) Haar random up to gauge
transformation [1]. For the same logic, the Polyakov line restricted to the
confined sector of the partially-deconfined phase should be slowly varying
SU($N-M$) Haar random up to gauge transformation.

‘Slowly varying’ is admittedly a vague expression. Below, we give an order
estimate for the rate of change. This estimate is used for a consistency check
in Sec. 3.2. We want to know how slowly a typical Polyakov line changes,
assuming the confined vacuum represented by a wave packet localized around a
pure-gauge configuration. Because a typical Polyakov line must lead to a large
enhancement factor, $V_{\vec{n}}^{-1}V_{\vec{n}+\hat{\mu}}$ has to be
sufficiently close to the identity so that the ground-state wave function does
not move too much (Fig. 3). Specifically, if
$V_{\vec{n}}^{-1}V_{\vec{n}+\hat{\mu}}$ is written as
$V_{\vec{n}}^{-1}V_{\vec{n}+\hat{\mu}}\sim e^{i\ell\tilde{X}}$ where
$\tilde{X}$’s entries are of order $1$, then $\ell$ must be smaller than the
radius of the wave packet. Link variable $U_{\vec{n},\mu}$ and continuum gauge
field $A_{\vec{n},\mu}$ are related by $U_{\vec{n},\mu}=e^{iag_{\rm
YM}A_{\vec{n},\mu}}$. Coupling $g_{\rm YM}$ and $A_{\vec{n},\mu}$ are of order
$a^{0}$ (there can be logarithmic dependence on $a$, but our argument here is
not precise enough to distinguish logarithm and constant), and hence, the size
of the wave packet is of order $a^{1}$ on the group manifold. Therefore,
$l\lesssim a$ is required for a large enhancement factor in the partition
function. A natural value of $v$ in (4) is of order one. Note that this is
just a crude estimate and should not be taken too precisely including order
one multiplicative factor or logarithmic corrections.

Let us also comment on the ‘acceleration’, i.e., the rate of change of $v$.
Except for $(1+1)$ dimensions, $v$ in (4) has to be smooth because
$P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{y}}$ depends only on the endpoints
$\vec{x}$ and $\vec{y}$ and does not change depending on the path connecting
these two points. For this reason, the Polyakov lines exhibit random walks
with slowly changing velocity, which is different from the most common random
walk in which the velocity at each step changes abruptly. In
$(1+1)$-dimensional theories, $v$ does not have to be continuous and the
Polyakov line can literally random walk. This has an important consequence as
we will discuss in Sec. 3.3.

##  3 Casimir scaling from random walk

In this section, we assume the exact Haar randomness and derive the linear
confinement potential with Casimir scaling. Hence, our task is to estimate
$\left\langle\chi_{\rm
r}(P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}})\right\rangle$ assuming
the random walk (4). To evaluate (4), let us introduce step size $a$, where
$Ka=L$, and take a limit of $a\to 0$. (We use $a$ because it is analogous to
the lattice spacing in lattice regularized theories.) We define $W_{j}$
($j=1,2,\cdots,K$) by

$\displaystyle
W_{j}\equiv(P^{\prime}_{\vec{x}+(j-1)a\hat{u}})^{-1}P^{\prime}_{\vec{x}+ja\hat{u}}\,.$
(24)

These can be combined to

$\displaystyle
P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}+L\hat{u}}=W_{1}W_{2}\cdots W_{K}\,,$
(25)

as illustrated in Fig. 4.

Figure 4: Contour corresponding to
$P_{\vec{x}}^{\prime-1}P^{\prime}_{\vec{x}^{\prime}}$ to illustrate the
definition of $W_{j}$ in (24).

Below, let us first consider the vanilla random walk as a solvable warm-up
exercise, and then consider a more realistic setup that takes into account the
correlation between $W_{i}$’s.

###  3.1 Vanilla random walk

As a crude approximation, let us treat $W_{1},\cdots,W_{K}$ totally random,
even neglecting the correlation of $W_{i}$ and $W_{i+1}$. This is the vanilla
random walk. Apparently, we cannot take the limit of $a\to 0$ in this case. To
keep this in mind, let us use $s$ instead of $a$. Physically, $s$ is a length
scale needed for $v$ in (4) to change significantly and the correlation
between $v(L)$ and $v(L+s)$ becomes small. With this approximation, the
average of a product becomes the product of the average:

$\displaystyle\left\langle W^{\rm(r)}_{1}W^{\rm(r)}_{2}\cdots
W^{\rm(r)}_{k}\right\rangle=\left\langle
W^{\rm(r)}_{1}\right\rangle\left\langle
W^{\rm(r)}_{2}\right\rangle\cdots\left\langle
W^{\rm(r)}_{k}\right\rangle=\left(\left\langle
W^{\rm(r)}\right\rangle\right)^{k}\,.$ (26)

This expression already tells us an exponential decay. To see the Casimir
scaling, we need a little bit more math. We assume $W=e^{i\Delta X}$, where
$X$ is random Gaussian with variance 1 and $\Delta$ is a small number. Then,

$\displaystyle\left\langle W^{\rm(r)}\right\rangle$
$\displaystyle=\left\langle\textbf{1}+i\Delta
x^{\alpha}T^{\rm(r)}_{\alpha}-\frac{\Delta^{2}x^{\alpha}x^{\beta}}{2}T^{\rm(r)}_{\alpha}T^{\rm(r)}_{\beta}+\cdots\right\rangle$
$\displaystyle=\textbf{1}-\frac{\Delta^{2}}{2}(T^{\rm(r)}_{\alpha})^{2}+\cdots$
$\displaystyle=\textbf{1}-2\Delta^{2}C_{\rm r}\textbf{1}+\cdots\,.$ (27)

Hence, for small $\Delta$,

$\displaystyle\left\langle W^{\rm(r)}\right\rangle\simeq e^{-2\Delta^{2}C_{\rm
r}}\textbf{1}$ (28)

and hence

$\displaystyle\left\langle W^{\rm(r)}_{1}\right\rangle\left\langle
W^{\rm(r)}_{2}\right\rangle\cdots\left\langle
W^{\rm(r)}_{K}\right\rangle\simeq e^{-2\Delta^{2}C_{\rm r}K}\textbf{1}\,.$
(29)

Therefore,

$\displaystyle\langle(\chi_{\rm r}(P_{\vec{x}}))^{\ast}\chi_{\rm
r^{\prime}}(P_{\vec{x}+L\hat{u}})\rangle\simeq\delta_{rr^{\prime}}e^{-2\Delta^{2}C_{\rm
r}K}$ (30)

up to $O(\Delta^{4})$ terms in the exponent. Therefore, by identifying
$\Delta$ with physical parameters as

$\displaystyle\Delta^{2}=\beta s\sigma_{\rm string}^{(0)}\,,$ (31)

string tension for representation r is $C_{\rm r}\sigma_{\rm string}^{(0)}$ if
$O(\Delta^{4})$ terms are negligible. The dependence on $s$ is natural because
after $K$ steps the distance traveled on the group manifold is $\sim
N\Delta\sqrt{K}\sim N\sqrt{L}$, which is a typical case of random walk.

That terms of order $\Delta^{4}$ or higher break Casimir scaling would be good
news because lattice simulations suggest Casimir scaling is not exact.

Several references [13, 21, 22] assumed the uncorrelated nature of color flux
and explained the Casimir scaling in a way similar to the vanilla random walk
in the limit of $\Delta\to 0$. Our discussion above is based on the
microscopic mechanism applicable to generic confining gauge theories, and
hence, we can give a more precise picture that takes into account a finite
amount of correlations, as we will show next. It would be interesting to see
if the hypothesis regarding the QCD vacuum proposed in the past can be
justified based on our approach.

###  3.2 Random walk with gradually changing velocity

In reality, we want the ‘velocity’ $W_{i}$ to change gradually. To incorporate
this feature, let us consider a model of random walk with $W_{i}$ defined by

$\displaystyle W_{i}=Z_{i}\cdots Z_{i+p-1}\,,$ (32)

where $Z_{j}=e^{i\epsilon X_{j}}$ and $X_{j}$ is Gaussian random with a unit
variance, and take the average of the product of $W$’s as before:

$\displaystyle\langle W_{1}W_{2}\cdots W_{K}\rangle\,.$ (33)

Then, $W_{i}$ changes gradually, and the correlation disappears completely
after $p$ steps. This corresponds to a physical length scale $ap$, and hence,
we will consider a continuum limit $a\to 0$ fixing $ap$. We are interested in
the behavior of the correlator as a function of physical distance $Ka$, or
equivalently, $\frac{K}{p}$.

As explained in Sec. 2.3, we expect $W_{i}\sim e^{ia\tilde{X}_{i}}$, where
$\tilde{X}_{i}$’s entries are of order $1$. Note that this is a crude order
estimate that does not even take possible logarithmic corrections into
account. If we further assume $\epsilon$ is very small, we obtain
$a\tilde{X}_{i}\sim\epsilon\sum_{j=0}^{p-1}X_{i+j}$. Each matrix entry on the
right-hand side is a sum of $p$ random numbers of order $\epsilon$, which is
of order $\epsilon\sqrt{p}$. Therefore, $a\sim\epsilon\sqrt{p}$ would be
expected. We will consider a continuum limit $a\to 0$, $\epsilon\to 0$ fixing
$ap\sim\epsilon p^{3/2}$.

Unfortunately, we were not able to solve this problem analytically. We
performed numerical experiments for SU(2) and SU(3). The results are very
clean, and we believe our findings apply to other groups as well.

Numerically, we can confirm the exponential decay with the exponent
proportional to $\frac{K}{p}C_{\rm r}$, by taking $\epsilon\to 0$ and
$p\to\infty$ fixing $\epsilon p^{3/2}$; see Fig. 5. Furthermore, the
dependence on $\epsilon p^{3/2}$ is small if $\epsilon p^{3/2}$ is not too
large. Indeed, as shown in Fig. 6, by taking the horizontal axis to be
$(K/p)\cdot(\epsilon p^{3/2})^{2}=K(\epsilon p)^{2}$, we can see that
exponential decays with different $p$ and $\epsilon$ line up. Therefore, by
identifying $(\epsilon p)^{2}$ with $\beta a\sigma_{\rm string}^{(0)}$, we can
reproduce $e^{-\beta L\sigma_{\rm string}^{(0)}C_{\rm r}}$. This $\epsilon p$
is a counterpart of $\Delta$ in the vanilla random walk, but now the finite
amount of correlation is taken into account.

As we have discussed, $a\sim\epsilon\sqrt{p}$ is a natural identification.
Then, $(\epsilon p)^{2}\sim\beta a\sigma_{\rm string}^{(0)}$ implies
$\beta\sigma_{\rm string}^{(0)}\sim\epsilon p^{3/2}$. If we took this estimate
very precisely including order one multiplicative factor, we could have a
problem that $\epsilon p^{3/2}$ cannot be taken very large without spoiling
the Casimir scaling. However, because this is just a crude order estimate, we
do not think this is a serious problem; a small correction to the estimate
could allow the counterpart of $\beta$ to be large. Therefore, we take it to
be good news that such a crude estimate gives $a\sim\epsilon p^{1/2}$, which
can be arbitrarily small, and $\beta\sigma_{\rm string}^{(0)}\sim\epsilon
p^{3/2}$, which can be of order one. More precise estimates can be achieved
via lattice gauge theory simulations.

We emphasize that we only provided a model of the random walk of Polyakov
lines satisfying essential features. We expect that different confining
theories can show different short-distance behaviors because of differences in
the details. Specifically, $X_{j}$ may not be exactly Gaussian, or $W_{i}$ may
change in a way different from (32).

Figure 5:  $\textrm{Tr}\langle W_{1}...W_{K}\rangle$ vs $(K/p)\times C_{\rm
r}$ for $p^{3/2}\epsilon=0.01\times 10^{3/2}\simeq 0.31622$, $p=10,20,40$.
[Top] SU(2), spin $j=\frac{1}{2},1,\frac{3}{2},2,\frac{5}{2}$. [Bottom] SU(3),
fundamental, rank-2 symmetric, adjoint, rank-3 symmetric. We normalized
$\textrm{Tr}\langle W_{1}...W_{K}\rangle$ by dividing $d_{\rm
r}=\textrm{Tr}\textbf{1}$. The number of samples for the averaging is 1
million for both SU(2) and SU(3).

Figure 6:  $\textrm{Tr}\langle W_{1}...W_{K}\rangle$ vs $(K/p)\cdot(\epsilon
p^{3/2})^{2}$ for several values of $p$ and $\epsilon p^{3/2}$. [Top] SU(2),
spin $j=\frac{1}{2}$ and $1$. [Bottom] SU(3), fundamental and rank-2
symmetric. We normalized $\textrm{Tr}\langle W_{1}...W_{K}\rangle$ by dividing
$d_{\rm r}=\textrm{Tr}\textbf{1}$. The number of samples for the averaging is
1 million for SU(2) and SU(3), $p=10,20,40$, and 500K for SU(2), $p=200$.

###  3.3 Simplification for $(1+1)$-dimensional theories

In Sec. 3.2, we have seen that $v$ in (4) gradually changes is responsible for
the absence of the Casimir scaling at short distance. In fact, the situation
is different for $(1+1)$-dimensional theories. As commented in Sec. 2.3,
$(1+1)$-dimensional theories are special in that $v$ in (4) does not have to
be continuous. On a lattice, the vanilla random walk is justified for any
lattice spacing, as long as the corrections to Haar randomness are neglected.
Therefore, quite generally, we expect the Casimir scaling also at short
distances.

It is known that the Casimir scaling is exact for $(1+1)$-dimensional pure
Yang-Mills theory [23, 24, 25]. With matter fields, infrared behaviors may
change drastically. However, we do not expect the short-distance behavior to
change, as long as the system remains confined. This provides us with a
nontrivial consistency check of the relationship between the Casimir scaling
and a random walk on the group manifold.

##  4 String breaking and loss of Casimir scaling

String breaking is almost trivial when the corrections to Haar-random
distribution are taken into account. If there is no particular global symmetry
that forces the Polyakov loop with representation r to become zero, we expect
to have a nonzero expectation value of the one-point function of the form

$\displaystyle\langle\chi_{\rm r}(P_{\vec{x}})\rangle=\langle\chi_{\rm
r}(P_{\vec{x}+L\hat{u}})\rangle\sim e^{-\beta m_{\rm r}}\,,$ (34)

where $m_{\rm r}$ is the mass of the lowest excitation in representation r.
Therefore, we need to take the disconnected part of the two-point function
into account:

$\displaystyle\langle(\chi_{\rm r}(P_{\vec{x}}))^{\ast}\chi_{\rm
r}(P_{\vec{x}+L\hat{u}})\rangle\sim e^{-\beta C_{\rm r}\sigma_{\rm
string}^{(0)}L}+e^{-2\beta m_{\rm r}}\,,$ (35)

where the second term is the disconnected part which is dominant at long
distances ($C_{\rm r}L\sigma_{\rm string}^{(0)}\gg m_{\rm r}$).

There must be corrections to the connected part as well. For QCD, the details
of the corrections to the connected part may not be important because the
disconnected part is nonzero for any representations and dominates long-
distance physics. The situation is somewhat different for theories with center
symmetry such as SU($N$) pure Yang-Mills. As long as the center symmetry is
not broken, the disconnected part vanishes. Furthermore, loops with different
center-symmetry charges cannot mix with each other, because the connected part
must vanish as well. In the random-walk approach, (11) receives some
corrections because $P^{\prime}_{\vec{x}}$ is not exactly Haar random, and we
need to consider a random walk on group manifold with some other weight than
the Haar measure. We do not find an immediate reason that forbids the mixing
in the same center-symmetry-charge sector, and hence, we expect that string
tension takes the same value in the same center-symmetry-charge sector
regardless of the representation. Precise dependence on charge requires a
better understanding and is out of the scope of this paper. See Refs. [14, 15,
16] for past attempts to build a random-walk model in the long-distance
regime.

##  5 Discussion

In the low-temperature regime of confining gauge theories such as QCD,
Polyakov lines are slowly varying Haar random modulo exponentially small
corrections with respect to the inverse temperature [1, 2, 3]. In this paper,
we suggested that linear confinement potential with Casimir scaling follows
naturally from such a random walk. Strictly speaking, we introduced a model of
random walk using a random matrix product. We expect the exponential decay
with approximate Casimir scaling to be universal and does not depend much on
the details of the random process, while short-distance behaviors depend on
details. With exponentially small corrections to Haar randomness, string
breaking and loss of Casimir scaling at long distance follow. So, what we have
seen is an approximate Casimir scaling in the intermediate distance scale,
which is precisely what we expect in confining gauge theories. Our picture can
be tested by calculating the Polyakov loops via lattice gauge theory
simulation.

The mechanism described for the Polyakov loop may generalize to the case of a
rectangular Wilson loop of the size $L_{t}$ for the time direction and $L$ for
the spatial direction. Specifically, we could take a product of Wilson lines
like in Fig. 7. In this case, $\epsilon p^{3/2}\sim L_{t}$ would be expected.
For this scenario to work, we need to take $L_{t}$ sufficiently large so that
the temporal edge inherits the Haar randomness of the Polyakov loop. Note also
that, if $L_{t}$ is too small, the conformal behavior can set in and slowly-
varying nature can be lost.444We thank Nadav Drukker for pointing this out for
us. It would be interesting to use lattice simulation to confirm this
mechanism.

Figure 7:  A Wilson line which can be used to build a rectangular Wilson loop
of the size $L_{t}$ for the time direction and $L$ for the spatial direction.

Note that we essentially assumed the mass gap to claim that the correction to
Haar randomness is exponentially small. Therefore, we are assuming one of the
important features of confinement, and hence, we do not claim we proved
confinement. Still, what we observed is a highly nontrivial consequence of the
mass gap, and we can find nontrivial connections between good old results and
deepen our understanding of confinement.

More than two decades ago, it was realized that large-$N$ theory in the weak-
coupling limit exhibits confinement/deconfinement transition that resembles
some aspects of the counterpart at strong coupling such as the $N$-dependence
of entropy [26, 27]. Most crucially, deconfinement can be understood as the
condensation of long strings both at weak coupling and at strong coupling.
However, it was not clear how linear confinement potential can follow from
this picture. In the meantime, a series of work on partial deconfinement [17,
18, 19, 20] led to the discovery of a microscopic mechanism of confinement in
the large-$N$ gauge theories [1, 2, 3] which is essentially the same as Bose-
Einstein condensation. (See Sec. 2.) This mechanism is based only on symmetry
and is applicable to both weak and strong coupling. Let us call it Confinement
$=$ BEC mechanism.555You may prefer $\simeq$ or $\sim$ depending on what you
mean by BEC, confinement, or $=$. The slowly varying Haar-randomness is a
consequence of the Confinement $=$ BEC mechanism. Therefore, we understood how
the Confinement $=$ BEC mechanism, which is based on a very basic property of
gauge theory, leads to linear confinement potential with Casimir scaling.

Confinement $=$ BEC mechanism leads to various phenomena including partial
deconfinement and intriguing pattern of global symmetry breaking [28, 29],
coexistence of parton-like and string-like degrees of freedom [30],666 In Ref.
[30], the partially-deconfined phase of the large-$N$ limit of strongly-
coupled lattice gauge theory was studied, and linear confinement potential in
the confined sector was observed. We expect the same pattern in the continuum
theory, by applying the same reasoning used in this paper to the confined
sector. Specifically, the Polyakov line restricted to the confined sector in
the partially-deconfined phase should be slowly varying Haar random.
instanton condensation [1, 2], and so on. Hopefully, this mechanism leads us
to a unified picture of an even richer list of nontrivial phenomena.

Acknowledgement

The authors would like to thank Pavel Buividovich, Nadav Drukker, Jeff
Greensite, Jack Holden, Enrico Rinaldi, Andreas Schäfer, Hidehiko Shimada,
Yuya Tanizaki, Masaki Tezuka, and Hiromasa Watanabe. V. G. thanks STFC for the
Doctoral Training Programme funding (ST/W507854-2021 Maths DTP). M. H. thanks
his STFC consolidated grant ST/X000656/1. G. B. is funded by the Deutsche
Forschungsgemeinschaft (DFG) under Grant No. 432299911 and 431842497.

## References

  * [1] M. Hanada and H. Watanabe, “On thermal transition in QCD,” arXiv:2310.07533 [hep-th].
  * [2] M. Hanada, H. Ohata, H. Shimada, and H. Watanabe, “A new perspective on thermal transition in QCD,” arXiv:2310.01940 [hep-th].
  * [3] M. Hanada, H. Shimada, and N. Wintergerst, “Color confinement and Bose-Einstein condensation,” JHEP 08 (2021) 039, arXiv:2001.10459 [hep-th].
  * [4] S. Gupta, K. Huebner, and O. Kaczmarek, “Renormalized Polyakov loops in many representations,” Phys. Rev. D 77 (2008) 034503, arXiv:0711.2251 [hep-lat].
  * [5] J. Ambjorn, P. Olesen, and C. Peterson, “Stochastic Confinement and Dimensional Reduction. 1. Four-Dimensional SU(2) Lattice Gauge Theory,” Nucl. Phys. B 240 (1984) 189–212.
  * [6] J. Ambjorn, P. Olesen, and C. Peterson, “Stochastic Confinement and Dimensional Reduction. 2. Three-dimensional SU(2) Lattice Gauge Theory,” Nucl. Phys. B 240 (1984) 533–542.
  * [7] L. Del Debbio, M. Faber, J. Greensite, and S. Olejnik, “Casimir scaling versus Abelian dominance in QCD string formation,” Phys. Rev. D 53 (1996) 5891–5897, arXiv:hep-lat/9510028.
  * [8] G. S. Bali, “Casimir scaling of SU(3) static potentials,” Phys. Rev. D 62 (2000) 114503, arXiv:hep-lat/0006022.
  * [9] S. Deldar, “Static SU(3) potentials for sources in various representations,” Phys. Rev. D 62 (2000) 034509, arXiv:hep-lat/9911008.
  * [10] S. Gupta, K. Huebner, and O. Kaczmarek, “Polyakov loop in different representations of SU(3) at finite temperature,” Nucl. Phys. A 785 (2007) 278–281, arXiv:hep-lat/0608014.
  * [11] A. Mykkanen, M. Panero, and K. Rummukainen, “Casimir scaling and renormalization of Polyakov loops in large-N gauge theories,” JHEP 05 (2012) 069, arXiv:1202.2762 [hep-lat].
  * [12] E. Poppitz and M. E. Shalchian T., “String tensions in deformed Yang-Mills theory,” JHEP 01 (2018) 029, arXiv:1708.08821 [hep-th].
  * [13] A. M. Brzoska, F. Lenz, J. W. Negele, and M. Thies, “Diffusion of Wilson loops,” Phys. Rev. D 71 (2005) 034008, arXiv:hep-th/0412003.
  * [14] G. Arcioni, S. de Haro, and P. Gao, “A Diffusion model for SU(N) QCD screening,” Phys. Rev. D 73 (2006) 074508, arXiv:hep-th/0511213.
  * [15] P. V. Buividovich and V. I. Kuvshinov, “Kramers-Moyall cumulant expansion for the probability distribution of parallel transporters in quantum gauge fields,” Phys. Rev. D 73 (2006) 094015, arXiv:hep-th/0605207.
  * [16] P. V. Buividovich and M. I. Polikarpov, “Random walks of Wilson loops in the screening regime,” Nucl. Phys. B 790 (2008) 28–41, arXiv:0704.3367 [hep-ph].
  * [17] M. Hanada and J. Maltz, “A proposal of the gauge theory description of the small Schwarzschild black hole in AdS${}_{5}\times$S5,” JHEP 02 (2017) 012, arXiv:1608.03276 [hep-th].
  * [18] D. Berenstein, “Submatrix deconfinement and small black holes in AdS,” JHEP 09 (2018) 054, arXiv:1806.05729 [hep-th].
  * [19] M. Hanada, G. Ishiki, and H. Watanabe, “Partial Deconfinement,” JHEP 03 (2019) 145, arXiv:1812.05494 [hep-th]. [Erratum: JHEP 10, 029 (2019)].
  * [20] M. Hanada, A. Jevicki, C. Peng, and N. Wintergerst, “Anatomy of Deconfinement,” JHEP 12 (2019) 167, arXiv:1909.09118 [hep-th].
  * [21] A. Di Giacomo, H. G. Dosch, V. I. Shevchenko, and Y. A. Simonov, “Field correlators in QCD: Theory and applications,” Phys. Rept. 372 (2002) 319–368, arXiv:hep-ph/0007223.
  * [22] J. Greensite, K. Langfeld, S. Olejnik, H. Reinhardt, and T. Tok, “Color Screening, Casimir Scaling, and Domain Structure in G(2) and SU(N) Gauge Theories,” Phys. Rev. D 75 (2007) 034501, arXiv:hep-lat/0609050.
  * [23] A. A. Migdal, “Recursion Equations in Gauge Theories,” Sov. Phys. JETP 42 (1975) 413.
  * [24] E. Witten, “On quantum gauge theories in two-dimensions,” Commun. Math. Phys. 141 (1991) 153–209.
  * [25] M. Nguyen, Y. Tanizaki, and M. Ünsal, “Noninvertible 1-form symmetry and Casimir scaling in 2D Yang-Mills theory,” Phys. Rev. D 104 no. 6, (2021) 065003, arXiv:2104.01824 [hep-th].
  * [26] B. Sundborg, “The Hagedorn transition, deconfinement and N=4 SYM theory,” Nucl. Phys. B 573 (2000) 349–363, arXiv:hep-th/9908001.
  * [27] O. Aharony, J. Marsano, S. Minwalla, K. Papadodimas, and M. Van Raamsdonk, “The Hagedorn - deconfinement phase transition in weakly coupled large N gauge theories,” Adv. Theor. Math. Phys. 8 (2004) 603–696, arXiv:hep-th/0310285.
  * [28] M. Hanada and B. Robinson, “Partial-Symmetry-Breaking Phase Transitions,” Phys. Rev. D 102 no. 9, (2020) 096013, arXiv:1911.06223 [hep-th].
  * [29] M. Hanada, J. Holden, M. Knaggs, and A. O’Bannon, “Global symmetries and partial confinement,” JHEP 03 (2022) 118, arXiv:2112.11398 [hep-th].
  * [30] V. Gautam, M. Hanada, J. Holden, and E. Rinaldi, “Linear confinement in the partially-deconfined phase,” JHEP 03 (2023) 195, arXiv:2208.14402 [hep-th].