# Topological Data Analysis for Particulate Gels

Alexander Smith${}^{\text{a}}$, Gavin J. Donley${}^{\text{b}}$, Emanuela Del
Gado${}^{\text{b,*}}$, and Victor M. Zavala${}^{\text{c,d,*}}$
${}^{\text{a}}$Department of Chemical Engineering and Material Science
University of Minnesota, Minneapolis, MN, USA
${}^{\text{b}}$Department of Physics
Georgetown University, Washington, DC, USA
${}^{\text{c}}$Department of Chemical and Biological Engineering
University of Wisconsin-Madison, Madison, WI, USA
${}^{\text{d}}$Mathematics and Computer Science Division
Argonne National Laboratory, Lemont, IL, USA
∗Corresponding Authors<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Soft gels, formed via the self-assembly of particulate organic materials,
exhibit intricate multi-scale structures that provides them with flexibility
and resilience when subjected to external stresses. This work combines
molecular simulations and topological data analysis (TDA) to characterize the
complex multi-scale structure of soft gels. Our TDA analysis focuses on the
use of the Euler characteristic, which is an interpretable and
computationally-scalable topological descriptor that is combined with
filtration operations to obtain information on the geometric (local) and
topological (global) structure of soft gels. We reduce the topological
information obtained with TDA using principal component analysis (PCA) and
show that this provides an informative low-dimensional representation of gel
structure. We use the proposed computational framework to investigate the
influence of gel preparation (e.g., quench rate, volume fraction) on soft gel
structure and to explore dynamic deformations that emerge under oscillatory
shear in various response regimes (linear, nonlinear, and flow). Our analysis
identifies specific scales and extents at which hierarchical structures in
soft gels are affected; moreover, correlations between structural deformations
and mechanical phenomena (such as shear stiffening) are explored. In summary,
we show that TDA facilitates the mathematical representation, quantification,
and analysis of soft gel structures, extending traditional network analysis
methods to capture both local and global organization.

##  1 Introduction

Particulate gel technology has emerged as a versatile and innovative solution
with widespread applications across diverse industries, including
pharmaceuticals, foods, and construction. In the pharmaceutical and medical
device industries, these gels play a pivotal role, serving as effective drug
delivery vehicles, encapsulation systems, and scaffolds. Their deformability
enables the seamless incorporation of various compounds, medications, or
diagnostic agents, contributing to enhanced bioavailability and therapeutic
efficacy in patient care [1, 2, 3]. Within the food industry, soft gel
encapsulation is employed for the delivery of functional ingredients,
vitamins, and supplements. The encapsulation of sensitive compounds not only
enhances their stability but also facilitates precise dosage control, meeting
the growing consumer demand for convenient and palatable nutritional solutions
[4, 5, 6, 7]. In the construction sector, the utility of soft gels extends to
the realm of binders, sealants, adhesives, and coatings, contributing to
improved performance and durability of building materials [8, 9, 10, 11].

Gel formation may proceed through phase separation, aggregation, and self-
assembly of polymers, colloids, or other soft matter components [12]. The
stochastic and multi-scale nature of these mechanisms produces hierarchical,
amorphous structures, and provides texture, elasticity, and stability [12, 13,
14, 15, 16]. The multi-scale structure of gels can be represented
mathematically as networks defined by the interactions of their building
blocks (e.g., particles, droplets, polymer aggregates). These network
representations provide insight into the deformation of soft gels, and aid in
understanding their mechanical failure [12, 13, 14, 17, 18]. For example,
while soft gels can accommodate large deformations, micro-scale cracks that
result from the accumulation of local tension can nucleate and trigger
macroscopic failure, in susceptible points determined by both the local
spatial distribution (geometry) and the large scale connectivity (topology) of
the gel network [12, 16]. The network structures are not static but highly
dynamically re-configurable; deformations, such as squeezing or stretching,
and deformation rates can modify the mechanical characteristics of a gel by
altering its network structure [16, 12, 13, 19, 20, 21]. In addition,
relaxation of the stresses accumulated in the gel through the gel formation
itself, under different environmental conditions, is another important source
of restructuring [12, 13, 19].

The adaptable nature of particulate gel materials, their flexibility and
tunability, provide exciting opportunities to precisely control and design
both geometry and topology of the particle networks, potentially opening up a
broad space for metamaterials design and discovery [22, 23, 24]. The geometry
of the structure (its local properties) can be controlled through the physical
chemistry of particles and solvent, particles surface properties, or depletion
interactions [25], and, in certain cases, is directly accessible through
confocal microscopy imaging[26, 27, 28]. However, getting access to the large-
scale (global) organization of the full network and its topology remain a
challenge. Computer simulations can provide this type of access, but many of
the methods used in understanding the structure of gels focus exclusively on
geometry or topology, but not both. From a geometric perspective, the
computation of the particle radial distribution functions provides a link to
scattering experiments and constitutes an excellent tool to quantify the local
structure [29, 30, 31] but does not capture the topology of the network [32,
33, 34]. Voronoi and Delaunay tesselation can capture properties such as
cavity size distribution [35, 36, 37], but provide limited quantification
tools of the topology of a gel structure.

Graph-theoretic methods use graph descriptors (e.g., modularity, average path
length, bond number, and minimal cycle basis) to quantify the structure of a
network [12, 38, 18]. These descriptors, however, are based on averages of
local topological structures and thus might fail to capture the multi-scale
nature of soft gels. Parametric methods based for example on distance
thresholds or $k$-nearest neighbors commonly used to identify a graph [39]
provide results which are often highly susceptible to the selection of
parameter (e.g., distance threshold, $k$-neighbors), which is problematic if
the parameter selection is not obvious from the data and rooted in the physics
of the system. It is also worth noting that some of these graph-theoretic
methods are computationally expensive. For example, the identification of a
minimal cycle basis for a graph is at least polynomial in the number of nodes
and edges present [40], requiring significant computational time for
simulations of soft gel structures from simulations which contain hundreds of
thousands of nodes and edges to be able to reproduce the heterogeneities
typical of real materials. It is also important to highlight that a
network/graph is inherently a 2-dimensional object; as such, these
representations might miss important information of the 3-dimensional
structure of soft gels.

In this work, we propose an approach, based on molecular simulations and
topological data analysis (TDA), that captures both the geometry and topology
of the 3-dimensional structure of soft particulate gels. Specifically, we
focus on the application of the Euler characteristic (EC) [32, 41], which is
an intuitive, interpretable, and computationally-scalable topological
descriptor. This approach quantifies the structure of the gel network using
basic topological invariants such as the number of connected components,
cycles (holes), and voids. We use a mathematical technique known as a
filtration to quantify how topological invariants emerge and disappear at
different length scales. This information is condensed in a topological
summary known as the EC curve. These approaches have been successfully applied
in the analysis of complex materials, simulations, and flow networks [42, 43,
44, 32, 41, 45, 46, 47]. The approaches allow us to bypass the need for
parameter selection, improving the robustness of our analysis in comparison to
other parametric methods (e.g., distance thresholding, neighbor selection).
Moreover, these approaches allow us to capture both topological (global) and
geometrical (local) characteristics of soft gel structures [48, 49].

We find that our proposed computational framework effectively represents,
quantifies, and summarizes the multi-scale nature of gel structures obtained
from molecular simulations. We reduce the topological information obtained
with TDA using principal component analysis (PCA) and show that this provides
an informative low-dimensional representation of gel structure. The
application of the proposed computational framework reveals the influence of
variations of the gelation kinetics on the local organization of the gels,
alongside with the impact of volume fraction changes on both local and global
structural aspects. Our analysis also identifies topological transitions
underpinning the onset of nonlinear response and flow when the gels are
subjected to large amplitude oscillatory tests, clearly disentangling
recoverable/unrecoverable changes from damage accumulation and revealing their
hierarchical nature; this ultimately reveals a connection between the
topological nature of the soft gel material properties and their rheological
response. The methods presented here are computationally scalable (analysis of
large simulations can be done in seconds-minutes and on a personal computer),
interpretable, and require minimal data. The proposed approach can also be
used for analysis of experimental data (such as confocal imaging), providing a
route for a more direct comparison of molecular simulation and experiments
[50, 51, 52]. All data and scripts needed for reproducing the results are
shared as open-source code.

##  2 Results

To understand and quantify the structure of particulate gel simulations, we
combine molecular simulations with the analysis of topological invariants of
the gel structures based on the EC, filtration operations, and PCA. These
methods are combined to represent and quantify the gel structures across
multiple length scales and capture the changes they undergo during shear. PCA
is used to visualize topological changes and provides insight into the impact
of diverse drivers on gel response (such as strain stiffening/hardening or
yielding). A summary of the computational workflow to quantify the
hierarchical topology of a soft gel is presented in Figure 1.

Figure 1: (left) Visualization of a soft gel for a single simulation frame. We
provide this visualization to illustrate how the various topological features:
components, holes, voids, could evolve during a filtration. (right) EC curve
for the gel simulation along with a visualization of the changes in the
simulation topology at various filtration diameters. The EC curve captures
various topological characteristics found at multiple scales within the soft
gel and show how the topological features such as holes ($\sim$2.3 particle
diameters) and large voids ($\sim$3.7 particle diameters) emerge. We also note
that there is minimal change in topology between $D=1$ and $D=2$ particle
diameters. This stable topology reflects the bonded structure of the gel, and
is often what is studied via network analysis tools. As we pass to $D>2$
particle diameters, we see large topological transitions in the material,
which now reflects the 3-dimensional structure of the gel and captures overlap
between particles that are not necessarily directly bonded.

###  2.1 Topological Analysis of Molecular Simulations

We use a particulate gel model that captures most aspects of microscopic
dynamics and rheology of this class of materials. Each gel is composed of $N$
particles, which interact via attractive short-range interactions of maximum
strength $\epsilon$, mediated by the solvent in which they are immersed and
through which their thermal motion is overdamped. Surface roughness, shape
irregularity and sintering processes can limit the relative motion of
particles as they aggregate in real materials. These effects are included in
the model through an angular modulation of the net attraction that introduces
a bending rigidity of the interparticle bonds [19, 53, 13, 54, 55]. We use a
cubic box of size $L$ with periodic boundary conditions and number density
$N/L^{3}$ which corresponds to an approximate solid volume fraction
$\phi=\frac{N\pi d^{3}/6}{(Ld)^{3}}$. For each value of $\phi$, various gel
microstructures are obtained by tuning the rate $\Gamma$ at which the relative
strength of the attractive interactions (with respect to $k_{B}T$) is
increased to induce gelation during the sample preparation (see Methods).
Finally, we subject a subset of the gels to large amplitude oscillatory shear
(LAOS) deformation and extract the rheological response from the virial stress
tensor (see Methods).

The 3-dimensional gel structures obtained from molecular simulation are
characterized in terms of three topological invariants: the number of
connected components $\beta_{0}$, of holes $\beta_{1}$, and of voids captured
by $\beta_{2}$ (these so-called Betti numbers are combined to obtain the EC
number) [32]. A void is an empty cavity within a shape that is surrounded on
all sides by a solid boundary. To provide an intuitive understanding of the
difference between voids and holes, we can think of an open cell versus closed
cell gel or foam. An open cell structure contains many holes, whereas a closed
cell structure contains many voids. To quantify the multi-scale structure of
our soft gel simulations we perform filtrations using increasing particle
diameter and quantify the topology of the resulting Čech complex (simplicial
complex) through the EC. The filtration diameter is scaled such that $D=1$
represents the true diameter of the simulation particles. For example, a
filtration radius $D=2$ represents a Euclidean ball that is twice as wide as
the original gel particle. We illustrate the filtration process on a gel
simulation snapshot in Figure 1, and we also visualize the topological changes
for the full simulation, which contains $\sim$16,000 particles, as particle
diameter is increased. We do this to illustrate how the various topological
invariants: components, cycles, voids, could evolve during a filtration.

Figure 1 shows how the multi-scale structure of soft gels can be directly
quantified through topology and filtrations. The topology of the gel goes
through various topological phase transitions as the filtration diameter is
increased. At $D=1$ to $D=2$ particle diameters we see a stable network
topology. This stable topology reflects the bonded structure of the gel, and
is often what is studied via network analysis tools. This network is dominated
by cycles, and yields a negative Euler characteristic. As we pass to $D>2$
particle diameters we see large topological phase transitions. Here we see a
large increase in the number of cycles, an increasingly negative EC value, and
the formation of an open cell structure. The measured structure now reflects
the 3-dimensional structure (e.g., packing) of the gel and captures overlap
between particles that are not necessarily directly bonded. Finally, at
$D>3.5$ particle diameters we see a transition from the open cell structure
dominated by cycles to the closed cell structure that is dominated by voids
driving the EC to a positive value.

We leverage these topological characterizations to understand the relationship
between soft gel structure and variance in the preparation parameters for the
soft gel. In Figure 2 we show a direct comparison of the EC curves for soft
gels with varying volume fraction (for a constant quench rate at which the gel
was formed in the simulations) and for soft gels formed with varying quench
rates (at a constant volume fraction of $0.10$). We see that there is a smooth
deformation of the EC curves with respect to changes in these two parameters,
but the way in which these parameters change the gel structure is different.
For the soft gels with varying quench rate shown in Figure 2(a), we see that
the changes in the soft gel are primarily impacting the local bonding
structure ($D=1-2.3$ particle diameters) and we see minimal changes in the
larger scale structure ($D>2.3$ particle diameters). This suggests that
changing the quench rate may modify the local structure of the gel, but less
so its mesoscopic 3-dimensional organization. Thus, quench rate could be used
to tune/control the local organization of a soft gel without necessarily
impacting its large-scale organization. In Figure 2(b) we see that the same is
not true when volume fraction is changed; specifically, here we see that there
is a large change in gel topology across multiple length scales.
Interestingly, the impact on local bonding topology ($D=1-2$ particle
diameters) for increasing volume fraction follows a similar trend to
decreasing quench rate, but diverges at larger scales ($D>2$ particle
diameters). This suggests that a combination of quench rate and volume
fraction tune/control could be used to create soft gels with target
structures.

(a)

(b)

Figure 2: (a) EC curves computed from gels at various quench rates with same
volume fraction (0.1). (b) EC curves computed from gels at various volume
fractions and same quench rate ($1\times 10^{-6}\epsilon/k_{B}\tau_{0}$, with
$\epsilon$ and $\tau_{0}$ respectively the unit energy and time in the
simulations). There is a continuous change in the topology/geometry of the
gels as quench rate and volume fraction are changed. There are significant
differences in the way the larger scale structure of the gel is changed when
either volume fraction or quench rate is changed. This suggests that these
parameters can be used to create bespoke gel structures.

Another important aspect in the application of TDA methods is computational
scalability. The methods employed here are able to process gel simulations
with 16,000 to 160,000 particles in a few minutes and on a single laptop
computer. Other methods that attempt to capture similar information, such as
minimal basis cycles of a graph, often run in polynomial time on the number of
edges and nodes in a graph [56]. Furthermore, these graph-theoretic methods
are focused at a single length scale; as such, attempting to perform these
computations across various length scales would compound computational costs
and affect overall scalability. Moreover these methods represent structure as
a 2-dimensional object and thus might miss important information of how the
network graph is embedded in 3-dimensional space.

###  2.2 Gel Dynamics - Shear Analysis

The ability to rapidly quantify the topology of these gel structures allows us
to perform high-throughput analysis of large temporal datasets, as those found
in the analysis of gel rheology [57, 58]. Here, we explore the topological and
geometric deformations of soft gels as they undergo oscillatory shear. The
time-averaged amplitude sweep response of the gel studied here is shown in
Figure 3(a). Due to the density of information contained in the time-resolved
rheological responses, we restrict our analyses to three key transition
amplitudes [58]: the transitions (1) between the linear and non-linear elastic
regimes, (2) the non-linear elastic and yielding regimes, and (3) the yielding
and flow regimes. These transitions occur in this specific gel system at
amplitudes of 0.1, 0.35, and 1.0 strain units, respectively (see Figure 3(a)).
At each of these amplitudes, we display time-resolved rheological and
structural data extracted at a frequency of 32 points per period of
oscillation.

Figure 3(b) illustrates the dynamic evolution of EC curves for a soft gel
undergoing various amplitudes of oscillatory shear that induce linear,
nonlinear, and flow regime responses. The same initialized gel is used (volume
fraction: 0.10, quench rate: $1\times 10^{-6}\epsilon/k_{B}\tau_{0}$, using
the simulations reduced units) but is exposed to different shear amplitudes.
This figure also shows the initial EC curve (initial topology) in the lightest
color and the darkest color representing the final EC curve (final topology).
We see that each simulation begins at the exact same EC curve but evolves
differently as the shear amplitude is changed, though the specific changes are
difficult to discern. To more clearly visualize the changes in the EC curves
as the system is sheared, we used PCA to project the EC curves to a low-
dimensional space (see Methods section). Figure 3(c) visualizes the projection
of all EC curves onto the first two principal components (orthogonal dynamic
modes) of the collective dataset. We found that the projection onto the first
two orthogonal dynamic modes accounts for approximately 98 % of the variance
in the EC data.

(a)

(b)

(c)

(d)

Figure 3: (a,b) EC curves computed at various temporal snapshots for a soft
gel that is undergoing linear, nonlinear, and flow response to shear. The
dynamic evolution is illustrated by difference in color, with the initial EC
curves in the lightest color and the darkest color representing the final EC
curve. We see that each simulation begins at the exact same EC curve but
evolve very different as the shear amplitude is changed. (b) PCA applied to
the EC curve data to create a low-dimensional representation of the dynamics
of the gel during three different shear amplitudes. PCA captures $98\%$ of the
variance in the data; we see two orthogonal dynamic modes that emerge, one
associated with the fast oscillations of the gel induced by oscillatory shear
and a slower mode capturing cumulative structure changes. (c) An optimal rigid
rotation matrix is applied to the dynamic modes so that we may observe the
cumulative and oscillatory components independently.

Figure 3(c) shows that PCA effectively captures the maximum amount of
variation in the data with the minimum number of dynamic modes. This analysis
reveals that there are two orthogonal dynamic behaviors captured in this low-
dimensional representation: one that contains much of the fast oscillatory
dynamics associated with the gel oscillation (we call this the oscillatory
component) and one that contains a slow cumulative structural dynamic change
within the soft gels, in particular for the nonlinear and flow regimes (we
call this the cumulative component). To directly analyze these dynamic modes
in the data, we can apply a rigid rotation matrix to the first two principal
component so that we may observe the cumulative and oscillatory components
independently. The resulting rotated visualization is found in Figure 3(d).
Details around the rigid rotation and criterion for optimal rotation can be
found in the Methods section.

We transformed the high-dimensional set of EC curves into a couple of simple
scalar values that describe the dominating dynamic modes emerging during
shearing. Figures 4(a) and 4(c) show how these values evolve during the three
simulations. Because PCA is a linear method, we can identify what portions of
the EC curve are associated with each dynamic mode by observing the weighting
of the EC curve by the oscillatory and cumulative components shown in Figures
4(b) and 4(d). Here, a larger gray area represents a higher weighting of that
specific portion of the EC curve. This information reveals that the
oscillatory component is capturing changes in the mesoscale structure of the
soft gel ($D>2.3$ particle diameters) while the cumulative component is
capturing local changes in the soft gel network ($D<2.3$ particle diameters).

Figure 4(a) captures much of the oscillatory behavior that is induced by the
oscillatory shearing of the gel. The linear response regime shows minimal
deformation and all deformation seems to be elastic as there are no continuous
offsets in the data. The same is true for the nonlinear response, but we see
larger oscillations which are expected as the shear amplitude has increased.
The flow response shows multiple interesting characteristics that
differentiate it from the linear and nonlinear regimes. First, there is a
large shift in topology in the first initial oscillations of the gel, this
suggests that there is immediate deformation in the mesoscale topology of the
gel during flow. Furthermore, we see that as shear continues the gel reaches a
new _topological_ steady state in which the oscillations are no longer clearly
defined, suggesting that there is a potential topological transition occurring
at the mesoscale in the material as it undergoes flow. Physically, these
changes in the EC suggest a near immediate formation of larger cavities and
holes within the gel without necessarily disrupting the local bonding
structure.

Moving to the cumulative component shown in Figure 4(c), we see that the
linear regime shows almost no change over time. The nonlinear and flow regimes
show similar dynamics but with a larger magnitude change for flow. We also
note there is no dramatic offset in the flow response which is seen in the
oscillatory component. This dynamic mode is associated with the local
network/bonding structure of the soft gel, suggesting that this cumulative
structural change is capturing the breaking and forming of bonds within the
soft gel, which is supported by the weighting shown in Figure 4(c). This mode
is of particular interest because the structural changes are similar to those
induced by reducing quench rate in gel preparation (Figure 2(a)), suggesting
that this mode is capturing potential hardening of the gel. We explore the
relationship between stiffening/hardening and the cumulative component in
Figures 4(e) and 4(f). Figure 4(e) shows the shear stress experienced by the
soft gel at a shear amplitude of 0.350. This amplitude induces a nonlinear
response in the gel and also stiffens the gel, this is evidenced by the
increase in peak shear stress over time. Interestingly, we see a similar
dynamic behavior in the cumulative component for the nonlinear response. We
test whether there is a correlation between this dynamic mode and stiffening
of the gel by comparing peak shear during oscillation and the corresponding
value of the cumulative component at these time points, this is shown in
Figure 4(f). We see there is an obvious correlation between the variables
($R^{2}=0.97$), suggesting that the cumulative component could be used as a
possible predictor of changes in peak shear stress. This also provides
evidence that the mechanism for shear stiffening here is based on local bond
reordering with minimal changes to the mesoscopic structure of the gel,
similar to what is observed when quench rates are reduced.

Combining the analysis of these components also provides insight into why this
same stiffening is not observed at the amplitude corresponding to the flow
transition. The cumulative component at this amplitude shows a much greater
change in the local bonding structure of the gel, which would suggest a
greater stiffening of the gel given the relationship found in Figure 4(f) for
the amplitude corresponding to the maximum non-linear elasticity. However,
when comparing the oscillatory component between these amplitudes we see
minimal non-elastic deformation in the mesoscale structure of the soft gel
($D>2$ particle diameters) in the nonlinear regime, whereas in the flow regime
there is a significant reorganization of the mesoscale structure of the soft
gel, which may be weakening the gel overall. Similar evidence was found in
[19], where significant bond breaking and forming was observed at high shear
amplitudes (but the gel was becoming weaker). The authors hypothesized that
there was a large scale reorganization of the gel, which is further supported
by our topological analysis.

We can further explore the physical connections by comparing the oscillatory
and cumulative dynamics identified with PCA and the stress and strain
undergone by the soft gel in Figure 5. Here, we qualitatively compare the
evolution of stress vs. strain over each cycle and compare this to the
dynamics of the oscillatory and cumulative component versus strain over each
cycle. We observe strong similarities between the evolution of oscillatory
component and stress when plotted against strain, illustrating similar changes
in dynamics as shear amplitude increases and the dramatic change in the
oscillatory component during flow that is not experienced during the linear
and nonlinear shear amplitudes. This reinforces the physical connection
between the topology of the soft gel and its rheological properties. We also
see that the cumulative component smoothly increases over each oscillation
during the nonlinear and flow regimes, while showing minimal change during the
linear regime.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 4: (a,c) Visualization of the oscillatory and cumulative component
identified through PCA for gels undergoing oscillatory shear. Each component
evolves over time, which illustrates the different dynamics expressed for
linear, nonlinear, and flow regime responses. (b,d) Visualization of the
oscillatory and cumulative component weightings for the EC curves. The area of
the weightings (gray) represent higher emphasis placed on that particular
portion of the EC. We see that the oscillatory component focuses on the
mesoscale structure of the soft gel ($D>2.3$ particle diameters), whereas the
cumulative component focuses on the local bonding/network structure of the gel
($D<2.3$ particle diameters). This suggests that these different scales of the
material undergo different dynamics during oscillatory shear. (e,f)
Demonstrates a clear correlation identified between the phenomena of shear
stiffening of a gel responding nonlinearly, evidence by the increase in peak
shear stress (e), and the dynamics of the cumulative component (c). We compute
an $R^{2}=0.97$ which suggests that the cumulative component is an excellent
predictor for the increase in peak shear stress and the phenomena of shear
stiffening.  Figure 5: Comparison of the dynamics of stress, the oscillatory
component, and the cumulative component versus strain in the material for the
linear (left column), nonlinear (middle column), and flow (right column)
regimes. We observe similar changes in oscillatory dynamics as shear amplitude
increases for both stress and the oscillatory component and a dramatic change
in the oscillatory component during flow that is not experienced during the
linear and nonlinear shear amplitudes. We also visualize the cumulative
component showing how the total cumulative damage is smoothly increased over
each oscillation of the gel in the nonlinear and flow regimes, while showing
minimal change in the linear regime.

##  3 Discussion

The topological and geometrical analysis of soft gels can be done effectively
and efficiently through the use of topological data analysis (TDA). We have
shown that the Euler characteristic (EC) of the Čech complex of the gel
network, coupled with filtration operations, captures the multi-scale
structure of soft gels that is missed by methods focusing independently on
topology (e.g., bonding networks) or geometry (e.g., radial particle
distributions in the physical space). We use EC curves as a concise summary of
the topology and geometry of soft gels and use this information to understand
the physical relationships between preparation parameters (e.g., quench rate,
volume fraction). We find that there are clear continuous deformations in the
structure of soft gels that correspond with changes in quench rate and volume
fraction, but that the changes in the structure of the gel are not equivalent.
In particular, we find that changing quench rate impacts the local
organization of the gel and changes in volume fraction impact both local and
global structure. This suggests that fine tuning of the quench rate and volume
fraction could produce bespoke soft gels with state-of-the-art physical
properties. It also indicates that applying our approach to gels formed
through a range of kinetic processes can gain novel insight into the
notoriously complex interplay between kinetics and gel structures. We explored
the dynamics of soft gels under oscillatory shear and identified significant
structural changes in soft gels that occur during linear, nonlinear, and flow
responses. We identified a couple of dominant fast/slow dynamic modes within
the sheared soft gels through the use of TDA and PCA. These modes represent
dynamics of the gel at various length scales. We show that these modes can be
directly correlated to physical phenomena such as shear stiffening or
hardening, and provide insight into the multi-scale dynamics of complex
nonlinear flow behaviors. Furthermore, these analysis methods are
computationally efficient and scalable. We are able to perform computations on
large simulations (in excess of 16,000 particles) in a few seconds on a single
laptop computer. This is in stark contrast with other methods, such as minimum
cycle basis algorithms, which scale nonlinearly in the number of
particles/bonds (nodes/edges) found in a structure and only focus on a single
length scale. The computational scalability and portability make the approach
proposed here particularly interesting to interface data-mining and data-
learning tools.

This work primarily explores the topology of Čech complexes and filtrations to
understand the structure of soft gels, but there are many other tools and
methods from topology and geometry that can be incorporated into these
analyses. The integration of topological concepts from knot theory, to start
with, can help to understand the intertwining of soft gel structures and how
the presence or absence of knots/links in the gel structure contribute to its
physical characteristics [59]. The use of other methods from integral
geometry, such as the full set of Minkowski functionals, can also help to
understand better the full geometry of the gel structures as they undergo
filtration. The Minkowski functionals (i.e., intrinsic volumes) represent
orthogonal measures of intrinsic geometry and for 3-dimensional objects such
as soft gels there are three independent measures: volume, surface area, and
the EC [33, 60, 61]. Future work could therefore explore how much information
those other measures can gain. Finally, the EC and more advanced topological
characteristics such as those obtained from persistence homology, can be used
to construct null-hypothesis and test for significant differences in soft gel
structures [62].

Overall, these methods can provide significant new insight into the behavior
and potential design of new soft materials that involve complex and disordered
network structures. Especially, there is often the need to identify where
these specific topological structures exist within a given soft gel or
material. While it is not possible to map directly from the EC curve back to
the soft gel, methods such as the Hodge Laplacian can be used to identify the
topological structures and understand what portions of a material are
contributing to their physical behavior [63]. The generalizability of the
mathematical foundations for these methods makes them applicable to a wide
range of data types, such as experimental images or videos. Recent
advancements in 3-dimensional imaging of colloidal gels provides data that is
directly amenable to the approaches proposed here and can be used to explore
physical experiments, providing a route for direct comparison of molecular
simulation and experiments [50, 51, 52].

##  4 Methods

###  4.1 Topological Data Analysis

Given a set of points $v\in\mathbb{R}^{m}$, in our case the positions of the
center of mass of the particles that constitute a gel, we construct a
_simplicial complex_ built from $k$-simplices of varying dimension $k$. A
$k$-simplex is a convex set spanned by $k$+1 affinely independent points,
denoted as:

$\displaystyle\sigma=\\{v_{1},v_{2},...,v_{n},v_{k+1}\\}$ (4.1)

Simplices of varying dimension are illustrated in Figure 6. A simplicial
complex (denoted as $\mathcal{K}$) is obtained by connecting simplicies of
varying dimension and can be used to describe the topology of complex shapes.
In Figure 7, we use a simple 2-dimensional soft gel structure to illustrate
the application of simplicial complexes to soft gel simulations. Here, the
soft gel particles are represented as points in 2-dimensional space
$x_{i}\in\mathbb{R}^{2}$, with a diameter $D=1$. When a couple of particles
overlap we add an edge (1-simplex) to our simplicial complex; three a triangle
(2-simplex), and so on. The resulting simplicial complex is known as the
_nerve_ of the overlapping particles, and more specifically the Čech complex
(see next section for more details). Simplicial complexes defined in this way
are exact representations of the topology of the particle system at a given
particle diameter (same number of connected component, cycles, voids, etc.)
while being able to be encoded and measured using simple algebraic operations.
In Figure 7(a), we show how at a $D=1$ particle diameter, many of the
particles overlap with their immediate neighbors. Here, we begin to form edges
or $1$-dimensional simplices that begin to capture a graph structure in our
data which results in four small holes and 5 different connected components,
reflecting a common network representation of a soft gel. However, we can see
that these holes and connected components are not randomly distributed in
space, there is a larger cyclic structure that is not captured when we
consider connectivity at a single particle diameter. As we move from $D=1$
particle diameters to $D=1.5$ particle diameters we see instances where $3$
particles begin to overlap, which yields the formation of $2$-dimensional
simplices (e.g., filled triangles). As we continue to increase $D$ we see
increasing levels of overlap and our resulting simplicial complex is
increasingly connected, revealing the larger cyclic structure captured within
the Čech complex. Eventually we will reach a point where every particle is
overlapped and we end up with a fully connected convex hull of our data, shown
at $D=3$ particle diameters. Through the Čech complex we gain an understanding
of the geometry and topology of the soft gel at varying length scales
(diameters).

(a)

(b)

(c)

Figure 6: Examples of $k$-dimensional simplexes. A simplex is a generalization
of a triangle to higher (or lower) dimensions. (a) $0$-simplices are vertices
(points), (b) $1$-simplices are edges, and (c) $2$-simplices are triangles.

(a)

(b)

Figure 7: (a) Illustration of the changing topology of a 2-dimensional system
as the diameter of the particles are increased. The overlap of the particles
defines the formation of edges and higher-dimensional simplicies in the Čech
complex (shown below). The Čech complex is a simplicial complex representation
of the connectivity of the particles shown above. The simplicial complex
representation allows us to compute and quantify the topology of the system.
(b) The Euler characteristic (EC) curve computed from the filtration of
structure shown in (a). We measure the EC $\chi(D)$ of the Čech complex at
increasing particle diameters $D_{i}\in\mathbb{R}$. The ordered pairing of
these values $\\{\chi(D_{i}),D_{i}\\}$ is the EC curve. The EC curve
summarizes the topological and geometrical changes in the structure at various
scales.

####  4.1.1 Čech Complex

To construct a Čech complex from our soft gel molecular simulations, we place
at the center of each point (i.e., particle) $x_{i}\in\mathbb{R}^{n}$ an
$n$-dimensional Euclidean ball
$B(x_{i},D):=\\{y\in\mathbb{R}^{n}:||x_{i}-y||_{2}<D/2\\}$ of diameter
$D\in\mathbb{R}$. The set of these Euclidean balls is known as a _cover_
$\mathcal{U}:=\cup\\{B_{i}\\}_{i\in I}$ of the points $x_{i}$, where
$I\in\mathbb{Z}$ represents the total number of particles in our simulation.
As we increase the diameter $D$ of the Euclidean balls $B(x_{i},D)$, we will
obtain some level of overlap between balls. This overlap is what defines the
connectivity, and subsequently the geometry and topology of our data. Thus,
from the cover of our points $\mathcal{U}$ we obtain a shape that describes
the structure of our soft gel at a length scale defined by our ball radius $D$
(see Figure 7). From the cover of balls, we can construct a _nerve_ which is a
simplicial complex that represents exactly the topology of the overlapping
balls.

The nerve of collection $\mathcal{U}:=\cup\\{B_{i\in I}\\}$ is the simplicial
complex with vertices $I$ and $k$-simplices built from
$\\{i_{0},i_{1},...,i_{k}\\}$ if and only if $B_{i_{0}}\cap
B_{i_{1}}\cap...\cap B_{i_{k}}\neq\\{\varnothing\\}$. The nerve of the cover
$\mathcal{U}$ is known as the Čech complex that is defined by the particle
positions and selected ball diameter. We then quantify the topology of the
Čech complex through the Euler characteristic.

####  4.1.2 The Euler Characteristic

The Euler characteristic (EC), introduced by Leonard Euler in 1758, is a
scalar value that quantifies the _topological invariants_ of a shape. A
topological invariant is a characteristic of a shape that is unaffected by
continuous deformation of the shape (e.g., stretching, bending), but is
affected by discontinuous deformations (e.g., cutting, glueing). In our work,
we are focused on 3-dimensional objects, which means that we are interested in
three particular topological invariants: connected components (1-dimensional),
holes (2-dimensional), and voids (3-dimensional). The total number of unique
$i$-dimensional invariants in a shape are referred to as the Betti numbers
$\beta_{i}$ where $i$ represents the dimension of the invariant being counted.
The EC $\chi$ is defined as the alternating sum of the Betti numbers:

$\displaystyle\chi=\sum_{i=1}^{n}(-1)^{i}\beta_{i}$ (4.2)

Computationally, there are many ways to quantify the Betti numbers and the
overall EC of a given Čech complex nerve [44, 64, 65]. These methods range
from the use of simplicial algebra to discrete Hodge theory [66]. One of the
simplest methods (we use this in our work) is an extension of Euler’s original
formula for polyhedra to simplicies given by [67]:

$\displaystyle\chi=\sum_{i=1}^{n}(-1)^{i}k_{i}$ (4.3)

where $k_{i}$ represents the number of $k$-dimensional simplicies in the Čech
complex.

####  4.1.3 Filtration and Euler Characteristic Curves

The EC is able to effectively quantify the topology of a shape, but it does
not account for all possible particle diameters and is thus sensitive to the
diameter selected. We can account for this issue by applying a _filtration_ to
our data where we measure the EC at multiple, increasing particle diameters.
For example, through the analysis of the simplified 2-dimensional system in
Figure 7, we notice how the topology of the resulting Čech complex is
sensitive to the selection of particle diameter $D$. This change in topology
with respect to diameter provides us an understanding of the multi-scale
nature of the system. For example, we see that at a radius $D=1$ particle
diameters, we capture multiple small cycles in the structure, while at a
larger radius $D=1.5$ particle diameters we see that there is a much larger
cycle captured. For a set of increasing diameters $D_{1}<D_{2}<...<D_{l}$, we
can compute the EC $\chi(D)$ of the resulting Čech complex at each diameter.
From this, we can create an ordered pairing of these values to construct an
Euler characteristic curve. An example curve for our simple 2-dimensional gel
is found in Figure 7. These curves summarize the evolution of the soft gel
topology as we vary the diameter of the particles (moving from local to global
structures). We can see in Figure 7 that we begin our filtration at $D=1$
particle diameter which captures a series of 4 small cycles and a total of 5
connected components giving $\chi(1)=5-1$, this structure reflects what might
be analyzed in a network representation of the system. We then reach a
threshold at $D=1.5$ where we see that the smaller cycles in our graph are
collapsed and we retain a single large cycle, $\chi(1.5)=1-1$. Eventually the
Čech complex becomes fully connected and we end with a single connected
component at $D=3$ particle diameters, $\chi(3)=1-0$. We can see from Figure 7
that the EC curve quantifies and summarizes these various topological changes.
Furthermore, the EC curve can be represented directly as a vector, which can
be integrated into common data analysis methods such as PCA and for use in
quantifying relationships between the structure of the system and its physical
behavior.

###  4.2 Principal Component Analysis

To apply PCA we construct vectors $X_{t}\in\mathbb{R}^{n}$ of EC values
$\chi(D)\in\mathbb{Z}$ of the Čech complex measured at $n$ diameters
$D\in\mathbb{R}$ for each simulation snapshot time point $t\in\mathbb{Z}$.
This vector $X_{t}:=[\chi(D_{1}),\chi(D_{2}),...,\chi(D_{n})]$ represents the
EC values measured at increasing diameter $D_{1}<D_{2}<...<D_{n}$. From this
we can construct a matrix $\mathbb{M}$ by stacking each of our
$t\in\mathbb{Z}$ simulation snapshot EC curves:
$\mathbf{M}:=[X_{1}^{T},X_{2}^{T},...X_{t}^{T}]^{T}\in\mathbb{R}^{t\times n}$.
We construct a matrix $\mathbf{M}_{i}$ for each shear amplitude
$i\in{0.100,0.350,1.00}$. In order to ensure the PCA projections are
consistent across all simulations we perform PCA on all $\mathbf{M}_{i}$
simultaneously so that the derived principal components are the same for each
simulation. Figure 3(c) visualizes the projection of all simulation snapshot
EC curves onto the first two principal components of the collective
$\mathbf{M}_{i}$.

The principal component analysis conducted provides an effective and
interpretable dimensionality reduction for the EC curve data measured from the
soft gel simulations. Here, we report the total explained variance captured in
the 10 leading principal components and justify the optimal rotation of the
first two principal components to separate the oscillatory and cumulative
components. Figure 8(a) shows the cumulative variance percentage contained in
the leading 10 principal components. We see that the first two components
capture over 98$\%$ of the total variance in the dataset supporting our
selection of the first two principal components for our analysis. Figure 8(b)
illustrates the optimal selection of rotation, in radians, for the data
projected on the first two principal components in order to identify the
oscillatory and cumulative components which are used in our analysis. Given
data projected onto the first two principal components
$\mathbf{v}=[v_{1},v_{2}]\in\mathbb{R}^{2\times 320}$, we apply a rotation
matrix $\mathbf{R}(\theta)$ given as:

$\mathbf{R}(\theta):=\begin{bmatrix}\cos(\theta)&-\sin(\theta)\\\
\sin(\theta)&\cos(\theta)\end{bmatrix}$

where $\theta$ represents the rotation angle in radians. We can identify the
optimal rotation $\hat{\theta}$ by finding the $\theta$ value that maximizes
the amplitude of oscillation in the oscillatory component and minimizes the
oscillatory amplitude in the cumulative component. We measure the oscillatory
amplitude through a simple Fourier transform. Figure 8(b) shows how these two
values evolve as $\theta$ is increased from 0, showing an optimal value
$\hat{\theta}=0.4$. Thus, we obtain and optimally rotated principal component
projection as $\tilde{\mathbf{v}}=\mathbf{R}(.4)\mathbf{v}$.

(a)

(b)

Figure 8: (a) Cumulative explained variance percentage captured with
consecutive leading principal components. It is shown that the first two
principal components explain over 98$\%$ of the total variance within the
dataset (b) Illustration of the absolute amplitude of oscillation contained
within the first and second component as they are rotated. We see that there
is an obvious minimum amplitude in the cumulative component and maximum
amplitude in the oscillatory component at a rotation of 0.40 radians which is
used in our analysis.

###  4.3 Numerical Gel Model and Simulation

We use an established numerical model for particulate gels [54, 55, 12, 19,
58] consisting of self-assembling spherical particles of diameter $d$ that
interact via a short-range attraction, $U_{2}$ and a three body term, $U_{3}$
which limits the bond angles and introduces a bending stiffness. This three-
body term is meant to model hindrance of the relative particle motion upon
aggregation, arising from sources such as particle surface roughness or the
irregularly shaped aggregates. The Molecular Dynamics (MD) simulations are
implemented in a system of $N$ particles with position vectors
{$\mathbf{r}_{1},...,\mathbf{r}_{N}$} and interacting with the potential
energy:

$U\left(\mathbf{r}_{1},...,\mathbf{r}_{N}\right)=\epsilon\left[\sum_{i>j}U_{2}\left(\frac{\textbf{
r}_{ij}}{d}\right)+\sum_{i}\sum_{j>k}^{j,k\neq
i}U_{3}\left(\frac{\textbf{r}_{ij}}{d},\frac{\textbf{r}_{ik}}{d}\right)\right]$
(4.4)

where $\textbf{r}_{ij}=\textbf{r}_{j}-\textbf{r}_{i}$, $\epsilon$ is the depth
of the attractive well $U_{2}$ and sets the energy scale, and $d$ is the
particle diameter, representing the unit length scale. In typical (colloidal)
systems, $d$ corresponds approximately to the range $d\simeq 10$ to $100$ nm
and $\epsilon\simeq 10$ to $100$ $k_{B}T$, where $k_{B}$ is the Boltzman
constant and T, typically room temperature. The functional forms of the two-
body ($U_{2}$) and three-body ($U_{3}$) terms in Eq. (4.4) are detailed
elsewhere [54, 12, 19, 58].

###  4.4 Gel Preparation

The preparation of the particulate gel configurations studied here has been
detailed in prior works [53, 12, 58]; here we cover the basics of the process.
The preparation protocol consists of two parts.

In a cubic simulation box of size $L$, we start from particles initially
equilibrated at $k_{B}T/\epsilon=0.5$ and brought into a spontaneously self-
assembled network at $k_{B}T/\epsilon=0.05$ through NVT (using a Nose-Hoover
thermostat). The network is composed of strands (particles that have
coordination number $z=2$) connected by branching points ($z=3$) and the
system is further equilibrated for additional $2\cdot 10^{4}$ MD steps.

In the second part of the gel preparation, a damped dynamics is used to bring
the configuration obtained as just described to a local minimum that more
likely corresponds to a mechanically stable state. This is achieved by
withdrawing the kinetic energy of the system to $\sim 10^{-10}$ of its initial
value with an overdamped dissipative dynamics:

$m\frac{d^{2}\textbf{r}_{i}}{dt^{2}}=-\nabla_{\textbf{r}_{i}}U-\zeta\frac{d\textbf{r}_{i}}{dt},$
(4.5)

where $m$ is the mass of each particle and $\zeta$ represents the drag
coefficient due to the surrounding solvent.

For the static structural analysis detailed in section 2.1, we study the
effect of varying the initial quench rate between $10^{-2}$ and
$10^{-6}\epsilon/k_{B}\tau_{0}$ at constant volume fraction $\phi=0.1$, and we
study the effect of varying the volume fraction from $0.05$ to $0.15$ at a
constant quench rate $\Gamma\approx 10^{-6}\epsilon/k_{B}\tau_{0}$. For the
rheological measurements in section 2.2, we exclusively study a system
prepared with volume fraction $\phi=0.1$ and quench rate $\Gamma\approx
10^{-5}\epsilon/k_{B}\tau_{0}$, which has been extensively studied in prior
work [19, 58].

###  4.5 Rheological Simulations

The rheological response is measured as in [58] by imposing an oscillatory
strain signal $\gamma(t)=\gamma_{0}\sin\omega t$ in the $xy$-plane of the
simulation box through

$m\frac{d^{2}{\textbf{r}}_{i}}{dt^{2}}=-\nabla_{{\textbf{r}}_{i}}U-\zeta\left(\frac{d{\textbf{r}}_{i}}{dt}-\dot{\gamma}(t)y_{i}{\hat{\textbf{x}}}\right)$
(4.6)

while updating the Lees-Edwards boundary conditions at every time step. Here
$\gamma_{0}$ is the strain amplitude and $\hat{\textbf{x}}$ denotes the unit
vector in the $x$-direction. In all the rheological measurements we used a
Stokes-like drag with $m/\zeta=0.5\tau_{0}$, having verified that the results
do not change qualitatively with further decreasing the $m/\zeta$ ratio.

To systematically study how the gels’ nonlinear rheological responses depend
on their microstructures, we perform large amplitude oscillatory shear (LAOS)
using a rate-controlled deformation:

$\dot{\gamma}(t)=\gamma_{0}\omega cos(\omega t).$ (4.7)

We specifically chose to perform an amplitude sweep by varying the amplitude
over the range $0.001<\gamma_{0}<10$ strain units and holding the frequency of
oscillation at $\omega$ = 0.0025$\tau_{0}^{-1}$. This range of strain
amplitudes spans from the linear viscoelastic regime through the point where
deformations are large enough to yield the material. It is worth noting that,
every amplitude in this amplitude sweep corresponds to a separate rheological
test that starts from the same unperturbed gel configuration [58]. This
differs from most experimental amplitude sweep tests, where the different
amplitudes are typically performed sequentially.

Under applied deformation, we compute the instantaneous shear stress
$\sigma_{xy}(t)$. The stresses are computed from the interaction part of the
global stress tensor using the standard virial equation [68] while neglecting
other contributions (kinetic and viscous terms) as in previous studies [69,
19, 58].

The complex viscoelastic modulus $G^{*}(\omega)$ is obtained from the Fourier
transforms of the stress output ${\tilde{\sigma}}(\omega)$ and the strain
input ${\tilde{\gamma}}(\omega)$ signals, as
$G^{*}(\omega)={\tilde{\sigma}}(\omega)/{\tilde{\gamma}}(\omega)$, from which
we compute the storage modulus $G^{\prime}(\omega)$ and the loss modulus
$G^{\prime\prime}(\omega)$ defined respectively as the real and imaginary part
of $G^{*}(\omega)$. The data used for these moduli is taken from cycles once
the response has reached steady alternance [58] (typically achieved after 4-6
cycles of deformation). Prior oscillatory shear and steady shear measurements
[12, 19, 58] indicate that the dynamic moduli are linear at small strain
amplitudes ($\gamma_{0}$) while the response becomes nonlinear at large
$\gamma_{0}$ values.

For three specific stress amplitudes identified from the amplitude sweep
($0.1$, $0.35$, and $1.0$ strain units) additional oscillatory tests were run
for ten cycles of deformation, with rheological data and structural
configurations being extracted from the system 32 times per period for use in
the topological analysis discussed in section 2.2.

## Acknowledgements

This work was supported by the U.S. Department of Energy, Office of Science,
Advanced Scientific Computing Research, under contract number DE-
AC02-06CH11357. Victor M. Zavala acknowledges partial support of the U.S.
National Science Foundation through the University of Wisconsin Materials
Research Science and Engineering Center (DMR-2309000). Gavin J. Donley and
Emanuela Del Gado acknowledge support from the U. S. National Science
Foundation (grant DMREF CBET—2118962).

## Code and Data Availability

All code for the topological data analysis is available via Zenodo at:
https://doi.org/10.5281/zenodo.10909693. Raw data for molecular simulations is
available upon request.

## References

  * [1] Paulomi Burey, BR Bhandari, Tony Howes, and MJ Gidley.  Hydrocolloid gel particles: formation, characterization, and application.  Critical reviews in food science and nutrition, 48(5):361–377, 2008\.
  * [2] Tharwat F Tadros.  Interfacial Phenomena and Colloid Stability: Industrial Applications.  Walter de Gruyter GmbH & Co KG, 2015.
  * [3] d J Livage and J Lemerle.  Transition metal oxide gels and colloids.  Annual Review of Materials Science, 12(1):103–122, 1982.
  * [4] Yiping Cao and Raffaele Mezzenga.  Design principles of food gels.  Nature Food, 1(2):106–118, 2020.
  * [5] Eric Dickinson.  Microgels—an alternative colloidal ingredient for stabilization of food emulsions.  Trends in Food Science & Technology, 43(2):178–188, 2015.
  * [6] Jakob M Townsend, Emily C Beck, Stevin H Gehrke, Cory J Berkland, and Michael S Detamore.  Flow behavior prior to crosslinking: The need for precursor rheology for placement of hydrogels in medical applications and for 3d bioprinting.  Progress in polymer science, 91:126–140, 2019.
  * [7] Huanan Wang, Otto C Boerman, Kemal Sariibrahimoglu, Yubao Li, John A Jansen, and Sander CG Leeuwenburgh.  Comparison of micro-vs. nanostructured colloidal gelatin gels for sustained delivery of osteogenic proteins: Bone morphogenetic protein-2 and alkaline phosphatase.  Biomaterials, 33(33):8695–8703, 2012.
  * [8] Haimei Zhang.  Building materials in civil engineering.  Elsevier, 2011.
  * [9] Libo Yan, Bohumil Kasal, and Liang Huang.  A review of recent research on the use of cellulosic fibres, their fibre fabric reinforced cementitious, geo-polymer and polymer composites in civil engineering.  Composites Part B: Engineering, 92:94–132, 2016.
  * [10] Katerina Ioannidou, Konrad J. Krakowiak, Mathieu Bauchy, Christian G. Hoover, Enrico Masoero, Sidney Yip, Franz-Josef Ulm, Pierre Levitz, Roland J.-M. Pellenq, and Emanuela Del Gado.  Mesoscale texture of cement hydrates.  Proceedings of the National Academy of Sciences, 113(8):2029–2034, 2016.
  * [11] Rongjin Liu, Zhihui Sun, Qingjun Ding, Ping Chen, and Kangli Chen.  Mitigation of early-age cracking of concrete based on a new gel-type superabsorbent polymer.  Journal of Materials in Civil Engineering, 29(10):04017151, 2017\.
  * [12] Mehdi Bouzid and Emanuela Del Gado.  Network topology in soft gels: Hardening and softening materials.  Langmuir, 34(3):773–781, 2018.
  * [13] Mehdi Bouzid and Emanuela Del Gado.  Mechanics of soft gels: Linear and nonlinear response.  Handbook of Materials Modeling: Applications: Current and Emerging Materials, pages 1719–1746, 2020.
  * [14] Oliver Lieleg, Jürgen Kayser, Gianfranco Brambilla, Luca Cipelletti, and Alexandra R Bausch.  Slow dynamics and internal stress relaxation in bundled cytoskeletal networks.  Nature materials, 10(3):236–242, 2011.
  * [15] JC Conrad, HM Wyss, Véronique Trappe, Suliana Manley, K Miyazaki, LJ Kaufman, Andrew B Schofield, DR Reichman, and DA Weitz.  Arrested fluid-fluid phase separation in depletion systems: Implications of the characteristic length on gel formation and rheology.  Journal of rheology, 54(2):421–438, 2010.
  * [16] Nick Koumakis, Esmaeel Moghimi, Rut Besseling, Wilson CK Poon, John F Brady, and George Petekidis.  Tuning colloidal gels by shear.  Soft Matter, 11(23):4640–4648, 2015.
  * [17] Federica Burla, Simone Dussi, Cristina Martinez-Torres, Justin Tauber, Jasper van der Gucht, and Gijsje H. Koenderink.  Connectivity and plasticity determine collagen network fracture.  Proceedings of the National Academy of Sciences, 117(15):8326–8334, 2020.
  * [18] Mohammad Nabizadeh, Farzaneh Nasirian, Xinzhi Li, Yug Saraswat, Rony Waheibi, Lilian C. Hsiao, Dapeng Bi, Babak Ravandi, and Safa Jamali.  Network physics of attractive colloidal gels: Resilience, rigidity, and phase diagram.  Proceedings of the National Academy of Sciences, 121(3):e2316394121, 2024.
  * [19] Jader Colombo and Emanuela Del Gado.  Stress localization, stiffening, and yielding in a model colloidal gel.  Journal of rheology, 58(5):1089–1116, 2014.
  * [20] Suliana Manley, Benny Davidovitch, Neil R Davies, L Cipelletti, AE Bailey, Rebecca J Christianson, U Gasser, V Prasad, PN Segre, MP Doherty, et al.  Time-dependent strength of colloidal gels.  Physical review letters, 95(4):048302, 2005.
  * [21] Iana Sudreau, Marion Servel, Eric Freyssingeas, Fran çois Liénard, Szilvia Karpati, Stéphane Parola, Xavier Jaurand, Pierre-Yves Dugas, Lauren Matthews, Thomas Gibaud, Thibaut Divoux, and Sébastien Manneville.  Shear-induced stiffening in boehmite gels: A rheo-x-ray-scattering study.  Phys. Rev. Mater., 7:115603, Nov 2023.
  * [22] Esmaeel Moghimi, Alan R Jacob, Nick Koumakis, and George Petekidis.  Colloidal gels tuned by oscillatory shear.  Soft Matter, 13(12):2371–2383, 2017.
  * [23] Thomas Gibaud, Noémie Dagès, Pierre Lidon, Guillaume Jung, L Christian Ahouré, Michael Sztucki, Arnaud Poulesquen, Nicolas Hengl, Frédéric Pignon, and Sébastien Manneville.  Rheoacoustic gels: tuning mechanical and flow properties of colloidal gels with ultrasonic vibrations.  Physical Review X, 10(1):011028, 2020.
  * [24] Daniel Bonn, Morton M Denn, Ludovic Berthier, Thibaut Divoux, and Sébastien Manneville.  Yield stress materials in soft condensed matter.  Reviews of Modern Physics, 89(3):035005, 2017.
  * [25] Jacob Israelachvili.  Intermolecular and Surface Forces.  Academic Press, third edition, 2011.
  * [26] A. D. Dinsmore, V. Prasad, I. Y. Wong, and D. A. Weitz.  Microscopic structure and elasticity of weakly aggregated colloidal gels.  Phys. Rev. Lett., 96:185502, May 2006.
  * [27] Kathryn A. Whitaker, Zsigmond Varga, Lilian C. Hsiao, Michael J. Solomon, James W. Swan, and Eric M. Furst.  Colloidal gel elasticity arises from the packing of locally glassy clusters.  Nature Communications, 10(1):2237, May 2019.
  * [28] C Patrick Royall, Malcolm A Faers, Sian L Fussell, and James E Hallett.  Real space analysis of colloidal gels: triumphs, challenges and future directions.  Journal of Physics: Condensed Matter, 33(45):453002, aug 2021.
  * [29] Lilian C Hsiao, Richmond S Newman, Sharon C Glotzer, and Michael J Solomon.  Role of isostaticity and load-bearing microstructure in the elasticity of yielded colloidal gels.  Proceedings of the National Academy of Sciences, 109(40):16029–16034, 2012.
  * [30] Takehiro Ohtsuka, C Patrick Royall, and Hajime Tanaka.  Local structure and dynamics in colloidal fluids and gels.  Europhysics Letters, 84(4):46002, 2008.
  * [31] Priya Varadan and Michael J Solomon.  Direct visualization of long-range heterogeneous structure in dense colloidal gels.  Langmuir, 19(3):509–512, 2003.
  * [32] Alexander Smith and Victor M Zavala.  The euler characteristic: A general topological descriptor for complex data.  Computers & Chemical Engineering, 154:107463, 2021.
  * [33] Hubert Mantz, Karin Jacobs, and Klaus Mecke.  Utilizing minkowski functionals for image analysis: a marching square algorithm.  Journal of Statistical Mechanics: Theory and Experiment, 2008(12):P12015, 2008.
  * [34] Yang Jiao, FH Stillinger, and S Torquato.  Modeling heterogeneous materials via two-point correlation functions: Basic principles.  Physical review E, 76(3):031110, 2007.
  * [35] S Arizzi, PH Mott, and UW Suter.  Space available to small diffusants in polymeric glasses: analysis of unoccupied space and its connectivity.  Journal of Polymer Science Part B: Polymer Physics, 30(5):415–426, 1992.
  * [36] Peter HF Hansen, Sandra Rödner, and Lennart Bergström.  Structural characterization of dense colloidal films using a modified pair distribution function and delaunay triangulation.  Langmuir, 17(16):4867–4875, 2001.
  * [37] Valeria Lotito and Tomaso Zambelli.  Pattern detection in colloidal assembly: A mosaic of analysis techniques.  Advances in colloid and interface science, 284:102252, 2020.
  * [38] John J Karnes, Todd H Weisgraber, James S Oakdale, Magi Mettry, Maxim Shusteff, and Juergen Biener.  On the network topology of cross-linked acrylate photopolymers: a molecular dynamics case study.  The Journal of Physical Chemistry B, 124(41):9204–9215, 2020.
  * [39] Aldo Glielmo, Brooke E Husic, Alex Rodriguez, Cecilia Clementi, Frank Noé, and Alessandro Laio.  Unsupervised learning methods for molecular simulation data.  Chemical Reviews, 121(16):9722–9758, 2021.
  * [40] A Kaveh and GR Roosta.  Revised greedy algorithm for formation of a minimal cycle basis of a graph.  Communications in numerical methods in engineering, 10(7):523–530, 1994.
  * [41] Alexander Smith, Spencer Runde, Alex Chew, Atharva Kelkar, Utkarsh Maheshwari, Reid Van Lehn, and Victor Zavala.  Topological analysis of molecular dynamics simulations using the euler characteristic.  2022\.
  * [42] Jason W Rocks, Andrea J Liu, and Eleni Katifori.  Hidden topological structure of flow network functionality.  Physical Review Letters, 126(2):028102, 2021.
  * [43] Jason W Rocks, Andrea J Liu, and Eleni Katifori.  Revealing structure-function relationships in functional flow networks via persistent homology.  Physical Review Research, 2(3):033234, 2020.
  * [44] Alexander D Smith, Paweł Dłotko, and Victor M Zavala.  Topological data analysis: concepts, computation, and applications in chemical engineering.  Computers & Chemical Engineering, 146:107202, 2021.
  * [45] Moo K Chung, Alexander Smith, and Gary Shiu.  Reviews: Topological distances and losses for brain networks.  arXiv preprint arXiv:2102.08623, 2021.
  * [46] Ippei Obayashi, Takenobu Nakamura, and Yasuaki Hiraoka.  Persistent homology analysis for materials research and persistent homology software: Homcloud.  Journal of the Physical Society of Japan, 91(9):091013, 2022.
  * [47] Takenobu Nakamura, Yasuaki Hiraoka, Akihiko Hirata, Emerson G Escolar, and Yasumasa Nishiura.  Persistent homology and many-body atomic structure for medium-range order in the glass.  Nanotechnology, 26(30):304001, 2015.
  * [48] KR Mecke and Herbert Wagner.  Euler characteristic and related measures for random geometric sets.  Journal of Statistical Physics, 64:843–850, 1991.
  * [49] Victor Klee.  The euler characteristic in combinatorial geometry.  The American Mathematical Monthly, 70(2):119–127, 1963.
  * [50] AD Dinsmore and DA Weitz.  Direct imaging of three-dimensional structure and topology of colloidal gels.  Journal of Physics: Condensed Matter, 14(33):7581, 2002.
  * [51] Jun Dong, Francesco Turci, Robert L Jack, Malcolm A Faers, and C Patrick Royall.  Direct imaging of contacts and forces in colloidal gels.  The Journal of Chemical Physics, 156(21):214907, 2022.
  * [52] Piotr Habdas and Eric R Weeks.  Video microscopy of colloidal suspensions and colloidal crystals.  Current opinion in colloid & interface science, 7(3-4):196–203, 2002.
  * [53] Mehdi Bouzid, Jader Colombo, Lucas Vieira Barbosa, and Emanuela Del Gado.  Elastically driven intermittent microscopic dynamics in soft solids.  Nature Communications, 8:15846, Jun 2017.
  * [54] M. Bantawa, W. A. Fontaine-Seiler, P. D. Olmsted, and E. Del Gado.  Microscopic interactions and emerging elasticity in model soft particulate gels.  Journal of Physics Condensed Matter, 33, 2021.
  * [55] M. Bantawa, B. Keshavarz, M. Geri, M. Bouzid, T. Divoux, G.H. McKinley, and E. Del Gado.  The hidden hierarchical nature of soft particulate gels.  Nature Phys., 19:1178–1184, 2023.
  * [56] Joseph Douglas Horton.  A polynomial-time algorithm to find the shortest cycle basis of a graph.  SIAM Journal on Computing, 16(2):358–366, 1987.
  * [57] José Ruiz-Franco, Nicoletta Gnan, and Emanuela Zaccarelli.  Rheological investigation of gels formed by competing interactions: A numerical study.  The Journal of chemical physics, 150(2), 2019.
  * [58] Gavin J Donley, Minaspi Bantawa, and Emanuela Del Gado.  Time-resolved microstructural changes in large amplitude oscillatory shear of model single and double component soft gels.  Journal of Rheology, 66(6):1287–1304, 2022.
  * [59] Kate E Horner, Mark A Miller, Jonathan W Steed, and Paul M Sutcliffe.  Knot theory in modern chemistry.  Chemical Society Reviews, 45(23):6432–6448, 2016.
  * [60] Markus Hütter.  Heterogeneity of colloidal particle networks analyzed by means of minkowski functionals.  Physical Review E, 68(3):031404, 2003.
  * [61] Patrice Koehl, Arseniy Akopyan, and Herbert Edelsbrunner.  Computing the volume, surface area, mean, and gaussian curvatures of molecules and their derivatives.  Journal of Chemical Information and Modeling, 63(3):973–985, 2023\.
  * [62] Paweł Dłotko, Niklas Hellmer, Łukasz Stettner, and Rafał Topolnicki.  Topology-driven goodness-of-fit tests in arbitrary dimensions.  Statistics and Computing, 34(1):34, 2024.
  * [63] Ronald Koh Joon Wei, Junjie Wee, Valerie Evangelin Laurent, and Kelin Xia.  Hodge theory-based biomolecular data analysis.  Scientific Reports, 12(1):9699, 2022.
  * [64] Herbert Edelsbrunner, John Harer, et al.  Persistent homology-a survey.  Contemporary mathematics, 453(26):257–282, 2008.
  * [65] Herbert Edelsbrunner and John L Harer.  Computational topology: an introduction.  American Mathematical Society, 2022.
  * [66] Oliver Knill.  One can hear the euler characteristic of a simplicial complex.  arXiv preprint arXiv:1711.09527, 2017.
  * [67] David S Richeson.  Euler’s gem: the polyhedron formula and the birth of topology.  Princeton University Press, 2019.
  * [68] A. Thompson, S. Plimpton, and W. Mattson.  General formulation of pressure and stress tensor for arbitrary many-body interaction potentials under periodic boundary conditions.  Journal of Chemical Physics, 131, 2009.
  * [69] Mehdi Bouzid and Emanuela Del Gado.  Mechanics of Soft Gels: Linear and Nonlinear Response, pages 1–29.  Springer International Publishing, Cham, 2018.