revtex4-1Repair the float

††thanks: Corresponding author.††thanks: Corresponding author.

# Local stability of spheres via the convex hull and the radical Voronoi
diagram

Peter K. Morse<EMAIL_ADDRESS>Department of Chemistry, Princeton
University, Princeton, NJ 08544  Department of Physics, Princeton University,
Princeton, NJ 08544  Princeton Institute of Materials, Princeton University,
Princeton, NJ 08544  Eric Corwin<EMAIL_ADDRESS>Department of Physics,
University of Oregon, Eugene, OR 97403  Materials Science Institute,
University of Oregon, Eugene, OR 97403

(August 28, 2024)

###### Abstract

Jamming is an emergent phenomenon wherein the local stability of individual
particles percolates to form a globally rigid structure. However, the onset of
rigidity does not imply that every particle becomes rigid, and indeed some
remain locally unstable. These particles, if they become unmoored from their
neighbors, are called rattlers, and their identification is critical to
understanding the rigid backbone of a packing, as these particles cannot bear
stress. The accurate identification of rattlers, however, can be a time-
consuming process, and the currently accepted method lacks a simple geometric
interpretation. In this manuscript, we propose two simpler classifications of
rattlers based on the convex hull of contacting neighbors and the maximum
inscribed sphere of the radical Voronoi cell, each of which provides geometric
insight into the source of their instability. Furthermore, the convex hull
formulation can be generalized to explore stability in hyperstatic soft sphere
packings, spring networks, non-spherical packings, and mean-field non-central-
force potentials.

##  I Introduction

A rigid structure is one which holds its shape when perturbed infinitesimally.
If this structure consists of particles, this rigid structure is said to be
jammed Stillinger  _et al._ (1964); Liu and Nagel (1998); Torquato  _et al._
(2000); O’Hern  _et al._ (2003); Torquato and Stillinger (2010); Liu and Nagel
(2010); van Hecke (2009). While the system as a whole may be rigid, local
regions of it may still be unconstrained. The particles—or clusters of
particles—making up these locally unconstrained regions are generally termed
“rattlers” Stillinger  _et al._ (1964); Speedy (1998); rat  and are removed
from the consideration of the structure for many analyses.

The rigorous rattler detection scheme in the literature Donev  _et al._ (2004)
relies on linear programming and is both computationally expensive and lacks a
simple geometric interpretation. Another, based on an event-driven packing
protocol, gives direct physical meaning to rattler detection by using a
stability analysis to systematically prune compressive forces, leaving
rattlers fully unconstrained Lerner  _et al._ (2013). However, this method
scales poorly with system size and dimension, as it requires matrix inversion.
These methods are, however, exact, and the resulting stable networks which
they find are identical. In light of the complexity of these algorithms, a
naive rattler detection scheme via constraint counting has proliferated and
been used widely as a proxy, despite its shortcomings. The naive algorithm
exploits the fact that the minimum number of constraints necessary to
stabilize a particle in $d$ dimensions is $d+1$. Thus, the number of contacts
on each particle is counted, and those with fewer than $d+1$ contacts are
deemed rattlers. Some (but not all) of these proxy methods apply this
criterion recursively, thus more closely approximating the true stable
network. However, this method cannot account for the presence of particles
with at least $d+1$ stable contacting neighbors which are nevertheless not
geometrically constrained.

Here, we present an alternative scheme for identifying rattlers that is
intuitive, efficient, and physically meaningful. In fact, we have been using
it for some time without realizing that it was not yet present in the
literature Corwin  _et al._ (2013); Morse and Corwin (2014); Charbonneau  _et
al._ (2015, 2016); Morse and Corwin (2016a, b, 2017); Charbonneau  _et al._
(2019); Hagh  _et al._ (2019); Sartor and Corwin (2020); Dennis and Corwin
(2020); Morse  _et al._ (2020); Richard  _et al._ (2020); Morse  _et al._
(2021); Rissone  _et al._ (2021); Sartor  _et al._ (2021); Charbonneau and
Morse (2021); Charbonneau  _et al._ (2021); Hagh  _et al._ (2022); Stanifer
and Lisa Manning (2022); Sartor and Corwin (2022); Dennis  _et al._ (2022);
Charbonneau and Morse (2023). Our method is based on a fundamental link
between local rigidity and the local geometry of force carrying contacts, and
implemented through the computation of the convex hull of the set of
contacting particles. The stable network obtained by this algorithm is
identical to that found in Refs. Donev  _et al._ (2004); Lerner  _et al._
(2013).

The central thrust of our algorithm is based on a comment within Ref. Donev
_et al._ (2004), namely that a sphere can only be locally rigid if it has
greater than $d+1$ non-cohemispheric contacts. While the authors of Ref. Donev
_et al._ (2004) note that simple constructions can be done in low spatial
dimensions (a method adopted in Refs. Wentworth-Nice  _et al._ (2020); Zhang
_et al._ (2022)), ours is a dimensionally independent construction: a particle
whose center is $\mathbf{r}_{0}$ is locally stable if the sum of all forces
acting on it is zero, and if the surface of the convex hull of the particle’s
center and the centers of all of its contacting neighbors
$\\{\mathbf{r}_{i}\\}$ does not include $\mathbf{r}_{0}$, i.e.
$\mathbf{r}_{0}\notin\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$,
where $\partial\mathrm{Conv}$ is the surface of the convex hull. We also prove
a related theorem, which can be shown to be equivalent to this, which states
that a particle is locally stable if the maximum inscribed sphere of its
radical Voronoi cell is unique and identical to the particle itself.

The rest of this article is structured as follows. In Sec. II, we provide
definitions for the generalized packing models that we can consider and a
series of mathematical definitions which will allow us to prove the two main
theorems. In Sec. III, we provide a formal proof that each construction finds
the correct stable network. In Sec. IV, we address computational complexity,
noting that even in the worst case scenario, the convex hull algorithm is
faster than the linear programming algorithm in $d<6$. We conclude in Sec. V
by discussing extensions of this construction to other models.

##  II Definitions

In the following, bold letters denote vectors in $\mathbb{R}^{d}$,
$\mathbf{0}$ represents the zero-vector, $\mathbf{a}\cdot\mathbf{b}$ denotes
the dot product between vectors $\mathbf{a}$ and $\mathbf{b}$,
$\\{\mathbf{r}_{i}\\}$ denotes a finite set of points, where each point is
represented by a vector from the origin, and
$\\{\mathbf{r}_{i}\\}\setminus\mathbf{r}_{0}$ denotes the set
$\\{\mathbf{r}_{i}\\}$ excluding the point $\mathbf{r}_{0}$. All definitions
assume the standard Euclidean distance metric on $\mathbb{R}^{d}$, where the
distance between points $\mathbf{a}$ and $\mathbf{b}$ is denoted
$|\mathbf{a}-\mathbf{b}|$. To define our packing, and to aid in later
definitions and theorems, we define both open and closed balls.

######  Definition 1.

An open ball of radius $\sigma$ around $\mathbf{s}$ is defined as the set of
points contained within a distance $\sigma$ of $\mathbf{s}$. The notation we
will use is
${B_{\sigma}(\mathbf{s})\equiv\\{\mathbf{y}:|\mathbf{s}-\mathbf{y}|<\sigma\\}}$.

######  Definition 2.

A closed ball of radius $\sigma$ around $\mathbf{s}$ is defined as the set of
points contained within and including a distance $\sigma$ of $\mathbf{s}$. The
notation we will use is
${\overline{B}_{\sigma}(\mathbf{s})\equiv\\{\mathbf{y}:|\mathbf{s}-\mathbf{y}|<\sigma\\}}$.

We thus consider particles defined by
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})$ with a non-dimensional overlap
between particles $i$ and $j$ defined as

$h_{ij}\equiv
1-\frac{\absolutevalue{\mathbf{r}_{i}-\mathbf{r}_{j}}}{\sigma_{i}+\sigma_{j}},$
(1)

subject to an additive potential $U=\sum_{ij}u(h_{ij})$ where contacts
($h_{ij}\geq 0)$ coincide with the potential cutoff, i.e. $u(h_{ij}\leq 0)=0$.
This form includes (but is not limited to) standard soft-sphere contact power
law potentials where $u(h_{ij}>0)\propto h_{ij}^{\gamma}$ for $\gamma>0$
($\gamma=2$ for Hookean spheres, and $\gamma=2.5$ for Hertzian spheres) and
hard spheres, where $u(h_{ij}>0)=\infty$.

From this, the force on particle $i$ from particle $j$ can be defined as

$\mathbf{f}_{ij}\equiv\nabla u(h_{ij})=\absolutevalue{\nabla
u(h_{ij})}\frac{\mathbf{r}_{j}-\mathbf{r}_{i}}{\absolutevalue{\mathbf{r}_{j}-\mathbf{r}_{i}}}.$
(2)

Here the only salient feature is that the force points towards the particle
center from the point of contact. Unless otherwise mentioned, we consider only
packings which are in a local energy minimum, such that the sum of forces
acting on each particle is zero. Extensions to non energy minimized packings
will be considered in Sec. V.

######  Definition 3 (Adapted from Ref. Donev  _et al._ (2004)).

A particle is locally stable if the sum of the forces acting on it is zero and
the forces acting on it span $\mathbb{R}^{d}$. Particles which are not locally
stable are called unstable.

In an effort to make this work as self contained as possible, we have compiled
a list of the mathematical definitions necessary to follow the theorems and
proofs of Sec. III such that only basic knowledge of set theory and linear
algebra will be prerequisite. The definitions are adapted from Refs. Ziegler
(1995); Munkres (2000); Grünbaum  _et al._ (2003).

Figure 1: Here we demonstrate the concept of an extreme point by examining
three red particles labelled (a-c). While this example is embedded in $d=2$,
the demonstration extends naturally to higher dimensions, replacing lines with
$(d-1)$-planes. (a) No line can be drawn which separates the particle from all
other particles, so (a) is not an extreme point. (b) A line can be drawn which
separates the particle from all other points, and it is thus an extreme point
and will be shown to be on the surface of the convex hull. (c) No line can be
drawn which separates the particle from all other particles, so it is not an
extreme point. However, a line exists which contains the particle and which
divides space such that all particles exist (inclusively) in one of its half
spaces, thus the point is on the surface of the convex hull.

######  Definition 4.

An extreme point $\mathbf{r}_{0}$ of the finite set $\\{\mathbf{r}_{i}\\}$ is
a point which can be separated from all other points by a $(d-1)$-plane. Thus
there exists a vector $\mathbf{a}\in\mathbb{R}^{d}$ with at least one non-zero
element and $b\in\mathbb{R}$ for which ${\mathbf{a}\cdot\mathbf{r}_{0}-b>0}$
while ${\mathbf{a}\cdot\mathbf{r}_{j}-b\leq 0}$ for all
${\mathbf{r}_{j}\in\\{\mathbf{r}_{i}\\}\setminus\mathbf{r}_{0}}$. An
illustration of both extreme and non-extreme points is given in Fig. 1.

Remark: In our proofs, we only need the extreme points of finite sets. The
concept of an extreme point can of course be generalized to infinite sets
Munkres (2000), but this makes several of the theorems unwieldy. The
definition used here is non-standard but reduces to the common definition in
the case of finite sets.

Figure 2: A simple illustration that a convex set containing points
$\mathbf{a}$ and $\mathbf{b}$ contains all points on a straight line between
them.

######  Definition 5.

A set $K\subset\mathbb{R}^{d}$ is convex if for all $\mathbf{a},\mathbf{b}\in
K$, $\mathbf{c}=(t-1)\mathbf{a}+t\mathbf{b}\in K$ for all $t\in[0,1]$. Put
simply, if $\mathbf{a}$ and $\mathbf{b}$ are in $K$, then $K$ is convex if
every point $\mathbf{c}$ along the straight line between $\mathbf{a}$ and
$\mathbf{b}$ is also in $K$. This is illustrated in Fig. 2.

######  Definition 6.

From Ref. Grünbaum  _et al._ (2003), a compact convex set
${K\subset\mathbb{R}^{d}}$ is a convex polytope if the extreme points of $K$
form a finite set. In this work, all instances of the word polytope are
implied to be convex.

######  Definition 7.

The surface $\partial K$ of a polytope $K$ is defined as the infinite set of
points $\mathbf{s}\in K$ for which there exists $\mathbf{s}_{\mathrm{out}}\in
B_{\sigma}(\mathbf{s})$ where $\mathbf{s}_{\mathrm{out}}\notin K$ for all
$\sigma$.

Figure 3: Here we demonstrate the convex hull (orange) of a set of points.
Points on the surface of the convex hull are colored black, while points not
on the surface of the convex hull are in teal.

######  Definition 8.

The convex hull of a set of points $\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ is
the unique closed $d$-dimensional polytope containing all points
$\\{\mathbf{r}_{i}\\}$ whose vertices are members of $\\{\mathbf{r}_{i}\\}$.
The surface of the convex hull is denoted
$\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ and is shown visually in Fig. 3.

######  Definition 9.

For a sphere given by $\overline{B}_{\sigma}(\mathbf{r})$, the points
$\\{\mathbf{b}_{i}\\}\subset\partial\overline{B}_{\sigma}(\mathbf{r})$ are
cohemispheric if there exists $\mathbf{a}\in\mathbb{R}^{d}$ with at least one
non-zero element, where $\mathbf{a}\cdot(\mathbf{b}_{i}-\mathbf{r})\geq 0$ for
all $i$. Similarly, forces $\\{\mathbf{f}_{i}\\}$ are cohemispheric if
$\mathbf{a}\cdot\mathbf{f}_{i}\geq 0$ for all $i$. If no such $\mathbf{a}$
exists, the points or forces are non-cohemispheric.

Figure 4: The maximum inscribed sphere (teal) of a convex polytope (purple) in
$d=2$. Contact points between the MIS and the polytope are shown with stars.
(a) The generic case with no symmetries has $d+1$ contact points between the
polytope and the MIS. (b) When two of the contacting surfaces are parallel, it
is possible to have an MIS with only 2 contacts in any dimension. (c) Highly
symmetric polytopes (regular ones, as shown here, or those near jamming), may
have MIS which have greater than $d+1$ contacts with the polytope.

######  Definition 10.

For a polytope $K$, the maximum inscribed sphere $M(K)$ is the largest closed
ball fully contained in $K$. That is,
$M(K)=\mathrm{max}_{\sigma}[\overline{B}_{\sigma}(\mathbf{r}):\overline{B}_{\sigma}(\mathbf{r})\subset
K]$. An illustration of the concept, including generic, degenerate, and highly
symmetric cases is given in Fig. 4. We use MIS as an abbreviation when not
referring to a specific $M(K)$.

######  Definition 11.

In a packing of particles with positions $\\{\mathbf{r}_{i}\\}$, the Voronoi
cell of particle 0 is the set
$V(\mathbf{r}_{0})=\\{\mathbf{y}:|\mathbf{y}-\mathbf{r}_{0}|\leq|\mathbf{y}-\mathbf{r}_{i}|\,\,\forall
i\\}$.

######  Definition 12.

The power of a point $\mathbf{c}\in\mathbb{R}^{d}$ with respect to a sphere
with center $\mathbf{r}$ and radius $\sigma$ is given by
$\Pi_{\mathbf{r},\sigma}(\mathbf{c})=|\mathbf{r}-\mathbf{c}|^{2}-\sigma^{2}$.
Points on the interior of the sphere have negative power, points on the
surface of the sphere have zero power, and points outside of the sphere have
positive power.

######  Definition 13.

In a packing of particles with positions $\\{\mathbf{r}_{i}\\}$ and radii
$\sigma_{i}$, the radical Voronoi cell of particle 0 is the set
$R(\mathbf{r}_{0})=\\{\mathbf{y}:\Pi_{\mathbf{r}_{0},\sigma_{0}}(\mathbf{y})\leq\Pi_{\mathbf{r}_{i},\sigma_{i}}(\mathbf{y})\,\,\forall
i\\}$.

Trivially, we see that if all particles are the same size (i.e.
$\sigma_{i}=\sigma$ for all $i$), then the radical Voronoi cell reduces to
that of the standard Voronoi cell. Both the radical Voronoi cell and, by
extension, the Voronoi cell are convex polytopes, and it is from the
definitions that these cells tessellate space, i.e. there is no point in space
which is not contained in the radical Voronoi cell of a particle, and the only
points which can be contained in multiple radical Voronoi cells are on the
shared surfaces of two or more cells.

##  III Proofs of the stability theorems

In this section, we provide proofs of the two main stability theorems,
labelled Theorem 6 (Sec. III.1) and Theorem 10 (Sec. III.2). While Sec. III.1
is entirely self contained, Sec. III.2 uses theorems from Sec. III.1. Some of
the theorems are elementary or have been proven by simpler means elsewhere,
but we formulate our own versions here, as we believe that they help to build
the physical intuition for the main theorems.

###  III.1 Stability via the convex hull

######  Theorem 1 (The Krein-Milman theorem Krein and Milman (1940)).

A compact convex subset of a Hausdorff locally convex topological vector space
is equal to the closed convex hull of its extreme points.

The proof of this theorem is given in Ref. Krein and Milman (1940). For the
purposes of this work, we will use the fact that the standard vector space on
$\mathbb{R}^{d}$ with the Euclidean distance metric and standard inner product
is a Hausdorff locally convex topological vector space. For clarification of
these terms, we suggest any standard textbook on topology (for example, Ref.
Munkres (2000)).

######  Corollary 1.1.

The convex hull of a set of points $\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ is
equal to the closed $d$-dimensional polytope whose vertices are the extreme
points of $\\{\mathbf{r}_{i}\\}$.

###### Proof.

Given that the standard vector space on $\mathbb{R}^{d}$ is a Hausdorff
locally convex topological vector space, the Krein-Milman theorem states that
a closed convex polytope is the convex hull of its extreme points, which for a
convex polytope are its vertices. ∎

######  Corollary 1.2.

If $\mathbf{r}_{0}$ is an extreme point of $\\{\mathbf{r}_{i}\\}$ then
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$.

###### Proof.

To prove that $\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$,
we must show that there exists a point $\mathbf{s}_{\mathrm{out}}\in
B_{\sigma}(\mathbf{r}_{0})$ such that
$\mathbf{s}_{\mathrm{out}}\notin\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$. Because
$\mathbf{r}_{0}$ is an extreme point, there exists
$\mathbf{a}\in\mathbb{R}^{d}$ with at least one non-zero element and
$b\in\mathbb{R}$ such that $\mathbf{a}\cdot\mathbf{r}_{0}-b>0$ while
$\mathbf{a}\cdot\mathbf{r}_{j}-b\leq 0$ for all
$\mathbf{r}_{j}\in\\{\mathbf{r}_{i}\\}\setminus\mathbf{r}_{0}$. We can thus
construct
$\mathbf{s}_{\mathrm{out}}=\mathbf{r}_{0}+\frac{\sigma\mathbf{a}}{2\absolutevalue{\mathbf{a}}}$,
for which $|\mathbf{r}_{0}-\mathbf{s}_{\mathrm{out}}|=\frac{\sigma}{2}$, and
thus $\mathbf{s}_{\mathrm{out}}\in B_{\sigma}(\mathbf{r}_{0})$. By
construction, $\mathbf{s}_{\mathrm{out}}$ is an extreme point of the set
$\\{\mathbf{s}_{\mathrm{out}},\mathbf{r}_{i}\\}$, and thus by Corollary 1.1,
$\mathbf{s}_{\mathrm{out}}\notin\mathrm{Conv}(\mathbf{r}_{i})$. This statement
is true for any value of $\sigma$, and thus
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$. ∎

######  Theorem 2.

A full dimensional convex polytope is equivalently defined by either its
vertices (V-Representation) or the intersection of half-planes representing
its surface (H-Representation).

The proof of this theorem is contained in standard texts on convex polytopes,
for example following the proofs of Theorems 3.1.1 and 3.1.2 of Ref. Grünbaum
_et al._ (2003) or Theorem 1.1 of Ref. Ziegler (1995). The theorem only
applies to full dimensional polytopes (i.e. ones which are $d$-dimensional
objects), but if the polytope is a $d^{\prime}$ dimensional object, where
$d^{\prime}\neq d$, it is sufficient for our purposes to consider the
V-Representation and the H-Representation in $\mathbb{R}^{d^{\prime}}$, in
which the polytope is full dimensional.

######  Corollary 2.1.

A point $\mathbf{r}_{0}$ which is contained on a $(d-1)$-plane, which defines
a halfspace containing all $\mathbf{r}_{i}$ is contained in
$\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$. That is, if there exists
$\mathbf{a}\in\mathbb{R}^{d}$ with at least one non-zero element and
$b\in\mathbb{R}$ such that $\mathbf{a}\cdot\mathbf{r}_{0}-b=0$ and
$\mathbf{a}\cdot\mathbf{r}_{j}-b\leq 0$ for all
$\mathbf{r}_{j}\in\\{\mathbf{r}_{i}\\}\setminus\mathbf{r}_{0}$, then
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$

###### Proof.

There are two cases here which need to be proven. If $\mathbf{r}_{0}$ is an
extreme point, then
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ by Corollary
1.2. If $\mathbf{r}_{0}$ is not an extreme point, then the half-plane
representation described here is equivalent to that defining the
H-Representation of a convex polytope, and thus
$\mathbf{r}_{0}\in\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ by Theorem 2. The
further statement that
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\\{\mathbf{r}_{i}\\})$ comes directly
from the definition of a halfspace. ∎

######  Theorem 3.

Any set of $d$ or fewer points on the surface of a sphere are cohemispheric.
That is, for a sphere centered at $\mathbf{r}_{0}$ with radius $\sigma_{0}$
and points $\\{\mathbf{c}_{i}\\}$ satisfying
$|\mathbf{c}_{i}-\mathbf{r}_{0}|=\sigma_{0}$, there exists a vector
$\mathbf{a}\in\mathbb{R}^{d}$ such that
$\mathbf{a}\cdot(\mathbf{c}_{i}-\mathbf{r}_{0})\geq 0$ for all $i$.

###### Proof.

Here, we can relax the condition $|\mathbf{c}_{i}-\mathbf{r}_{0}|=\sigma_{0}$
and prove a more general theorem. A hyperplane in $\mathbb{R}^{d}$ can always
be formed which passes through the $d$ contact points. That is, there exist
$\mathbf{a}^{\prime}\in\mathbb{R}^{d}$ and $b\in\mathbb{R}^{d}$ such that
$\mathbf{a}^{\prime}\cdot\mathbf{c}_{i}=b$ for all $\mathbf{c}_{i}$. Note that
if we construct a matrix $C$ with rows $\mathbf{c}_{i}$, then this hyperplane
is not unique if $\mathrm{det}(C)=0$, but any of the infinitely many solutions
will suffice.

We can define $b^{\prime}\in\mathbb{R}^{d}$ by
$\mathbf{a}^{\prime}\cdot\mathbf{r}_{0}=b^{\prime}$, then
$\mathbf{a}^{\prime}\cdot(\mathbf{c}_{i}-\mathbf{r}_{0})=b-b^{\prime}$. If
$b^{\prime}\leq b$, then $b-b^{\prime}\geq 0$, and we can take
$\mathbf{a}=\mathbf{a}^{\prime}$, whereupon the theorem is proven. If
$b^{\prime}>b$, then we can take $\mathbf{a}=-\mathbf{a}^{\prime}$, whereupon
the theorem is proven. ∎

From this, we note that the minimal number of points $\mathbf{c}_{i}$ for
which this theorem no longer holds is $d+1$. This is not to say that any $d+1$
points on the surface are non-cohemispheric (see for example, Fig. 5a), but to
state that the minimal number of points on a sphere which are non-
cohemispheric is $d+1$.

######  Theorem 4.

Given $\mathbf{f}_{i1}\in\\{\mathbf{f}_{i}\\}$ and
$\mathbf{a}\in\mathbb{R}^{d}$ where $\mathbf{f}_{i1}\neq\mathbf{0}$,
$\mathbf{a}\neq\mathbf{0}$, and $\mathbf{a}\cdot\mathbf{f}_{i1}\neq 0$,
$\sum_{i}\mathbf{f}_{i}=\mathbf{0}$ only if there exists
$\mathbf{f}_{i2}\in\\{\mathbf{f}_{i}\\}\setminus\mathbf{f}_{i1}$ such that
$\mathrm{sign}(\mathbf{a}\cdot\mathbf{f}_{i1})=-\mathrm{sign}(\mathbf{a}\cdot\mathbf{f}_{i2})$.

###### Proof.

Here we project $\mathbf{a}$ onto the sum yielding
${\sum_{i}(\mathbf{a}\cdot\mathbf{f}_{i})=\mathbf{a}\cdot\mathbf{f}_{i1}+\sum_{i\neq
i1}(\mathbf{a}\cdot\mathbf{f}_{i})=0}$. This last equality can only be true if
there is at least one element of the sum which is of the opposite sign of
$\mathbf{a}\cdot\mathbf{f}_{i1}$, implying that there exists
$\mathbf{f}_{i2}\in\\{\mathbf{f}_{i}\\}\setminus\mathbf{f}_{i1}$ such that
$\mathrm{sign}(\mathbf{a}\cdot\mathbf{f}_{i1})=-\mathrm{sign}(\mathbf{a}\cdot\mathbf{f}_{i2})$.
∎

This theorem is meant to be a vector extension of the trivial theorem that a
sum of numbers can only be zero if either all elements are zero, or if it
contains both positive and negative elements. Setting aside the null case,
this theorem simply states that a sum of vectors with at least one non-zero
element can only be zero if it contains positive and negative elements when
projected onto (almost) any axis. A mild caveat must be added, namely that the
projection is not onto a vector normal to a chosen non-zero vector in the set.
This caveat is only a formality as the projecting vector $\mathbf{a}$ is
arbitrary.

######  Corollary 4.1.

A particle with zero net force and at least $d+1$ non-cohemispheric non-zero
forces is locally stable.

Remark: This theorem applies more generally to both point particles and any
shape of particle with forces pointing towards its center of mass. Such a
particle will be stable to translations, but not to rotations.

###### Proof.

We label the set of non-zero forces $\\{\mathbf{f}_{i}\\}$ and the particle
center by $\mathbf{r}$. Note that the minimum number of vectors needed to span
$\mathbb{R}^{d}$ is $d$, so a particle is unstable with fewer than $d$ forces
acting upon it. Furthermore, a particle with $d$ forces acting upon it is
unstable by Theorem 3, as these forces are necessarily cohemispheric, and thus
there exists $\mathbf{a}\in\mathbb{R}^{d}$ such that
$\mathbf{a}\cdot\mathbf{f_{i}}\geq 0$ for all $i$. By Theorem 4,
$\sum\mathbf{f}_{i}\neq 0$ unless $\mathbf{f}_{i}=0$ for all $i$, and thus a
particle with $d$ non-zero forces acting upon it is unstable.

By definition, if there are $d+1$ non-cohemispheric non-zero forces, then no
$\mathbf{a}$ exists for which $\mathbf{a}\cdot\mathbf{f}_{i}\geq 0$ for all
$i$. Thus, for all $\mathbf{a}\in\mathbb{R}^{d}$ with at least one non-zero
element, Theorem 4 states that there will be positive and negative
projections, and thus the net force can sum to zero without all forces being
trivially zero, and thus the particle is locally stable. ∎

######  Theorem 5.

A particle with center $\mathbf{r}_{0}$ and with contacting particles centered
at $\\{\mathbf{r}_{i}\\}$ is unstable if
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$.

###### Proof.

We have two instances to prove. If $\mathbf{r}_{0}$ is an extreme point, then
by Definition 4, there exist $\mathbf{a}\in\mathbb{R}^{d}$ and
$b\in\mathbb{R}$ where $\mathbf{a}\cdot\mathbf{r}_{0}-b>0$ while
$\mathbf{a}\cdot\mathbf{r}_{j}-b\leq 0$ for all
$\mathbf{r}_{j}\in\\{\mathbf{r}_{i}\\}\setminus\mathbf{r}_{0}$. The contact
forces on $\mathbf{r}_{0}$ are all of the form
$\mathbf{f}_{j}=c_{j}(\mathbf{r}_{0}-\mathbf{r}_{j})$ with
$c_{j}\in\mathbb{R}$ and $c_{j}\geq 0$. Thus
$\sum_{j}\mathbf{f}_{j}=\sum_{j}c_{j}(\mathbf{r}_{0}-\mathbf{r}_{j})$. Taking
the projection on $\mathbf{a}$, we have
$\sum_{j}\mathbf{a}\cdot\mathbf{f}_{j}=\sum_{j}c_{j}(\mathbf{a}\cdot\mathbf{r}_{0}-\mathbf{a}\cdot\mathbf{r}_{j})$.
Depending on the sign of $b$, the non-zero terms are either all positive or
all negative, meaning that the sum cannot be zero unless all $c_{j}$ are zero.
Thus, by Theorem 4, either $\sum_{j}\mathbf{f}_{j}\neq\mathbf{0}$, or
$\mathbf{f}_{j}=0$ for all $j$. Either condition means that the particle is
unstable.

If $\mathbf{r}_{0}$ were not an extreme point, then the sum
$\sum_{j}\mathbf{a}\cdot\mathbf{f}_{j}$ could only be $0$ if
$\mathbf{a}\cdot\mathbf{r}_{j}-b=0$ for all $j$. These forces would then all
be co-hemispheric, and thus by Theorem 4, either
$\sum_{j}\mathbf{f}_{j}\neq\mathbf{0}$, or $\mathbf{f}_{j}=0$ for all $j$, and
thus the particle is unstable. ∎

Here we note that this is a sufficient condition for instability, and not a
necessary one. If $\mathbf{r}_{0}$ is out of force balance with neighboring
contacts $\\{\mathbf{r}_{i}\\}$, but
$\mathbf{r}_{0}\notin\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$, then
$\mathbf{r}_{0}$ is still unstable.

Figure 5: We test whether the blue particle is stable by looking at the convex
hull (red) of its own center and the centers of its stable neighboring
particles (in black). Note that there may be other contacts with the blue
particle which have been determined to be unstable and are thus not shown. In
(a) the blue particle is unstable, because its center lies on the surface of
the convex hull. In (b) the blue particle is stable, because its center is not
on the surface of the convex hull.

######  Theorem 6.

A particle with center $\mathbf{r}_{0}$ and with a non-empty set of stable
contacting particles centered at $\\{\mathbf{r}_{i}\\}$ is locally stable if
and only if
$\mathbf{r}_{0}\notin\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$
and the sum of forces acting on the particle is zero.

###### Proof.

The statement that $\mathbf{r}_{0}$ is locally stable if
$\mathbf{r}_{0}\notin\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$,
and the sum of all forces acting on the particle from $\\{\mathbf{r}_{i}\\}$
is zero follows a recursive application of Definition 3 and Theorem 5.

Next, we must prove that
$\mathbf{r}_{0}\notin\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$
with stable contacts $\\{\mathbf{r}_{i}\\}$ and zero net force implies that
$\mathbf{r}_{0}$ is locally stable and thus has a set of stable forces acting
on the particle centered at $\mathbf{r}_{0}$ which both span $\mathbb{R}^{d}$
and sum to zero. Because
$\mathbf{r}_{0}\notin\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{i}\\})$,
we know that $\mathbf{r}_{0}$ is neither an extreme point of the convex hull,
nor is it on the surface. Thus no $\mathbf{a}$ exists for which the contact
forces, labelled $\\{\mathbf{f}_{i}\\}$ have the property
$\mathbf{a}\cdot\mathbf{f}_{i}\geq 0$ for all $i$. These forces are thus non-
cohemispheric, and so from Theorem 3, there must be $d+1$ of them. And because
this particle has zero net force acting upon it, by Corollary 4.1, the
particle is locally stable.

An illustration of this theorem is given in Fig. 5. ∎

###  III.2 Stability via the radical Voronoi cell

######  Theorem 7.

If $i$ and $j$ are hard particles with centers $\mathbf{r}_{i}$ and
$\mathbf{r}_{j}$ and radii $\sigma_{i}$ and $\sigma_{j}$ and $h_{ij}=0$, then
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})$
contains exactly one point $\mathbf{c}_{ij}$ where $\mathbf{c}_{ij}\in\partial
R(\mathbf{r}_{i})$ and $\mathbf{c}_{ij}\in\partial R(\mathbf{r}_{j})$.

Figure 6: The radical Voronoi diagram is shown between contacting particles
$i$ and $j$ with contact point $\mathbf{c}_{ij}$. A second contact point
between the two particles $\mathbf{c}_{ij}^{\prime}$ is assumed, so that we
can show $\mathbf{c}_{ij}^{\prime}=\mathbf{c}_{ij}$ via the triangle
inequality.

###### Proof.

We define

$\mathbf{c}_{ij}=\mathbf{r}_{i}+\sigma_{i}\frac{\mathbf{r}_{j}-\mathbf{r}_{i}}{|\mathbf{r}_{j}-\mathbf{r}_{i}|}$
(3)

and note that $|\mathbf{c}_{ij}-\mathbf{r}_{i}|=\sigma_{i}$ so that
$\mathbf{c}_{ij}\in\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})$ and
$|\mathbf{c}_{ij}-\mathbf{r}_{j}|=|(\mathbf{r}_{j}-\mathbf{r}_{i})|-\sigma_{i}$.
We then note that $h_{ij}=0$ implies
$\sigma_{j}=|(\mathbf{r}_{j}-\mathbf{r}_{i})|-\sigma_{i}$, and thus
$|\mathbf{c}_{ij}-\mathbf{r}_{j}|=\sigma_{j}$ and so
$\mathbf{c}_{ij}\in\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})$.

To show that the intersection contains only one point, we assume that
$\mathbf{c}_{ij}^{\prime}\in\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})$
so that $|\mathbf{c}^{\prime}_{ij}-\mathbf{r}_{i}|\leq\sigma_{i}$ and
$|\mathbf{c}^{\prime}_{ij}-\mathbf{r}_{j}|\leq\sigma_{j}$, but
$\mathbf{c}_{ij}^{\prime}\neq\mathbf{c}_{ij}$ (as in Fig. 6). By the triangle
inequality,
$|\mathbf{r}_{j}-\mathbf{r}_{i}|\leq|\mathbf{r}_{i}-\mathbf{c}_{ij}^{\prime}|+|\mathbf{r}_{j}-\mathbf{c}_{ij}^{\prime}|$,
which becomes the degenerate statement
$\sigma_{i}+\sigma_{j}\leq\sigma_{i}+\sigma_{j}$. The degeneracy implies a
triangle of zero area, so that $\mathbf{c}_{ij}^{\prime}$ lies on the line
between $\mathbf{r}_{i}$ and $\mathbf{r}_{j}$, and by simple algebra, we find
that $\mathbf{c}_{ij}^{\prime}=\mathbf{c}_{ij}$. This is a contradiction, and
thus the intersection
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})$
contains only one point.

To show that $\mathbf{c}_{ij}\in R(\mathbf{r}_{i})$ and $\mathbf{c}_{ij}\in
R(\mathbf{r}_{j})$, we calculate the power of $\mathbf{c}_{ij}$ with respect
to each sphere. Here we find that
$\Pi_{\mathbf{r}_{i},\sigma_{i}}(\mathbf{c}_{ij})=\Pi_{\mathbf{r}_{j},\sigma_{j}}(\mathbf{c}_{ij})=0$.
The only lower power would be negative (interior of a sphere), and because
these are hard spheres, that is not possible. Thus, $\mathbf{c}_{ij}\in
R(\mathbf{r}_{i})$ and $\mathbf{c}_{ij}\in R(\mathbf{r}_{j})$. ∎

######  Corollary 7.1.

In a hard particle system,
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial R(\mathbf{r}_{i})$
contains only the contact points between particle $i$ and its contacting
neighbors, centered at $\\{\mathbf{r}_{j}\\}$.

###### Proof.

We know from Theorem 7 that
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial R(\mathbf{r}_{i})$
contains the contact points between particle $i$ and its contacting neighbors,
so we need now only show that it contains no other points. Suppose
$\mathbf{b}\in\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial
R(\mathbf{r}_{i})$ and that $\mathbf{b}\neq\mathbf{c}_{ij}$ from Eq. (3) for
any $j$. Points on $\partial R(\mathbf{r}_{i})$ have equal power with respect
to at least one other sphere, which we will generically call $\mathbf{r}_{j}$.
We have so far covered the case of zero power, and now consider points with
negative power. As per Definition 12, points of negative power are on the
interior of both spheres, i.e. $\mathbf{b}\in
B_{\sigma_{i}}(\mathbf{r}_{i})\cap B_{\sigma_{j}}(\mathbf{r}_{j})$, but
because $i$ and $j$ are hard spheres $B_{\sigma_{i}}(\mathbf{r}_{i})\cap
B_{\sigma_{j}}(\mathbf{r}_{j})=\emptyset$. Thus points of negative power are
not in the intersection $\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial
R(\mathbf{r}_{i})$. Points of positive power are not contained within
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})$ and are thus not in the
intersection $\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial
R(\mathbf{r}_{i})$. Therefore,
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\partial R(\mathbf{r}_{i})$
contains only the contact points between particle $i$ and its contacting
neighbors, centered at $\\{\mathbf{r}_{j}\\}$. ∎

######  Theorem 8.

In a convex region $K$, if $\overline{B}_{\sigma}(\mathbf{a})\subset K$ and
$\overline{B}_{\sigma}(\mathbf{b})\subset K$, then
$\overline{B}_{\sigma}(\mathbf{c})\subset K$ for all
$\mathbf{c}=(t-1)\mathbf{a}+t\mathbf{b}$ where $t\in[0,1]$.

Figure 7: An illustration of the fact that if two closed balls exist within a
convex region (centered at $\mathbf{a}$ and $\mathbf{b}$ respectively), every
closed ball on the line between the two is also contained in the region.

###### Proof.

From Definition 5, this property is true for every individual point within the
closed ball, so it is true for the closed ball itself. An illustration of the
concept is given in Fig. 7, where every ball contained on the line between
$\mathbf{a}$ and $\mathbf{b}$ is contained in the convex region if the closed
balls centered at $\mathbf{a}$ and $\mathbf{b}$ are contained in the region. ∎

Remark: We note that a further generalization of Theorem 8 is true when we
have different radii balls at the endpoints
$\overline{B}_{\sigma_{a}}(\mathbf{a})$ and
$\overline{B}_{\sigma_{b}}(\mathbf{a})$, where then the interpolated ball has
radius $\sigma_{c}=(t-1)\sigma_{a}+t\sigma_{b}$. This generalization is,
however, not necessary for our purposes and would potentially obscure the
results.

######  Theorem 9.

If $M(R(\mathbf{r}_{0}))$ is not unique in a hard particle system, then the
particle centered at $\mathbf{r}_{0}$ is not locally stable.

Figure 8: An example of a radical Voronoi diagram (orange lines) for a set of
particles (blue) which is highly degenerate. Here, because the MIS of each
particle is not unique despite having radii equal to that of the particles,
none of the particles are locally stable.

###### Proof.

We assume $M(R(\mathbf{r}_{0}))$ is not unique, such that
$\overline{B}_{\sigma}(\mathbf{r}_{1})\subset R(\mathbf{r}_{0})$ and
$\overline{B}_{\sigma}(\mathbf{r}_{2})\subset R(\mathbf{r}_{0})$ with
$\mathbf{r}_{1}\neq\mathbf{r}_{2}$ and there is no solution to
$\overline{B}_{\sigma}^{\prime}(\mathbf{r}_{3})\subset R(\mathbf{r}_{0})$
where $\sigma^{\prime}>\sigma$. We then assume that the particle centered at
$\mathbf{r}_{0}$ is locally stable and try to find a contradiction. If the
particle is stable, there exist at least $d+1$ non-cohemispheric points
$\mathbf{c}_{ij}$ given by Eq. (3) which, by Theorem 7, have the property
$\mathbf{c}_{ij}\in\overline{B}_{\sigma_{0}}(\mathbf{r}_{0})\cap\partial
R(\mathbf{r}_{0})$. Because the particle centered at $\mathbf{r}_{0}$ is fully
locally constrained, there exist no dilations or translations which maintain
the hard sphere condition. We now have two scenarios to consider, which each
contain a contradiction: $\sigma<\sigma_{0}$ and $\sigma\geq\sigma_{0}$.

If $\sigma<\sigma_{0}$, then neither $\overline{B}_{\sigma}(\mathbf{r}_{1})$
nor $\overline{B}_{\sigma}(\mathbf{r}_{1})$ represent the MIS, because
$\overline{B}_{\sigma_{0}}(\mathbf{r}_{0})\subset R(\mathbf{r}_{0})$ has a
larger radius. If $\sigma\geq\sigma_{0}$, then
$\overline{B}_{\sigma_{0}}(\mathbf{r}_{1})\subset\overline{B}_{\sigma}(\mathbf{r}_{1})\subset
R(\mathbf{\mathbf{r}_{0}})$. Theorem 8 states that all closed balls of radius
$\sigma_{0}$ on the straight line between $\mathbf{r}_{0}$ and
$\mathbf{r}_{1}$ are also contained in $R(\mathbf{r}_{0})$. However, because
the particle centered at $\mathbf{r}_{0}$ with radius $\sigma_{0}$ is stable,
no translations $T$ exists such that
$T(\overline{B}_{\sigma_{0}}(\mathbf{r}_{0}))\subset R(\mathbf{r}_{0})$.
Because no case relating $\sigma$ and $\sigma_{0}$ exists without a
contradiction, this implies that if $M(R(\mathbf{r}_{0}))$ is not unique in a
hard particle system, then the particle centered at $\mathbf{r}_{0}$ is not
locally stable. ∎

A packing with highly degenerate (non-unique) maximum inscribed spheres is
illustrated in Fig. 8, where clearly the particles are not stable.

######  Corollary 9.1.

If $M(K)$ is unique for a polytope $K$, then $M(K)\cap\partial K$ contains at
least $d+1$ non-cohemispheric points.

###### Proof.

If $M(K)$ is unique, then there are no translations represented by the
transformation $T$ which can be done such that $T(M(K))\subset K$. Thus $M(K)$
is fully constrained by the boundary $\partial K$. By Corollary 4.1, if we
impose a fictive force on $M(K)$ from each point of contact
$\\{\mathbf{c}_{i}\\}$ between $M(K)$ and $\partial K$, then there must be at
least $d+1$ non-cohemispheric $\mathbf{c}_{i}$ for $M(K)$ to be stable. Thus
$M(K)\cap\partial K$ contains at least $d+1$ non-cohemispheric points. ∎

######  Theorem 10.

In a packing of hard particles, a particle with center $\mathbf{r}_{0}$ and
radius $\sigma_{0}$ is locally stable if and only if $M(R(\mathbf{r}_{0}))$ is
unique and has center $\mathbf{r}_{0}$ and radius $\sigma_{0}$.

###### Proof.

First, we must prove that in a hard sphere system, a particle with center
$\mathbf{r}_{0}$ and radius $\sigma_{0}$ being locally stable implies that
$M(R(\mathbf{r}_{0}))$ is unique and has center $\mathbf{r}_{0}$ and radius
$\sigma_{0}$. Following the logic of the proof of Theorem 9, we assume that
$M(R(\mathbf{r}_{0}))$ has center $\mathbf{r}_{1}$ and radius $\sigma_{1}$
with $\mathbf{r}_{1}\neq\mathbf{r}_{0}$ and $\sigma_{1}\neq\sigma_{0}$ and
find a contradiction to show that $\mathbf{r}_{1}\neq\mathbf{r}_{0}$ and
$\sigma_{1}=\sigma_{0}$. If $\sigma_{1}<\sigma_{0}$, then this does not
correspond to the maximum inscribed sphere. If $\sigma_{1}\geq\sigma_{0}$,
then
$\overline{B}_{\sigma_{0}}(\mathbf{r}_{1})\subset\overline{B}_{\sigma_{1}}(\mathbf{r}_{1})\subset
R(\mathbf{r}_{0})$ and thus by Theorem 8
$\overline{B}_{\sigma_{0}}(\mathbf{r}_{2})\subset R(\mathbf{r}_{0})$ for all
$\mathbf{r}_{2}$ on a straight line between $\mathbf{r}_{0}$ and
$\mathbf{r}_{1}$. But because $\overline{B}_{\sigma_{0}}(\mathbf{r_{0}})$ is
locally stable, no translations or dilations exist which remain in
$R(\mathbf{r}_{0})$, so $\mathbf{r}_{1}=\mathbf{r}_{0}$ and
$\sigma_{1}=\sigma_{0}$.

Second, we must prove that in a hard sphere system, for a particle centered at
$\mathbf{r}_{0}$ with radius $\sigma_{0}$, $M(R(\mathbf{r}_{0}))$ being unique
and having center $\mathbf{r}_{0}$ and $\sigma_{0}$ implies that the particle
is stable. This follows immediately from Corollary 9.1, as the particle has
$d+1$ non-cohemispheric points of contact with $R(\mathbf{r}_{0})$, which by
Corollary 7.1, correspond to contacts with neighboring particles. Thus, by
Corollary 4.1, the particle centered at $\mathbf{r}_{0}$ with radius
$\sigma_{0}$ is stable. ∎

Figure 9: The radical Voronoi diagram (black lines) is computed for a set of
particles of different radii with contacts displayed as blue lines. All blue
particles, labelled $i$, have
$M(R(\mathbf{r}_{i}))=\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})$ and are thus
stable. The red particle, which we will call $0$ is a rattler, and its MIS is
shown as a dashed magenta line. We see clearly that
$M(R(\mathbf{r}_{0}))\neq\overline{B}_{\sigma_{0}}(\mathbf{r}_{0})$

##  IV Algorithmic complexity

Theorems 6 and 10 provide a natural recursive algorithm for determining the
stable set of particles in a packing, and, through its complement, the set of
rattlers. The algorithm begins with a tentative statement that all particles
are stable, and it loops over each particle testing for stability, taking the
function $\mathrm{isStable}(i)$ from either Theorem 6, Theorem 10, or Eq. 12
of Ref. Donev  _et al._ (2004), considering only the stable set of particles.
The algorithm ends when no changes are made to the stable list in a full loop.
The structure of the algorithm is similar to that of Ref. Donev  _et al._
(2004), and as expected, it produces an identical stable list.

Algorithm 1 Global Stability Algorithm

1:$i\in\mathrm{stableList}\,\,\forall i$

2:$\mathrm{unstableList}=\emptyset$

3:$\mathrm{flip}\leftarrow true$

4:while flip do

5: $\mathrm{flip}\leftarrow false$

6: for $i\in\mathrm{stableList}$ do

7: if isUnstable(i) then

8: $\mathrm{flip}\leftarrow true$

9: Move $i$ from stableList to unstableList

10: end if

11: end for

12:end while

13:Return stableList

The worst-case scenario for this algorithm is a packing in which only a single
particle is initially unstable, but its removal destabilizes one of its
neighbors, and so on. Such a situation will require $N$ iterations through the
algorithm, each of which takes $\mathcal{O}(N)$ time, yielding a total worst
case runtime of $\mathcal{O}(N^{2})$. We note, however, that no typical case
approaches this complexity. The method of Ref. Lerner  _et al._ (2013),
meanwhile, scales as at least $\mathcal{O}(d^{3}N^{3})$ Charbonneau and Morse
(2023).

The only difference between the methods of Ref. Donev  _et al._ (2004),
Theorem 6, and Theorem 10 is the speed of the function $\mathrm{isStable}(i)$.
For a particle with $n$ contacting particles (where $n\sim\mathcal{O}(d)$),
the linear programming method scales as $\mathcal{O}(n^{2+a})$ where
$a=\frac{1}{18}$ Jiang  _et al._ (2021) while the convex hull scales as
$\mathcal{O}(n^{\lfloor d/2\rfloor})$ in the worst case scenario, where
$\lfloor\cdot\rfloor$ is the floor function Chazelle (1993). The radical
Voronoi diagram for an individual cell can be computed in
$\mathcal{O}(n^{\lceil d/2\rceil})$ where $\lceil\cdot\rceil$ is the ceiling
function. Thus while the radical Voronoi method is slower than the convex hull
method in odd dimensions, it is of the same order in even dimensions. The
calculation of the MIS is then either a linear programming problem Gritzmann
and Klee (1992, 1993) or a minimization problem Morse and Corwin (2014) whose
complexity has not yet been interrogated. The worst case scenario then makes
this calculation the rate determining step, and it is thus no faster than the
linear programming methods of Ref. Donev  _et al._ (2004). By comparison, we
see that the convex hull algorithm is faster than the linear programming
algorithm for at least $d<6$.

##  V Further extensions

We have shown that the convex hull and the radical Voronoi cell can be used to
quickly determine the stability of individual spheres in a packing with only
minimal requirements on the interparticle potential. It is straightforward to
show that the construction can be applied more generally in a variety of
cases. Here, we list several:

  1. 1.

In a spring network under compression, an individual node is unstable if it is
on the surface of the convex hull of its connecting nodes.

  2. 2.

A particle of any shape is unstable if the only forces acting on it are point
forces directed towards its center of mass, and the center of mass is on the
surface of the convex hull of the contact points and the center of mass.

  3. 3.

In Mari-Kurchan (MK) interactions Mari  _et al._ (2009); Mari and Kurchan
(2011), where the distance between particles is given by
${h_{ij}^{\mathrm{MK}}=\frac{\absolutevalue{\mathbf{r}_{i}-\mathbf{r}_{j}+\mathbf{\Lambda}_{ij}}}{\sigma_{i}+\sigma_{j}}}$
where $\mathbf{\Lambda}_{ij}$ is a random vector with
$\mathbf{\Lambda}_{ij}=-\mathbf{\Lambda}_{ji}$, a particle $\mathbf{r}_{0}$
with contacts $\\{\mathbf{r}_{j}\\}$ is unstable if
$\mathbf{r}_{0}\in\partial\mathrm{Conv}(\mathbf{r}_{0},\\{\mathbf{r}_{j}+\mathbf{\Lambda}_{ij}\\})$.
This method was used in Ref. Charbonneau and Morse (2023). Note that this is
true despite not technically being a central force potential.

  4. 4.

Several recent studies have analyzed soft sphere systems during energy
minimization Charbonneau and Morse (2023); Stanifer and Lisa Manning (2022);
Nishikawa  _et al._ (2021, 2022); Manacorda and Zamponi (2022), wherein it may
be important to study the evolution of rattlers and stable subsystems. Here,
the convex hull theorem may be used, with the additional caveat that a
particle is only locally stable if the sum of all forces acting on it is zero,
and if the forces acting on it span $\mathbb{R}^{d}$.

  5. 5.

Following the logic of Sec. III.2, we conjecture that Theorem 10 also holds
for additively-weighted Voronoi cells and any generalization of Voronoi cells
$G$ for which the contact point of two hard spheres ($i$ and $j$) is contained
on the surface of the generalized Voronoi cell, i.e.
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})\in\partial
G(\mathbf{r}_{i})$ and
$\overline{B}_{\sigma_{i}}(\mathbf{r}_{i})\cap\overline{B}_{\sigma_{j}}(\mathbf{r}_{j})\in\partial
G(\mathbf{r}_{j})$. However, these cells are generically non-convex, and so
some of the tools we have used do not suffice.

These extensions show the utility of our methods, which extend beyond simple
sphere packings. It is our hope that this work not only provides a simple
computational tool, but helps to illuminate the interplay between geometry and
mechanical rigidity.

## References

  * Stillinger  _et al._ (1964) F. H. Stillinger, Jr., E. A. DiMarzio, and R. L. Kornegay, J. Chem. Phys. 40, 1564 (1964).
  * Liu and Nagel (1998) A. J. Liu and S. R. Nagel, Nature 396, 21 (1998).
  * Torquato  _et al._ (2000) S. Torquato, T. M. Truskett, and P. G. Debenedetti, Phys. Rev. Lett. 84, 2064 (2000).
  * O’Hern  _et al._ (2003) C. S. O’Hern, L. E. Silbert, A. J. Liu, and S. R. Nagel, Phys. Rev. E 68, 011306 (2003).
  * Torquato and Stillinger (2010) S. Torquato and F. H. Stillinger, Rev. Mod. Phys. 82, 2633 (2010).
  * Liu and Nagel (2010) A. J. Liu and S. R. Nagel, Annu. Rev. Condens. Matter Phys. 1, 347 (2010).
  * van Hecke (2009) M. van Hecke, J. Phys.: Condens. Matter 22, 033101 (2009).
  * Speedy (1998) R. J. Speedy, J. Phys. Condens. Matter 10, 4185 (1998).
  * (9) In the original paper describing these unstable particles, they were called inclusions Stillinger  _et al._ (1964), though the terminology eventually switched to rattlers Speedy (1998), while the concept remained fixed.
  * Donev  _et al._ (2004) A. Donev, S. Torquato, F. H. Stillinger, and R. Connelly, Journal of Computational Physics 197, 139 (2004).
  * Lerner  _et al._ (2013) E. Lerner, G. Düring, and M. Wyart, Comput. Phys. Commun. 184, 628 (2013).
  * Corwin  _et al._ (2013) E. I. Corwin, R. Stinchcombe, and M. F. Thorpe, Phys. Rev. E 88, 014102 (2013).
  * Morse and Corwin (2014) P. K. Morse and E. I. Corwin, Phys. Rev. Lett. 112, 115701 (2014).
  * Charbonneau  _et al._ (2015) P. Charbonneau, E. I. Corwin, G. Parisi, and F. Zamponi, Phys. Rev. Lett. 114, 125504 (2015).
  * Charbonneau  _et al._ (2016) P. Charbonneau, E. I. Corwin, G. Parisi, A. Poncet, and F. Zamponi, Phys. Rev. Lett. 117, 045503 (2016).
  * Morse and Corwin (2016a) P. K. Morse and E. I. Corwin, Soft Matter 12, 1248 (2016a).
  * Morse and Corwin (2016b) P. K. Morse and E. I. Corwin, J. Stat. Mech. 2016, 074009 (2016b).
  * Morse and Corwin (2017) P. K. Morse and E. I. Corwin, Phys. Rev. Lett. 119, 118003 (2017).
  * Charbonneau  _et al._ (2019) P. Charbonneau, E. I. Corwin, L. Fu, G. Tsekenis, and M. van der Naald, Phys. Rev. E 99, 020901 (2019).
  * Hagh  _et al._ (2019) V. F. Hagh, E. I. Corwin, K. Stephenson, and M. F. Thorpe, Soft Matter 15, 3076 (2019).
  * Sartor and Corwin (2020) J. D. Sartor and E. I. Corwin, Phys. Rev. E 101, 050902 (2020).
  * Dennis and Corwin (2020) R. C. Dennis and E. I. Corwin, Phys. Rev. Lett. 124, 078002 (2020).
  * Morse  _et al._ (2020) P. Morse, S. Wijtmans, M. van Deen, M. van Hecke, and M. L. Manning, Phys. Rev. Res. 2, 023179 (2020).
  * Richard  _et al._ (2020) D. Richard, M. Ozawa, S. Patinet, E. Stanifer, B. Shang, S. A. Ridout, B. Xu, G. Zhang, P. K. Morse, J.-L. Barrat, L. Berthier, M. L. Falk, P. Guan, A. J. Liu, K. Martens, S. Sastry, D. Vandembroucq, E. Lerner, and M. L. Manning, Phys. Rev. Mater. 4, 113609 (2020).
  * Morse  _et al._ (2021) P. K. Morse, S. Roy, E. Agoritsas, E. Stanifer, E. I. Corwin, and M. L. Manning, Proc. Natl. Acad. Sci. U.S.A. 118, e2019909118 (2021).
  * Rissone  _et al._ (2021) P. Rissone, E. I. Corwin, and G. Parisi, Phys. Rev. Lett. 127, 038001 (2021).
  * Sartor  _et al._ (2021) J. D. Sartor, S. A. Ridout, and E. I. Corwin, Phys. Rev. Lett. 126, 048001 (2021).
  * Charbonneau and Morse (2021) P. Charbonneau and P. K. Morse, Phys. Rev. Lett. 126, 088001 (2021).
  * Charbonneau  _et al._ (2021) P. Charbonneau, E. I. Corwin, R. C. Dennis, R. Díaz Hernández Rojas, H. Ikeda, G. Parisi, and F. Ricci-Tersenghi, Phys. Rev. E 104, 014102 (2021).
  * Hagh  _et al._ (2022) V. F. Hagh, S. R. Nagel, A. J. Liu, M. L. Manning, and E. I. Corwin, Proc. Natl. Acad. Sci. U.S.A. 119, e2117622119 (2022).
  * Stanifer and Lisa Manning (2022) E. Stanifer and M. Lisa Manning, Soft Matter 18, 2394 (2022).
  * Sartor and Corwin (2022) J. D. Sartor and E. I. Corwin, Phys. Rev. Lett. 129, 188001 (2022).
  * Dennis  _et al._ (2022) R. C. Dennis, V. F. Hagh, and E. I. Corwin, Phys. Rev. E 106, 044901 (2022).
  * Charbonneau and Morse (2023) P. Charbonneau and P. K. Morse, “Jamming, relaxation, and memory in a structureless glass former,”  (2023), arxiv:2304.11100 [cond-mat] .
  * Wentworth-Nice  _et al._ (2020) P. Wentworth-Nice, S. A. Ridout, B. Jenike, A. Liloia, and A. L. Graves, Soft Matter 16, 5305 (2020).
  * Zhang  _et al._ (2022) A. L. Zhang, S. A. Ridout, C. Parts, A. Sachdeva, C. S. Bester, K. Vollmayr-Lee, B. C. Utter, T. Brzinski, and A. L. Graves, Phys. Rev. E 106, 034902 (2022).
  * Ziegler (1995) G. M. Ziegler,  _Lectures on Polytopes_ , Graduate Texts in Mathematics No. 152 (Springer-Verlag, New York, 1995).
  * Munkres (2000) J. R. Munkres,  _Topology_ , 2nd ed. (Prentice Hall, Inc, Upper Saddle River, NJ, 2000).
  * Grünbaum  _et al._ (2003) B. Grünbaum, V. Kaibel, V. Klee, and G. M. Ziegler,  _Convex Polytopes_ , 2nd ed., Graduate Texts in Mathematics No. 221 (Springer, New York, 2003).
  * Krein and Milman (1940) M. Krein and D. Milman, Stud. Math. 9, 133 (1940).
  * Jiang  _et al._ (2021) S. Jiang, Z. Song, O. Weinstein, and H. Zhang, in  _Proceedings of the 53rd Annual ACM SIGACT Symposium on Theory of Computing_, STOC 2021 (Association for Computing Machinery, New York, NY, USA, 2021) pp. 823–832.
  * Chazelle (1993) B. Chazelle, Discrete Comput Geom 10, 377 (1993).
  * Gritzmann and Klee (1992) P. Gritzmann and V. Klee, Discrete Comput. Geom. 7, 255 (1992).
  * Gritzmann and Klee (1993) P. Gritzmann and V. Klee, Math. Program. 59, 163 (1993).
  * Mari  _et al._ (2009) R. Mari, F. Krzakala, and J. Kurchan, Phys. Rev. Lett. 103, 025701 (2009).
  * Mari and Kurchan (2011) R. Mari and J. Kurchan, J. Chem. Phys. 135, 124504 (2011).
  * Nishikawa  _et al._ (2021) Y. Nishikawa, A. Ikeda, and L. Berthier, J. Stat. Phys. 182, 37 (2021).
  * Nishikawa  _et al._ (2022) Y. Nishikawa, M. Ozawa, A. Ikeda, P. Chaudhuri, and L. Berthier, Phys. Rev. X 12, 021001 (2022).
  * Manacorda and Zamponi (2022) A. Manacorda and F. Zamponi, J. Phys. A: Math. Theor. 55, 334001 (2022).