unitsize(1cm); pair[] points = (0, 0), (0, 1), (1, 1), (1, 0); int[][] conns =
0, 1, 1, 2, 2, 3, 3, 0, 0, 2; int[][] connsT = 0, 1, 0, 2, 0, 3, 1, 2, 2, 3;
pair[] dirs = SW, NW, NE, SE; pair[] es = (1/4, 3/4), (3/4, 1/4), (1/2, 5/3);
path[] edges; for(int i=0; i¡conns.length; ++i) edges.push(points[conns[i][0]]
– points[conns[i][1]]); path[] edgesT; for(int i=0; i¡connsT.length; ++i)
edgesT.push(points[connsT[i][0]] – points[connsT[i][1]]); void drawgraph(pen
pe=defaultpen) for(int i=0; i¡edges.length; ++i) draw(edges[i], p=pe); for(int
i=0; i¡points.length; ++i) dot(points[i], p=(i¿0?pe:pe+red)+6); void
drawgrapharrows(pen pe=defaultpen) for(int i=0; i¡edgesT.length; ++i)
draw(edgesT[i], p=pe, MidArrow); for(int i=0; i¡points.length; ++i)
dot(points[i], p=(i¿0?pe:pe+red)+6); void labelgraph() for(int i=0;
i¡points.length; ++i) label(points[i], ”$v_{"}+(string)i+"$”, dirs[i]); void
labelgraph(string[] values) for(int i=0; i¡points.length; ++i)
label(points[i], ”$v_{"}+(string)i+","+values[i]+"$”, dirs[i]); void
drawdualverts(pen pe, real size=6) for(int i=0; i¡es.length; ++i) dot(es[i],
p=pe+size); void drawdualedges(pen pe) draw(es[0] – es[1], p=pe); draw(es[0]
.. midpoint(edges[0]) .. (-.5, 1.3) .. es[2], p=pe); draw(es[0] ..
midpoint(edges[1]) .. es[2], p=pe); draw(es[1] .. midpoint(edges[2]) .. (1.5,
1) .. es[2], p=pe); draw(es[1] .. midpoint(edges[3]) .. (-.5, -.5) .. (-1, 1)
.. (0, 2) .. es[2], p=pe); void drawdual(pen pe=defaultpen) drawdualedges(pe);
drawdualverts(pe); pair cis(real angle) return (cos(angle), sin(angle)); pair
a=cis(0.5*pi); pair b=cis(0.1*pi); pair c=cis(1.7*pi); pair d=cis(1.3*pi);
pair e=cis(0.9*pi); void drawfive(real size=3) dot(a, black+size); dot(b,
black+size); dot(c, black+size); dot(d, black+size); dot(e, black+size); pair
u=(0, 0); pair v=(0, 1); pair w=(0, 2); pair x=(1, 0); pair y=(1, 1); pair
z=(1, 2); void drawsix(real size=3) dot(u, black+size); dot(v, black+size);
dot(w, black+size); dot(x, black+size); dot(y, black+size); dot(z,
black+size);

# Cycle and circuit chip-firing on graphs

Anton Dochtermann, Eli Meyers, Rahgav Samavedan, Alex Yi

###### Abstract.

Motivated by the notion of chip-firing on the dual graph of a planar graph, we
consider ‘cycle chip-firing’ on an arbitrary graph $G$. The chip-firing rule
is governed by $\operatorname{\mathcal{L}}^{*}(G)$, the dual Laplacian of $G$
determined by choosing a basis for the lattice of integral flows on $G$. We
show that any graph admits such a basis so that
$\operatorname{\mathcal{L}}^{*}(G)$ is an $M$-matrix, leading to a firing rule
on these basis elements that is avalanche finite. This follows from a more
general result on bases of integral lattices that may be of independent
interest. Our results provide a notion of $z$-superstable cycle configurations
that are in bijection with the set of spanning trees of $G$. We show that for
planar graphs, as well as for the graphs $K_{5}$ and $K_{3,3}$, one can find
such an M-basis consisting of circuits of the underlying graph. We consider
the question for arbitrary graphs and address some open questions.

##  1\. Introduction

The classical theory of _chip-firing_ involves a simple game played on the
vertices of a graph. Versions of these games play a role in physics in the
context of self-organized criticality [2, 10], and also have been studied for
their combinatorial properties and underlying algebraic structure [6, 5]. More
recently chip-firing has seen connections to various other disciplines
including algebraic geometry [3]. We refer to [9] and [14] for good
introductions to the subject.

In one variant of the game, a vertex is designated to be the _root_ , and a
nonnegative configuration of chips ${\bf
c}=(c_{1},c_{2},\dots,c_{n})\in{\mathbb{N}}^{n}$ is placed on the non-root
vertices. If the number of chips on a non-root vertex $i$ is greater than its
degree, the vertex $i$ can _fire_ , in which case the vertex passes chips to
each of its neighbors (one for each edge connected it to $i$).

A configuration is _stable_ if no non-root vertex can fire, and is _critical_
if it is stable and recurrent (see below for proper definitions). The set of
critical configurations on a graph $G$ forms an abelian group $\kappa(G)$
under coordinate-wise addition and _stabilization_ , known as the _critical
group_ of $G$. One can show that $\kappa(G)$ does not depend on the choice of
root vertex in the graph $G$. Although easy to define, $\kappa(G)$ is a subtle
invariant of the graph $G$ and for instance there has been recent interest in
determining the distribution of critical groups among random graphs [18]. It
is also well-known that the set of critical configurations on $G$ has a simple
duality with the so-called _superstable_ configurations of that graph (known
in some contexts as _$G$ -parking functions_).

The dynamics of chip-firing on a graph $G$ can be encoded by the reduced
Laplacian matrix
$\tilde{\operatorname{\mathcal{L}}}=\tilde{\operatorname{\mathcal{L}}}(G)$,
where the effect of firing a vertex $i$ can be seen to correspond to
subtracting the $i$th row of $\tilde{\operatorname{\mathcal{L}}}$. If we think
of
$\tilde{\operatorname{\mathcal{L}}}:\mathbb{Z}^{n-1}\rightarrow\mathbb{Z}^{n-1}$
as a linear operator, one can prove that the critical group $\kappa(G)$ is
recovered as the cokernel of the reduced Laplacian:

$\kappa(G)\cong\mathbb{Z}^{n-1}/\operatorname{im}\tilde{\operatorname{\mathcal{L}}}.$

By the well-known matrix tree theorem, we have that $|\kappa(G)|=|\tau(G)|$,
where $\tau(G)$ is the set of spanning trees of $G$ (recall that $G$ is
assumed to be connected). In particular the set of superstable configurations
is in bijection with $\tau(G)$. Merino [15] has shown that this bijection can
be chosen to preserve degree/activity, which connects the theory of chip-
firing on a graph $G$ to the properties of the Tutte polynomial $T_{G}(x,y)$.

A more general notion of ‘chip-firing’ can be encoded in any $n\times n$
integer matrix $L$ that has positive entries on the diagonal and non-positive
entries on the off-diagonal. Chip-firing at a site $i$ is then defined by
subtracting $L^{T}{\bf e}_{i}$, the $i$th row of $L$, from the current
configuration. One is interested in such matrices $L$ with the property that
all initial configurations eventually stabilize by repeated firings determined
by the matrix $L$. This ‘avalanche finite’ property is equivalent to $L$ being
a so-called $M$-matrix. Building on the work of Baker and Shokrieh from [4],
Guzmán and Klivans have shown [12] that a certain ‘energy minimizing’
perspective on $M$-matrices leads to good notions of so-called $z$-superstable
configurations, where a notion of firing multisets of vertices simultaneously
is relevant. In particular, the set of $z$-superstable configurations are in
bijection with the equivalence classes defined by the matrix $L$.

Returning to the world of graphs, we see that if a graph $G$ is embedded in
the plane with dual graph $G^{*}$, then a configuration on the vertices of
$G^{*}$ can naturally be thought of as a ‘circuit configuration’ on $G$. After
fixing an orientation on $G$ and an orientation of the plane, firing a vertex
in $G^{*}$ then corresponds to firing a directed circuit in $G$. A circuit $C$
sends chips to any circuit $C^{\prime}$ that it shares at least one edge with,
with the total number of chips sent determined by the number of edges that are
oriented consistently versus inconsistently. In addition, we see that the
number of ‘superstable circuit configurations’ (i.e. the superstable
configurations of $G^{*}$) are in bijection with the number of spanning trees
of $G$, where the number of trees with $k$ _internally_ passive edges is equal
to the number of superstable circuit configurations of degree $k$.

What happens if $G$ is not planar? In this paper, we study the notion of
circuit (and, more generally, cycle) chip-firing for general (not necessarily
planar) graphs. We adopt the perspective that traditional chip-firing on $G$
can be viewed as choosing a particular basis for the lattice of integral cuts
$\operatorname{\mathcal{C}}(G)$, given by all vertex cuts determined by non-
root vertices. Choosing this basis for $\operatorname{\mathcal{C}}(G)$ then
defines the reduced incidence matrix $\tilde{\delta}$ and
$\tilde{\delta}\circ\tilde{\delta}^{T}=\tilde{\operatorname{\mathcal{L}}}$
recovers the reduced Laplacian, which we have seen governs the usual rules of
chip-firing.

In the dual case we want to choose a basis for
$\operatorname{\mathcal{F}}(G)=\ker\delta\cap\mathbb{Z}^{m}$, the integral
kernel of the (signed) incidence matrix $\delta$ (sometimes called the
_lattice of integral flows_ for the graph $G$ [1]). Choosing such a basis
$\iota=\\{{\bf f}_{1}\dots,{\bf f}_{g}\\}$ gives rise to a dual Laplacian
$\iota^{T}\iota=\operatorname{\mathcal{L}}^{*}:\mathbb{Z}^{g}\rightarrow\mathbb{Z}^{g}$,
where $g=|E|-|V|+1$ is the _genus_ of the connected graph. We first observe
This dual Laplacian can be used to recover the critical group of the graph as
$\kappa(G)\cong\mathbb{Z}^{g}/\operatorname{im}\operatorname{\mathcal{L}}^{*}$
(see Proposition 3.3).

This result follows from established properties of integer lattices for
example from [1], and is also used in [13] (a fact we learned after preparing
this note). In Section 3 we give a mostly self-contained proof.

Although any basis for $\operatorname{\mathcal{F}}(G)$ gives rise to a dual
Laplacian $\operatorname{\mathcal{L}}^{*}$, this matrix will not in general
define a satisfactory avalanche-finite firing rule. We would like to interpret
$\operatorname{\mathcal{L}}^{*}$ as governing a chip-firing rule on the chosen
basis elements, with associated notions of $z$-superstable configurations and
‘good’ representatives for each chip-firing equivalence class. For this we
need the matrix $\operatorname{\mathcal{L}}^{*}$ to be an $M$-matrix, and our
next result shows that this can always be achieved.

######  Corollary 3.13.

Suppose $G$ is any graph. Then there exists a basis $\iota=\\{{\bf
f}_{1}\dots,{\bf f}_{g}\\}$ for its lattice of integral flows such that the
associated dual Laplacian $\operatorname{\mathcal{L}}^{*}$ is an $M$-matrix.

A basis for $\operatorname{\mathcal{F}}(G)$ satisfying the conditions in
Corollary 3.13 will be called an _cycle $M$-basis_ for $G$. Corollary 3.13
follows from a more general result (Theorem 3.11) regarding bases of any
integral lattice that may be of independent interest. Our proof of Theorem
3.11 is constructive, but our algorithm can lead to complicated linear
combinations of integer vectors with many large entries. To interpret our
basis in terms of ‘circuit chip-firing’, we would like the elements in our
cycle $M$-basis to encode directed circuits, that is to consist of $0,-1,1$
entries with the nonzero entries corresponding to the edges of some circuit in
the graph, with the sign determined by their orientation with respect to the
fixed orientation on $G$.

For the case of planar graphs we see that such a _circuit $M$-basis_ is given
by the bounded faces (appropriately oriented) in some embedding of $G$, see
Proposition 4.2. It is then natural to ask if any graph $G$ admits an
$M$-basis consisting of only circuits. This question is wide open but in
Section 4 we show that such bases do exist for the case of $G=K_{5}$ and
$G=K_{3,3}$ (found via a computer search).

######  Proposition 4.3.

The graphs $K_{5}$ and $K_{3,3}$ both admit cycle $M$-bases consisting of
circuits.

We also prove some results concerning the structure that circuit bases must
satisfy. For instance we show that if $K_{n}$ has a circuit $M$-basis then
some circuit must contain at least 4 edges, see Proposition 4.4. Similarly,
any circuit $M$-basis for $K_{m,n}$ with $m\geq n\geq 3$ must contain a
$6$-cycle.

A choice of cycle $M$-basis for a graph $G$ defines a dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ that governs an avalanche finite chip-firing
process. By results of [12] this then leads to notions of $z$-superstable
configurations, unique energy-minimizing elements in the equivalence classes
defined by $\operatorname{\mathcal{L}}^{*}$. By Proposition 3.3 we know the
size of this set of equivalence classes and hence we get the following.

######  Proposition 3.14.

Let $G$ be a connected graph, and let $\iota=\\{{\bf f}_{1}\dots,{\bf
f}_{g}\\}$ be a cycle $M$-basis with associated dual Laplacian
$\operatorname{\mathcal{L}}^{*}$. Then the number of $z$-superstable cycle
configurations of $G$ (with respect to $\operatorname{\mathcal{L}}^{*}$) is
given by $|\tau(G)|$, the number of spanning trees of $G$.

In the case of traditional chip-firing the degree sequence of superstable
configurations counts the number of spanning trees via the number of
externally passive edges (and in particular is independent of the choice of
sink vertex). If $G$ is a planar graph then a natural choice of circuit
$M$-basis given by a planar embedding leads to a dual Laplacian
${\mathcal{L}}^{*}$ whose $z$-superstable configurations count spanning trees
by internal passivity. In the general case the interpretation is not so clear.
In particular the degree sequence of $z$-superstables depends on the choice of
cycle $M$-basis for $\operatorname{\mathcal{F}}(G)$. One wonders if there is a
notion of ‘activity’ of a spanning tree that depends on the choice of basis
for $\operatorname{\mathcal{F}}(G)$ and that is reflected in the degree
sequence of the associated $z$-superstable circuit configurations . We leave
this for a future project.

The rest of the paper is organized as follows. In Section 1 we review some
basic notions from the classical theory of chip-firing on graphs as well as
the theory of generalized chip-firing defined by $M$-matrices. In Section 3 we
discuss our approach to circuit chip-firing in more detail and prove that any
graph admits a cycle $M$-basis. In Section 4 we show that $K_{5}$ and
$K_{3,3}$ both admit cycle $M$-bases consisting of circuits and discuss some
other examples. Here we also show that any circuit bases for nonplanar $K_{n}$
and $K_{m,n}$ must have large circuits. We end in Section 5 with some further
discussion and open questions.

Acknowledgements. We thank Johnny Guzmán for helpful comments. Most examples
were computed using SageMath [17], and we thank Dave Perkinson for sharing his
code regarding $z$-superstables. The work presented here was initiated as part
of Summer Honors Camp within the Mathworks Program at Texas State University.
We thank Mathworks and Texas State for the support and productive working
environment.

##  2\. Basics of classical chip-firing

Here we fix some notation and recall some results from the classical theory of
chip-firing on graphs, for the most part following conventions from [9].
Throughout the paper we let $G=(V,E)$ denote a connected undirected graph on
vertex set $V=\\{0,1,\dots,n\\}$, with the vertex $0$ designated as the root.
For $u\in V$ we let $\operatorname{deg}(u)$ denote the degree of the vertex,
so that

$deg(u)=|\\{w\in V:uw\in E\\}|.$

A configuration of chips on the (non-root) vertices is a vector ${\bf
c}=(c_{1},c_{2},\dots,c_{n})\in{\mathbb{N}}^{n}$. A non-root vertex $i$ can
fire when $c_{i}\geq\operatorname{deg}(u)$, so that the number of chips is
great than or equal to the number of vertices that are adjacent to $i$. In
this case the vertex passes chips to each of its neighbors (one for each edge
connecting it to $i$), resulting in a new configuration ${\bf c}^{\prime}$.

A first question to ask is which configurations can be reached from a given
configuration via these firing rules. A configuration ${\bf c}$ is stable if
no non-root vertex can fire, i.e. $c_{i}<\operatorname{deg}(i)$ for all
vertices $i=1,2,\dots,n$. In this case the root vertex can fire, passing a
chip to each vertex adjacent to $0$. A configuration ${\bf c}$ is recurrent if
it is stable and reappears in this process; that is after firing the root
vertex there exists a sequence of firings that returns to ${\bf c}$. The set
of stable and recurrent configurations form group $\kappa(G)$ under
coordinatewise addition (and stabilization), called the critical group of $G$.
It is well known that $\kappa(G)$ is a finite abelian group whose order is
given by the number of spanning trees of $G$, although determining the
structure of $\kappa(G)$ is not easy.

Dually, a configuration ${\bf c}$ is said to be superstable if no set of (non-
root) vertices can fire, so that for all $\sigma\subset\\{1,2,\dots,n\\}$ we
have $c_{i}<\operatorname{deg}_{\sigma}(i)$, where

$\operatorname{deg}_{\sigma}(i)=|\\{j\in[n]\backslash\sigma:ij\in E(G)\\}$

is the number of vertices outside $\sigma$ that are adjacent to $i$. The set
of superstable configurations are also called ‘$G$-parking functions’ in some
contexts. It turns out that there is a simple duality between superstable
configurations and critical configurations. If we let ${\bf
k}=(\operatorname{deg}(1)-1,\operatorname{deg}(2)-1,\dots,\operatorname{deg}(n)-1)$
denote the _canonical configuration_ of the graph $G$, one can show that a
configuration ${\bf c}$ is superstable if and only of ${\bf k}-{\bf c}$ is
critical. In particular the maximal superstables are in bijection with the
minimal critical configurations of the graph $G$.

###  2.1. Linear algebra

We can rephrase many of the constructions from above in terms of some
underlying linear algebra. If $G=(V,E)$ is our graph with
$V=\\{0,1,\dots,n\\}$ we let $\partial$ denote its (signed) incidence matrix.
Here we orient the edges of $G$ so that an edge $e=ij$ has $i\leftarrow j$ if
$i<j$; hence the column of $\partial$ corresponding to $e$ has a $-1$ in the
row corresponding to $i$ and a $1$ in the row corresponding to $j$ (and all
other entries $0$). Removing the row of $\partial$ corresponding to the sink
vertex gives the reduced incidence matrix $\tilde{\partial}$. By definition
the reduced Laplacian matrix $\tilde{\operatorname{\mathcal{L}}}$ of $G$ is
given by

$\tilde{\operatorname{\mathcal{L}}}\vcentcolon=\tilde{\partial}\tilde{\partial}^{t}.$

One can see that the this matrix plays the role of chip-firing in the
dynamical system described above. In particular, given a configuration of
chips ${\bf c}$ we have that firing a vertex $i$ corresponds to subtracting
the $i$th column of $\operatorname{\mathcal{L}}$ from ${\bf c}$ to obtain a
new configuration ${\bf c}^{\prime}$:

${\bf c}^{\prime}={\bf c}-\operatorname{\mathcal{L}}{\bf e}_{i},$

where ${\bf e}_{i}$ is the standard basis vector corresponding to the vertex
$i$. In this sense the set of configurations that are obtainable from another
by a sequence of firings can be rephrased as ${\bf c}\sim{\bf c}^{\prime}$ if
${\bf c}-{\bf c}^{\prime}\in\operatorname{im}\operatorname{\mathcal{L}}$. One
can prove that there is an isomorphism of groups

$\kappa(G)\cong\mathbb{Z}^{|V|-1}/\operatorname{im}\tilde{\operatorname{\mathcal{L}}}.$

We see that $\kappa(G)$ is a finite abelian group and (by the matrix tree
theorem) has cardinality given by the number of spanning trees of $G$.

######  Example 2.1.

As an example, one can take the graph $G$ on vertex set $\\{0,1,2,3\\}$
depicted in Figure 1, with edge set $\\{01,02,03,12,23\\}$. According to the
orientation described above the relevant matrices are

$\tilde{\partial}=\begin{bmatrix}1&0&0&-1&0\\\ 0&1&0&1&-1\\\
0&0&1&0&1\end{bmatrix}\hskip
14.45377pt\tilde{\operatorname{\mathcal{L}}}=\begin{bmatrix}2&-1&0\\\
-1&3&-1\\\ 0&-1&2\end{bmatrix}$

Figure 1. The Diamond graph.

###  2.2. Energy minimizers and M-matrices

The superstable configurations on a graph correspond to those configurations
in which no set of vertices can fire simultaneously (without the number of
chips at some vertex becoming negative). More recently it was shown that
superstable configurations can equivalently be seen as the unique energy
minimizers among the elements in the equivalence class determined by the
Laplacian ([4], [12]). Here one defines a norm on chip configurations defined
by the (inverse of) the reduced Laplacian. In [4] Baker and Shokrieh show that
the energy minimizers correspond to superstable configurations, and hence are
unique per equivalence class.

In an attempt to expand the notion of chip-firing to more general contexts,
Gabrielov ([11]) studied a class of _abelian avalanche models_ , extending the
work of Dhar ([10]). Under this model we fix a finite set $V$ (here taken to
be $[n]=\\{1,\dots,n\\}$) and a _redistribution matrix_ $\Delta$ with indices
in $V$ satisfying

$\Delta_{ii}>0,\text{ for all $i$};\;\;\Delta_{ij}\leq 0,\text{ for $j\neq
i$}.$

Inspired by considerations in physics, this captures the idea that firing a
state involves losing chips at that site and increasing chips at neighboring
sites. A vector ${\bf h}=(h_{1},\dots,h_{n})\in\mathbb{Z}^{n}$ defines a
configuration, and a site $i$ is allowed to _fire_ if $h_{i}\geq\Delta_{ii}$.
In this case firing the site $i$ is defined by replacing ${\bf h}$ with the
vector ${\bf h}-\Delta^{T}{\bf e}_{i}$, subtracting the $i$th row of $\Delta$
from the configuration. Two configurations ${\bf c}$ and ${\bf d}$ are
equivalent if their difference ${\bf c}-{\bf d}$ is in
$\operatorname{im}\Delta$, the $\mathbb{Z}$-image of the matrix $\Delta$.

Mimicking the setup for graphs we say that a configuration is stable if
$c_{i}<\Delta_{ii}$. If $\Delta=\tilde{\operatorname{\mathcal{L}}}(G)$ is the
reduced Laplacian of a graph $G$ we know that if ${\bf c}$ is any initial chip
configuration, then any sequence of chip-firing moves will terminate in a
stable configuration. It is then natural question to ask which matrices
$\Delta$ have the property that the analogous chip-firing process will always
eventually stabilize.

Suppose $L$ is any $n\times n$ matrix. We say that $L$ is a Z-matrix if
$L_{ij}\leq 0$ for all $i\neq j$. Then from [12] we have the following.

######  Proposition 2.2.

[12, Definition 2.2] Suppose $L$ is an $n\times n$ Z-matrix. Then the
following are equivalent.

  1. (1)

$L$ is avalanche finite.

  2. (2)

The real part of the eigenvalues of $L$ are all positive.

  3. (3)

The inverse $L^{-1}$ exists and all entries of $L^{-1}$ are non-negative.

  4. (4)

There exists a vector ${\bf x}\in{\mathbb{R}}^{n}$ with ${\bf x}\geq 0$ such
that $L{\bf x}$ has all positive entries.

The equivalence of the first three conditions can be found in Plemmons [16],
and the connection to the first condition is due to Gabrielov [11]. If any
(and hence all) of the conditions above hold we say that $L$ is a non-singular
$M$-matrix. Such matrices appear in disparate fields, including economics,
operations research, finite element analysis ([16]).

In [12] it is shown that if $L$ is an $M$-matrix then energy-minimizing
configurations exists and are unique per equivalence class. Given a
configuration ${\bf c}\in\mathbb{Z}^{n}$ one defines the energy (or norm) to
be

$E({\bf c})={\bf c}^{T}L{\bf c}.$

Recall that in this set up two configurations ${\bf f},{\bf
g}\in\mathbb{Z}^{n}$ are considered equivalent if ${\bf g}-{\bf
f}\in\operatorname{im}L$, in which case we write ${\bf f}\sim{\bf g}$,. A
configuration ${\bf f}$ is effective, denoted ${\bf f}\geq 0$, if we have
$f_{i}\geq 0$ for all $i=1,\dots,n$. Given any ${\bf f}\in\mathbb{Z}^{n}$ with
${\bf f}\geq 0$ we consider the minimization problem

(1) $\min_{{\bf g}\sim{\bf f},{\bf g}\geq 0}E({\bf g}).$

######  Definition 2.3.

A configuration ${\bf f}\in\mathbb{Z}^{n}$ with ${\bf f}\geq 0$ is
$z$-superstable (with respect to an $M$-matrix $L$ if for every ${\bf
z}\in\mathbb{Z}^{n}$ with ${\bf z}\geq 0$ and ${\bf z}\neq{\bf 0}$ there
exists $i=1,\dots,n$ such that

${\bf f}_{i}-(L{\bf z})_{i}<0.$

In other words there is no ‘multiset-firing’ that can be performed on the
given configuration ${\bf f}$ without resulting in a negative value on some
site. In [12], Guzmán and Klivans show that $z$-superstables are energy
minimizers in the following sense.

######  Theorem 2.4.

[12] Let $L$ be an $M$-matrix. A vector ${\bf f}\in\mathbb{Z}^{n}$ with ${\bf
f}\geq 0$ is $z$-superstable if and only if it is the (unique) minimizer of

$\min_{{\bf g}\sim{\bf f},{\bf g}\geq 0}E({\bf g}).$

As a corollary we see that in each equivalence class defined by $L$ there
exists a unique $z$-superstable configuration. We mention that a natural
notion of critical configurations for $M$-matrices is also discussed in [12],
but we will not need that here.

##  3\. Cycle firing and M-bases

Here we discuss in more detail our approach to cycle and circuit chip-firing
on a graph $G$. As mentioned in Section 1 the basic idea is to view
traditional chip-firing on $G$ as determined by the incidence matrix of $G$
and to apply matroidal/Gale duality.

We first recall some notions from lattice theory from [1]. For this let
${\mathbb{R}}^{m}$ denote Euclidean space with the usual inner product
$\langle\cdot,\cdot\rangle$. Suppose $F$ and $C$ are orthogonal subspaces of
${\mathbb{R}}^{m}$ and further suppose that both are rational (so that each
has a basis consisting of vectors from ${\mathbb{Q}}^{m}$). Let
$\operatorname{\mathcal{F}}$ and $\operatorname{\mathcal{C}}$ denote the
lattices given by $\operatorname{\mathcal{F}}=\mathbb{Z}^{m}\cap F$ and
$\operatorname{\mathcal{C}}=\mathbb{Z}^{m}\cap C$.

For a lattice $\Lambda\subset{\mathbb{R}}^{m}$ we define its dual lattice as

$\Lambda^{\sharp}=\\{x\in{\mathbb{R}}^{m}:\langle
x,\lambda\rangle\in\mathbb{Z}\;\text{for all $\lambda\in\Lambda$}\\}.$

From [1] we have the following.

######  Lemma 3.1.

[1, Lemma 1] With notation as above,

$\operatorname{\mathcal{C}}^{\sharp}/\operatorname{\mathcal{C}}\cong\operatorname{\mathcal{F}}^{\sharp}/\operatorname{\mathcal{F}}.$

In fact both are isomorphic to the group
$\mathbb{Z}^{m}/(\operatorname{\mathcal{C}}\oplus\operatorname{\mathcal{F}})$.

The group $\Lambda^{\sharp}/\Lambda$ is called the determinant group of the
lattice $\Lambda$. As was pointed out in [13], we have another way to compute
determinant groups in terms of a basis for the lattice.

######  Lemma 3.2.

[13, Theorem 12] Suppose $\Lambda$ is a rank $r$ sublattice of
$\mathbb{Z}^{m}\subset{\mathbb{R}}^{m}$, and $\\{z_{1},\dots,z_{r}\\}$ is an
integer basis for $\Lambda$. Then the determinant group of $\Lambda$ is given
by

$\Lambda^{\sharp}/\Lambda\cong\mathbb{Z}^{r}/\operatorname{im}(\langle
z_{i},z_{j}\rangle)_{i,j=1,\dots,r}.$

To apply these results in our context we let $G=(V,E)$ be a finite simple
graph on vertex set $V=\\{0,1,\dots,n\\}$ and edge set $E$. Let
$\partial=\partial_{G}:{\mathbb{R}}^{E}\rightarrow{\mathbb{R}}^{V}$ denote the
signed incidence matrix of $G$ (after specifying some orientation). We let

$\operatorname{\mathcal{F}}(G)=\mathbb{Z}^{E}\cap\ker(\partial)$

denote the lattice of integral flows and

$\operatorname{\mathcal{C}}(G)=\mathbb{Z}^{E}\cap\operatorname{im}(\partial^{T})$

denote the lattice of integral cuts. We use $\operatorname{\mathcal{F}}$ and
$\operatorname{\mathcal{C}}$ to denote these lattices if the context is clear.
Note that (assuming $G$ is connected) an integer basis for
$\operatorname{\mathcal{C}}$ is obtained by removing a single row
corresponding to any (say the sink) vertex.

Dually, suppose we have an integer basis $\\{{\bf f}_{1},\dots,{\bf f}_{g}\\}$
for $\operatorname{\mathcal{F}}$, where $g=|E|-|V|+1$ is the so-called genus
of the graph. Let $\iota$ be the $|E|\times g$ matrix with columns ${\bf
f}_{i}$ and define $\operatorname{\mathcal{L}}^{*}\vcentcolon=\iota^{T}\iota$
to be the dual Laplacian (with respect to this choice of basis). We recover
the critical group as the cokernel of $\operatorname{\mathcal{L}}^{*}$, as
follows.

######  Proposition 3.3.

Suppose $G$ is any graph. Let $\iota$ be the matrix associated to any basis
for the lattice of integral flows $\operatorname{\mathcal{F}}(G)$, and let
$\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$ denote the associated dual
Laplacian. We then have an isomorphism of groups

$\kappa(G)\cong\mathbb{Z}^{g}/\operatorname{im}\operatorname{\mathcal{L}}^{*}.$

###### Proof.

If $\partial$ denotes the incidence matrix of $G$ (under some orientation)
then we see that $F=\ker\partial$ and $C=\operatorname{im}\partial^{T}$
satisfy the conditions of Lemma 3.1. Recall that we have an integral basis for
$\operatorname{\mathcal{C}}(G)=\mathbb{Z}^{E}\cap C$ given by all rows of
$\partial$ corresponding to nonroot vertices. From Lemma 3.2 we see that

$\operatorname{\mathcal{F}}(G)^{\sharp}/\operatorname{\mathcal{F}}(G)\cong\mathbb{Z}^{n}/\operatorname{im}\tilde{\operatorname{\mathcal{L}}}\cong\kappa(G).$

But also from Lemma 3.2 we have that

$\operatorname{\mathcal{C}}(G)^{\sharp}/\operatorname{\mathcal{C}}(G)\cong\mathbb{Z}^{g}/\operatorname{im}\operatorname{\mathcal{L}}^{*},$

and so the result follows from Lemma 3.1. ∎

If $G$ is a planar graph, we can take an embedding of $G$ in the plane to
define a dual graph $G^{*}$. By definition $G$ has vertices given by regions
in the plane determined by the embedding, and adjacency given by regions
sharing an edge (see Figure 2). After choosing an orientation of the graph
$G$, as well as an orientation of the plane, we then have a basis for
$\operatorname{\mathcal{F}}$ given by the bounded regions of the embedding.
Proposition 3.3 in particular tells us that $\tau(G^{*})\cong\tau(G)$, as was
first established in [8].

Figure 2. The embedded Diamond graph and its dual.

######  Example 3.4.

We continue with the graph $G$ from Example 2.1. For the choice of basis for
$\operatorname{\mathcal{F}}(G)$ described above we have the reduced Laplacian
$\tilde{\operatorname{\mathcal{L}}}$ given by

$\tilde{\operatorname{\mathcal{L}}}=\begin{bmatrix}2&-1&0\\\ -1&3&-1\\\
0&-1&2\\\ \end{bmatrix}.$

The critical group of $G$ is given by
$\mathbb{Z}^{3}/\operatorname{im}\tilde{\operatorname{\mathcal{L}}}\cong\mathbb{Z}/8\mathbb{Z}$,
which has representatives (superstable configurations) given by
$\\{121,111,120,021,101,110,011,020\\}$. In the dual graph depicted in Figure
2, we take the vertex corresponding to the unbounded face as the (red) sink
vertex and order the other vertices arbitrarily. In this case we have

$\operatorname{\mathcal{L}}^{*}=\begin{bmatrix}3&-1\\\ -1&3\end{bmatrix}.$

The ‘dual’ critical group is given by
$\mathbb{Z}^{2}/\operatorname{\mathcal{L}}^{*}\cong\mathbb{Z}/{8\mathbb{Z}}$,
and the superstable configurations of $G^{*}$ are
$\\{22,21,12,02,11,20,01,10\\}$.

###  3.1. Finding an M-basis

As we have seen the (usual, reduced) Laplacian
$\tilde{\operatorname{\mathcal{L}}}$ plays a prominent role in the traditional
theory of chip-firing on a graph $G$. We would like to use the dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ to define a similar rule for ‘cycle chip-
firing’ on $G$, where elements in our basis for ${\mathcal{F}}(G)$ lend chips
to other ‘neighboring elements’. However, an arbitrary choice of basis for
$\operatorname{\mathcal{F}}(G)=\ker\tilde{\partial}$ may lead to a dual
Laplacian that is not a $Z$-matrix, and in particular may have positive
entries off the diagonal.

######  Example 3.5.

For the graph in our running Example 2.1 if we take the set
$\left\\{\begin{bmatrix}1\\\ -1\\\ 0\\\ 1\\\
0\end{bmatrix},\begin{bmatrix}0\\\ -1\\\ 1\\\ 0\\\ -1\end{bmatrix}\right\\}$
as a basis for $\ker\tilde{\partial}$ we obtain a dual Laplacian
$\operatorname{\mathcal{L}}^{*}=\begin{bmatrix}3&1\\\ 1&3\end{bmatrix}$.

To simplify notation we make the following definition.

######  Definition 3.6.

Suppose $G$ is a connected graph and let $\operatorname{\mathcal{B}}=\\{{\bf
f}_{1},\dots,{\bf f}_{g}\\}$ be a basis for the lattice of integral flows
$\operatorname{\mathcal{F}}(G)$. Let $\iota$ be the $|E|\times g$ matrix with
columns given by $\operatorname{\mathcal{B}}$. We say that
$\operatorname{\mathcal{B}}$ is an cycle $M$-basis if the matrix
$\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$ is an $M$-matrix.

As expected the situation is fairly straightforward for the case of planar
graphs.

######  Proposition 3.7.

Any planar graph admits a cycle $M$-basis.

###### Proof.

Suppose $G$ is a planar graph embedded in the plane and let
$R_{1},\dots,R_{g}$ denote the set of bounded regions. Orient $R_{1}$
arbitrarily, and for any $R_{i}$ that shares an edge with $R_{1}$ orient
$R_{i}$ so that the edges are oriented in opposite directions. Continue this
process until all regions have been oriented. If some $R_{j}$ does not share
an edge with the regions oriented so far, then orient $R_{j}$ and proceed as
above. A choice of compatible orientations is guaranteed since
${\mathbb{R}}^{2}$ itself is orientable.

Now for each region $R_{i}$ we construct a vector ${\bf
f}_{i}\in\mathbb{Z}^{E}$ according to the orientation describe above relative
to the fixed orientation on $G$ determined by the (signs of) entries of
$\partial$. Namely, the vector ${\bf f}_{i}$ has a $1$ in the entry
corresponding to an edge $e_{j}$ if the orientation in $R_{i}$ agrees with
that of $\delta$, $-1$ if the orientation is opposite, and $0$ if the edge
$e_{j}$ does not appear in $R_{i}$. The collection $\\{{\bf f}_{1},\dots,{\bf
f}_{g}\\}$ forms a basis for $\ker\tilde{\delta}$, and we let $\iota$ denote
the matrix with these vectors as columns.

It follows that the matrix $\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$
equals $\tilde{\operatorname{\mathcal{L}}}(G^{*})$, the reduced Laplacian of
the dual graph of $G$ with respect to the given embedding. Here the designated
root of $G^{*}$ corresponds to the unbounded face of $G$ in this embedding.
Reduced Laplacians of graphs are known to be $M$-matrices, so the result
follows. ∎

In the case of planar graphs we can import the usual notions of superstable
configurations, etc. from the theory of chip-firing on $G^{*}$, where $G^{*}$
is the dual graph defined by some embedding of $G$. The bounded faces of $G$
(thought of as circuits) will then correspond to the nonsink vertices of
$G^{*}$, so firing a vertex of $G^{*}$ corresponds to firing a circuit of $G$,
and we have a reasonable notion of ‘cycle chip-firing’ (see Figure 2).

If $G$ is not planar we still have a good notion of
$\operatorname{\mathcal{L}}^{*}$ (after choosing an integral basis for $F$)
and we would like to interpret this matrix as describing some rule for chip-
firing. We first observe that the only obstruction to being an $M$-matrix is
the sign pattern of the entries in $\operatorname{\mathcal{L}}^{*}$.

######  Lemma 3.8.

Suppose $G$ is a graph and $\\{{\bf f}_{1},\dots,{\bf f}_{g}\\}$ is a basis
for $\operatorname{\mathcal{F}}(G)$. Then the resulting dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ is a cycle $M$-basis if and only if the
matrix $\operatorname{\mathcal{L}}^{*}$ is a $Z$-matrix.

###### Proof.

Let $\iota$ denote the matrix with columns ${\bf f}_{1},\dots\bf f_{g}$. By
definition we have that $\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$ so
that $\operatorname{\mathcal{L}}^{*}$ is a symmetric positive definite matrix
(note that $\operatorname{\mathcal{L}}^{*}$ is invertible). From this it
follows that $\operatorname{\mathcal{L}}^{*}$ will have positive real
eigenvalues. Hence as long as $\operatorname{\mathcal{L}}^{*}$ is a
$Z$-matrix, it will satisfy the conditions of being an $M$-matrix listed in
Proposition 2.2. ∎

######  Example 3.9.

Consider $K_{5}$, the complete graph on vertex set $\\{1,2,3,4,5\\}$ with an
orientation on edges given by $(i,j)$ whenever $i<j$. With this orientation
the matrix $\partial$ is given below (with columns ordered in lexicographical
order):

$\partial=\begin{bmatrix}1&1&1&1&0&0&0&0&0&0\\\ -1&0&0&0&1&1&1&0&1&0\\\
0&-1&0&0&-1&0&0&1&1&0\\\ 0&0&-1&0&0&-1&0&-1&0&1\\\
0&0&0&-1&0&0&-1&0&-1&-1\end{bmatrix}$

A basis for $\operatorname{\mathcal{C}}=\ker\partial$ and resulting dual
Laplacian $\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$ are given by

$\iota=\begin{bmatrix}1&0&-1&0&1&0\\\ 0&-1&0&0&-1&1\\\ 0&0&0&1&0&-1\\\
-1&1&1&-1&0&0\\\ 0&1&-1&1&0&0\\\ 1&0&0&-1&0&0\\\ 0&-1&0&0&1&0\\\
0&0&0&0&-1&1\\\ 0&0&-1&1&0&0\\\ 1&0&0&0&-1&0\end{bmatrix},\hskip
14.45377pt\operatorname{\mathcal{L}}^{*}=\begin{bmatrix}4&-1&-2&0&0&0\\\
-1&4&0&0&0&-1\\\ -2&0&4&-3&-1&0\\\ 0&0&-3&5&0&-1\\\ 0&0&-1&0&5&-2\\\
0&-1&0&-1&-2&3\end{bmatrix}$

Figure 3. The cycles represented by the basis for $\ker\partial$ in Example
3.9.

According to Lemma 3.8 we see that $\operatorname{\mathcal{L}}^{*}$ is an
$M$-matrix. Note that the sum of the entries in the third row is negative.
Hence there will exist circuit configurations that are stable under set
firings but not $z$-superstable.

For example, one can check that the configuration $\\{1,1,0,0,0,1\\}$ is
$\chi$-superstable. However, firing the multiset $\\{5,3,8,6,4,6\\}$ leads to
the zero configuration.

For the case of general graphs we will need a slightly modified version of the
Gram-Schmidt algorithm which we record below.

######  Lemma 3.10.

Suppose $\\{{\bf v}_{1},\dots\bf v_{n}\\}$ is a linearly independent
collection of vectors in ${\mathbb{R}}^{n}$. Then there exists an orthogonal
set of vectors $\\{{\bf q}_{1},\dots,{\bf q}_{n}\\}$ with the property that
for all $i=1,\dots,n$ the vector ${\bf q}_{i}$ is in the integer span of
$\\{{\bf v}_{1},\dots\bf v_{i}\\}$.

###### Proof.

We apply the Gram-Schimdt algorithm but scale at each step so that the
constructed vectors are in the desired integer span. In particular we set
${\bf q}_{1}\vcentcolon={\bf v}_{1}$ and for $i>1$ define the vectors ${\bf
q}_{i}$ as

${\bf q}_{i}\vcentcolon=\left(\prod^{i-1}_{j=1}\|{\bf
q}_{j}\|^{2}\right)\left({\bf v}_{i}-\sum^{i-1}_{j=1}\dfrac{{\bf
v}_{i}\cdot{\bf q}_{j}}{\|{\bf q}_{j}\|^{2}}{\bf q}_{j}\right).$

We see that ${\bf q}_{i}$ is an integer linear combination of the set $\\{{\bf
v}_{i},{\bf q}_{1},\dots,{\bf q}_{i-1}\\}$. By induction, for all
$j=1,\dots,i-1$ each ${\bf q}_{j}$ is an integer linear combination of
$\\{{\bf v}_{1},\dots,{\bf v}_{j}\\}$. The vector ${\bf q}_{i}$ is nonzero
(since the original set was linearly independent) and a calculation shows that
it is orthogonal to ${\bf q}_{j}$ for all $j<i$. ∎

With this we prove our main algebraic result.

######  Theorem 3.11.

Suppose $\\{{\bf v}_{1}\dots,{\bf v}_{g}\\}$ is an integral basis for a
lattice $\Lambda\subset{\mathbb{R}}^{d}$. Then there exists an integral basis
$\\{{\bf f}_{1},\dots,{\bf f}_{g}\\}$ for $\Lambda$ with the property that
${\bf f}_{i}\cdot{\bf f}_{j}\leq 0$ for all $i\neq j$.

###### Proof.

Let $\\{{\bf q}_{1},\dots\bf q_{g}\\}$ be the set of orthogonal vectors
obtained from $\\{{\bf v}_{1}\dots,{\bf v}_{g}\\}$ as in the proof of Lemma
3.10. We set ${\bf f}_{1}\vcentcolon={\bf v}_{1}$ and construct the rest of
our vectors inductively. For ${\bf f}_{2}$ we want ${\bf f}_{2}={\bf
v}_{2}+t{\bf q}_{1}$ for some $t\in\mathbb{Z}$ satisfying ${\bf
f}_{1}\cdot({\bf v}_{2}+t{\bf q}_{1})\leq 0$. Since ${\bf q}_{1}={\bf f}_{1}$
we have that ${\bf q}_{1}\cdot{\bf f}_{1}>0$ and hence we need
$t\leq\dfrac{-{\bf f}_{1}\cdot{\bf v}_{2}}{{\bf f}_{1}\cdot{\bf q}_{1}}$.
Therefore we set

$t=\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{2}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor,$

so that ${\bf f}_{2}\vcentcolon={\bf v}_{2}+\left\lfloor\dfrac{-{\bf
f}_{1}\cdot{\bf v}_{2}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf q}_{1}$.
Note that ${\bf f}_{2}$ is obtained by performing integer column operations to
the matrix with columns $\\{{\bf v}_{1}{\bf v}_{2}\\}$ and hence
$\operatorname{span}_{\mathbb{Z}}\\{{\bf f}_{1},{\bf
f}_{2}\\}=\operatorname{span}_{\mathbb{Z}}\\{{\bf v}_{1},{\bf v}_{2}\\}$. By
construction we have ${\bf f}_{1}\cdot{\bf f}_{2}\leq 0$.

To illustrate one more step of the process, we define the vector ${\bf f}_{3}$
as

${\bf f}_{3}\vcentcolon={\bf v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf
v}_{3}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf
q}_{1}+\left\lfloor\dfrac{-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf f}_{2}\cdot{\bf
q}_{2}}\right\rfloor{\bf q}_{2}.$

Since ${\bf q}_{2}$ is orthogonal to ${\bf v}_{1}$ we have that

${\bf f}_{1}\cdot{\bf f}_{3}={\bf f}_{1}\cdot{\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf f}_{1}\cdot{\bf q}_{1}\leq 0.$

To see that ${\bf f}_{2}\cdot{\bf f}_{3}\leq 0$ note that

${\bf f}_{2}\cdot{\bf f}_{3}={\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)+{\bf
f}_{2}\cdot\left\lfloor\dfrac{-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf f}_{2}\cdot{\bf
q}_{2}}\right\rfloor{\bf q}_{2},$

so that

${\bf f}_{2}\cdot{\bf f}_{3}-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)=\left\lfloor\dfrac{-{\bf
f}_{2}\cdot\left({\bf v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf
v}_{3}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf
f}_{2}\cdot{\bf q}_{2}}\right\rfloor{\bf f}_{2}\cdot{\bf q}_{2}.$

But ${\bf f}_{2}\cdot{\bf q}_{2}>0$ so that

$\left\lfloor\dfrac{-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf f}_{2}\cdot{\bf
q}_{2}}\right\rfloor{\bf f}_{2}\cdot{\bf q}_{2}\leq\left(\dfrac{-{\bf
f}_{2}\cdot\left({\bf v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf
v}_{3}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf
f}_{2}\cdot{\bf q}_{2}}\right){\bf f}_{2}\cdot{\bf q}_{2}$

Hence we have

${\bf f}_{2}\cdot{\bf f}_{3}-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)\leq\left(\dfrac{-{\bf
f}_{2}\cdot\left({\bf v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf
v}_{3}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf q}_{1}\right)}{{\bf
f}_{2}\cdot{\bf q}_{2}}\right){\bf f}_{2}\cdot{\bf q}_{2}={-{\bf
f}_{2}\cdot\left({\bf v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf
v}_{3}}{{\bf f}_{1}\cdot{\bf q}_{1}}\right\rfloor{\bf q}_{1}\right)}$

and we conclude

${\bf f}_{2}\cdot{\bf f}_{3}\leq-{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)+{\bf f}_{2}\cdot\left({\bf
v}_{3}+\left\lfloor\dfrac{-{\bf f}_{1}\cdot{\bf v}_{3}}{{\bf f}_{1}\cdot{\bf
q}_{1}}\right\rfloor{\bf q}_{1}\right)=0.$

For the general case we let $\\{{\bf q}_{1},\dots\bf q_{g}\\}$ be the set of
vectors constructed from Lemma 3.10. Assuming we have constructed vectors
$\\{{\bf f}_{1},\dots,{\bf f}_{n}\\}$ with the desired properties we define

${\bf f}_{n+1}\vcentcolon={\bf v}_{n+1}+\sum^{n}_{j=1}\left\lfloor\dfrac{-{\bf
f}_{j}\cdot{\bf a}_{n,j}}{{\bf f}_{j}\cdot{\bf q}_{j}}\right\rfloor{\bf
q}_{j},$

where

${\bf a}_{n,j}={\bf v}_{n+1}+\sum^{j-1}_{k=1}\left\lfloor\dfrac{-{\bf
f}_{k}\cdot{\bf a}_{n,k}}{{\bf f}_{k}\cdot{\bf q}_{k}}\right\rfloor{\bf
q}_{k},$

with ${\bf a}_{n,1}=v_{n+1}$.

First we note that ${\bf f}_{i}\cdot{\bf q}_{i}>0$ for all $i=1,\dots,n$, so
that in particular these formulas make sense. To see this let
$\alpha_{j}=\left\lfloor\dfrac{-{\bf f}_{j}\cdot{\bf a}_{i-1,j}}{{\bf
f}_{j}\cdot{\bf q}_{j}}\right\rfloor$ and since the set $\\{{\bf
q}_{1},\dots,{\bf q}_{i-1}\\}$ is orthogonal observe that

${\bf f}_{i}\cdot{\bf q}_{i}=({\bf v}_{i}+\alpha_{1}{\bf
q}_{1}+\cdots+\alpha_{i-1}{\bf q}_{i-1})\cdot{\bf q}_{i}={\bf v}_{i}\cdot{\bf
q}_{i}.$

But ${\bf v}_{i}\cdot{\bf q}_{i}>0$ since ${\bf v}_{i}$ and ${\bf q}_{i}$ are
on the same side of the hyperplane $\operatorname{span}\\{{\bf
q}_{1},\dots,{\bf q}_{i-1}\\}$ in the space $\operatorname{span}\\{{\bf
q}_{1},\dots,{\bf q}_{i}\\}$. Note that ${\bf
q}_{i}\notin\operatorname{span}\\{{\bf q}_{1},\dots,{\bf q}_{i-1}\\}$ since
the set $\\{{\bf q}_{1},\dots,{\bf q}_{i}\\}$ is orthogonal. Hence for any
$\beta\in{\mathbb{R}}$ whenever $t=\left\lfloor\dfrac{\beta}{{\bf
f}_{i}\cdot{\bf q}_{i}}\right\rfloor$ we have that $t({\bf f}_{i}\cdot{\bf
q}_{i})\leq\beta$.

We claim that ${\bf f}_{k}\cdot{\bf f}_{\ell}\leq 0$ for all $1\leq k<\ell\leq
n+1$. To see this first note that

${\bf f}_{k}\cdot{\bf f}_{\ell}={\bf f}_{k}\cdot{\bf
v}_{\ell}+\sum^{k}_{j=1}\left\lfloor\dfrac{-{\bf f}_{j}\cdot{\bf
a}_{\ell-1,j}}{{\bf f}_{j}\cdot{\bf q}_{j}}\right\rfloor{\bf f}_{k}\cdot{\bf
q}_{j},$

since ${\bf q}_{m}$ is orthogonal to ${\bf f}_{k}$ for all $m>k$ (recall ${\bf
f}_{k}$ is constructed as a linear combination of $\\{{\bf v}_{1},\dots,{\bf
v}_{k}\\}$).

But we see that

$\left\lfloor\dfrac{{-\bf f}_{k}\cdot{\bf a}_{\ell-1,k}}{{\bf f}_{k}\cdot{\bf
q}_{k}}\right\rfloor{\bf f}_{k}\cdot{\bf q}_{k}\leq-{\bf f}_{k}\cdot{\bf
a}_{\ell-1,k}=-{\bf f}_{k}\cdot\left({\bf
v}_{\ell}+\sum^{k-1}_{j=1}\left\lfloor\dfrac{-{\bf f}_{j}\cdot{\bf
a}_{\ell-1,j}}{{\bf f}_{j}\cdot{\bf q}_{j}}\right\rfloor{{\bf q}_{j}}\right).$

Hence ${\bf f}_{k}\cdot{\bf f}_{\ell}\leq 0$, as desired.

Finally it is clear that $\operatorname{span}_{\mathbb{Z}}\\{{\bf
f}_{1},\dots,{\bf f}_{n+1}\\}=\operatorname{span}_{\mathbb{Z}}\\{{\bf
v}_{1},\dots\bf v_{n+1}\\}$ since ${\bf v}_{n+1}$ was constructed by adding
integer multiples of vectors from the set $\\{{\bf v}_{1},\dots,{\bf
v}_{n}\\}$ to the vector ${\bf v}_{n+1}$. The result follows. ∎

######  Example 3.12.

As an small illustration of the algorithm, suppose

$\\{{\bf v}_{1},{\bf v}_{2},{\bf v}_{3}\\}=\left\\{\begin{bmatrix}1\\\ 0\\\
2\\\ 1\end{bmatrix},\begin{bmatrix}1\\\ 1\\\ 0\\\
1\end{bmatrix},\begin{bmatrix}2\\\ 0\\\ 1\\\ 1\end{bmatrix}\right\\}$

is a given basis for a lattice $\Lambda\in{\mathbb{R}}^{4}$. Then Lemma 3.10
produces the orthogonal set

$\\{{\bf q}_{1},{\bf q}_{2},{\bf q}_{3}\\}=\left\\{\begin{bmatrix}1\\\ 0\\\
2\\\ 1\end{bmatrix},\begin{bmatrix}4\\\ 6\\\ -4\\\
4\end{bmatrix},\begin{bmatrix}396\\\ -288\\\ -144\\\
-108\end{bmatrix}\right\\}.$

We then apply the algorithm in the proof of Theorem 3.11 to get

$\\{{\bf f}_{1},{\bf f}_{2},{\bf f}_{3}\\}=\left\\{\begin{bmatrix}1\\\ 0\\\
2\\\ 1\end{bmatrix},\begin{bmatrix}0\\\ 1\\\ -2\\\
0\end{bmatrix},\begin{bmatrix}-3\\\ -6\\\ 3\\\ -4\end{bmatrix}\right\\},$

a basis for $\Lambda$ with the desired properties.

######  Corollary 3.13.

Any graph $G$ admits a cycle $M$-basis. That is, there exists a basis
$\iota=\\{{\bf f}_{1}\dots,{\bf f}_{g}\\}$ for its flow space
$\operatorname{\mathcal{F}}(G)$ such that the associated dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ is an $M$-matrix.

###### Proof.

From Lemma 3.8 we have that $\operatorname{\mathcal{L}}^{*}$ is an $M$-matrix
if and only if it is a $Z$-matrix. From Theorem 3.11 we have that a basis
$\\{{\bf f}_{1},\dots,{\bf f}_{g}\\}$ for the lattice
$\operatorname{\mathcal{F}}(G)$ can always be chosen so that ${\bf
f}_{i}\cdot{\bf f}_{j}\leq 0$ for all $i\neq j$. The result follows.

∎

###  3.2. z-superstable cycle configurations

Suppose $G$ is a graph and $\iota=\\{{\bf f}_{1}\dots,{\bf f}_{g}\\}$ is a
cycle $M$-basis. Then results from [12] imply that the set of $z$-superstable
configurations represent a unique element from each equivalence class of
$\mathbb{Z}^{g}/\operatorname{im}({\operatorname{\mathcal{L}}^{*}})$. Recall
that a configuration is $z$-superstable if no multiset of sites can be fired
without resulting in a negative number of chips at some site. In the context
of dual Laplacians we will call these $z$-superstable cycle configurations.

######  Proposition 3.14.

Let $G$ be a connected graph, and let $\iota=\\{{\bf f}_{1}\dots,{\bf
f}_{g}\\}$ be a cycle $M$-basis with associated dual Laplacian
$\operatorname{\mathcal{L}}^{*}$. Then the number of $z$-superstable cycle
configurations of $G$ is given by $|\tau(G)|$, the number of spanning trees of
$G$.

###### Proof.

By results of [12] discussed above, there exists a unique $z$-superstable
configuration in each equivalence class determined by
$\operatorname{\mathcal{L}}^{*}$. Hence the number of such $z$-superstable
configurations is given by

$|\mathbb{Z}^{g}/\operatorname{im}\operatorname{\mathcal{L}}^{*}|.$

From Proposition 3.3, we have

$\kappa(G)\cong\mathbb{Z}^{g}/\operatorname{im}\operatorname{\mathcal{L}}^{*},$

and since $|\kappa(G)|=|\tau(G)|$ the result follows. ∎

Recall for classical chip-firing on a graph $G$, the number of superstable
configurations of degree $d$ is given by the number of spanning trees of $G$
with $d$ externally passive edges. It is not clear if our $z$-superstable
configurations have some other combinatorial interpretation in terms of some
statistic on the set of spanning trees. See the last section for further
discussion.

######  Example 3.15.

For the graph $K_{5}$ we use the dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ from Example 3.9 which produces $125$
$z$-superstable circuit configurations, corresponding to the $5^{5-2}=125$
spanning trees of $K_{5}$. This collection of multisets is clearly closed
under taking subsets, so it suffices to describe the maximal elements. Here
they are given by

$\displaystyle\\{000112,000211,010022,010210,010300,020021,020040,020111,021020,021110,030101,$
$\displaystyle
100102,101020,101110,130020,130110,200021,200111,210020,210110,300020,310000\\}.$

The degree sequence is given by ${\bf c}=(1,6,19,38,39,19,3)$, where the
number of $z$-superstable configurations of degree $d-1$ is given by the entry
$c_{d}$.

##  4\. Circuit M-bases

From the above results we see any graph $G$ admits an $M$-basis, a basis for
the flow space $\operatorname{\mathcal{F}}(G)$ with the property that the
associated dual Laplacian $\operatorname{\mathcal{L}}^{*}$ is an $M$-matrix.
In the proof of Theorem 3.11, however, we see that the $M$-basis elements are
constructed via linear combinations of a given basis. A natural question to
ask if one can find a $M$-basis for $\operatorname{\mathcal{F}}(G)$ consisting
of _circuits_. More precisely we want the entries of each basis vector to be
$0$ or $\pm 1$, with the the nonzero entries corresponding to some circuit of
in the graph $G$ (with some chosen orientation). Note that to be an element of
$\ker\tilde{\partial}$ we need those cycles to be oriented, so the sign of any
entry will be determined by how this orientation compares to the fixed
orientation on $G$.

######  Definition 4.1.

Suppose $G$ is a connected graph. A circuit $M$-basis for $G$ is a cycle
$M$-basis such that each basis element corresponds to a circuit as above.

In the proof of Proposition 3.7 we see that if $G$ is a planar graph then any
embedding of $G$ gives rise to a basis for $\operatorname{\mathcal{F}}(G)$
that consists of circuits. Hence we in fact have the following stronger
result.

######  Proposition 4.2.

Any planar graph admits a circuit $M$-basis.

The next natural question to ask is whether any graph $G$ admit a circuit
$M$-basis. Although we suspect that the answer is no, we have not found any
counterexamples. However, via a computer search we have been able to find
circuit $M$-bases for the two ‘smallest’ nonplanar graphs.

######  Proposition 4.3.

The complete graph $K_{5}$ and the complete bipartite graph $K_{3,3}$ both
admit circuit $M$-bases.

###### Proof.

For $K_{5}$ we take the basis for $\operatorname{\mathcal{F}}(K_{5})$
described in Example 3.9. There we saw that the associated dual Laplacian
$\operatorname{\mathcal{L}}^{*}$ was an $M$-matrix.

For $K_{3,3}$ we can choose a basis for $\operatorname{\mathcal{F}}(K_{3,3})$
given by the columns of the matrix $\iota$ below. The corresponding dual
Laplacian $\operatorname{\mathcal{L}}^{*}=\iota^{T}\iota$ is also indicated.

$\iota=\begin{bmatrix}-1&0&1&0\\\ 0&1&-1&0\\\ 1&-1&0&0\\\ 1&0&0&-1\\\
0&-1&1&1\\\ -1&1&-1&0\\\ 0&0&-1&1\\\ 0&0&0&-1\\\ 0&0&1&0\end{bmatrix},\hskip
14.45377pt\operatorname{\mathcal{L}}^{*}=\begin{bmatrix}4&-2&0&-1\\\
-2&4&-3&-1\\\ 0&-3&6&0\\\ -1&-1&0&4\end{bmatrix}$

Here the vertex set of $K_{3,3}$ is given by bipartition
$\\{0,1,2\\}\cup\\{3,4,5\\}$, and the rows in $\iota$ are labeled
lexicographically $(03,04,\dots,24,25)$. See Figure 4 for an illustration of
the columns of $\iota$. From Lemma 3.8 we see that
$\operatorname{\mathcal{L}}^{*}$ is an $M$-matrix. ∎

Figure 4. A circuit $M$-basis for the graph $K_{3,3}$ from the proof of
Proposition 4.3.

###  4.1. Obstructions to a circuit M-basis

In this section we study restrictions on possible circuit $M$-bases for
nonplanar graphs. For instance, having found a circuit $M$-basis for the
complete graph $K_{5}$ consisting of circuits (see Example 3.9) it is natural
to ask whether it is possible to find a circuit $M$-basis consisting of all
3-cycles. In fact this is not the case.

######  Proposition 4.4.

Let $n\geq 5$, and suppose $\operatorname{\mathcal{B}}$ is a circuit $M$-basis
for $K_{n}$. Then at least one element of $\operatorname{\mathcal{B}}$ has
more than three nonzero entries (so that at least one of the circuits is not a
triangle).

###### Proof.

For a contradiction suppose $\operatorname{\mathcal{B}}=\\{{\bf
f}_{1},\dots,{\bf f}_{g}\\}$ is a circuit $M$-basis for $K_{n}$ consisting of
all triangles. Let $\iota$ be the matrix with columns given by the elements of
$\operatorname{\mathcal{B}}$. In what follows we will also think of elements
of $\operatorname{\mathcal{B}}$ as oriented triangles.

If an edge $e$ is contained in three or more elements of
$\operatorname{\mathcal{B}}$ then since any pair of triangles shares at most
edge, we have that some pair $\\{{\bf f}_{i},{\bf f}_{j}\\}$ of these elements
would have this edge oriented in the same direction. Hence the corresponding
$(i,j)$ entry in $\iota\iota^{T}$ would be positive, a contradiction. We
conclude that every row of $\iota$ must have either one nonzero entry or two
nonzero entries with different signs.

Recall that the genus of $K_{n}$ is $\frac{(n-2)(n-1)}{2}$. Hence since
$\operatorname{\mathcal{B}}$ is a circuit basis consisting only of triangles,
we have that $\iota$ has a total of $\frac{3(n-2)(n-1)}{2}$ nonzero entries.
However, because there is a maximum of two nonzero entries per row and there
are $\frac{n(n-1)}{2}$ edges, then there is a maximum of $n(n-1)$ nonzero
values. If $n>6$, then $\frac{3(n-2)(n-1)}{2}>n(n-1)$ and hence no such basis
exists.

If $n=6$, then $\frac{3(n-2)(n-1)}{2}=n(n-1)=30$. In this case every edge must
appear in exactly two triangles, with an opposite sign in each. The sum of
these cycles is the zero cycle, and hence the set ${}_{1},\dots,{\bf
f}_{g}\\}$ is linearly dependent. This is a contradiction to the assumption
that $\operatorname{\mathcal{B}}$ is a basis.

If $n=5$, then we need 18 nonzero entries in 10 rows. Hence we must have two
rows in $\iota$ that have a single nonzero entry. But in this case the sum
${\bf f}_{1}+\cdots+{\bf f}_{g}$ of the basis elements would be a vector with
two nonzero entries, which cannot be an element of the flow space
$\operatorname{\mathcal{F}}(G)=\ker\tilde{\partial}$. The result follows. ∎

Note that since $K_{4}$ is planar with an embedding with bounded regions given
by $3$-cycles, we see that $n\geq 5$ is tight.

######  Proposition 4.5.

Let $3\leq m\leq n$, and suppose $\operatorname{\mathcal{B}}$ is a circuit
$M$-basis for $K_{m,n}$. Then at least one element of
$\operatorname{\mathcal{B}}$ has more than four nonzero entries (so that at
least one of the circuits is not a $4$-cycle).

###### Proof.

For the sake of contradiction, suppose $\operatorname{\mathcal{B}}=\\{{\bf
f}_{1},\dots,{\bf f}_{g}\\}$ is a circuit $M$-basis for $K_{m,n}$ consisting
of all $4$-cycles. Let $\iota$ be the matrix with columns given by the
elements of $\operatorname{\mathcal{B}}$. In what follows we will also think
of elements of $\operatorname{\mathcal{B}}$ as oriented $4$-cycles.

As in the proof of Proposition 4.4 we claim there cannot exist an edge that
appears in three or more element of $\operatorname{\mathcal{B}}$. Here a pair
of distinct $4$-cycles can share two edges but these must be adjacent and
furthermore must be oriented in opposite directions in the respective
$4$-cycles. Hence three or more $4$-cycles involving these edges would lead to
a positive entry in the corresponding entry of
$\operatorname{\mathcal{L}}^{*}$.

Suppose $m=3$. For this case the matrix $\iota$ must have $8(n-1)$ nonzero
entries, with at most two of these entries in each of the $3n$ rows. Hence we
need $8n-8<6n$, so that $n<4$. Note that if $8n-8=6n$ then we get a linear
dependence among the elements of $\operatorname{\mathcal{B}}$.

For the case of $K_{3,3}$ the total number of nonzero entries in $\iota$ is
$16$, and there are $9$ edges in $K_{3,3}$. We’ve seen that the maximum number
of possible entries in any row of $\iota$ row is two (a $-1$ and a $1$), and
in order for $\operatorname{\mathcal{B}}$ to span the space
$\operatorname{\mathcal{F}}(G)$ we need $\iota$ to have a nonzero entry in
every row. Hence the only way to distribute the entries would be to have 2
entries in 7 rows and only one in the last 2 rows. Similar to the $K_{5}$ case
in 4.4, this can not create a cycle basis.

For the general case of $K_{m,n}$ with $4\leq m,n$ we need $4(mn-
m-n+1)=4mn-4m-4n+1$ nonzero entries in $\iota$, with at most two of these
appearing in any of the $mn$ rows. Hence we need $2mn<4n+4m-4$, so that
$(m-2)(n-2)<2$. Since $m,n\geq 4$ it is evident that this is never true. The
result follows. ∎

##  5\. Further thoughts

In this section we discuss some open questions and further directions. Some of
these have been mentioned above but we collect them here for convenience. The
most natural open question concerns the existence of a cycle $M$-basis
consisting of only circuits of the underlying graph.

######  Question 5.1.

Does any graph $G$ admit a circuit $M$-basis?

We suspect that the answer is no but have not found any counterexamples. We
have seen that all planar graphs, as well as $K_{5}$ and $K_{3,3}$, admit
circuit $M$-bases but we know of no other nontrivial examples. Is there a
naturally defined nontrivial class of graphs that admit a circuit $M$-basis?
Do complete graphs have this property?

######  Question 5.2.

Can one find an explicit bijection between the set of $z$-superstable and the
set $\tau(G)$ of spanning trees?

After fixing a cycle $M$-basis for a graph $G$, we have seen that the set of
$z$-superstable configurations represent a unique member of each equivalence
class determined by the underlying dual Laplacian
$\operatorname{\mathcal{L}}^{*}$. These configurations can be interpreted as
‘minimal energy’ representatives, or as cycle configurations for which no
superset can be fired. From Proposition 3.14 we see that the set of such
$z$-superstable configurations is in bijection with the set of spanning trees.
One wonders if an explicit bijection exists in the spirit of those constructed
in (for instance) [7].

######  Question 5.3.

What do the $z$-superstable configurations count?

Since the set of $z$-superstable configurations are in bijection with the set
of spanning trees, a natural question to ask is whether the collection of such
configurations is related to some statistic on the spanning trees.

In the case of $z$-superstable configurations arising from the dual graph $G*$
of a planar graph $G$ the number of $z$-superstable configurations of degree
$d$ are given by the number of spanning trees of $G$ with $d$ internally
passive edges. This is not the case for general graphs, and for instance the
degree of a $z$-superstable configuration can be larger than $|V(G)|-1$ (see
Example 3.15). However, it would be interesting to see if there is some notion
of activity of a spanning tree, dependent on the choice of $M$-basis, that
does correspond to the degree sequence of $z$-superstable configurations.

Regardless of the status of Question 5.1, we can ask if there are any
restrictions on the sizes of elements in any circuit $M$-basis for certain
families of graphs. We have seen that the graphs $K_{n}$ and $K_{m,n}$ cannot
have circuit bases with all ‘small’ circuits. Experimental evidence suggests
the following.

######  Question 5.4.

Is it true that any circuit $M$-basis for $K_{n}$ contains a $n$-cycle? Does
any circuit basis for $K_{m,n}$ with $m\leq n$ contain a $2m$-cycle?

Finally, much of the theory and applications discussed above make sense in the
more general of setting of (totally) unimodular matrices and unimodular
matroids. It would be interesting to pursue these constructions in those
settings.

## References

  * [1] R. Bacher, P. De La Harpe, and T. Nagnibeda, _The lattice of integral flows and the lattice of integral cuts on a finite graph_ , Bull. Soc. Math. France 125, no. 2 (1997), pp. 167–198.
  * [2] P. Bak, C. Tang, K. Wiesenfeld, _Self-organized criticality_ , Phys. Rev. A, 38 (1988), pp. 364–374.
  * [3] M. Baker and S. Norine, _Riemann-Roch and Abel-Jacobi theory on a finite graph_ , Adv. Math., 215(2007), pp. 766–788.
  * [4] M. Baker, F. Shokrieh, _Chip-firing games, potential theory on graphs, and spanning trees_ , J. Combin. Theory Ser. A, 120, Issue 1 (2013), pp. 164–182.
  * [5] N. Biggs. _Chip-firing and the critical group of a graph_ , J. Algebraic Combin., 9 (1999), pp. 25–45.
  * [6] A. Björner, L. Lovász, P. W. Shor, _Chip-firing games on graphs_ , European J. Combin., 12 (1991), pp. 283–291.
  * [7] D. Chebikin and P. Pylyavskyy, _A family of bijections between G-parking functions and spanning trees_ , J. Combin. Theory Ser. A, 110, Issue 1 (2005), pp. 31–41.
  * [8] R. Cori and D. Rossin, _On the Sandpile Group of Dual Graphs_ , Europ. J. Combinatorics, 21 (2000), pp. 447–459.
  * [9] S. Corry and D. Perkinson, _Divisors and Sandpiles: An Introduction to Chip-Firing_ , American Mathematical Society, 2018.
  * [10] D. Dhar _Self-organized critical state of sandpile automaton models_ , Phys. Rev. Lett., 64 (14) (1990), pp. 1613–1616.
  * [11] Andrei Gabrielov, _Abelian avalanches and Tutte polynomials_ , Phys. A 195 (1993), no. 1-2, pp. 253–274.
  * [12] J. Guzmán, C. Klivans, _Chip-firing and energy minimization on M-matrices_ , J. Combin. Theory Ser. A, 132, (2015), pp. 14–31.
  * [13] B. Jacobson, _Critical groups of graphs_ , unpublished thesis.
  * [14] C. Klivans, _The Mathematics of Chip-Firing_ , Chapman & Hall / CRC Press, 2018.
  * [15] C. Merino, _Chip firing and the tutte polynomial_ , Ann. Comb. 1 (1997) pp. 253–259.
  * [16] R.J. Plemmons, _M-matrix characterizations. I. Nonsingular M-matrices_ , Linear Algebra Appl., 18 (2) (1977), pp. 175–188.
  * [17] W. A. Stein et al. Sage Mathematics Software (Version 9.0), The Sage Development Team, 2020, http://www.sagemath.org.
  * [18] M. Wood, _The distribution of sandpile groups of random graphs_ , J. Amer. Math. Soc. 30 (2017), pp. 915–958.