# Explicit formulas for matrices associated to ladder, circular ladder, and
Möbius ladder graphs

Ali Azimi  Department of Mathematics and Applied Mathematics, Xiamen
University Malaysia, 43900, Sepang, Selangor Darul Ehsan, Malaysia
<EMAIL_ADDRESS><EMAIL_ADDRESS>and  M. Farrokhi D. G  Research
Center for Basic Sciences and Modern Technologies (RBST), Institute for
Advanced Studies in Basic Sciences (IASBS), Zanjan 45137-66731, Iran
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

We give explicit formulas for resistance distance matrices and Moore-Penrose
inverses of incidence and Laplacian matrices of ladder, circular ladder, and
Möbius ladder graphs. As a result, we compute the Kirchhoff index of these
graphs and give new combinatorial formulas for the number of their spanning
trees.

###### Key words and phrases:

Moore-Penrose inverse, resistance distance, Kirchhoff index, ladder graph,
circular ladder graph, Möbius ladder graph

###### 2000 Mathematics Subject Classification:

Primary 05C20, 05C50, Secondary 15A09.

##  1\. introduction

The incidence matrix $Q(\Gamma)$ of an oriented graph $\Gamma$ with vertex-set
$V(\Gamma)=\\{v_{1},\ldots,v_{n}\\}$ and edge-set
$E(\Gamma)=\\{e_{1},\ldots,e_{m}\\}$ is an $n\times m$ matrix defined as
follows: The rows and columns of $Q$ are indexed by vertices and edges of
$\Gamma$, respectively. The $(i,j)$ entry of $Q$ is $0$ if the vertex $v_{i}$
and the edge $e_{j}$ are not incident, otherwise it is $1$ or $-1$ depending
on whether $e_{j}$ starts from or ends at $v_{i}$. The Laplacian matrix
$L(\Gamma)$ of $\Gamma$, that is actually equal to
$Q(\Gamma)Q^{\prime}(\Gamma)$, is an $n\times n$ matrix whose rows and columns
are indexed by vertices of $\Gamma$. The $(i,j)$ entry of $L(\Gamma)$ is equal
to the degree of the vertex $v_{i}$ if $i=j$, and it is $-1$ or $0$ if the
vertices $v_{i}$ and $v_{j}$ are adjacent or non-adjacent, respectively.

The Moore-Penrose inverse of an $m\times n$ real matrix $A$, denoted by
$A^{+}$, is the (unique) $n\times m$ real matrix that satisfies the following
equations

$AA^{+}A=A,\hskip 11.38092ptA^{+}AA^{+}=A^{+},\hskip
11.38092pt(AA^{+})^{\prime}=AA^{+},\hskip 11.38092pt(A^{+}A)^{\prime}=A^{+}A$

(see [9, 10]). The Moore-Penrose inverse of the incidence matrix of a graph
was first studied by Ijira in [17]. We refer the interested reader to [1, 2,
3] for previous studies on the Moore-Penrose inverses of matrices associated
to various classes of graphs.

Incidence matrices and their inverses, such as the Moore-Penrose inverses, are
important tools used in analysis of graphs. Let $\Gamma$ be a simple connected
graph. The resistance distance between two vertices $v_{i}$ and $v_{j}$ can be
computed via the formula

$r(v_{i},v_{j})=L^{+}_{v_{i},v_{i}}+L^{+}_{v_{j},v_{j}}-L^{+}_{v_{i},v_{j}}-L^{+}_{v_{j},v_{i}},$

where $L^{+}_{v_{i},v_{j}}$ denotes the $(i,j)$-entry of the Moore-penrose
inverse $L^{+}$ of $L=L(\Gamma)$. The notion of resistance distance between
vertices of a graph is introduced by Klein and Randic in 1993 [19] as the
effective resistance between the given vertices (computed by Ohm’s law) when a
battery is attached between them assuming that all edges are unit resistors.
Beside many applications in physics, chemistry, and network sciences (see [13,
18, 22]), resistance distances can also be studied from a purely mathematical
point of view as they represent new metrics on graph as well as making
connections to graph theoretical objects like the Moore-Penrose inverse of
Laplacian matrices above, see [5, 7, 8, 12, 14, 15, 16]. A strictly related
parameter is the so-called Kirchhoff index of the graph, namely

$Kf(\Gamma)=\frac{1}{2}\sum_{i=1}^{n}\sum_{j=1}^{n}r(v_{i},v_{j}).$

In this paper, we consider three familiar classes of graphs, namely ladder,
circular ladder, and Möbius graphs, and study their resistance distances and
the Moore-Penrose inverse of their Laplacian matrices. Recall that the ladder
graph $L_{n}$ on $2n$ vertices is composed of two paths of length $n-1$ with
spokes/rungs connecting corresponding vertices. The circular ladder graph
$CL_{n}$ defines analogous to ladder graphs with paths replaced by cycles.
Also, the Möbius ladder graph $M_{n}$ obtains from the circular ladder graph
by switching two non-spoke edges of a square.

##  2\. Preliminaries

In this section, we introduce some recursive sequences arising from the study
of ladder, circular ladder, and Möbius ladder graphs.

######  Definition.

By a $(p,q)$-sequence $\\{x_{n}\\}$ we mean a sequence satisfying
$x_{n}=px_{n-1}-qx_{n-2}$ for all $n$ provided that $x_{n-1}$ and $x_{n-2}$
are defined.

  * •

$\\{a_{n}\\}_{n=0}^{\infty}$ is the $(4,1)$-sequence A001835 with initial
values $(a_{0},a_{1})=(1,1)$.

  * •

$\\{s_{n}\\}_{n=0}^{\infty}$ is the $(4,1)$-sequence A001353 with initial
values $(s_{0},s_{1})=(0,1)$.

Lemma 2.1 presents a list of identities satisfied by these sequences and will
be used frequently in the rest of the paper.

######  Lemma 2.1.

The following identities hold among sequences $\\{a_{n}\\}$ and $\\{s_{n}\\}$:

(1) $\displaystyle s_{n}$ $\displaystyle=\sum_{i=1}^{n}a_{i},$ (2)
$\displaystyle 2s_{n}$ $\displaystyle=a_{n+1}-a_{n},$ (3) $\displaystyle
a_{n+1}-a_{n}$ $\displaystyle=a_{k+1}a_{n-k+1}-a_{k}a_{n-k},$ (4)
$\displaystyle s_{n}$ $\displaystyle=a_{k}s_{n-k}+a_{n-k+1}s_{k},$ (5)
$\displaystyle(a_{n+1}+a_{n})^{2}$ $\displaystyle=3(a_{n+1}-a_{n})^{2}+4,$ (6)
$\displaystyle\sum_{i+j=n+1}^{n}a_{i}a_{j}$
$\displaystyle=\frac{n(a_{n+1}+a_{n})+(a_{n+1}-a_{n})}{6}=\frac{(n+1)s_{n}+na_{n}}{3}$

for all $n\geq k\geq 0$.

###### Proof.

The equations (1) and (2) follow by induction on $n$ and the equation (3)
follows by induction on $k$. Let $\lambda:=2+\sqrt{3}$,
$\alpha_{a}:=(3-\sqrt{3})/6$, $\beta_{a}:=(3+\sqrt{3})/6$, and
$\alpha_{s}=-\beta_{s}=\sqrt{3}/6$. Then $a_{n}$ and $s_{n}$ are given by the
Binet formulas

$a_{n}=\alpha_{a}\lambda^{n}+\beta_{b}\lambda^{-n}\quad\text{and}\quad
s_{n}=\alpha_{s}\lambda^{n}+\beta_{s}\lambda^{-n}.$

The equations (4) and (5) can be obtained simply from the Binets formulas.
Now, we prove equation (6). Let $b_{n}:=\sum_{i=1}^{n}a_{i}a_{n+1-i}$ for all
$n\geq 1$. First observe that

$\displaystyle b_{n}=\sum_{i=1}^{n}a_{i}a_{n+1-i}$
$\displaystyle=\sum_{i=1}^{n}(4a_{i}a_{n-i}-a_{i}a_{n-1-i}),$

which implies that

(7) $b_{n}-4b_{n-1}+b_{n-2}=a_{n}-a_{n-1}$

On the other hand, by equation (3),

$\displaystyle b_{n}-b_{n-2}$
$\displaystyle=a_{1}a_{n}+\left((a_{2}a_{n-1}-a_{1}a_{n-2})+\cdots+(a_{n-1}a_{2}-a_{n-2}a_{1})\right)+a_{n}a_{1}$
$\displaystyle=2a_{n}+(n-2)(a_{n}-a_{n-1})$

so that

(8) $b_{n}-b_{n-2}=na_{n}-(n-2)a_{n-1}.$

Summing up the equations (8) when $n$ is replaced by $3,\ldots,n$, it follows
that

(9) $b_{n}+b_{n-1}=na_{n}+s_{n-1}.$

Now, by solving the system of linear equations (7), (8), and (9) for $b_{n}$,
$b_{n-1}$, and $b_{n-2}$, the result follows.

∎

##  3\. ladder graph

Let $L_{n}$ denote the ladder graph of order $2n$. Suppose vertices along the
top and bottom of $L_{n}$, as depicted in Fig. 1, are labeled
$u_{1}^{+},\ldots,u_{n}^{+}$ and $u_{1}^{-},\ldots,u_{n}^{-}$, respectively.
Also, suppose the edges of $L_{n}$ are oriented as
$e_{i}^{\varepsilon}=(u_{i}^{\varepsilon},u_{i+1}^{\varepsilon})$ for
$i=1,\ldots,n-1$ with $\varepsilon=\pm 1$, and $f_{i}=(u_{i}^{+},u_{i}^{-})$
for $i=1,\ldots,n$, where $f_{i}$ are the spokes.

$u_{1}^{-}$$u_{2}^{-}$$u_{3}^{-}$$u_{n-2}^{-}$$u_{n-1}^{-}$$u_{n}^{-}$$u_{1}^{+}$$u_{2}^{+}$$u_{3}^{+}$$u_{n-2}^{+}$$u_{n-1}^{+}$$u_{n}^{+}$
Figure 1. Ladder graph $L_{n}$

Let $S(\Gamma)$ denote the set of all spanning trees of the underlying graph
of $\Gamma$ and $s(\Gamma):=\\#S(\Gamma)$. The graph obtained from $\Gamma$ by
contracting edge $e$ is denoted by $\Gamma/e$. It is known from [20, p.35]
that $s(L_{n})=4s(L_{n-1})-s(L_{n-2})$ for all $n\geq 2$. Since, $s(L_{0})=0$
and $s(L_{1})=1$, we observe that $s(L_{n})=s_{n}$ for all $n\geq 0$.

Let $Q=[q_{u,e}]$ be the incidence matrix of $L_{n}$ with the orientation
given in Fig. 1 whose rows and columns are indexed by ordered sets

$V(L_{n})=\\{u_{1}^{+},\ldots,u_{n}^{+},u_{1}^{-},\ldots,u_{n}^{-}\\}$

and

$E(L_{n})=\\{f_{1},\ldots,f_{n},e_{1}^{+}\ldots,e_{n-1}^{+},e_{1}^{-}\ldots,e_{n-1}^{-}\\},$

respectively. Also, let

(10) $H=\begin{pmatrix}B&-B\\\ C&D\\\ D&C\end{pmatrix},$

be the $(3n-2)\times 2n$ matrix defined as follows: The rows and columns of
$H$ are indexed by ordered sets $E(L_{n})$ and $V(L_{n})$, respectively. Let
$B=[b_{i,j}]$ be the $n\times n$ symmetric matrix given by
$b_{i,j}=a_{i}a_{n-j+1}/2s_{n}$, and $C=[c_{i,j}]$ and $D=[d_{i,j}]$ be the
$(n-1)\times n$ matrices given by

$c_{i,j}=\begin{cases}\delta_{i,j}-\frac{i}{2n}-\frac{a_{n-j+1}s_{i}}{2s_{n}},&i\leq
j,\\\ \frac{1}{2}-\frac{i}{2n}+\frac{a_{j}s_{n-i}}{2s_{n}},&i>j,\end{cases}$

and

$d_{i,j}=\begin{cases}-\frac{i}{2n}+\frac{a_{n-j+1}s_{i}}{2s_{n}},&i\leq j,\\\
\frac{1}{2}-\frac{i}{2n}-\frac{a_{j}s_{n-i}}{2s_{n}},&i>j,\end{cases}$

respectively. In the following series of lemmas, we provide the machinery to
prove that $H$ is the Moore-Penrose inverse of $Q$.

######  Lemma 3.1.

$H\mathbf{1}=0$.

###### Proof.

It is obvious from the definition that $(C+D)\mathbf{1}=0$. Thus

$H\mathbf{1}=\begin{bmatrix}(B-B)\mathbf{1}\\\ (C+D)\mathbf{1}\\\
(D+C)\mathbf{1}\end{bmatrix}=0,$

as required ∎

######  Lemma 3.2.

$QH=I-\frac{1}{2n}J$.

###### Proof.

Let $\delta^{\prime}_{i,j}:=1-\delta_{i,j}$. For $1\leq i\leq j\leq n$, we
have

$\displaystyle(QH)_{u_{i}^{+},u_{j}^{+}}$
$\displaystyle=b_{i,j}-\delta^{\prime}_{i,1}c_{i-1,j}+\delta^{\prime}_{i,n}c_{i,j}$
$\displaystyle=\frac{a_{i}a_{n-j+1}}{2s_{n}}-\delta^{\prime}_{i,1}\left(-\frac{i-1}{2n}-\frac{a_{n-j+1}s_{i-1}}{2s_{n}}\right)+\delta^{\prime}_{i,n}\left(\delta_{i,j}-\frac{i}{2n}-\frac{a_{n-j+1}s_{k}}{2s_{n}}\right)$
$\displaystyle=\delta_{i,j}-\frac{1}{2n},$

by using equation (1) when $i=n$. Also, an analogous argument yields

$(QH)_{u_{j}^{+},u_{i}^{+}}=(QH)_{u_{i}^{-},u_{j}^{-}}=(QH)_{u_{j}^{-},u_{i}^{-}}=\delta_{i,j}-\frac{1}{2n}.$

Similarly, by discussing the cases where $i=1$, $i=n$, and $i\neq 1,n$ and
either $i\leq j$, $i=j+1$, and $i>j+1$, we obtain

$(QH)_{u_{i}^{+},u_{j}^{-}}=(QH)_{u_{i}^{-},u_{j}^{+}}=-b_{i,j}-\delta^{\prime}_{1,i}d_{i-1,j}+\delta^{\prime}_{i,n}d_{i,j}=-\frac{1}{2n}$

for all $1\leq i,j\leq n$. The proof is complete. ∎

######  Lemma 3.3.

$HQ$ is a symmetric matrix.

###### Proof.

Let $1\leq i,j\leq n$. Since $B$ is a symmetric matrix,

$(HQ)_{f_{i},f_{j}}=2b_{i,j}=2b_{j,i}=(HQ)_{f_{j},f_{i}}.$

From the definition of $H$, we get

$(HQ)_{e_{i}^{+},f_{j}}=-(HQ)_{e_{i}^{-},f_{j}}=c_{i,j}-d_{i,j}=\begin{cases}\delta_{i,j}-\frac{a_{n-j+1}s_{i}}{s_{n}},&i\leq
j,\\\ \frac{a_{j}s_{n-i}}{s_{n}},&i>j,\end{cases}$

and

$(HQ)_{f_{j},e_{i}^{+}}=-(HQ)_{f_{j},e_{i}^{-}}=b_{j,i}-b_{j,i+1}=\begin{cases}\frac{a_{i}a_{n-j+1}-a_{i+1}a_{n-j+1}}{2s_{n}},&i<j,\\\
\frac{a_{j}a_{n-i+1}-a_{j}a_{n-i}}{2s_{n}},&i\geq j.\end{cases}$

Utilizing equations (1), (2), and (4) and a simple case-by-case analysis, it
yields $(HQ)_{e_{i}^{+},f_{j}}=(HQ)_{f_{j},e_{i}}$ and
$(HQ)_{e_{i}^{-},f_{j}}=(HQ)_{f_{j},e_{i}^{-}}$. On the other hand,

$(HQ)_{e_{i}^{+},e_{j}^{+}}=c_{i,j}-c_{i,j+1}=\begin{cases}\delta_{i,j}+\frac{(a_{n-j}-a_{n-j+1})s_{i}}{2s_{n}},&i\leq
j,\\\ \frac{(a_{j}-a_{j+1})s_{n-i}}{2s_{n}},&i>j+1,\\\
-\frac{1}{2}+\frac{a_{j}s_{n-i}+a_{n-j}s_{i}}{2s_{n}},&i=j+1,\end{cases}$

$(HQ)_{e_{i}^{+},e_{j}^{-}}=d_{i,j}-d_{i,j+1}=\begin{cases}\frac{(a_{n-j+1}-a_{n-j})s_{i}}{2s_{n}},&i\leq
j,\\\ \frac{1}{2}-\frac{a_{j}s_{n-i}+a_{n-j}s_{i}}{2s_{n}},&i=j+1,\\\
\frac{(a_{j+1}-a_{j})s_{n-i}}{2s_{n}},&i>j+1,\end{cases}$

$(HQ)_{e_{j}^{-},e_{i}^{+}}=c_{j,i}-c_{j,i+1}=\begin{cases}\frac{(a_{i}-a_{i+1})s_{n-j}}{2s_{n}},&j>i+1,\\\
-\frac{1}{2}+\frac{a_{i}s_{n-j}+a_{n-i}s_{j}}{2s_{n}},&j=i+1,\\\
\frac{(a_{n-i}-a_{n-i+1})s_{j}}{2s_{n}},&i\geq j\end{cases}$

and

$(HQ)_{e_{i}^{-},e_{j}^{-}}=d_{i,j}-d_{i,j+1}=\begin{cases}\frac{(a_{n-j+1}-a_{n-j})s_{i}}{2s_{n}},&i\leq
j,\\\ \frac{1}{2}-\frac{a_{j}s_{n-i}+a_{n-j}s_{i}}{2s_{n}},&i=j+1,\\\
\frac{(a_{j+1}-a_{j})s_{n-i}}{2s_{n}},&i>j+1,\end{cases}$

from which the result follows. ∎

######  Theorem 3.4.

Let $L_{n}$ be a ladder graph with incident matrix $Q$. Then $Q^{+}=H$.

###### Proof.

By Lemmas 3.2 and 3.3, $QH$ and $HQ$ are symmetric. On the other hand, $QHQ=Q$
for $\mathbf{1}^{\prime}Q=0$. Since $HQH=H$ by Lemmas 3.1 and 3.2, it follows
that $H=Q^{+}$. ∎

The following lemma will be used frequently in the reminder of this section.

######  Lemma 3.5 ([4, Corollary 3.2]).

Let $\Gamma$ be a connected graph and $P$ be a path of length $d$ between
vertices $u,v\in V(\Gamma)$. If all edges in $P$ have the same direction from
$u$ to $v$ then,

$r(u,v)=\sum_{e\in E(P)}(q^{+}_{e,u}-q^{+}_{e,v}),$

where $r(u,v)$ is the resistance distance between vertices $u$ and $v$.

######  Corollary 3.6.

If $\mathcal{B}=2s(L_{n})B$. Then $\mathcal{B}$ has eigenvalue $s(L_{n})$ with
associated eigenvector $\mathbf{1}$. Also, the diagonal entries of
$\mathcal{B}$ satisfy the equations

$\mathcal{B}_{i,i}=s(L_{n}/f_{i})$

for $i=1,\ldots,n$.

###### Proof.

By [3, Lemma 2.1]

$\sum_{i=1}^{n}q^{+}_{f_{i},u_{j}^{+}}=1-\frac{n}{2n}=\frac{1}{2},$

On the other hand, by Theorem 3.4,

$\sum_{i=1}^{n}q^{+}_{f_{i},u_{j}^{+}}=\frac{a_{i}\sum_{i\leq
j}a_{n-j+1}+a_{n-i+1}\sum_{i>j}a_{j}}{2s(L_{n})}$

so that

$s(L_{n})=a_{i}\sum_{i\leq
j}a_{n-j+1}+a_{n-i+1}\sum_{i>j}a_{j}=\sum_{j=1}^{n}b_{i,j},$

that is $\mathcal{B}\mathbf{1}=s(L_{n})\mathbf{1}$. Now, from Lemma 3.5 and
Theorem 3.4, and the definition of the resistance distance between two
vertices, we get

$\frac{s(L_{n}/f_{i})}{s(L_{n})}=r(u_{i}^{+},u_{i}^{-})=2q^{+}_{f_{i},u_{i}^{+}}=\frac{a_{i}a_{n-i+1}}{s(L_{n})}$

(see [6, p. 133]). Therefore

(11) $s(L_{n}/f_{i})=a_{i}a_{n-i+1}=\mathcal{B}_{i,i},$

as required. ∎

In [11] the author uses circus reduction to obtain the resistance distance and
Kirchhoff index of ladder graphs. Here we use the Moore-Penrose inverses of
Laplacian matrices of ladder graphs to drive new formulas for resistance
distance and Kirchhoff index of ladder graphs.

######  Theorem 3.7.

Let $R=[r(u,v)]$ be the resistance distance matrix of the ladder graph $L_{n}$
of order $2n$. Then

$r(u_{i}^{\varepsilon},u_{j}^{\varepsilon^{\prime}})=\frac{d}{2}-\varepsilon\varepsilon^{\prime}\frac{a_{i}a_{n-j+1}}{2s_{n}}+\alpha(i,j),$

where $d=d(u_{i}^{+},u_{j}^{+})$, $\varepsilon,\varepsilon^{\prime}=\pm 1$,
and $\alpha(i,j)=(s(L_{n}/f_{i})+s(L_{n}/f_{j}))/4s_{n}$.

###### Proof.

Let $d=j-i=d(u_{i}^{+},u_{j}^{+})$ and $P:e_{i}^{+},\ldots,e_{j-1}^{+}$ be the
shortest path between $u_{i}^{+}$ and $u_{j}^{+}$, where $1\leq i<j\leq n$.
Then

$\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{i}^{+}}=\sum_{k=i}^{j-1}c_{k,i}=c_{i,i}+\sum_{k=i+1}^{j-1}c_{k,i},$

Thus

$\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{i}^{+}}=\left(1-\frac{i}{2n}-\frac{a_{n-i+1}s_{i}}{2s_{n}}\right)+\sum_{k=i+1}^{j-1}\left(\frac{1}{2}-\frac{k}{2n}+\frac{a_{i}s_{n-k}}{2s_{n}}\right).$

By equations (2) and (3),

$\displaystyle\sum_{k=i}^{j}q^{+}_{e_{k}^{+},u_{i}^{+}}$
$\displaystyle=\frac{d+1}{2}-\frac{id}{2n}-\frac{d(d-1)}{4n}+\frac{a_{i}(a_{n-i}-a_{n-j+1})-a_{n-i+1}(a_{i+1}-a_{i})}{4s_{n}}$
$\displaystyle=\frac{d+1}{2}-\frac{id}{2n}-\frac{d(d-1)}{4n}+\frac{a_{n}-a_{n+1}+a_{i}a_{n-i+1}-a_{i}a_{n-j+1}}{4s_{n}}$
$\displaystyle=\frac{d+1}{2}-\frac{id}{2n}-\frac{d(d-1)}{4n}+\frac{-2s_{n}+s(L_{n}/f_{i})-a_{i}a_{n-j+1}}{4s_{n}}$
$\displaystyle=\frac{d}{2}-\frac{id}{2n}-\frac{d(d-1)}{4n}-\frac{a_{i}a_{n-j+1}}{4s_{n}}+\frac{s(L_{n}/f_{i})}{4s_{n}}.$

The same argument for the vertex $u_{j}^{+}$ yields

$\displaystyle\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{j}^{+}}=\sum_{k=i}^{j-1}c_{k,j}$
$\displaystyle=-\frac{id}{2n}-\frac{d(d-1)}{4n}+\frac{a_{n-j+1}(a_{i}-a_{d+i})}{4s_{n}}$
$\displaystyle=-\frac{id}{2n}-\frac{d(d-1)}{4n}+\frac{a_{i}a_{n-j+1}}{4s_{n}}-\frac{s(L_{n}/f_{j})}{4s_{n}}.$

Therefore,

$\displaystyle
r(u_{i}^{+},u_{j}^{+})=\frac{d}{2}-\frac{a_{i}a_{n-j+1}}{2s_{n}}+\frac{s(L_{n}/f_{i})+s(L_{n}/f_{j})}{4s_{n}}.$

From the block structure of $H$ it follows that
$r(u_{i}^{-},u_{j}^{-})=r(u_{i}^{+},u_{j}^{+})$. Finally, from the equality

$\displaystyle\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{j}^{-}}=\sum_{k=i}^{j-1}d_{k,j}$
$\displaystyle=\sum_{k=i}^{j-1}\left(-\frac{k}{2n}+\frac{a_{n-j+1}s_{k}}{2s_{n}}\right)$
$\displaystyle=-\frac{id}{2n}-\frac{d(d-1)}{4n}-\frac{a_{i}a_{n-j+1}}{4s_{n}}+\frac{s(L_{n}/f_{j})}{4s_{n}}$

we obtain

$\displaystyle r(u_{j}^{-},u_{i}^{+})=r(u_{i}^{+},u_{j}^{-})$
$\displaystyle=\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{i}^{+}}+b_{j,i}-\sum_{k=i}^{j-1}q^{+}_{e_{k}^{+},u_{j}^{-}}+b_{j,j}$
$\displaystyle=\frac{d}{2}+\frac{a_{i}a_{n-j+1}}{2s_{n}}+\frac{s(L_{n}/f_{i})+s(L_{n}/f_{j})}{4s_{n}},$

as required. ∎

Utilizing Theorem 3.7, we can compute the Kirchhoff index of ladder graphs.

######  Corollary 3.8.

For any $n\geq 1$,

$Kf(L_{n})=\frac{n^{2}}{3}\left(n+1+\frac{a_{n}}{s_{n}}\right).$

###### Proof.

We have

$\displaystyle Kf(L_{n})$ $\displaystyle=\sum_{\\{u,v\\}\subseteq
V(L_{n})}r(u,v)$
$\displaystyle=\frac{4n\sum_{i=1}^{n}s(L_{n}/f_{i})}{4s_{n}}+2Kf(P_{n})$
$\displaystyle=\frac{n\sum_{i=1}^{n}a_{i}a_{n+1-i}}{s_{n}}+\frac{n(n^{2}-1)}{3}$

from which the result follows by equation (6). ∎

##  4\. Circular ladder graph

In this section, we compute the Moore-Penrose inverse of incidence matrices of
circular ladder graphs. Accordingly, we obtain the resistance distance matrix
and Kirchhoff index of circular ladders. In what follows, $CL_{n}$ denotes the
circular ladder graph of order $2n$ and $u_{1}^{+},\ldots,u_{n}^{+}$ and
$u_{1}^{-},\ldots,u_{n}^{-}$ stand for the vertices of the internal and
external cycles of $CL_{n}$, respectively, as in Fig. 2. Also,
$u_{i+kn}^{\varepsilon}$ denotes the vertex $u_{i}^{\varepsilon}$ for all
$1\leq i\leq n$, integer $k$, and $\varepsilon=\pm 1$.

$u_{1}^{+}$$u_{2}^{+}$$u_{3}^{+}$$u_{4}^{+}$$u_{n-2}^{+}$$u_{n-1}^{+}$$u_{n}^{+}$$u_{1}^{-}$$u_{2}^{-}$$u_{3}^{-}$$u_{4}^{-}$$u_{n-2}^{-}$$u_{n-1}^{-}$$u_{n}^{-}$
Figure 2. Circular ladder $CL_{n}$

Let $\Gamma$ be a planar graph embedded on a shaded cylinder as in Fig. 3. A
Jordan path is a smooth simple path starting from outer border and ending at
the inner border such that it meets every face and edge at most once and never
paths through vertices (see [20]). Removing edges crossed by a Jordan path
$\mathcal{J}$ results in a graph whose spanning trees depend only on
$\mathcal{J}$. On the other hand, every spanning tree of $\Gamma$ determines a
Jordan path. Therefore, the set of all spanning trees of $\Gamma$ can be
partitioned into sets $S_{\mathcal{J}}(\Gamma)$ where each set
$S_{\mathcal{J}}(\Gamma)$ corresponds to spanning trees associated to a Jordan
path $\mathcal{J}$. Let $s_{\mathcal{J}}(\Gamma):=\\#S_{\mathcal{J}}(\Gamma)$.

Figure 3. Jordan path

######  Theorem 4.1.

Let $CL_{n}$ be the circular ladder graph of order $2n$. If $f$ is a spoke,
then

$s(CL_{n}/f)=ns(L_{n}).$

###### Proof.

Consider a natural embedding of circular ladder $CL_{n}$ on a cylinder as in
Fig. 3 and assume without loss of generality that $f=u_{1}^{+}u_{1}^{-}$. We
have two (equivalent up to reflection) classes of Jordan paths:

  * (1)

Jordan paths $\mathcal{J}_{i,j}^{1}$ through edges
$u_{i}^{+}u_{i+1}^{+},u_{i+1}^{+}u_{i+1}^{-},\ldots,u_{j-1}^{+}u_{j-1}^{-},u_{j-1}^{-}u_{j}^{-}$
for $1\leq i<j\leq n+1$;

  * (2)

Jordan paths $\mathcal{J}_{i,j}^{2}$ through edges
$u_{j}^{+}u_{j-1}^{+},u_{j-1}^{+}u_{j-1}^{-},\ldots,u_{i+1}^{+}u_{i+1}^{-},u_{i}^{-}u_{i+1}^{-}$
for $1\leq i<j\leq n+1$.

Note that the Jordan paths in parts (1) and (2) coincide when $j=i+1$. Also
$s_{\mathcal{J}_{i,j}^{1}}(CL_{n}/f)=s_{\mathcal{J}_{i,j}^{2}}(CL_{n}/f)$.
Utilizing the above Jordan paths in conjunction with equation (11), yields

$\displaystyle s(CL_{n}/f)=$ $\displaystyle 2\sum_{1\leq i<j\leq
n+1}s_{\mathcal{J}_{i,j}^{1}}(CL_{n}/f)-\sum_{i=1}^{n}s_{\mathcal{J}_{i,i+1}^{1}}(CL_{n}/f)$
$\displaystyle=$ $\displaystyle 2\sum_{1\leq i<j\leq
n+1}s(L_{n+1+i-j}/f_{i})-\sum_{i=1}^{n}s(L_{n}/f_{i})$ $\displaystyle=$
$\displaystyle 2\sum_{1\leq i<j\leq
n+1}a_{i}a_{n+2-j}-\sum_{i=1}^{n}a_{i}a_{n+1-i}$ $\displaystyle=$
$\displaystyle 2\sum_{2\leq j\leq
n+1}s_{j-1}a_{n+2-j}-\sum_{i=1}^{n}a_{i}a_{n+1-i}$ $\displaystyle=$
$\displaystyle 2\sum_{1\leq i\leq
n}a_{i}s_{n+1-i}-\sum_{i=1}^{n}a_{i}a_{n+1-i},$

where $L_{s}/f_{t}$ denotes the contraction of the $t$’s spoke of $L_{s}$ from
right or left. Let $x_{n}:=\sum_{i=1}^{n}a_{i}a_{n+1-i}$ and
$y_{n}:=\sum_{i=1}^{n}a_{i}s_{n+1-i}$ for all $n\geq 1$. Clearly,
$y_{n}-y_{n-1}=x_{n}$. Also, from equation (4), it follows that
$y_{n}+y_{n-1}=ns_{n}$. Thus $2y_{n}=ns_{n}+x_{n}$, which implies that

$s(CL_{n}/f)=2y_{n}-x_{n}=ns_{n},$

as required. ∎

######  Theorem 4.2.

Let $CL_{n}$ be the circular ladder graph of order $2n$. Then

$s(CL_{n})=\frac{1}{2}n(a_{n+1}+a_{n}-2).$

###### Proof.

Analogous to the proof of Theorem 4.1, we consider two families of (equivalent
up to reflection) Jordan paths given by $\mathcal{J}_{i,j}^{1}$ and
$\mathcal{J}_{i,j}^{2}$ as in Theorem 4.1 where $i$ is behind $j$ regarding
clockwise rotations. From the rotational symmetry of $CL_{n}$ and the equality
of $\mathcal{J}_{i,j}^{1}$ and $\mathcal{J}_{i,j}^{2}$ for consecutive values
of $i$ and $j$, it follows that

$\displaystyle s(CL_{n})$
$\displaystyle=2n\sum_{i=1}^{n}s(L_{i})-ns(L_{n})=\frac{1}{2}n(a_{n+1}-a_{n}-2)$

by equation (2). ∎

In what follows, $e_{i}^{+}$, $e_{i}^{-}$, and $f_{i}$ stand for directed
edges $(u_{i}^{+},u_{i+1}^{+})$, $(u_{i}^{-},u_{i+1}^{-})$, and
$(u_{i}^{+},u_{i}^{-})$ for $i=1,\ldots,n$. Accordingly, we consider the
orientation of $CL_{n}$ whose directed edges are $e_{i}^{+}$, $e_{i}^{-}$, and
$f_{i}$ for $i=1,\ldots,n$.

######  Theorem 4.3.

Let $CL_{n}$ be the oriented circular ladder graph of order $2n$ with
incidence matrix $Q$. The Moore-Penrose inverse $Q^{+}=[q_{e,u}^{+}]$ of $Q$
is given by

$\displaystyle q^{+}_{e_{i+t}^{\varepsilon^{\prime}},u_{i}^{\varepsilon}}$
$\displaystyle=-\varepsilon\varepsilon^{\prime}\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\varepsilon\varepsilon^{\prime}\frac{a_{t+1}+\varepsilon\varepsilon^{\prime}1}{4}-\frac{2t+1}{4n},$
$\displaystyle\varepsilon q^{+}_{f_{i+t},u_{i}^{\varepsilon}}$
$\displaystyle=\frac{a_{t+1}+a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\frac{a_{t+1}-a_{t}}{4}$

for all $1\leq i\leq n$, $t=0,\ldots,n-1$, and
$\varepsilon,\varepsilon^{\prime}=\pm 1$.

###### Proof.

Let $d=d(f_{i},u_{j})$. From the symmetry of graph, it follows that
$q^{+}_{f_{i},u_{i}^{-}}=-q^{+}_{f_{i},u_{i}^{+}}$. Hence

$2q^{+}_{f_{i},u_{i}^{+}}=r(u_{i}^{+},u_{i}^{-})=\frac{s(CL_{n}/f_{i})}{s(CL_{n})}=\frac{s_{n}}{a_{n+1}+a_{n}-2}$

by Lemma 3.5 and Theorems 4.1 and 4.2.

Consider the edge cut-set $\\{e_{j}^{+},f_{j+1},e_{j+1}^{+}\\}$, where $j=i+t$
with $0\leq t<n$. From [3, Lemma 2.1], it follows that

(12)
$q^{+}_{e_{j+1}^{+},u_{i}^{+}}=\delta_{i+1,j}-\frac{1}{2n}-q^{+}_{f_{j+1},u_{i}^{+}}+q^{+}_{e_{j}^{+},u_{i}^{+}}.$

Analogously, by considering the edge cut-set
$\\{e_{j}^{-},f_{j+1},e_{j+1}^{-}\\}$, we obtain

(13)
$q^{+}_{e_{j+1}^{-},u_{i}^{+}}=-\frac{1}{2n}+q^{+}_{f_{j+1},u_{i}^{+}}+q^{+}_{e_{j}^{-},u_{i}^{+}}.$

for all $j\geq i$. If we replace $j$ by $i-1$, and use the equalities

$q^{+}_{e_{i-1}^{+},u_{i}^{+}}=-q^{+}_{e_{i}^{+},u_{i}^{+}}\quad\text{and}\quad
q^{+}_{e_{i-1}^{-},u_{i}^{+}}=-q^{+}_{e_{i}^{-},u_{i}^{+}}$

obtained from the symmetry of the graph, then we get the initial values

$q^{+}_{e_{i}^{+},u_{i}^{+}}=\frac{1}{2}\left(1-\frac{1}{2n}-q^{+}_{f_{i},u_{i}^{+}}\right)\quad\text{and}\quad
q^{+}_{e_{i}^{-},u_{i}^{+}}=\frac{1}{2}\left(\frac{-1}{2n}+q^{+}_{f_{i},u_{i}^{+}}\right).$

On the other hand, the incidence vector of the cycle induced by
$\\{u_{j}^{+},u_{j+1}^{+},u_{j}^{-},u_{j+1}^{-}\\}$ lies in the left null
space of $Q^{+}$, from which we obtain

$q^{+}_{f_{j+1},u_{i}^{+}}-q^{+}_{e_{j}^{-},u_{i}^{+}}-q^{+}_{f_{j},u_{i}^{+}}+q^{+}_{e_{j}^{+},u_{i}^{+}}=0$

so that

(14)
$q^{+}_{f_{j+1},u_{i}^{+}}=q^{+}_{e_{j}^{-},u_{i}^{+}}+q^{+}_{f_{j},u_{i}^{+}}-q^{+}_{e_{j}^{+},u_{i}^{+}}$

for all $j\geq i$. Let

$C:=\begin{pmatrix}2&-1&-1\\\ -1&2&1\\\
-1&1&1\end{pmatrix}\quad\text{and}\quad D:=-\frac{1}{2n}\begin{pmatrix}1\\\
1\\\ 0\end{pmatrix}.$

If $X_{j}$ denotes the column matrix $(q^{+}_{e_{j}^{+},u_{i}^{+}}\
q^{+}_{e_{j}^{-},u_{i}^{+}}\ q^{+}_{f_{j},u_{i}^{+}})^{T}$ for $i\leq j<i+n$,
then from the equalities (13) and (14), it follows that

$X_{j+1}=CX_{j}+D$

for all $i\leq j<i+n$. Using induction on $t\geq 0$, one can prove that

$C^{t}=\begin{pmatrix}\frac{a_{t+1}+1}{2}&-\frac{a_{t+1}-1}{2}&-s_{t}\\\
-\frac{a_{t+1}-1}{2}&\frac{a_{t+1}+1}{2}&s_{t}\\\
-s_{t}&s_{t}&a_{t}\end{pmatrix}$

and

$I+C+\cdots+C^{t-1}=\begin{pmatrix}\frac{s_{t}+t}{2}&-\frac{s_{t}-t}{2}&-\frac{a_{t}-1}{2}\\\
-\frac{s_{t}-t}{2}&\frac{s_{t}+t}{2}&\frac{a_{t}-1}{2}\\\
-\frac{a_{t}-1}{2}&\frac{a_{t}-1}{2}&s_{t-1}+1\end{pmatrix}.$

Now, from the equality

$X_{i+t}=C^{t}X_{i}+(I+C+\cdots+C^{t-1})D,$

it follows that

$\displaystyle q^{+}_{e_{i+t}^{+},u_{i}^{+}}$
$\displaystyle=-\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\frac{a_{t+1}+1}{4}-\frac{2t+1}{4n},$
$\displaystyle q^{+}_{e_{i+t}^{-},u_{i}^{+}}$
$\displaystyle=\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\frac{a_{t+1}-1}{4}-\frac{2t+1}{4n},$
$\displaystyle q^{+}_{f_{i+t},u_{i}^{+}}$
$\displaystyle=\frac{a_{t+1}+a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\frac{a_{t+1}-a_{t}}{4}$

for all $t=0,\ldots,n-1$. Finally, the result follows by applying the same
argument or using the symmetry of graph for vertices
$u_{1}^{-},\ldots,u_{n}^{-}$. ∎

######  Corollary 4.4.

If $n=2k+1$, then

$s(CL_{n})=ns(L_{n})\left(\frac{s(L_{k+1})+s(L_{k})}{s(L_{k+1})-s(L_{k})}\right)$

and if $n=2k$, then

$s(CL_{n})=ns(L_{n})\left(\frac{s(L_{k+1})+2s(L_{k})+s(L_{k-1})}{s(L_{k+1})-s(L_{k-1})}\right).$

###### Proof.

We know that
$s(CL_{n}/f_{i})/s(CL_{n})=r(u_{i}^{+},u_{i}^{-})=2q^{+}_{f_{i},u_{i}^{+}}$.
Thus

$s(CL_{n})=ns(L_{n})\cdot\frac{a_{n+1}+a_{n}-2}{a_{n+1}-a_{n}}.$

Using Binet’s formulas, one can easily show that

$\frac{a_{n+1}+a_{n}-2}{a_{n+1}-a_{n}}=\frac{s_{k+1}+s_{k}}{s_{k+1}-s_{k}}\quad\text{or}\quad\frac{s_{k+1}+2s_{k}+s_{k-1}}{s_{k+1}-s_{k}}$

according to $n=2k+1$ or $n=2k$, respectively. The result follows. ∎

######  Corollary 4.5.

Let $1\leq i\leq n$, $0\leq t<n$, and $\varepsilon,\varepsilon^{\prime}=\pm
1$. Then

$r(u_{i}^{\varepsilon},u_{i+t}^{\varepsilon^{\prime}})=-\varepsilon\varepsilon^{\prime}\frac{a_{t+1}+a_{t}-\varepsilon\varepsilon^{\prime}2}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\varepsilon\varepsilon^{\prime}\frac{a_{t+1}-a_{t}}{4}+\frac{t}{2}-\frac{t^{2}}{2n}.$

###### Proof.

First observe that
$r(u_{i}^{-},u_{i+t}^{-})=r(u_{i}^{+},u_{i+t}^{+})=r(u_{0}^{+},u_{t}^{+})$ by
Theorem 4.3 and Lemma 3.5. We have

$\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{+}}=-\frac{a_{t+1}+a_{t}-2}{8}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\frac{a_{t+1}-a_{t}}{8}+\frac{t}{4}-\frac{t^{2}}{4n}$

and

$\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{-}}=\frac{a_{t+1}+a_{t}-2}{8}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\frac{a_{t+1}-a_{t}}{8}+\frac{t}{4}-\frac{t^{2}}{4n}.$

Also,

$\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{t}^{+}}=\sum_{s=0}^{t-1}q^{+}_{e_{n+s-t}^{+},u_{0}^{+}}=-\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{+}}$

and

$\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{t}^{-}}=\sum_{s=0}^{t-1}q^{+}_{e_{n+s-t}^{+},u_{0}^{-}}=-\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{-}}$

by symmetry of the graph. Let $1\leq i\leq n$ and $0\leq t<n$. Then, by Lemma
3.5,

$\displaystyle r(u_{i}^{+},u_{i+t}^{+})$
$\displaystyle=r(u_{0}^{+},u_{t}^{+})=2\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{+}}$
$\displaystyle=-\frac{a_{t+1}+a_{t}-2}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\frac{a_{t+1}-a_{t}}{4}+\frac{t}{2}-\frac{t^{2}}{2n}$

and

$\displaystyle r(u_{i}^{+},u_{i+t}^{-})$
$\displaystyle=r(u_{0}^{+},u_{t}^{-})=\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{+}}-\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{t}^{-}}+q^{+}(f_{t},u_{0}^{+})-q^{+}(f_{t},u_{t}^{-})$
$\displaystyle=\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{+}}+\sum_{s=0}^{t-1}q^{+}_{e_{s}^{+},u_{0}^{-}}+q^{+}(f_{t},u_{0}^{+})+q^{+}(f_{0},u_{0}^{+})$
$\displaystyle=\frac{a_{t+1}+a_{t}+2}{4}\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\frac{a_{t+1}-a_{t}}{4}+\frac{t}{2}-\frac{t^{2}}{2n},$

as required. ∎

######  Corollary 4.6.

For every $n\geq 3$, we have

$Kf(CL_{n})=n^{2}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}+\frac{1}{6}n(n^{2}-1)$

It is known from [19, Theorem F] that $Kf(\Gamma)=n\mathrm{tr}(L^{+}(\Gamma))$
for every graph $\Gamma$ of order $n$, where $L(\Gamma)$ denotes the Laplacian
matrix of $\Gamma$. Also, from the definition, we know that

$r(u,v)=L^{+}_{u,u}+L^{+}_{v,v}-L^{+}_{u,v}-L^{+}_{v,u}$

for all $u,v\in V(\Gamma)$, where $L^{+}$ denotes the Moore-Penrose inverse of
$L=L(\Gamma)$ and $r(u,v)$ is the resistance distance between vertices $u$ and
$v$ of $\Gamma$. Utilizing the fact that circular ladders are vertex
transitive graphs, we obtain the following result immediately.

######  Corollary 4.7.

Let $L:=L(CL_{n})$ be the Laplacian matrix of $CL_{n}$ and
$L^{+}=[l_{uv}^{+}]$ be its Moore-Penrose inverse. Then

$l_{u_{i}^{\varepsilon},u_{i+t}^{\varepsilon^{\prime}}}=\varepsilon\varepsilon^{\prime}\frac{a_{t+1}+a_{t}}{8}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}-2}-\varepsilon\varepsilon^{\prime}\frac{a_{t+1}-a_{t}}{8}-\frac{t}{4}+\frac{t^{2}}{4n}+\frac{1}{24}\left(n-\frac{1}{n}\right)$

for all $1\leq i\leq n$, $0\leq t<n$, and
$\varepsilon,\varepsilon^{\prime}=\pm 1$.

##  5\. Möbius ladder graph

Following the same techniques used in the previous section, we shall compute
the Moore-Penrose inverse of incidence matrices of Möbius ladder graphs.
Likewise, we obtain the resistance distance matrix and Kirchhoff index of
Möbius ladders. In what follows, $M_{n}$ denotes the Möbius ladder graph of
order $2n$ and $u_{1},\ldots,u_{n}$ and $v_{1},\ldots,v_{n}$ stand for the
vertices of the external and internal cycles of $M_{n}$, respectively, as in
Fig. 4(a). Also, we may set $u_{i+n}:=v_{i}$ for all $i=1,\ldots,n$, and that
$u_{i+2nk}:=u_{i}$ for all $1\leq i\leq 2n$ and integers $k$. In what follows,
$f_{i}:=u_{i}v_{i}$ ($1\leq i\leq n$) is the $i$-th spoke and
$e_{i}=u_{i}u_{i+1}$ is the $i$-th edge of the $2n$-cycle
$u_{1},\ldots,u_{2n}$ ($1\leq i\leq 2n$), see the second drawing in Fig. 4(b).

$u_{1}$$u_{2}$$u_{3}$$u_{4}$$u_{n-3}$$u_{n-2}$$u_{n-1}$$u_{n}$$v_{1}$$v_{2}$$v_{3}$$v_{4}$$v_{n-3}$$v_{n-2}$$v_{n-1}$$v_{n}$ |  | $u_{1}$$u_{2}$$u_{3}$$u_{n-2}$$u_{n-1}$$u_{n}$$u_{n+1}$$u_{n+2}$$u_{n+3}$$u_{2n-2}$$u_{2n-1}$$u_{2n}$
---|---|---
(a) |  | (b)
Figure 4. Möbius ladder graph $M_{n}$

######  Theorem 5.1.

Let $M_{n}$ be the Möbius ladder graph of order $2n$ and $Q$ be its incidence
matrix. Then the Moore-Penrose inverse of $Q=[q^{+}_{e,v}]$ is given by

$\displaystyle q^{+}_{e_{i+t},u_{i}}$
$\displaystyle=-\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}+\frac{a_{t+1}+1}{4}-\frac{2t+1}{4n},$
$\displaystyle q^{+}_{e_{n+i+t},u_{i}}$
$\displaystyle=\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}-\frac{a_{t+1}-1}{4}-\frac{2t+1}{4n},$
$\displaystyle\varepsilon_{i+t}q^{+}_{f_{i+t},u_{i}}$
$\displaystyle=\frac{a_{t+1}+a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}-\frac{a_{t+1}-a_{t}}{4}$

for all $t=0,\ldots,n-1$.

###### Proof.

First observe that $s(M_{n}/f_{i})=s(CL_{n}/f_{i})$ and that by [20, p. 41]

$s(M_{n})=s(CL_{n})+2n=\frac{n}{2}(a_{n+1}+a_{n}+2).$

Therefore, by Lemma 3.5 and Theorem 4.1,

$q^{+}_{f_{i},u_{i}}=-q^{+}_{f_{i},v_{i}}=\frac{r(u_{i},v_{i})}{2}=\frac{s(M_{n}/f_{i})}{2s(M_{n})}=\frac{s_{n}}{a_{n+1}+a_{n}+2}.$

Consider the edge cut-set $\\{e_{j},f_{j+1},e_{j+1}\\}$. From [3, Lemma 2.1],
it follows that

(15)
$q^{+}_{e_{j+1},u_{i}}=\delta_{i,j+1}-\frac{1}{2n}+q^{+}_{e_{j},u_{i}}-\varepsilon_{j+1}q^{+}_{f_{j+1},u_{i}},$

where $\varepsilon_{j}=(-1)^{[(j-1)/n]}$. Since
$q^{+}_{e_{i},u_{i}}=-q^{+}_{e_{i-1},u_{i}}$ by the symmetry of the graph,
equation (15) yields the initial values

$q^{+}_{e_{i},u_{i}}=\frac{1}{2}\left(1-\frac{1}{2n}-\varepsilon_{i}q^{+}_{f_{i},u_{i}}\right)$

for all $i=1,\ldots,2n$. Likewise, from $f_{n+i}=f_{i}$ and
$q^{+}_{e_{n+i},u_{i}}=-q^{+}_{e_{n+i-1},u_{i}}$ it follows that

$q^{+}_{e_{n+i},u_{i}}=\frac{1}{2}\left(-\frac{1}{2n}+\varepsilon_{i}q^{+}_{f_{i},u_{i}}\right)$

for all $i=1,\ldots,2n$. Notice that the values of $q^{+}_{e_{i+t},u_{i}}$ are
independent of $i$ by the symmetry of the graph, so we may restrict ourselves
to $i$ with $1\leq i\leq n$ and assume $\varepsilon_{i}=1$.

On the other hand, the incidence vector of the cycle
$u_{j},u_{j+1},u_{n+j},u_{n+j+1}$ lies in the left null space of $Q^{+}$, from
which it follows that

$q^{+}_{e_{j},u_{i}}+\varepsilon_{j+1}q^{+}_{f_{j+1},u_{i}}-q^{+}_{e_{n+j},u_{i}}-\varepsilon_{j}q^{+}_{f_{j},u_{i}}=0$

or

(16)
$\varepsilon_{j+1}q^{+}_{f_{j+1},u_{i}}=q^{+}_{e_{n+j},u_{i}}-q^{+}_{e_{j},u_{i}}+\varepsilon_{j}q^{+}_{f_{j},u_{i}}.$

Let

$C:=\begin{pmatrix}2&-1&-1&0\\\ -1&2&0&-1\\\ -1&1&1&0\\\ 1&-1&0&1\\\
\end{pmatrix}\quad\text{and}\quad D:=-\frac{1}{2n}\begin{pmatrix}1\\\ 1\\\
0\\\ 0\end{pmatrix}.$

If $X_{j}$ denotes the column matrix $(q^{+}_{e_{j},u_{i}}\
q^{+}_{e_{n+j},u_{i}}\ \varepsilon_{j}q^{+}_{f_{j},u_{i}}\
\varepsilon_{n+j}q^{+}_{f_{n+j},u_{i}})^{T}$ for $i\leq j<i+n$, then from the
equalities (15) and (16), it follows that

$X_{j+1}=CX_{j}+D$

for all $i\leq j<n+i-1$. Using induction on $t\geq 0$, one can prove that

$C^{t}=\begin{pmatrix}\frac{a_{t+1}+1}{2}&-\frac{a_{t+1}-1}{2}&-\frac{s_{t}+t}{2}&\frac{s_{t}-t}{2}\\\
-\frac{a_{t+1}-1}{2}&\frac{a_{t+1}+1}{2}&\frac{s_{t}-t}{2}&-\frac{s_{t}+t}{2}\\\
-s_{t}&s_{t}&\frac{a_{t}+1}{2}&-\frac{a_{t}-1}{2}\\\
s_{t}&-s_{t}&-\frac{a_{t}-1}{2}&\frac{a_{t}+1}{2}\\\ \end{pmatrix}$

and

$I+C+\cdots+C^{t-1}=\begin{pmatrix}\frac{s_{t}+t}{2}&-\frac{s_{t}-t}{2}&-\frac{a_{t}+t^{2}-t-1}{4}&\frac{a_{t}-t^{2}+t-1}{4}\\\
-\frac{s_{t}-t}{2}&\frac{s_{t}+t}{2}&\frac{a_{t}-t^{2}+t-1}{4}&-\frac{a_{t}+t^{2}-t-1}{4}\\\
-\frac{a_{t}-1}{2}&\frac{a_{t}-1}{2}&\frac{s_{t-1}+t+1}{2}&-\frac{s_{t-1}-(t-1)}{2}\\\
\frac{a_{t}-1}{2}&-\frac{a_{t}-1}{2}&-\frac{s_{t-1}-(t-1)}{2}&\frac{s_{t-1}+t+1}{2}\\\
\end{pmatrix}.$

Now, from the equality

$X_{i+t}=C^{t}X_{i}+(I+C+\cdots+C^{t-1})D,$

it follows that

$\displaystyle q^{+}_{e_{i+t},u_{i}}$
$\displaystyle=-\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}+\frac{a_{t+1}+1}{4}-\frac{2t+1}{4n},$
$\displaystyle q^{+}_{e_{n+i+t},u_{i}}$
$\displaystyle=\frac{2a_{t+1}-a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}-\frac{a_{t+1}-1}{4}-\frac{2t+1}{4n},$
$\displaystyle\varepsilon_{i+t}q^{+}_{f_{i+t},u_{i}}$
$\displaystyle=\frac{a_{t+1}+a_{t}}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}-\frac{a_{t+1}-a_{t}}{4}$

for all $t=0,\ldots,n-1$, as required. ∎

######  Corollary 5.2.

Let $1\leq i\leq 2n$ and $0\leq t\leq n$. Then

$r(u_{i},u_{i+t})=-\frac{a_{t+1}+a_{t}-2}{4}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}+\frac{a_{t+1}-a_{t}}{4}+\frac{t}{2}-\frac{t^{2}}{2n}.$

###### Proof.

The result follows from Lemma 3.5 and Theorem 5.1 in conjunction with the fact
that

$\sum_{j=0}^{t-1}q^{+}_{e_{i+j},u_{i+t}}=-\sum_{j=0}^{t-1}q^{+}_{e_{i+j},u_{i}}$

by the symmetry of the graph. ∎

Theorem 5.2 in conjunction with equation (5) yields

######  Corollary 5.3.

For every $n\geq 3$, we have

$Kf(CL_{n})=n^{2}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}+\frac{1}{6}n(n^{2}-1).$

Finally, the Moore-Penrose inverse of the Laplacian matrix of $M_{n}$ can be
obtained analogous to Corollary 4.7.

######  Corollary 5.4.

Let $L:=L(M_{n})$ be the Laplacian matrix of $M_{n}$ and $L^{+}=[l_{uv}^{+}]$
be its Moore-Penrose inverse. Then

$l^{+}_{u_{i},u_{i+t}}=\frac{a_{t+1}+a_{t}}{8}\cdot\frac{a_{n+1}-a_{n}}{a_{n+1}+a_{n}+2}-\frac{a_{t+1}-a_{t}}{8}-\frac{t}{4}+\frac{t^{2}}{4n}+\frac{1}{24}\left(n-\frac{1}{n}\right)$

for all $1\leq i\leq 2n$ and $0\leq t<n$.

## References

  * [1] A. Azimi and R. B. Bapat, Moore-Penrose inverse of the incidence matrix of a distance regular graph, Linear Algebra Appl. 551 (2018), 92–103.
  * [2] A. Azimi and R. B. Bapat, The Moore-Penrose inverse of the incidence matrix of complete multipartite and bi-block graphs, Discrete Math. 342(8) (2019), 2393–2401.
  * [3] A. Azimi, R. B. Bapat, and E. Estaji, Moore-Penrose inverse of incidence matrix of graphs with complete and cyclic blocks, Discrete Math. 342 (2019), 10–17.
  * [4] A. Azimi, R. B. Bapat, and M. Farrokhi D. G., Resistance distance of blow-ups of trees, Discrete Appl. Math. 320 (2022), 311–323.
  * [5] P. Ali, F. Atik, and R. B. Bapat, Identities for minors of the Laplacian, resistance and distance matrices of graphs with arbitrary weights, Linear Multilinear Algebra 68(2) (2020), 323–336.
  * [6] R. B. Bapat, Graphs and Matrices, Second Edition, Springer, 2014.
  * [7] R. B. Bapat and S. Gupta, Resistance distance in wheels and fans, Indian J. Pure Appl. Math. 41(1) (2010), 1–13.
  * [8] N. L. Biggs, Potential theory on distance-regular graphs, Combin. Probab. Comput. 2(3) (1993), 243–255.
  * [9] A. Ben-Israel and T. N. E. Greville, Generalized Inverses: Theory and Applications, Second Edition, Springer, New York, 2003.
  * [10] S. L. Campbell and C. D. Meyer, Generalized Inverses of Linear Transformations, Pitman, London, 1979.
  * [11] Z. Cinkir, Effective resistances and Kirchhoff index of ladder graphs, J. Math. Chem. 54(4) (2016), 955–966.
  * [12] K. C. Das and Y. Yang, Eigenvalues of the resistance-distance matrix of complete multipartite graphs, J. Inequal. Appl. (2017), No. 296, 11 pp.
  * [13] E. Ernesto and M. Delio, Hubs-biased resistance distances on graphs and networks, J. Math. Anal. Appl. 507, no. 1, (2022), Paper No. 125728, 22 pp.
  * [14] X. Gao, Y. Luo, and W. Liu, Resistance distances and the Kirchhoff index in Cayley graphs, Discrete Appl. Math. 159(17) (2011), 2050–2057.
  * [15] S. V. Gervacio, Resistance distance in complete $n$-partite graphs, Discrete Appl. Math. 203 (2016), 53–61.
  * [16] Q. Huang, H. Chen, and Q. Deng, Resistance distances and the Kirchhoff index in double graphs, J. Appl. Math. Comput. 50 (2016), 1–14.
  * [17] Y. Ijiri, On the generalized inverse of an incidence matrix, J. Soc. Indust. Appl. Math. 13 (1965), 827–836.
  * [18] M. Kagan and B. Mata, A physics perspective on the resistance distance for graphs, Math. Comput. Sci. 13(1-2) (2019), 105–115.
  * [19] D. J. Klein and M. Randić, Resistance distance, J. Math. Chem. 12 (1993), 81–95.
  * [20] M. Rubey, Counting Spanning Trees, Diplomarbeit, Universitat Wein, Mai, 2000.
  * [21] N. J. A. Sloane, The On-Line Encyclopedia of Integer Sequences, https://oeis.org.
  * [22] Y. Yi and S. Yue, Resistance-distance betweenness: Definition, properties, and application in complex linguistic networks, Phys. A, 610, (2023), Paper No. 128408.


# Depleting states dictate the ideal glass and physics of glass transition

Anshul D. S. Parmar<EMAIL_ADDRESS>Institute of Physical Chemistry,
Westfälische Wilhelms-Universität Münster, Corrensstraße 28/30, 48149 Münster,
Germany Helmholtz-Institute Münster (IEK-12) Forschungszentrum Jülich GmbH,
Corrensstraße 46, 48149 Münster, Germany Andreas Heuer  andheuer@uni-
muenster.de Institute of Physical Chemistry, Westfälische Wilhelms-Universität
Münster, Corrensstraße 28/30, 48149 Münster, Germany

###### Abstract

Understanding the properties of supercooled fluids in equilibrium, even below
the calorimetric glass transition $T_{g}$, is an elusive challenge. This is
even more true for the properties of the ideal glass, defined as the minimum
of the potential energy landscape (PEL) in the non-crystalline regime.
Although its existence is a mathematical necessity due to the finite range of
energies, its properties, and physical relevance are still undecided. Here we
combine the million-fold acceleration of Swap Monte Carlo with PEL analysis to
study a non-network 2D model glass former in equilibrium for a wide range of
system sizes and temperatures. We observe the transition from fragile to
strong behavior and provide a generic perspective to such observations in many
experimental non-network glass formers, or the saturation of structural
disorder upon cooling. Furthermore, we can identify a particular system size
that shows quantitative agreement with the macroscopic limit, while allowing a
complete characterization of the potential energy landscape (PEL) down to its
global minimum. This implies the availability of configurational entropy down
to the zero temperature limit without the detour via liquid entropy, thereby
quantifying the putative entropy crisis. Through appropriate system size
analysis, we can identify a low-energy depletion regime, including the ideal
glass, and reveal its physical relevance for equilibrium properties over a
wide temperature range, including below $T_{g}$.

Introduction— Glass, a disordered material, has remained ubiquitous in human
civilization for millennia. Furthermore, the intriguing and elusive problem of
the nature of glasses, emerging upon cooling a supercooled liquid, persists as
an unsolved and notorious enigma in the realm of condensed matter physics
Anderson (1995). With its broad range of applicability, glass finds relevance
across a spectrum of materials, from polymers and alloys to colloids, thereby
marking its significance in a diverse array of applications Debenedetti and
Stillinger (2001). The equilibrium window of supercooling, metastability below
melting, ranges from 100ps to 100s of structural relaxation in the laboratory,
marking the onset ($T_{onset}$), and experimental glass transition temperature
($T_{g}$), respectively Schmidtke  _et al._ (2012). Glass-forming liquids fall
into two classes based on their temperature dependence of dynamical behavior:
“strong” liquids with exponential dependence on inverse temperature (known as
the Arrhenius relation) and “fragile” liquids with a super-Arrhenius (SA)
behavior, exhibiting a stronger-than-exponential dependence Angell (1991a).
Network formers such as silica and water display both types of dynamics,
separated by a characteristic Fragile to Strong Crossover (FSC) in the
supercooled regime Saika-Voivod  _et al._ (2001); Saksaengwijit  _et al._
(2004).

Thermodynamics-wise, a relation between dynamics and entropy has been
formulated in well-known approaches such as Adam-Gibbs Adam and Gibbs (1965)
or the Random first-order transition (RFOT) theory Lubchenko and Wolynes
(2007); Kirkpatrick  _et al._ (1989). For network formers, the emergence of
the FSC indeed goes along with an inflection of the entropy Saika-Voivod  _et
al._ (2001). The distinct characteristics intrinsic to the system present a
formidable obstacle in constructing a comprehensive theoretical framework to
decipher the nature of glass.

The potential energy landscape (PEL), as a general approach, one studies the
distribution of the inherent structures (IS), i.e., local energy minima
(denoted IS-density), and their topology. This information helps to elucidate
the underlying nature of thermodynamic and dynamical properties Goldstein
(1969); Angell (1991b); Sastry (2001); Sciortino (2005); La Nave  _et al._
(2006); Heuer (2008); Baity-Jesi  _et al._ (2021). If the complete density of
states, including its ground state, were known, important thermodynamic
properties such as entropy or specific heat could be directly deduced for all
temperatures. Thus, the PEL approach has the potential to provide answers to
fundamental questions of glass physics, including the nature of fragility, the
ideal glass corresponding to the ground state of the IS, and the nature of
approaching the zero-temperature entropy, which so far remains unanswered due
to the astronomical time scales of relaxation Berthier and Biroli (2011).

A PEL study requires the use of relatively small system sizes because for
large systems, the information content is reduced due to the presence of many
uncorrelated subsystems Doliwa and Heuer (2003); Heuer (2008). Thus, one looks
for system sizes $N$ that are as small as possible, while ensuring minimal
finite-size effects for thermodynamic, structural, and dynamical observables.
Denoting $N_{c}$ as the corresponding system size, reasonable choices were
$N_{c}=65$ for binary 3D Lennard-JonesKob and Andersen (1994); Doliwa and
Heuer (2003), $N_{c}=99$ for 3D BKS-silica Saksaengwijit  _et al._ (2004), and
$N_{c}=80$ for 2D silica Roy and Heuer (2022). For the network formers, it was
even possible to extract the complete IS-density for $N=N_{c}$, displaying a
Gaussian energy dependence with a sharp decay at low energies Saksaengwijit
_et al._ (2004); Roy and Heuer (2022). An alternative structural approach to
silica can also be found Yu  _et al._ (2022). The emergence of this cutoff has
been related to the presence of well-defined low-energy states void of defects
such as dangling bonds Saksaengwijit  _et al._ (2004), and a subtle structural
change Yu  _et al._ (2022); Roy and Heuer (2022). The FSC emerged as a
consequence of that cutoff. Beyond the special case of network formers, the
complete IS-density is unknown for system sizes around $N_{c}$ due to the
necessity of exponentially demanding numerical burden at low energies.

Here we study a 2D non-network glass former via the swap Monte Carlo technique
Ninarello  _et al._ (2017). With the removal of previous sampling limitations,
we generate a large ensemble of equilibrium configurations for system sizes
$N\in[33,1056]$ and for the temperature range $T\gtrsim T_{g}/2$, i.e., even
exceeding the experimentally relevant regime $T>T_{g}$, while no
crystallization was observed. This allows the reconstruction of the IS-density
relevant for the range $T\gtrsim T_{g}/2$. Importantly, after identifying an
appropriate system size $N_{c}$, we can show that for $N\leq N_{c}$ even the
complete IS-density, including the bottom state, can be determined.
Consequently, for that range of system sizes, thermodynamic and structural
properties can be completely described for all temperatures $T>0$ (neglecting
quantum-mechanical effects), and the configurational entropy can be obtained
without starting from the liquid entropy and using techniques like
thermodynamic integration. For $N>N_{c}$, it will turn out that this physical
insight is at least relevant for $T\gtrsim T_{g}/2$, i.e., for all
experimentally relevant and exceeding supercooling. Furthermore, also an FSC
can be identified at low $T$, providing unprecedented information about the
presence and origin of the FSC in non-network formers. Furthermore, a
comparison with corresponding experiments and conclusions about Kauzmann’s
paradox is possible.

Specifically, we analyze a discrete-polydisperse mixture of purely repulsive
spheres with $10$ particle types in two dimensions. The composition is kept
strictly the same for all system sizes, ensuring consistency in the study; see
the Methods section for details Ninarello  _et al._ (2017); Parmar  _et al._
(2020). Equilibrated configurations are generated for a wide range of system
sizes $N\in\\{33,66,132,264,528,1056\\}$ and temperatures $T\in[0.0263,0.5]$,
where the calorimetric glass transition temperature is $T_{g}=0.0555$,
associated to time scales of 100 s. The lowest equilibrated sample from the
swap moves corresponds to astronomical time scales of $10^{21}$s. The well-
equilibrated independent samples from the swap MC were used as the initial
sample for the NVT-MD simulations with LAMMPS Plimpton (1995). The dynamics of
the system explored for $T\in[0.069,0.5]$, spanning 7(10) decades of
diffusion(relaxation) below the onset temperature($T_{onset}=0.25$). The
simulation at the lowest T is pushed to $t=6\cdot 10^{8}$ MD time, roughly
corresponding to the millisecond time scale. In what follows,
$e_{IS}=E_{IS}/N$ is the IS energy per particle.

Figure 1: Finite-size effects and dynamics. a The energy fluctuation per
particle for two different temperatures. It decreases proportionally to
$N^{-1/2}$. b The finite-size effect of the average IS energy relative to the
fluctuations, for the same two temperatures. Included is an empirical $1/N$
asymptotic behavior (dashed line). c The inverse diffusion constant $1/D$
(scaled with $D_{o}$, diffusion at the onset temperature) covers 7 decades of
glassy dynamics. It reveals the crossover from fragile to strong behavior at
low T($\lesssim 0.1$). The inset presents the low and high-T activation energy
ratio for all system sizes.

Quantifying finite-size effects— We start by evaluating the size dependence of
the standard deviation and the average energy of the IS-density at two
different temperatures. In the large system size limit, one expects that a
system of size $2N$ can be considered as an independent superposition of two
systems with size $N$. This superposition hypothesis implies that the standard
deviation of the energy distribution (per particle) scales like $N^{-1/2}$ as
indeed observed for all N, see Fig. 1(a). Generally, the energy $E_{N}$ will
deviate by a constant $c/2$ from the large-size limit $Ne_{\infty}$. Then the
relation $E_{2N}=2E_{N}$, resulting from the superposition hypothesis, yields
$e_{2N}=E_{N}/N=e_{\infty}+c/(2N)$, i.e. $e_{N}-e_{\infty}\propto 1/N$. Also,
this relation holds very well (Fig. 1b) for the whole range of $N$-values and
for both temperatures. Thus no deviations from the superposition hypothesis
are visible for $N\geq 33$.

The diffusion constant is an appropriate observable for the identification of
finite-size effects of the transport behavior Rehwald  _et al._ (2010). The
results are shown in fig. 1c. The emergence of an FSC at low temperatures
($T\lesssim 0.1$) can be observed for all system sizes. The ratio of the low-T
($A_{L}$) and high-T ($A_{H}$) activation energies are displayed in fig. 1c(
inset). Thus, a system as small as $N=33$ reflects the dynamics in the
thermodynamic limit on a semi-quantitative level, whereas for $N=66$, even a
nearly quantitative agreement can be observed. Detailed analysis shows, also
for structural relaxation (SI (VI)), that the low-T dynamics is incompatible
with alternative standard descriptions for supercooled liquids. The resulting
choice $N_{c}=66$ is also well supported by the weak $N$-dependence of the
participation ratio for $N\geq N_{c}$, quantifying the number of mobile
particles in relaxation events during an IS-IS transition (see SI (V)) Vogel
_et al._ (2004)).

Figure 2: From depletion to the bottom of the PEL. The per-particle potential
energy of IS with system sizes (a) N=33 and (b) N=66, respectively. The color
bar for each temperature represents the equilibrium probability distribution.
The deviation from the observed Gaussian behavior in the initial supercooled
regime is stressed with the dashed line fitted for the temperature range
$T\in[0.1429,0.25]$. The horizontal lines represent the lowest energy of IS
for both system sizes. c &d The Boltzmann energy distributions, for
$N\in[66,264]$ compared with the distribution of subsystems of size $N/2$,
suggest the subsystem does satisfactorily characterize the larger system.
Since the $N=66$ system accesses much lower energies than the $N=33$ system,
deviations for this comparison have to occur in the low-energy region (as in
c).

The energy distribution— A meticulous analysis of the energy distribution
across various system sizes, temperatures, and energy levels yields further
insights and understanding. First, we report the Boltzmann distribution
$P(E_{IS},T)$ of IS energies, determined at different temperatures, as well as
its first moment for $N=33$ and $N=66$, see Fig. 2(a&b). For $T>0.1$ we
observe a linear behavior of the average IS energy with respect to the inverse
temperature. This is known to reflect a Gaussian distribution of IS Heuer and
Büchner (2000). For all systems sizes, around $T\lesssim 0.1$, the slope
starts to decrease for decreasing temperature (also for $N>66$, see SI (IV)).
This is caused by the depletion of low-energy states, as seen from the
emerging asymmetric distribution. Importantly, the occurrence of the FSC
manifests at exactly the temperature where the non-Gaussian regime of the PEL
becomes relevant, independent of system size. This suggests a very strong
relationship between both observations.

At still significantly lower temperatures, the bottom of the landscape is
found, namely for $N=33$ at $T\approx 0.06$ (close to the calorimetric
$T_{g}$), and for $N=66$ at much lower temperatures (around $T_{g}/2$). Its
identification required these extensive simulations. The claim that the bottom
of the PEL is indeed found is further supported by statistical arguments (see
SI (VII) and SI (VIII)), stating that for $N=66$, all the 16 lowest energies
are found for both the lowest two temperatures. The IS at the bottom has lower
energy for $N=66$. Most importantly, for the first time, the whole IS-
distribution is now accessible for system size $N_{c}$.

One might conclude that the different temperatures, required for the
identification of the bottom of the PEL, is related to the different energies
(per particle) at the bottom. However, as argued in SI (VI), this effect would
also be present if the superposition hypothesis is strictly valid, implying
identical bottom energies. This is also the reason why it is not possible to
find the bottom of the PEL for $N=132$ (estimated required temperature: 0.02).

The Boltzmann distributions $P(E_{IS},T)$ can be further explored to check the
superposition hypothesis, see Fig. 2(c&d). Here we show the actual
distributions at different sizes $N$ and their prediction via the
superposition hypothesis, based on the distribution at size $N/2$ (see
Methods). The temperatures are chosen far below $T_{g}$. For the transition
from $N=33$ to $N=66$, we see a very good agreement for higher energies,
whereas deviations are visible on the low-energy side. This was expected due
to the shifted bottom of the PEL, as observed in Fig. 2(a&b). Remarkably, we
see a perfect superposition for the transitions $N=66\rightarrow N=132$ and
$N=132\rightarrow N=264$ for the whole energy range, supporting the choice
$N_{c}=66$ again.

Figure 3: From the energy distribution to the configurational entropy. a The
energy distribution of IS for various system sizes across the ultra-stable
limit. The dashed-colored lines show a Gaussian fit from the higher energy
supercooled regime, revealing faster depletion of states than predicted for a
pure Gaussian. The dashed-black line depicts the distribution constructed from
the respective half-sized subsystem. The vertical-colored lines represent the
lowest energy for IS for N=33 and 66. b For $N=33,66$, the configurational
entropy is calculated from the complete IS energy distribution, including the
information of the basin’s volumes from the Hessian (see SI (IX)). For $N\geq
132$, the entropy curves are shifted to agree with N=66 at $T_{MCT}=0.105$.
The inset compares the ‘shifted’ entropy for $N=33$ to $N=66$ at $T_{MCT}$ and
its range of (nearly) perfect agreement. The vertical dashed line represents
$T_{g}/2$.

A more profound picture of the PEL can be obtained from knowledge of the IS-
density, written as $G(E_{IS})$. It is defined such that the equilibrium
probability is given by $P(E_{IS},T)\propto G(E_{IS})\exp(-\beta E_{IS})$ with
$\beta=1/{k_{B}T}$. Thus, $G(E_{IS})$ can be obtained by reweighting
Saksaengwijit  _et al._ (2004). Strictly speaking, the IS-density additionally
contains the information about the effective volume of the basin, i.e.,
$G(E_{IS})=\sum_{i}Y_{i}\delta(E-E_{i,IS})$ where the factor $Y_{i}$ expresses
the volume contribution. Thus, one should rather speak of a volume-weighted
IS-density. In harmonic approximation, $Y_{i}$ is related to its curvature at
the minimum and is estimated from the eigenvalues of the Hessian Sciortino
(2005) (for more details and the negligible impact of additional anharmonic
corrections, see SI (IX)). In Fig. 3(a) $G(E_{IS})$ is displayed in a
normalized continuous representation. One can explicitly observe the depletion
of IS at low energies, giving rise to a transition from Gaussian to non-
Gaussian density, as already reflected in Fig. 2. Consistent with Fig. 2(c&d),
the superposition hypothesis is also valid for $G(E_{IS})$ for the whole range
of energies. The only deviations are seen for the transition $N=33\rightarrow
N=66$ close to the bottom of the PEL. It cannot be excluded that the
superposition hypothesis is also violated for larger $N$ in the range of very
low energies. However, these differences would, if at all, only be visible for
$T<T_{g}/2$, and thus not relevant for experiments.

Note that the non-Gaussian energy regime becomes much more dominant for large
system sizes simply because of the superposition effect, i.e., the underlying
convolution of the energy distributions. In contrast, for $N=33$, the non-
Gaussian regime is quite small. Note that even the $N=33$ system can be viewed
as a superposition of even smaller entities, as evidenced by much smaller
values of the participation ratio (see SI (V)). Thus, the non-Gaussian regime
even for the $N=33$ system can be seen as a consequence of the finite energy
range of the PEL. We have reported that the onset temperature of the FSC
appears when the depletion of low-energy states and thus non-Gaussian effects
become relevant. Based on our observations, this suggests that the FSC is
influenced by the presence of the bottom of the PEL of small-scale units. As
already discussed exclusively for network formers, the presence of a low-
energy bottom generally leads to low-T Arrhenius behavior Saksaengwijit  _et
al._ (2004).

Of utmost relevance for the discussion of the properties of glass-forming
systems is the configurational entropy, defined via the Shannon description as
$S_{c}(T)=-\sum_{i}P(E_{i,IS},T)\ln P(E_{i,IS},T)$ where the sum is over all
the IS. Since $P(E_{i,IS},T)\propto Y_{i}\exp(-\beta E_{i,IS})$, the
configurational entropy can be rewritten as (see SI (IX)) $S_{c}=-\langle
logY_{i}\rangle+\beta\langle E_{i,IS}\rangle+\ln(\sum_{i}Y_{i}\exp(-\beta
E_{i,IS}))$. Whereas the first two terms are thermodynamic averages, the last
term requires knowledge about the absolute number of IS, which is available
for $N=33$ and $N=N_{c}=66$. Thus, the complete knowledge of the PEL for both
systems unveils a unique prospect for evaluating the configurational entropy
for all temperature ranges without resorting to the liquid entropy. For larger
$N$, where the ground state is unknown, the entropy can be directly obtained
apart from a temperature-independent shift factor and can be calculated for
the temperature range, for which equilibrium simulations are available.

Fig. 3b shows the configurational entropy for various system sizes. For
$N=N_{c}=66$ its temperature dependence, involving a transition from positive
to negative curvature at low but finite temperature (see also SI (IX)),
reflects the low-energy state depletion, including the cutoff. For the present
model glass, it stems a unique footing to quantify Kauzmann’s picture of
entropy crisis Kauzmann (1948) and provides a direct view of the thermodynamic
behavior Angell (2008) such as the nature of the strong decrease of entropy at
low temperatures Berthier  _et al._ (2019) due to the direct relation to
individual states.

One might have expected that the entropy for $N=33$ only differs for low
temperatures because of the lack of low-energy states as seen, e.g., in
Fig.3a. Indeed, for $T>0.04$, the shape of the entropy curve is identical.
However, the absolute value differs significantly as seen in Fig. 3b, i.e.,
the number of IS for all energies is by a factor
$exp[(s_{c}(N=66)-s_{c}(N=33))\cdot 66]\approx 10^{7}$ higher for $N=66$ as
expected from the factorization properties. This reveals a surprising entropy
variation, although other thermodynamic properties (i.e., average energy) or
the diffusivity are qualitatively similar for $N=33$ and $N=66$ in the whole
temperature range. A simple suggestion for this thermodynamic finite-size
effect of the absolute number of states is provided in the SI (IX), suggesting
a dominant $\ln(N)/N$-dependence of the configurational entropy.

Figure 4: Structural features across the PEL. a & b The lowest energy sample
was explored during a wide range of temperatures. For N=33(66), the lowest
energy configuration is shown for T=0.0385(0.0270) and T=0.0357(0.0263) with
filled and dashed styles, respectively. The color scheme represents the
diameter of the particle. c The disorder depicts a similar temperature
dependence as the energy. d The parametric plot between the structural
disorder and the potential energy provides a comprehensive approach to
studying the disorder. e The mean diameter for the neighbors of the particle
$i$ with diameter $\sigma_{i}\geq 1.4$ for a range of temperatures across
$T_{g}$. The distinct nature of N=33 suggests that the variety of possible
configurations is not large enough for the smallest system.

Lowest state and structure— To quantitatively characterize the low-energy
state, Fig. 4a&b presents the ground states observed at various temperatures
for $N=33$ and $N=66$, respectively. Evidently, the packaging suggests a
strong anti-correlation between larger and smaller particles.

Inspired by the network-forming liquids, the structural features of the IS are
expected to characterize the system’s disorder, depending on their energy Roy
and Heuer (2019). One can quantify the disorder ($\Theta$) as the deviation
from steric arrangements with the Voronoi neighbors (see methods) Tong and
Tanaka (2018). Upon cooling, already above $T_{g}$ one can see the approach
towards a constant disorder, see Fig. 4c. The temperature dependence is very
similar to that of the energy (Fig. 2). Indeed, this is a direct consequence
of the linear relationship between disorder and energy, see Fig. 4d,
complemented by the presence of the low-energy cutoff. Since structural
features are also accessible experimentally, a signature of the bottom of the
PEL might thus be observable. This structural feature does not display any
finite-size effects, even not for $N=33$.

Another structural property is a possible correlation of the diameters of
adjacent particles. For this purpose, we report the mean diameter of the
Voronoi neighbors for particles with $\sigma\geq 1.4$ for system sizes between
$N=33$ and $N=1056$, see Fig.4e. For lower energies, a stronger anti-
correlation between adjacent diameters is observed. The deviations for $N=33$
in the low-energy limit possibly reflect the smaller number of possible
arrangements. The stronger anti-correlations enable a more efficient packing
of the particles and thus the accessibility of lower energies. The analog of
this statistical feature is the disappearance of well-defined defects in the
case of network formers when reaching the low-energy limit. Such an order
parameter does not exist for the present non-network system, yielding a
contrasting perspective on both types of glass formers.

Experimental system— The study by Hecksher et al. Hecksher  _et al._ (2008)
suggests a possibility of VFT-like divergence not being a complete description
of the glassy dynamics, which is further corroborated by Mallamace et al.
Mallamace  _et al._ (2010) for a variety of experimental systems suggesting
strong behavior in the low-T limit in supercooled liquids. In the context of
FSC, we analyze a range of fragile experimental glassy systems and observe
that the low-T dynamics(below $T_{onset}$) possess 4 to 12 decades to be
described with the Arrhenius relation (see SI (XI)). The exact reasoning
behind such crossover is not clearly understood. Whereas, e.g., Mallamace et
al. Mallamace  _et al._ (2010) proposed dynamics heterogeneity’s role with a
crossover around a viscosity $10^{3}$ Poise, in other approaches the crossover
is predicted either in the vicinity of $T_{MCT}$ Casalini  _et al._ (2003a,
b), or via excluded volume based theory close to $T_{g}$ Rah (2007). The
state-of-the-art measurements for the 20 million-old amber Zhao  _et al._
(2013), metallic glass Zhang  _et al._ (2010), and viscosities beyond
laboratory glass transition Pogna  _et al._ (2015); Jadhao and Robbins (2017)
highlight the lack of understanding of the temperature dependence of dynamics.
Based on the arguments for such a wide range of systems, our study proposes a
general framework to comprehend the emergence of Arrhenius behavior at low-T
within the context of the PEL.

Discussion— The study accentuates a central issue regarding the relaxation
dynamics of supercooled liquids near the experimental glass transition and
beyond. Combining the idea of the study of relatively small systems with the
swap simulation, we can investigate the energy landscape for an equilibrium
relaxation time span of astronomical scale, thus achieving a view of the
complete energy landscape. The non-Gaussian nature of PEL appears to be a
general characteristic of glasses and not a special feature of network
formers. Consequently, the crossover to Arrhenius behavior can be seen as a
general glassy feature, avoiding any finite-T divergence. It also
contextualizes the recent observations of low-T Arrhenius behavior, reported
in other independent simulation studies Berthier and Ediger (2020); Guiselin
_et al._ (2022); Das and Sastry (2022). The results presented in this work
provide a broad perspective on glass physics, extending beyond experimental
time scales.

The novel finding of crossover for the PEL and dynamics may represent a
mechanism that can follow Kauzmann’s scenario for the glassy state, providing
valuable insights into the fundamental nature of glass. The corresponding
depletion of low-energy states is mainly a consequence of the energy cutoff,
seen for the smallest studied system. Additionally, this supports the idea
that systems cannot be aged without crystallizing to arbitrarily deep
energies; such observation could be crucial for modeling energy landscapes and
mesoscale modeling of glass. Finally, our investigation paves the way for new
studies for experimental and model glass formers within the context of PEL and
understanding the dynamics around and below $T_{g}$ and possibly also about
the impact of the dimensionality on the low-energy and thus low-temperature
properties.

Glass forming in-silico model— We study a discrete-polydisperse mixture of $N$
purely repulsive particles, consisting of $M$ species of equal mass $m$,
pairwise interacting via power-law potential. The purely repulsive interaction
potential for particles $i$ and $j$, at the distance $r_{ij}$ is defined as

$e(r_{ij})=\epsilon\left(\frac{\sigma_{ij}}{r_{ij}}\right)^{12}+\sum_{k=\\{0,2\\}}c_{2k}\left(\frac{r_{ij}}{\sigma_{ij}}\right)^{2k},$
(1)

with interaction range $r_{ij}/\sigma_{ij}(=r_{c,ij})\leq 1.25$. The fourth-
order polynomial ensures that the potential, force, and second-order
derivatives are continuous at the cutoff. We use $\epsilon$ as the unit of
energy and temperature and average particle diameter $\left<\sigma\right>$ as
the unit of length. The molecular dynamics time is measured in reduced units
of $\sqrt{m\left<\sigma\right>^{2}/\epsilon}$, while the Monte Carlo time is
reported as the MC steps performed. We study the system sizes ranging from
N=33 to 1056, consists of $M=10$ species, following the diameter-distribution
$P(\sigma)\propto\sigma^{-3}$ Ninarello  _et al._ (2017), with
$\sigma\in[0.73,1.62]$. Such a strategy for modeling allows for rigorous
comparison of various system sizes while keeping the size-disorder identical
Parmar  _et al._ (2020)(see SI(I)). In addition to polydispersity, to ensure
good glass-forming ability and avoid crystallization, we utilize non-addictive
interactions in diameter:
$\sigma_{ij}=0.5{(\sigma_{i}+\sigma_{j})}(1-0.2|\sigma_{i}-\sigma_{j}|)$. A
system with N particles inside a square cell of area V is simulated under
periodic boundary conditions with a number density of $\rho=N/V=1$.

Sample preparation and dynamics— We perform swap Monte Carlo simulations to
achieve equilibrated samples in a deeply supercooled regime. One swap Monte
Carlo step is consistent with 80% translation moves(local), and the rest are
swap (nonlocal) Ninarello  _et al._ (2017), while the Metropolis rule ensures
a detailed balance for equilibrium sampling. Detailed examinations of the
aging and relaxation ensure that samples reside in perfect equilibrium (see
SI(II)).

For local dynamics, we perform molecular dynamics simulations. Well-
equilibrated samples from swap Monte Carlo simulations are subjected to
constant temperature (NVT) simulations using the Nośe-Hoover thermostat in the
temperature range of $[0.5000,0.0690]$. Simulations have been carried out for
the system sizes $N\in[33,1056]$, with the sampling size $n_{s}\in[1056,100]$,
for averaging.

Dynamical observable and time scales— The long-wavelength Mermin-Wagner
fluctuations can mask the particle motions relevant to the relaxation,
nonetheless irrelevant to the structure and cage-escape process Illing  _et
al._ (2017). We characterize the dynamics using cage relative positions
${\bf{\bar{r}}}_{i}(t)={\bf{r}}_{i}(t)-\sum_{j}{\bf{r}}_{j}(t)/\mathcal{N}_{i}$,
where the summation goes over all the neighbors, $N_{i}$, of the particle $i$,
defined at the initial time t=0. The neighbors are defined by the first
coordination shell of the pair-correlation function, a pair $\left<ij\right>$
is considered a neighbor if $r_{ij}/\sigma_{ij}\leq 1.33$. To characterize the
dynamics, we estimate the diffusion coefficient ($D$) from the mean squared
displacement ($\Delta(t)$) of the system.

$\Delta(t)=\frac{1}{N}\left<\sum_{i=1,N}({\bf{\bar{r}}}_{i}(t)-{\bf{\bar{r}}}_{i}(0))^{2}\right>,$
(2)

where the angular bracket represents the averaging over the independent
samples. For dynamics, we use diffusivity($D$), computed from the long time
limit of the MSD $\Delta(t)=\lim_{t\rightarrow\infty}4Dt$.

At the high temperature, the system acts as a simple liquid, and the dynamics
expressed in the Arrhenius manner,

$D^{-1}(T)=D_{\infty}^{-1}\exp(-\beta A),$ (3)

where $D_{\infty}$ and $A$ are infinite temperature diffusion and the
activation energy, respectively. For the high temperature, a departure from
Arrhenius nature marks the onset temperature, $T_{onset}$=0.25.

The mode coupling temperature, $T_{MCT}=0.105$, is marked by fitting the
power-law form- $D^{-1}(T)\propto|T-T_{MCT}|^{-\gamma}$, with exponent
$\gamma=2.5$, in the first two decades of supercooled dynamics. The system’s
dynamical nature is also studied with the parabolic and Vogel-Fulcher-
Tammann’s depiction in the first three decades of the glassy dynamics. Such
projection describes the growth of the dynamics only within the first three
decades of the supercooling; which systemically overestimates at the low
temperature. While the initial supercooled liquid is fragile, the
four(five)-decades of diffusion(relaxation) at extreme super-cooling
(T$\lesssim$0.1) adequately describes the strong behavior(further details see
SI(VI)).

The experimental glass transition temperature, $T_{g}=0.0555$, corresponds to
the first 12-decade of glassy relaxation. The present study advances the
equilibrium description of dynamics nearer to the glass transition till
$1.24T_{g}$. The seven decades of glassy relaxation provide the mili-second
narration of the dynamics approaching the glass transition. With the advantage
of swap-MC, studying the potential energy landscape extends to much lower
temperatures exceeding the billions of years for relaxation in the ultrastable
limit.

The potential energy landscape(PEL)— To accurately assess the PEL, we analyze
the energy of the inherent structures ($E_{IS}$). An equilibrated sample is
subjected to the conjugate gradient minimization up to machine precision using
LAMMPS Plimpton (1995). The complete energy distribution is constructed for
the system sizes $N\in[33,1056]$ with a sample size of $1\times 10^{7}$ to
$5\times 10^{5}$.

The probability distribution for inherent structure $P(E_{IS},T)$ can be re-
weighted to estimate the energy distribution of the PEL from

$P_{W}(E_{IS},T)\propto P(E_{IS},T)\exp(\beta E_{IS}).$ (4)

For reliable construction of the probability distribution($G(E_{IS})$) of PEL,
the reweighted distributions glued with the average value for the range
$\left<E_{IS}\right>\pm 1.2\sigma$, with $\sigma$ is the variance of the
$E_{IS}$. We focus on the supercooled regime below the onset temperature, and
the Gaussian nature of the PEL can be marked in $T\in[0.105,0.25]$. The
equilibrium sampling is ensured with stringent assessments- (i) the particle
size autocorrelation similar to the structural correlation, (ii) the system is
not aging, (iii) and the energy-wise ratio of the weighted distribution for
different temperatures ($T_{i}$,$T_{j}$) remains constant, i.e.,
$P_{W}(E_{IS},T_{i})/P_{W}(E_{IS},T_{j})$=constant, see also SI, section II.

The study of the smaller system is also extended to describe the bulk system.
A bulk system of N particles and energy $E_{IS}$, can be considered as two
coupled subsystems i(and j) of a microcanonical ensemble of size
$N_{i}$($N_{j}$) and energy $E_{i,IS}$($E_{j,IS}$). Therefore, the total
probability for the compound system can be expressed as

$P(E_{IS};N)=\int dE_{i}P(E_{i,IS},N_{i})\cdot P(E-E_{i,IS},N-N_{i}),$ (5)

where the integration is performed over all the energies $E_{i,IS}$. We
strictly compare small and bulk systems past the glass transition temperature
with such a strategy.

The disorder— The structural disorder, $\Theta$, measures the deviation of
compactness from the steric arrangements Tong and Tanaka (2018). For the
particle ‘i’ with the Voronoi neighbors $\langle j,k\rangle$, the disorder
highlights the angular deviation with the central particle $i$ in the
following manner

$\Theta_{i}=\frac{1}{\mathcal{N}_{i}}\sum_{\left<j,k\right>}|\theta_{j,k}-\theta^{o}_{j,k}|,$
(6)

where $\theta_{j,k}$ and $\theta^{o}_{j,k}$ represents the angle between i and
$\langle j,k\rangle$ while the position of particle for inherent structure and
all the particles in contact with each other, respectively. The
$\mathcal{N}_{i}$ is the number of pairs of neighbors next to each other. For
our model, the system’s mean disorder captures the pairs-wise energetic
contribution and an excellent structural descriptor similar to energy.

## References

  * Anderson (1995) P. W. Anderson, Science 267, 1615 (1995).
  * Debenedetti and Stillinger (2001) P. G. Debenedetti and F. H. Stillinger, Nature 410, 259 (2001).
  * Schmidtke  _et al._ (2012) B. Schmidtke, N. Petzold, R. Kahlau, M. Hofmann, and E. Rössler, Physical Review E 86, 041507 (2012).
  * Angell (1991a) C. Angell, Journal of Non-Crystalline Solids 131, 13 (1991a).
  * Saika-Voivod  _et al._ (2001) I. Saika-Voivod, P. H. Poole, and F. Sciortino, Nature 412, 514 (2001).
  * Saksaengwijit  _et al._ (2004) A. Saksaengwijit, J. Reinisch, and A. Heuer, Physical review letters 93, 235701 (2004).
  * Adam and Gibbs (1965) G. Adam and J. H. Gibbs, The journal of chemical physics 43, 139 (1965).
  * Lubchenko and Wolynes (2007) V. Lubchenko and P. G. Wolynes, Annu. Rev. Phys. Chem. 58, 235 (2007).
  * Kirkpatrick  _et al._ (1989) T. R. Kirkpatrick, D. Thirumalai, and P. G. Wolynes, Physical Review A 40, 1045 (1989).
  * Goldstein (1969) M. Goldstein, The Journal of Chemical Physics 51, 3728 (1969).
  * Angell (1991b) C. Angell, in  _Hydrogen-Bonded Liquids_ (Springer, 1991) pp. 59–79.
  * Sastry (2001) S. Sastry, Nature 409, 164 (2001).
  * Sciortino (2005) F. Sciortino, Journal of Statistical Mechanics: Theory and Experiment 2005, P05015 (2005).
  * La Nave  _et al._ (2006) E. La Nave, S. Sastry, and F. Sciortino, Physical Review E 74, 050501 (2006).
  * Heuer (2008) A. Heuer, Journal of Physics: Condensed Matter 20, 373101 (2008).
  * Baity-Jesi  _et al._ (2021) M. Baity-Jesi, G. Biroli, and D. R. Reichman, The European Physical Journal E 44, 1 (2021).
  * Berthier and Biroli (2011) L. Berthier and G. Biroli, Reviews of modern physics 83, 587 (2011).
  * Doliwa and Heuer (2003) B. Doliwa and A. Heuer, Journal of Physics: Condensed Matter 15, S849 (2003).
  * Kob and Andersen (1994) W. Kob and H. C. Andersen, Physical review letters 73, 1376 (1994).
  * Roy and Heuer (2022) P. K. Roy and A. Heuer, The Journal of Chemical Physics 157, 174506 (2022).
  * Yu  _et al._ (2022) Z. Yu, D. Morgan, M. Ediger, and B. Wang, Physical Review Letters 129, 018003 (2022).
  * Ninarello  _et al._ (2017) A. Ninarello, L. Berthier, and D. Coslovich, Physical Review X 7, 021039 (2017).
  * Parmar  _et al._ (2020) A. D. Parmar, M. Ozawa, and L. Berthier, Physical Review Letters 125, 085505 (2020).
  * Plimpton (1995) S. Plimpton, Journal of computational physics 117, 1 (1995).
  * Rehwald  _et al._ (2010) C. Rehwald, O. Rubner, and A. Heuer, Physical review letters 105, 117801 (2010).
  * Vogel  _et al._ (2004) M. Vogel, B. Doliwa, A. Heuer, and S. C. Glotzer, The Journal of chemical physics 120, 4404 (2004).
  * Heuer and Büchner (2000) A. Heuer and S. Büchner, Journal of Physics: Condensed Matter 12, 6535 (2000).
  * Kauzmann (1948) W. Kauzmann, Chemical reviews 43, 219 (1948).
  * Angell (2008) C. Angell, MRS bulletin 33, 544 (2008).
  * Berthier  _et al._ (2019) L. Berthier, P. Charbonneau, A. Ninarello, M. Ozawa, and S. Yaida, Nature communications 10, 1508 (2019).
  * Roy and Heuer (2019) P. K. Roy and A. Heuer, Physical Review Letters 122, 016104 (2019).
  * Tong and Tanaka (2018) H. Tong and H. Tanaka, Physical Review X 8, 011041 (2018).
  * Hecksher  _et al._ (2008) T. Hecksher, A. I. Nielsen, N. B. Olsen, and J. C. Dyre, Nature Physics 4, 737 (2008).
  * Mallamace  _et al._ (2010) F. Mallamace, C. Branca, C. Corsaro, N. Leone, J. Spooren, S.-H. Chen, and H. E. Stanley, Proceedings of the National Academy of Sciences 107, 22457 (2010).
  * Casalini  _et al._ (2003a) R. Casalini, M. Paluch, and C. M. Roland, Journal of Physics: Condensed Matter 15, S859 (2003a).
  * Casalini  _et al._ (2003b) R. Casalini, M. Paluch, and C. M. Roland, The Journal of chemical physics 118, 5701 (2003b).
  * Rah (2007) K. Rah, Physica A: Statistical Mechanics and its Applications 378, 167 (2007).
  * Zhao  _et al._ (2013) J. Zhao, S. L. Simon, and G. B. McKenna, Nature communications 4, 1 (2013).
  * Zhang  _et al._ (2010) C. Zhang, L. Hu, Y. Yue, and J. C. Mauro, The Journal of chemical physics 133, 014508 (2010).
  * Pogna  _et al._ (2015) E. A. A. Pogna, C. Rodríguez-Tinoco, G. Cerullo, C. Ferrante, J. Rodríguez-Viejo, and T. Scopigno, Proceedings of the National Academy of Sciences 112, 2331 (2015).
  * Jadhao and Robbins (2017) V. Jadhao and M. O. Robbins, Proceedings of the National Academy of Sciences 114, 7952 (2017).
  * Berthier and Ediger (2020) L. Berthier and M. D. Ediger, The Journal of Chemical Physics 153, 044501 (2020).
  * Guiselin  _et al._ (2022) B. Guiselin, C. Scalliet, and L. Berthier, Nature Physics 18, 468 (2022).
  * Das and Sastry (2022) P. Das and S. Sastry, Journal of Non-Crystalline Solids: X 14, 100098 (2022).
  * Illing  _et al._ (2017) B. Illing, S. Fritschi, H. Kaiser, C. L. Klix, G. Maret, and P. Keim, Proceedings of the National Academy of Sciences 114, 1856 (2017).