Network Analysis of Particles and Grains
Lia Papadopoulos Mason A. Porter Karen E. Daniels Danielle S. Bassett






Granular Networks



Network Analysis of Particles and Grains
    
Lia Papadopoulos

Department of Physics & Astronomy, University of Pennsylvania, USA

epapad@sas.upenn.edu

Mason A. Porter

Department of Mathematics, University of California Los Angeles, USA
 
Mathematical Institute, University of Oxford, UK

CABDyN Complexity Centre, University of Oxford, UK

mason@math.ucla.edu

Karen E. Daniels

Department of Physics, North Carolina State University, USA

kdaniel@ncsu.edu

Danielle S. Bassett^*

Departments of Bioengineering and Electrical & Systems Engineering, University of Pennsylvania, USA

^*Corresponding author: dsb@seas.upenn.edu
    December 30, 2023
===========================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================




The arrangements of particles and forces in granular materials have a complex organization on multiple spatial scales that ranges from local structures to mesoscale and system-wide ones. This multiscale organization can affect how a material responds or reconfigures when exposed to external perturbations or loading. The theoretical study of particle-level, force-chain, domain, and bulk properties requires the development and application of appropriate physical, mathematical, statistical, and computational frameworks. Traditionally, granular materials have been investigated using particulate or continuum models, each of which tends to be implicitly agnostic to multiscale organization. Recently, tools from network science have emerged as powerful approaches for probing and characterizing heterogeneous architectures across different scales in complex systems, and a diverse set of methods have yielded fascinating insights into granular materials. In this paper, we review work on network-based approaches to studying granular matter and explore the potential of such frameworks to provide a useful description of these systems and to enhance understanding of their underlying physics. We also outline a few open questions and highlight particularly promising future directions in the analysis and design of granular matter and other kinds of material networks. 
Granular materials, particulate systems, networks, network science, multiscale organization










§ GLOSSARY OF TERMS







  Granular materials: Granular materials are collections of discrete, macroscopic particles that interact with each other through contact (rather than long-range) forces. Importantly, these systems are non-equilibrium: the particles are large enough to avoid rearrangement under thermal fluctuations, and they lose energy through frictional and inelastic interactions
with neighboring particles. 


  Particulate materials: Like granular materials, particulate materials are collections of discrete, macroscopic elements. However, the elements making up the system may be entities — such as bubbles, foams, colloids, or suspensions — that include multiple phases of matter.
The term “particulate material” is a more general one than “granular material”. 


  Packing fraction: The fraction of a granular material that consists of the particles. One calculates the packing fraction as the ratio of the total volume of all particles to the volume of the region that they occupy. The packing fraction is also sometimes called the “packing density” or “volume fraction”. 


  Force chain: Force chains are typically described as the subset of inter-particle contacts in a granular material that carry the largest forces in the system. They often form filamentary networks that align preferentially with the principal stress axes under which a material is loaded.


  Jamming: As certain system parameters change, disordered, particulate materials typically undergo a transition from an underconstrained, liquid-like state to a rigid, solid-like state characterized by the onset of mechanical stability. The transition to/from a jammed state may arise through increases/decreases in quantities like packing fraction or contact number, which can occur due to an applied load. A formal theory of jamming exists for idealized situations (with soft, frictionless, and spherical particles).


  Isostatic: A jammed packing is isostatic when it has exactly the minimum number of contacts that are required for maintaining mechanical stability through force balance and torque balance. One typically examines isostaticity by calculating the mean number of contacts per particle. “Hyperstatic” and “hypostatic” packings have more and fewer contacts, respectively. 


  Mono/bi/poly-disperse: A particulate material is monodisperse if it is composed of particles of a single species (with the same size, shape, and material properties). Bidisperse materials have particles of two species, and polydisperse materials can either have particles from three or more discrete species or have particles with a continuum of properties.


  Structural rigidity theory: For a specified structure composed of fixed-length rods connected to one another by hinges, structural rigidity theory studies the conditions under which the associated structural graph is able to resist deformations and support applied loads. 


  Stress: A stress is a force applied to an object's surfaces. (The units measure a force per unit area.) Shear stress arises from the component of the applied force that acts in a direction parallel to the object's cross-section, and normal stress arises from the perpendicular component.


  Strain: Strain is the fractional (unitless) deformation of a material that arises due to an applied stress. One calculates strain from the relative displacement of particles in a material, excluding rigid-body motion such as translation or rotation. Like stress, strain has both shear and normal components. 


  Pure shear: One can use the term “pure shear” to describe either stresses or strains in which an object is elongated along one axis and shortened in the perpendicular direction without inducing a net rotation.



  Axial compression: In axial compression, one applies inward forces to an object in one direction (uniaxial), two directions (biaxial), or all directions (isotropic compression). These forces result in uniaxial strain, biaxial strain, or isotropic strain, respectively. 


  Cyclic shear/compression: These consist of repeated cycles of shear or compression applied to the same system.


  Shear band: A shear band is a narrow region of a particulate material in which most of the strain is localized, whereas other regions remain largely undeformed. A shear band is also sometimes called a region of “strain localization”. 


  Strain softening/hardening: As a material is loaded and undergoes deformation, continuing deformation can become either easier (strain softening) or harder (strain hardening). Eventually, after much deformation, the material can reach a critical state in which there are no further changes in the resistance to deformations. 


  Stress ratio: The stress ratio, which is analogous to Coulomb's Law, is the ratio of shear to normal stresses. Frictional failure occurs when the shear force exceeds the product of the normal force and the coefficient of friction.


  Photoelasticity/birefringence: Photoelasticity is an optical technique for quantifying internal stresses based on the transmission of polarized light through “birefringent” materials, which have preferentially fast and slow directions for the propagation of light.


  DEM or MD simulations: The Discrete (or Distinct) Element Method and Molecular Dynamics simulations are related numerical techniques that compute the motions of all particles in a system (such as a granular material). In each method, a computer algorithm treats each particle as an object subject to Newton's laws of motion, where forces consist of body forces (e.g., gravity) and those that arise from interactions with the object's neighbors.






§ INTRODUCTION



Granular materials comprise a subset of the larger set of particulate matter <cit.>. People engage with such materials — which include sands, beans, grains, powders such as cornstarch, and more — often in their daily lives. One can define a granular material as a large collection of discrete, macroscopic particles that interact only when in contact. Granular materials are inherently non-equilibrium in two distinct ways, characterized by (1) the lack of rearrangement under thermal fluctuations and (2) the loss of energy through frictional and inelastic dissipation during contact between grains. Nonetheless, they phenomenologically reproduce equilibrium states of matter, exhibiting characteristics of solids (forming rigid materials), liquids (flowing out of a container), or gases (infrequent contacts between grains), depending on the type and amount of driving. In this review, we focus mainly on granular solids and slow (non-inertial) flows <cit.>; these are dense materials in which sustained inter-particle contacts provide the dominant contribution to material properties.

The functional properties of granular materials are related in a nontrivial way to the complex manner in which particles interact with one another and to the spatial scales (particle, chain, domain, and bulk) and time scales over which those interactions occur. For example, pairs of particles can exert force on one another in a local neighborhood. However, as particles push on adjacent particles, the combined effect can transmit forces over long distances via important mesoscale structures commonly called force chains <cit.>.
The idea of networks has been invoked for many years to help provide a quantitative understanding and explanation of force-chain organization <cit.>.
Broadly speaking, force chains form a network of filamentary-like structures that are visually apparent in images from experiments, like the one shown in Fig. <ref>. In such images, the brighter particles carry larger forces <cit.>. Furthermore, force chains tend to align preferentially along the principal stress axes <cit.>. It can be helpful to think of a force-chain network as the backbone of strong forces that span a system, providing support for both static <cit.> and dynamic <cit.> loading. However, weaker forces can also play a stabilizing role, much as guy-wires do on an aerial tower <cit.>.

It is also possible for sets of particles to cluster together into larger geographical domains, with potentially distinct properties, that can have weak structural boundaries between them <cit.>. At the largest scale, granular materials as a whole exhibit bulk properties, such as mechanical stability or instability in response to sheer or compression <cit.>. All of the aforementioned spatial scales are potentially relevant for understanding phenomena such as transmission of acoustic waves <cit.>, thermal conductivity and heat transfer <cit.>, electrical properties <cit.>, and more. The time scales of interactions in granular materials are also important, and they can vary over many orders of magnitude. For example, in systems under compression, statistical fluctuations of grain displacements depend fundamentally on the length of the strain step (i.e., “increment”) over which one makes measurements, as fluctuations over short windows are consistent with anomalous diffusion and those over longer windows are consistent with Brownian behavior <cit.>.



The principled study of such diverse characteristics and organization in a single system can be very challenging, and the development of appropriate physical, mathematical, statistical, and computational models is important to attain a mechanistic understanding of granular materials. Traditionally, it has been common to model granular materials using either particulate-based or continuum-based frameworks <cit.>. However, both of these approaches are often implicitly agnostic to intermediate-scale organization, which is important for understanding both static granular packings <cit.> as well as granular dynamics  <cit.>. Recently, tools from network science <cit.> and related mathematical subjects — which include approaches that can account explicitly for mesoscale structures <cit.> — have been used successfully to study properties of granular materials across multiple spatial and temporal scales. The most common representation of a network, an important idea for the study of complex systems of interacting entities <cit.>, is as a graph <cit.>. A graph consists of a set of nodes (to represent the entities) and a set of edges, each of which represents an interaction between a pair of entities (or between an entity and itself). Increasingly, more complicated network representations (such as multilayer networks <cit.>) are also employed. Moreover, there is also an increasing recognition that it is important to consider the impact of other features, such as spatial embedding and other spatial effects <cit.>, on network structure and dynamics, rather than taking an approach that promises that “one size fits all.” Network science offers methods for quantitatively probing and analyzing large, interacting systems whose associated networks have heterogeneous patterns that defy explanations attained by considering exclusively all-to-all, regular, or lattice-like interactions <cit.>.

There are several open problems in granular physics that may benefit from network-science approaches. In particular, because granular materials have multiple relevant length and time scales <cit.>, it can be challenging to model and quantify their structural organization, material properties, and responses to external loads <cit.>. However, although complex, the pairwise inter-particle interactions that underlie and govern the structure and behavior of granular systems (and other particulate matter) render them amenable to various network representations and network-based analyses. <cit.> were among the first to explicitly suggest and formalize the use of ideas from network science to begin to study some of the difficult questions in granular physics. In their paper, they highlighted the ability of a network-based perspective to complement traditional methods for studying granular materials and to open new doors for the analysis of these complex systems. One place in which network analysis may be especially useful is in quantifying how local, pairwise interactions between particles in a granular packing yield organization on larger spatial scales (both mesoscale and system-level). For example, in sheared or compressed granular packings, such organization can manifest as force chains or other intermediate-sized sets of particles that together comprise a collective structure. Network science provides approaches to extract and quantitatively characterize heterogeneous architectures at microscale, mesoscale, and macroscale sizes, and one can use these methods to understand important physical phenomena, including how such multiscale organization relates to bulk material properties or to spatial and temporal patterns of force transmission through a material.

Network-based approaches should also be able to provide new insights into the mechanisms that govern the dynamics of granular materials. For example, as we will discuss, network analysis can helpfully describe certain aspects of complex dynamics (such as granular flows), and can provide quantitative descriptions of how the structure of a dense granular material evolves as a system deforms under external loads (such as those induced by compression, shear, tapping, or impact). It seems sensible to use a network-based approach when a network is changing on temporal scales slower than the time that it takes for information to propagate along it. We also expect ideas from temporal networks <cit.> or adaptive networks <cit.> to be fruitful for studying faster dynamics, and investigation of granular dynamics in general should benefit from the development of both novel network representations and methods of network analysis that are designed specifically to understand temporally evolving systems. Another important problem in the study of granular materials is to predict when and where a granular system will fail. There has been some progress made in this area using network-based approaches, but it is important to continue to develop and apply tools from network analysis and related areas to gain a deeper understanding of which network features regulate or are most indicative of eventual failure. Another exciting direction for future work is to combine network-based approaches with questions about material design. In particular, can one use network-based approaches to help engineer granular systems — or other materials that are amenable to network representations — with desired and specialized properties?

It is also important to note that network-based representations and methods of analysis can provide insightful descriptions of granular materials with various additional complexities, such as systems that are composed of differently-shaped particles, 3-dimensional (3D) materials, and so on. This flexibility makes the application of tools from network science a powerful approach for studying the structural properties and dynamics of granular networks. Such a framework also allows one to compare network architectures in diverse situations, such as between simulations and experiments, across systems that are composed of different types of particles or exposed to different loading conditions, and more. Exploiting these capabilities will yield improved understanding of which properties and behaviors of granular materials are general, versus which are specific to various details of a system. With the continued development of physically-informed network-analysis tools, network-based approaches show considerable promise for further development of both qualitative and quantitative descriptions of the organization and complex behavior of granular materials.

The purpose of our paper is to review the nascent application of network theory (and related topics) to the study of granular materials. We begin in Sec. <ref> with a mathematical description of networks. In Sec. <ref>, we briefly review a set of measures that one can calculate on graphs and which have been useful in past investigations of granular materials. In Sec. <ref>, we review several different ways in which granular materials have been represented as networks, and we discuss investigations of such networks to quantify heterogeneous, multiscale organization in granular materials and to understand how these systems evolve when exposed to external perturbations. We also point out insights into the underlying physics that have resulted from network-based investigations of granular matter. We close in Sec. <ref> with some thoughts on the many remaining open questions, and we describe a few specific future directions that we feel are important to pursue. We hope that our review will be helpful for those interested in using tools from network science to better understand the physics of granular systems, and that it will spur interest in using these techniques to inform material design.  



§ NETWORK CONSTRUCTION AND CHARACTERIZATION




 §.§ What is a network?


 


It is often useful to model a complex system as a network, the simplest type of which is a graph <cit.>. For most of our article, we will use the terms network and graph synonymously, but the former concept is more general than the latter.[Indeed, it is increasingly important to examine network representations that are more complicated than graphs (see Secs. <ref> and <ref>) — such as multilayer networks <cit.>, simplicial complexes <cit.>, and others — and it is also essential to study dynamical processes on networks, rather than focusing exclusively on structural characteristics <cit.>.]  A graph G consists of nodes (i.e., vertices), where pairs of nodes are adjacent to each other via edges (i.e., links). 
We denote the set of nodes by 𝒱 and the set of edges by ℰ. 
A node can also be adjacent to itself via a self-edge (which is also sometimes called a self-loop), and a multi-edge describes the presence of two or more edges that are attached to (i.e., incident to) the same pair of nodes. (Unless we state otherwise, we henceforth assume that our networks have neither self-edges nor multi-edges.) The number of nodes in a graph is the size of the graph, and we also use the word “size” in the same way for other sets of nodes. A subgraph of a graph G is a graph constructed using a subset of G's nodes and edges.

An edge between two nodes represents some sort of relationship between them. For example, edges can represent information flow between different parts of the internet <cit.>, friendship or other social interactions between people <cit.>, trading between banks <cit.>, anatomical or functional connections between large-scale brain regions <cit.>, physical connections between particles in contact <cit.>, and so on. Edges can be either unweighted or weighted, and they can be either undirected or directed <cit.>. In an unweighted (i.e., binary) network, an edge between two nodes is assigned a binary value (traditionally 0 or 1) to encode the absence or presence of a connection. In a weighted network, edges can take a variety of different values to convey varying strengths of relationships between nodes. In an undirected network, all edges are bidirectional, so one assumes that all relationships are reciprocal. In a directed network, however, edges have a direction that encodes a connection from one node to another.

An adjacency matrix is a useful way to represent the information in a graph. For an unweighted and undirected graph, an element of the adjacency matrix 𝐀 of an N-node network is 

    A_ij = {[ 1 ,  if there is an edge between nodes i and j ,;                                 0 ,  otherwise , ].

where i,j ∈{1, …, N}. For a network in which nodes do not have labels, one can apply a permutation to A's rows and columns — one uses the same permutation for each — to obtain another adjacency matrix that represents the same network. In this paper, when we refer to “the” adjacency matrix A of a graph with unlabeled nodes, we mean any one of these matrices, which have the same properties (spectra, etc.).

For a weighted graph, if nodes i and j are adjacent via an edge, we denote the corresponding edge weight by w_ij (which is usually given by a nonnegative real number (e.g., the value of the normal or tangential component of the force between two contacting particles)[We do not consider edges with negative weights, although it may be interesting to do so in future work if there is an appropriate physical reason.]). An element in the associated weighted adjacency matrix (which is sometimes called a weight matrix) 𝐖 is

    W_ij = {[ w_ij ,  if there is an edge between nodes i and j ,;                                    0 ,  otherwise . ].


For the more general case of a weighted, directed graph, if there is an edge from node j to node i, then we let w_ij represent the weight of that edge <cit.>. The associated weighted and directed adjacency matrix 𝐖 is 

    W_ij = {[ w_ij ,  if there is an edge from node j to node i ,;                                    0 ,  otherwise . ].


An adjacency matrix associated with an undirected network is symmetric, but an adjacency matrix for a directed network need not be (and is symmetric if and only if all directed edges are reciprocated). In the present review, we primarily consider undirected networks, although we will occasionally make remarks about directed situations.

For weighted graphs, it is often also important to consider a binary adjacency matrix 𝐀 associated with a weight matrix 𝐖. Note that 𝐀 captures only the connectivity of nodes (i.e., their adjacencies), irrespective of how strongly they interact with each other. In terms of 𝐖, the corresponding binary network (which can be either directed or undirected) is

    A_ij = {[ 1 ,  if  W_ij≠ 0 ,;   0 ,  otherwise . ].

It is common to use terms like network topology when discussing structural properties of 𝐀, and sometimes one uses terms like network geometry when discussing properties that also depend on edge weights. Because we will also employ ideas from subjects like algebraic topology (see Sec. <ref>), we will need to be very careful with such terminology.

The network representations that have been used to study granular matter (and other kinds of materials) employ diverse definitions of edges (both weighted and unweighted, and both directed and undirected), and some generalizations of graphs have also been considered. See Fig. <ref> for a schematic showing possible choices of nodes and edges for a network representation of a granular packing. A variety of tools and measures from network analysis have been used to study granular networks. We discuss some of these ideas in Sec. <ref>.



 §.§ Some tools for characterizing granular networks



Network theory <cit.> provides myriad ways to characterize and quantify the topological and geometrical organization of complex networks. Thus, different network methods can reveal different important features of the underlying system, and these features, in turn, can help explain how a system behaves in certain situations. In the context of granular matter, for example, it is often desirable to understand the stability of a material, mechanical responses to external stresses, or wave propagation through a system. Recent investigations have demonstrated that network analysis can inform understanding of the mechanisms that underlie these phenomena. In this section, we discuss several network concepts; and in Sec. <ref>, we describe how they have been used for the study of granular materials. We are, of course, not presenting anything close to an exhaustive list of tools from network science. See <cit.> and other books and reviews (and references therein) for discussions of other tools from network science. For simplicity, we primarily give definitions for undirected networks, though many of the ideas that we present also have directed counterparts. We start with basic network diagnostics, and also discuss some more complicated methods.



  §.§.§ Degree.



One local property of a network is node degree. In an undirected network, a node's degree is equal to the number of edges that are attached to it (see Fig. <ref>d). We denote the degree of node i as k_i, and we recall that N denotes the total number of nodes. For an unweighted graph with adjacency matrix A, one can calculate k_i with the formula

    k_i = ∑_j=1^N A_ij .
 

One can generalize the idea of degree to strength (i.e., weighted degree) using a weight matrix 𝐖 <cit.>. The strength s_i of node i is equal to the sum of the weights of the edges that are attached to that node: 

    s_i = ∑_j=1^N W_ij .
 
Its associated degree k_i is still given by Eq. (<ref>).

A common network representation of granular materials is to treat particles as nodes and physical contacts between particles as either unweighted or weighted edges (see Fig. <ref>a–c). In this representation, node degree and node strength quantify information at the scale of single particles. 

One can compute the mean degree (a global property) of a network by calculating

    ⟨k⟩ = 1/N∑_i k_i ,

and one can similarly compute the mean strength with the formula

    ⟨s⟩ = 1/N∑_i s_i .

In an undirected network, mean degree is related to N and the total number m of edges through the relation ⟨k⟩ = 2m/N. It is sometimes useful to characterize a network using its degree distribution P(k) <cit.>, which gives the probability that a node has degree k.

When one represents a granular packing as a contact network (see Fig. <ref> and Sec. <ref>), which is a binary network (i.e., unweighted network), the degree k_i of node i is known more commonly in the physics literature as its contact number or coordination number Z_i. If every node has the same degree (or even if most nodes have this degree), such as in a regular lattice, one often refers to the mean coordination number Z of the lattice or packing. It is well-known that coordination number is related to the stability of granular packings <cit.> and plays a critical role in the jamming transition <cit.>, a change of phase from an underconstrained state to a rigid state that is characterized by the onset of mechanical stability.[Unless we note otherwise, we use the phrase jamming in the formal sense of the jamming transition as defined by <cit.>. Packings of particles above the jamming point (a critical point related to the jamming transition) are rigid and overconstrained (i.e., “hyperstatic”), those at this point are marginally stable and exactly constrained (i.e., “isostatic”), and those below this point are underconstrained (i.e., “hypostatic”). Additionally, packings below the jamming point are sometimes called “unjammed”, and those above the jamming point are called “jammed”.] We discuss these ideas further throughout the review. 



  §.§.§ Walks and paths.
 


In a network, a walk is an alternating sequence of nodes and edges that starts and ends at a node, such that consecutive edges are both incident to a common node. A walk thus describes a traversal from one node to another node (or to itself) along edges of a network. A path is a walk that does not intersect itself or visit the same node (or the same edge) more than once, except for a closed path, which starts and ends at the same node (see Sec. <ref>). One can compute the number of (unweighted) walks of a given length from a binary, undirected adjacency matrix A <cit.>. The length l of an unweighted walk is defined as the number of edges in the associated sequence (counting repeated edges the number of times that they appear). Letting Ξ^l_ij denote the number of walks of length l between nodes i and j, one calculates

    Ξ^l_ij = [𝐀^l]_ij .

Various types of random walks yield short paths between nodes in a network, and such ideas (and their relation to topics such as spectral graph theory) are very insightful for studying networks <cit.>.

In an undirected network, a path from node i to node j is necessarily also a path from node j to node i. However, this is not typically true in directed networks, which have sometimes been utilized in studies of granular force chains (see, e.g., <cit.>). Depending on a network's structure, there may be several or no paths between a given pair of nodes. An undirected network is called connected if there exists a path between each pair of nodes, and a directed network is called strongly connected if there is a path between each pair of nodes <cit.>. (A directed network is called weakly connected if its associated undirected network is connected, and a strongly connected network is necessarily also weakly connected.) In networks of granular packings, both the existence and lengths of paths can impact system behavior. A connected network consists of a single component. When a network has multiple components, it is common to study one component at a time (e.g., focusing on a largest connected component (LCC), which is one that has the largest number of nodes).

The length of an unweighted path is the number of edges in the associated sequence, and it is sometimes also called the hop distance (or occasionally, unfortunately, the topological distance). Paths in a network can also be weighted by defining some (possibly abstract) notion of distance associated with the edges of the network. For example, in a spatially-embedded network <cit.>, distance may refer to actual physical distance along an edge, in which case the length of a weighted path in the network is given by the sum of the physical distances along the sequence of edges in the path. However, one can also consider “distance" more abstractly. For example, in a transportation or flow network, one can define a distance between two adjacent nodes to be some measure of resistance between those nodes, and then the length of a weighted path in such a network is given by the sum of the resistances along the sequence of edges in the path.[One can also calculate distances between nodes if they occupy positions in a metric space (such as a latent one that determines the probability that a pair of nodes is adjacent to each other) <cit.>, and the properties of that metric space can influence the distances (namely, the ones along a network) that concern us.]

We use the term network distance to indicate a distance between two nodes (which can be either unweighted or weighted) that is computed by summing along edges in a path. A geodesic path — i.e., a shortest path (which need not be unique) — between two nodes can be particularly relevant (though other short paths are often also relevant), and a breadth-first search (BFS) algorithm <cit.> is commonly employed to find geodesic paths in a network. The diameter of a graph is the maximum geodesic distance between any pair of nodes. 

Denoting the shortest, unweighted network distance (i.e., shortest-path distance) between nodes i and j as d_ij, the mean shortest-path distance L between pairs of nodes in a graph is <cit.>

    L = 1/N(N-1)∑_i,j    (i ≠ j) d_ij .

Note that one must be cautious when computing the mean shortest-path distance on disconnected networks (i.e., on networks that have more than one component), because the usual convention is to set the distance between two nodes in different components to be infinite <cit.>. Therefore, in a network with multiple components, the distance L from Eq. (<ref>) is infinite. One solution to this problem is to compute L for each component separately. Another network notion that relies on paths is network efficiency <cit.>

    E = 1/N(N-1)∑_i,j   (i ≠ j)1/d_ij .


One can generalize measures based on walks and paths to incorporate edge weights <cit.>. Letting d^w_ij denote the shortest, weighted network distance between nodes i and j, one can define a weighted mean shortest-path distance L^w and weighted efficiency E^w as in Eqs. (<ref>) and  (<ref>), respectively, but now one uses d_ij^w instead of d_ij. The network efficiency E is a normalized version of the Harary index of a graph <cit.>. Additionally, the convention d_ij = ∞ (or d^w_ij = ∞) if there is no path from i to j allows one to use Eq. (<ref>) (or its weighted counterpart E^w) on connected graphs or on graphs with more than one component. For both unweighted and weighted scenarios, large values of network efficiency tend to correspond to small values of mean shortest-path length, and vice versa. One can also readily generalize notions of paths, distances, and efficiency to directed networks <cit.>.

In later sections, we will describe the use of paths, walks, and related ideas for investigating the structure of granular materials and their response to perturbations — including, but not limited to, how these quantities change as a granular packing is compressed and goes through the jamming transition <cit.> — and we will also describe their use in specific applications, such as in understanding heat transfer through a granular material <cit.>.




  §.§.§ Cycles.
 


A cycle (i.e., a closed walk) in a network is a walk that begins and ends at the same node <cit.>. As with other walks, one way to characterize a cycle is by calculating its length. An l-cycle is a cycle in which l edges are traversed (counting repeated edges the number of times that they appear in the cycle). A simple cycle is a cycle that does not include repeated nodes or edges, aside from one repetition of the origin node at the termination of a closed cycle. Thus, for example, a simple 3-cycle in an undirected network is a triangle. For the remainder of this review, we assume that cycles are simple cycles, unless we explicitly state otherwise. In the context of a granular packing, one can directly map particle contact loops — sets of physically-connected grains arranged in a circuit — to cycles in a corresponding graphical representation (see Fig. <ref>d). The length l is odd for an odd cycle and even for an even cycle.

We briefly note a few related concepts that are used to examine cycles in graphs because of their relevance to several network-based studies of granular materials. These are the notions of cycle space, cycle basis, and minimum cycle basis <cit.>. The cycle space of an undirected graph is the set of all simple cycles in a graph along with all subgraphs that consist of unions of edge-disjoint simple cycles (i.e., they can share nodes but not edges) <cit.>. A cycle basis is a minimal set of simple cycles such that any element of the cycle space can be written as a symmetric difference of cycles in the cycle basis <cit.>. Finally, for unweighted networks, a minimum cycle basis is a basis in which the total length of all cycles in the basis is minimal. For weighted networks, it is a basis in which the sum of the weights of all cycles in the basis is minimal.

Minimum cycle bases can provide useful information about the structure and organization of cycles in a network, so several algorithms have been developed to extract them (see, for example, <cit.>). Once one has determined a minimum cycle basis, one can examine the distribution of cycle lengths or define measures to quantify the participation of different nodes in cycles of different lengths. For example, <cit.> defined the concept of a cycle-participation vector X_i^cycle = [x_i^0,x_i^3,…,x_i^l] for each node i. The elements of this vector count the number of cycles of each length in which node i participates. In this definition, x_i^3 is the number of 3-cycles in which node i participates, x_i^4 is the number of 4-cycles in which node i participates, and so on (up to cycles of length l). If a node is not part of any cycle, then x_i^0 = 1 and x_i^j = 0 for all j ≥ 3; otherwise, x_i^0 = 0.

One reason to examine cycles in granular networks <cit.> is that they can help characterize mesoscale structural features of a network. Cycles (that are nontrivial) involve more than a single node, but they do not typically embody global structures of a large network. This makes them appealing for studying network representations of granular materials, because mesoscale features seem to play a role in the behavior of these systems <cit.>. Perhaps the most important motivation, however, is that cycles appear to be relevant for stability and rigidity of a system. Specifically, in the context of structural rigidity theory, 3-cycles tend to be stabilizing structures that can maintain rigidity under applied forces <cit.>, whereas 4-cycles can bend or deform (see Sec. <ref>). In Sec. <ref>, we discuss in more detail how cycles can help characterize granular systems.



  §.§.§ Clustering coefficients.



Clustering coefficients are commonly-used diagnostics to measure the density of triangles either locally or globally in a network <cit.>. For an unweighted, undirected network, the local clustering coefficient C_i is usually defined as the number of triangles involving node i divided by the number of triples centered at node i <cit.>. A triple is a set of three nodes that can include either three edges (to form a 3-cycle) or just two of them. In terms of the adjacency matrix and node degree, the local clustering coefficient is

    C_i = ∑_hj A_hj A_ih A_ij/k_i(k_i - 1)

for k_i ≥ 2 (and C_i = 0 if k_i ∈{0,1}). One can then calculate a global clustering coefficient of a network as the mean of C_i over all nodes:

    C = 1/N∑_i^N C_i .

There is also another (and simpler) common way of defining a global clustering coefficient in a network that is particularly useful when trying to determine analytical approximations of expectations over ensembles of random graphs <cit.>.



The notion of a local clustering coefficient has also been extended to weighted networks in several ways <cit.>. In one formulation <cit.>, a local, weighted clustering coefficient C^w_i is defined as

    C^w_i = 1/s_i(k_i - 1)∑_j,h(W_ij + W_ih)/2 A_ij A_ih A_jh

for strength s_i > 0 and degree k_i ≥ 2. The quantity C^w_i = 0 if either s_i = 0 (so that k_i = 0) or k_i = 1. 
Recall that 𝐖 and 𝐀 are, respectively, associated weighted and unweighted adjacency matrices. The mean of C^w_i over all nodes gives a weighted clustering coefficient C^w of a network. As we will discuss later (see Secs. <ref> and <ref>), clustering coefficients have been employed in several studies of granular materials. For example, they have been used to examine stability in granular packings <cit.>. See Fig. <ref>a for an example of the spatial distribution of a clustering coefficient in a granular packing. 




  §.§.§ Centrality measures.



In network analysis, one calculates centrality measures to attempt to quantify the importance of particular nodes, edges, or other structures in a network <cit.>. Different types of centralities characterize importance in different ways. The degree centrality (i.e., degree) of a node, for example, is simply the number of edges attached to it (see Sec. <ref>). A few other types of centrality that have been used to study granular materials are closeness centrality, node betweenness centrality, edge betweenness centrality, and subgraph centrality.

Notions of closeness centrality of a node measure how close that node is to other nodes in a network <cit.>. For a given node i, the most standard notion of closeness is defined as the inverse of the sum over the shortest-path lengths from node i to all other nodes j in a network. That is, node i's closeness centrality is

    H_i = N-1/∑_j ≠ i d_ij .

Note that if we use the convention that the distance between two nodes in different components is infinite, then Eq. (<ref>) only makes sense for connected networks. For any network with more than one component, Eq. (<ref>) yields a closeness centrality of 0.

The geodesic node betweenness centrality of node i is the fraction of geodesic paths (either unweighted or weighted) between distinct nodes (not including i) that traverse node i <cit.>. Let ψ_gh(i) denote the number of geodesic paths from node g to node h that traverse node i (with i ∉{g,h}), and let ψ_gh denote the total number of geodesic paths from node g to node h. The geodesic node betweenness centrality of node i is then

    B_i = ∑_g,h; g ≠ hψ_gh(i)/ψ_gh ,    i ∉{g,h} .


Geodesic node betweenness can probe the heterogeneity of force patterns in granular networks. See Fig. <ref>b for an example spatial distribution of a geodesic node betweenness centrality in an experimental granular packing. One can also compute a geodesic edge betweenness centrality of an edge by calculating the fraction of shortest paths (either unweighted or weighted) that traverse it <cit.>. Let ψ_gh(i,j) denote the number of geodesic paths from node g to node h that traverse the edge that is attached to nodes i and j, and let ψ_gh denote the total number of geodesic paths from node g to node h. The geodesic edge betweenness centrality of this edge is then

    B^e_ij = ∑_g,h;  g ≠ hψ_gh(i,j)/ψ_gh .


Another measure of node importance is subgraph centrality Y <cit.>, which quantifies a node's participation in closed walks of all lengths. Recall from Sec. <ref> that one can write the number of length-l walks from node i to node j in terms of powers of the adjacency matrix 𝐀. To calculate closed walks of length l that begin and end at node i, we take i = j in Eq. (<ref>). The subgraph centrality of node i, with a specific choice for how much we downweight longer paths, is then given by

    Y_i = ∑_l=0^∞[𝐀^l]_ii/l! .


Because shorter walks are weighted more strongly than longer walks in Eq. (<ref>), they contribute more to the value of subgraph centrality. (In other contexts, centrality measures based on walks have also been used to compare the spatial efficiencies of different networks, and such ideas are worth exploring in granular materials <cit.>.) One can also express subgraph centrality in terms of the eigenvalues and eigenvectors of the adjacency matrix <cit.>. Let v^i_α denote the ith component of the αth eigenvector v_α of 𝐀, and let λ_α denote the corresponding αth eigenvalue. 

One can then write

    Y_i = ∑_α=1^n (v^i_α)^2 e^λ_α .


One can then calculate a mean subgraph centrality Y by averaging Y_i over the nodes in a network. In one study of granular materials <cit.>, a subgraph centrality was examined for weighted networks by considering the eigenvalues and eigenvectors of the weight matrix 𝐖 in Eq. (<ref>). One can also compute network bipartivity R <cit.> to quantify the contribution to mean subgraph centrality Y from closed walks of even length. In particular, the network bipartivity R_i of node i is 

    R_i = Y^even_i/Y_i ,
 
where Y^even_i is the contribution to the sum in Eq. (<ref>) from even values of l (i.e., even-length closed walks). As with other node diagnostics, one can average bipartivity over all nodes in a network to obtain a global measure, which we denote by R. 

In Sec. <ref>, we will discuss calculations of closeness, betweenness, and subgraph centralities in granular packings. Obviously, our discussion above does not give an exhaustive presentation of centrality measures, and other types of centralities have also been used in studies of granular materials (see, for example, <cit.>).
	


  §.§.§ Subgraphs, motifs, and superfamilies.
 


One can interpret the local clustering coefficient in Eq. (<ref>) as a relationship between two small subgraphs: a triangle and a connected triple. Recall that a subgraph of a graph G is a graph constructed using a subset of G's nodes and edges. Conceptually, one can interpret small subgraphs as building blocks or subunits that together can be used to construct a network. For example, in a directed network, there exist three possible 2-node subgraphs (i.e., dyads): the dyad in which node i is adjacent to node j by a directed edge, the dyad in which node j is adjacent to node i by a directed edge, and the dyad in which both of these adjacencies exist. In a directed, unweighted graph, there are 13 different connected 3-node subgraphs <cit.> (see Fig. <ref>). 



The term motif is sometimes used for a small subgraph that occurs often in a particular network or set of networks (typically relative to some null model, such as a randomly rewired network that preserves the original degree distribution) <cit.>. Borrowing terminology from genetics, these motifs appear to be overexpressed in a network (or set of networks). Unsurprisingly, the number of n-node subgraphs increases very steeply with n, so identifying subgraphs in large networks is computationally expensive, and many algorithms have been developed to estimate the number of subgraphs in an efficient (though approximate) way. See, for example, <cit.>. In applying algorithms for motif counting to data, one seeks to identify subgraphs that are present more often than expected in some appropriate random-network null model. 

The over-representation of a motif in a network is often interpreted as indicative of its playing a role in the function of that network (though one has to be cautious about drawing such conclusions). For example, 3-node motifs can form feedforward loops in which there are directed edges from node i_1 to node i_2, from node i_2 to node i_3, and from node i_3 to node i_1. The identification and characterization of motifs has yielded insights into the structure and function of a variety of systems, including food webs <cit.>, gene-regulation networks of yeast <cit.>, neuronal networks of the macaque monkey <cit.>, and others. For different types of networks, one can also identify so-called superfamilies, which are sets of networks that have similar motif-frequency distributions <cit.>. There also exists a less-stringent definition of a superfamily in which one disregards whether a subgraph is a motif in the sense of it being more abundant than expected from some random-graph null model and instead considers a superfamily to be a set of networks that have the same rank-ordering of the number of n-node subgraphs for some fixed value of n <cit.>. In either case, one can examine different superfamilies to help understand the role that specific motifs (or subgraphs) or sets of motifs (or subgraphs) may have in potentially similar functions of networks in a given superfamily. 

Subgraphs, motifs, and superfamilies have been examined in several studies that applied network analysis to granular materials <cit.>. They have revealed interesting insights into the deformation and reconfiguration that occurs in granular systems for different types of loading conditions and external perturbations. We discuss these ideas further in Secs. <ref> and <ref>.



  §.§.§ Community structure.
 


Many real-world networks also have structure on intermediate scales (mesoscales) that can arise from particular organizations of nodes and edges <cit.>. The most commonly-studied mesoscale network property is community structure <cit.>, which describes sets of nodes, called communities, that are densely (or strongly) interconnected to each other but only weakly connected to other dense sets of nodes. In other words, a community has many edges (or large total edge weight, in the case of weighted networks) between its own nodes, but the number and/or weight of edges between nodes in different communities is supposed to be small. Once one has detected communities in a network, one way to quantify their organization is to compute and/or average various network quantities over the nodes (or edges) within each community separately, rather than over an entire network. For example, one can compute the size (i.e., number of nodes) of a community, mean path lengths between nodes in a given community, or some other quantity to help characterize the architecture of different communities in a network. Studying community structure can reveal useful insights about granular systems, whose behavior appears to be influenced by mesoscale network features <cit.>. 

Community structure and methods for detecting communities have been studied very extensively <cit.>. We will briefly discuss the method of modularity maximization <cit.>, in which one optimizes an (occasionally infamous) objective function known as modularity, as this approach has been employed previously in several studies of granular materials (see, e.g., Sec. <ref>). However, myriad other approaches exist for studying community structure in networks. These include stochastic block models (SBMs) and other methods for statistical inference (which are increasingly favored by many scholars) <cit.>, approaches based on random walks (e.g., InfoMap <cit.>), various methods for detecting local community structure (see, e.g., <cit.>), edge-based communities <cit.>, and many others.

The goal of modularity maximization is to identify communities of nodes that are more densely (or more strongly) interconnected with other nodes in the same community than expected with respect to some null model. To do this, one maximizes a modularity objective function

    Q = ∑_i,j [W_ij - γ P_ij] δ(g_i,g_j) ,

where g_i is the community assignment of node i and g_j is the community assignment of node j, and where the Kronecker delta δ(g_i,g_j)=1 if g_i = g_j and δ(g_i,g_j)=0 otherwise. The quantity γ is a resolution parameter that adjusts the relative average sizes of communities <cit.>, where smaller values of γ favor larger communities and larger values of γ favor smaller communities <cit.>. The element P_ij is the expected weight of the edge between node i and node j under a specified null model. In many contexts, the most common choice is to determine the null-model matrix elements P_ij from the Newman–Girvan (NG) null model <cit.>, for which

    P^NG_ij = s_i s_j/2m ,
	
where s_i=∑_j W_ij is the strength (and, for unweighted networks, the degree k_i) of node i and m=1/2∑_i,j W_ij is the total edge weight (and, for unweighted networks, the total number of edges) in the network. There are several other null models, which are usually based on a random-graph model, and they can incorporate system features (such as spatial information) in various ways <cit.>. In the next part of this subsubsection, we discuss a physically-motivated null model that is useful for studying granular force networks.



Maximizing Q is NP-hard <cit.>, so it is necessary to use computational heuristics to identify near-optimal partitions of a network into communities of nodes <cit.>. Two well-known choices are the Louvain <cit.> and Louvain-like <cit.> locally greedy algorithms, which begin by placing all nodes in their own community, and they then iteratively agglomerate nodes when the resulting partition increases modularity Q. Because of the extreme near degeneracy of the modularity landscape (a very large number of different partitions can have rather similar values of the scalar Q), it is often useful to apply such an algorithm many times to construct an ensemble of partitions, over which one can average various properties to yield a consensus partition <cit.>.

 

Physical considerations. Community-detection tools, such as modularity maximization, have often been applied to social, biological, and other networks <cit.>. In applying these techniques to granular materials, however, it is important to keep in mind that the organization of particulate systems (such as the arrangements of particles and forces in a material) is subject to significant spatial and physical constraints, which can severely impact the types of organization that can arise in a corresponding network representation of the material. When studying networks that are embedded in real space or constructed via some kind of physical relationship between elements, it is often crucial to consider the spatial constraints — and, more generally, a system's underlying physics — and their effects on network architecture <cit.>. Such considerations also impact how one should interpret network diagnostics such as path lengths and centrality measures, the null models that one uses in procedures such as modularity maximization, and so on. The NG null model was constructed to be appropriate for networks in which a connection between any pair of nodes is possible. Clearly, in granular materials — as in other spatially-embedded systems <cit.> — this assumption is unphysical and therefore problematic. 

Bassett et al. <cit.> defined a null model that accounts explicitly for geographical (and hence spatial) constraints in granular materials, in which each particle can contact only its nearest neighbors <cit.>. In the context of granular networks with nodes representing particles and edges representing forces between those particles, the geographical null model P in <cit.> has matrix elements

    P_ij = ρ A_ij ,

where ρ is the mean edge weight in the network and A is the binary adjacency matrix of the network. In this particular application, ρ = f := ⟨ f_ij⟩ is the mean inter-particle force. As we illustrate in Fig. <ref>, modularity maximization with the geographical null model [Eq. (<ref>)] produces different communities than modularity maximization with the NG null model [Eq. (<ref>)] <cit.>.  




 

Generalization of modularity maximization to multilayer networks. Although studying community structure in a given granular packing can provide important insights, one is also typically interested in how such mesoscale structures reconfigure as a material experiences external perturbations, such as those from applied compression or sheer. To examine these types of questions, one can optimize a multilayer generalization of modularity to study multilayer granular force networks in which each layer represents a network at a different step in the evolution of the system (for example, at different time steps or for different packing fractions) <cit.>. In Fig. <ref>, we show a schematic of a multilayer construction that has been employed in such investigations. See <cit.> for reviews of multilayer networks (including generalizations of this construction).

One way to detect multilayer communities in a network is to use a generalization of modularity maximization <cit.>, which was derived for multilayer networks with interlayer edges between counterpart nodes in different layers. For simplicity, we suppose that all edges are bidirectional. One maximizes

    Q_multi = 1/2 η∑_ijqr [(𝒲_ijq - γ_q𝒫_ijq) δ_qr + ω_jqrδ_ij] δ(g_iq,g_jr) ,


where 𝒲_ijq is the (i,j)th component of the qth layer of the adjacency tensor 𝒲 <cit.> associated with the multilayer network, 𝒫_ijq is the (i,j)th component of the qth layer of the null-model tensor, γ_q is a resolution parameter (sometimes called a structural resolution parameter) for layer q, and ω_jqr is the interlayer coupling between layers q and r. (In the context of multilayer representations of temporal networks, if ω_jqr = ω for all j, q, and r, one can interpret ω as a temporal resolution parameter.)
More specifically, ω_jqr is the strength of the coupling that links node j in layer q to itself in layer r. (This type of interlayer edge, which occurs between counterpart nodes in different layers, is called a diagonal edge <cit.>.) The quantities g_iq and g_jr, respectively, are the community assignments of node i in layer q and node j in layer r. The intralayer strength of node j in layer q is s_jq = ∑_i𝒲_ijq, and the interlayer strength of node j in layer q is ζ_jq = ∑_rω_jqr, so the multilayer strength of node j in layer q is given by κ_jq = s_jq + ζ_jq. Finally, the normalization factor η = 1/2∑_jqκ_jq is the total strength of the adjacency tensor.[In the study of multilayer networks, it is common to use the term “tensor” to refer to a multidimensional array <cit.> (as is common in some disciplines <cit.>), and proper tensorial structures have been explored briefly in adjacency tensors <cit.>.] 

Maximizing multilayer modularity [Eq. (<ref>)] allows one to examine phenomena such as evolving communities in time-dependent networks or communities that evolve with respect to some other parameter, and communities in networks with multiple types of edges. Capturing such behavior has been useful in many applications, including financial markets <cit.>, voting patterns <cit.>, international relations <cit.>, international migration <cit.>, disease spreading <cit.>, human brain dynamics <cit.>, and more. In the context of granular matter, multilayer community detection allows one to examine changes in community structure of a force network, in which communities can either persist or reconfigure, with respect to both particle content and the mean strength of nodes inside a community, due to applied loads on a system.  



  §.§.§ Flow networks.
 


One can examine many natural and engineered systems — such as animal and plant vasculature, fungi, and urban transportation networks <cit.> — from the perspective of flow networks (which are often directed) that transport a load (of fluids, vehicles, and so on) along their edges. It is of considerable interest to examine how to optimize flow through a network <cit.>. A well-known result from optimization theory is the maximum-flow–minimum-cut theorem <cit.>: for a suitable notion of flow and under suitable assumptions, the maximum flow that can pass from a source node to a sink node is given by the total weight of the edges in the minimum cut, which is the set of edges with smallest total weight that, when removed, disconnect the source and the sink. A related notion, which applies to networks in which there is some cost associated with transport along network edges, is that of maximum-flow–minimum-cost. In this context, one attempts to find a route through a network that maximizes flow transmission from source to sink, while minimizing the cost of flow along network edges <cit.>. The maximum-flow–minimum-cut and maximum-flow–minimum-cost problems are usually examined under certain constraints, such as flow conservation at each node and an upper bound (e.g., limited by a capacitance) on flow through any edge. One can examine granular networks from such a perspective by considering a flow of force along a network formed by contacting grains. We discuss relevant studies in Sec. <ref>. 
	


  §.§.§ Connected components and percolation.



Sometimes it is possible to break a network into connected subgraphs called components (which we introduced briefly in Sec. <ref>). A component, which is sometimes called a cluster, is a subgraph G_C of a graph G such that at least one path exists between each pair of nodes in G_C <cit.>. Components are maximal subsets in the sense that the addition of another node of G to it destroys the property of connectedness. An undirected graph is connected when it consists of a single component. Networks with more than one component often have one component that has many more nodes than the other components, so there can be one large component and many small components. One can find the components of a graph using a breadth-first search (BFS) algorithm <cit.>, and one can determine the number of components by counting the number of 0 eigenvalues of a graph's combinatorial Laplacian matrix <cit.>. To study graph components, one can also use methods from computational algebraic topology. Specifically, the zeroth Betti number β_0 indicates the number of connected components in a graph <cit.> (see Sec. <ref>). 

Percolation theory <cit.>, which builds on ideas from subjects such as statistical physics and probability theory, is often used to understand the emergence and behavior of connected components in a graph <cit.>. For example, in the traditional version of what is known as bond percolation (which is also traditionally studied on a lattice rather than on a more general network) <cit.>, edges are occupied with probability p, and one examines quantities such as the size distributions of connected components as a function of the parameter p, called the bond occupation probability. It is especially interesting to determine a critical value p_c, called the percolation threshold, at which there is a phase transition: below p_c, there is no percolating component (or cluster), which spans the system and connects opposite sides; above p_c, there is such a cluster <cit.>. In the latter scenario, it is common to say that there is a “percolating network". In percolation on more general networks, one can study how the size of the largest component, as a fraction of the size of a network, changes with p. Related ideas also arise in the study of components in Erdős–Rényi random graphs G(N,p), in which one considers an ensemble of N-node graphs and p is the independent probability that an edge exists between any two nodes <cit.>. In the limit N →∞, the size of the largest connected component (LCC) undergoes a phase transition at a critical probability p_c = 1/N. When p < p_c, the ER graph in expectation does not have a giant connected component (GCC); at p = p_c, a GCC emerges whose size scales linearly with N for p > p_c. Similarly, for bond percolation on networks, a transition occurs at a critical threshold p_c, such that for p > p_c, there is a GCC (sometimes also called a “giant cluster” or “percolating cluster”) whose size is a finite fraction of the total number N of nodes as N →∞ <cit.>. When studying percolation on networks, quantities of interest include the fraction of nodes in the LCC, the mean component size, the component-size distribution, and critical exponents that govern how these quantities behave just above the percolation threshold <cit.>.

We will see in Sec. <ref> that it can be informative to use ideas from percolation theory to study the organization of granular networks. For example, it is particularly interesting to examine how quantities such as the number and size of connected components evolve as a function of packing density (or another experimental parameter). <cit.>. Some studies have considered connectivity percolation transitions, which are characterized by the appearance of a connected component that spans a system (i.e., a percolating cluster, as reflected by an associated GCC in the infinite-size limit of a network); or rigidity percolation transitions, which can be used to examine the transition to jamming <cit.>. Rigidity percolation is similar to ordinary bond percolation (which is sometimes used to study connectivity percolation), except that edges represent the presence of rigid bonds between network nodes <cit.> and one examines the emergence of rigid clusters in the system as a function of the fraction of occupied bonds.  One can also study percolation in force networks by investigating the formation of connected components and the emergence of a percolating cluster of contacts as a function of a force threshold, which is a threshold applied to a force-weighted adjacency matrix (representing contact forces between particles) to convert it to a binary adjacency matrix  <cit.>. (See Sec. <ref> for additional discussion.) However, it is important to note that when studying networks of finite size, one needs to be careful with claims about GCCs and percolation phase transitions, which are defined mathematically only in the limit of infinite system size.




  §.§.§ Methods from algebraic topology and computational topology.
 


The tools that we have described thus far rely on the notion of a dyad (i.e., a 2-node subgraph) as the fundamental unit of interest (see Fig. <ref>a). However, recent work in algebraic topology and computational topology <cit.> offers a complementary view, in which the fundamental building blocks that encode relationships between elements of a system are k-simplices (each composed of k+1 nodes), rather than simply nodes and dyadic relations between them (see Fig. <ref>b). These structures can encode “higher-order” interactions and can be very useful for understanding the architecture and function of real-world networks (e.g., they yield a complementary way to examine mesoscale network features), and they have been insightful in studies of sensor networks <cit.>, contagion spreading <cit.>, protein interactions <cit.>, neuronal networks <cit.>, and many other problems. See <cit.> for further discussion and pointers to additional applications. The discussion in <cit.> is also useful.



A collection of simplices that are joined in a compatible way is called a simplicial complex, which is a generalization of a graph that can encode non-dyadic relations <cit.>. More precisely, and following <cit.>, we define an (abstract) simplicial complex 𝒳 as a pair of sets: V_𝒳, called the vertices (or nodes); and S_𝒳, called the simplices, each of which is a finite subset of V_𝒳, subject to the requirement that if σ∈ S_𝒳, then every subset τ of σ is also an element of S_𝒳. A simplex with k elements is called a (k-1)-simplex, and subsets τ⊂σ are called faces of σ. Using this notation, a 0-simplex is a node, a 1-simplex is an edge and its two incident nodes (i.e., a dyad), a 2-simplex is a filled triangle, and so on (see Fig. <ref>b). One type of simplicial complex that can be used to encode the information in a graph is a clique complex (sometimes also called a flag complex); we show an example in Fig. <ref>. To construct the clique complex of a graph G, one associates every k-clique (a complete — i.e., fully connected — subgraph of k nodes) in G with a (k-1)-simplex. One can thus think of building the clique complex of a graph G as “filling in” all of the k-cliques in G (see Fig. <ref>c). Note that we use the terms k-simplex and k-clique because they are standard, but it is important not to confuse the use of k in this context with the use of k as the (also standard) notation for node degree.

One important feature of a simplicial complex is the potential presence of cycles.[Although we use the term cycle, which is standard in algebraic topology, note that this concept of a cycle is distinct from (though related to) the standard network-science use of the word “cycle” (see Sec. <ref>). The latter is sometimes called a circuit, a term that we will use occasionally for clarity (especially given our focus on connected graphs).] A cycle can consist of any number of nodes, and a k-dimensional cycle is defined as a closed arrangement of k-simplices, such that a cycle has an empty boundary[The precise mathematical definition of a cycle requires a more detailed presentation than what we include in our present discussion. For more information and further details from a variety of perspectives, see <cit.>.].
For example, Fig. <ref>d illustrates a closed arrangement of 1-simplices (i.e., edges) that forms a 1-dimensional cycle. It important to distinguish between cycles that encircle a region that is filled by simplices with cycles that enclose a void (which is often called a “hole” for the case of 1-dimensional cycles). For example, the set of purple edges in the object in the upper portion of Fig. <ref>d constitute a 1-dimensional cycle that surrounds a region filled by 2-simplices (i.e., filled triangles), whereas the purple edges in the object in the lower portion of Fig. <ref>d constitute a 1-dimensional cycle that encloses a hole.
 
Characterizing the location and prevalence of void-enclosing cycles in the clique complex of a network representation of a granular packing can offer fascinating insights into the packing's structure <cit.>. One way to do this is by computing topological invariants such as Betti numbers <cit.>. 
The kth Betti number β_k counts the number of inequivalent k-dimensional cycles that enclose a void, where two k-dimensional cycles are equivalent if they differ by a boundary of a collection of (k+1)-simplices. In other words, the kth Betti number β_k counts the number of nontrivial equivalence classes of k-dimensional cycles and can thus also be interpreted as counting the number of voids (i.e., “holes" of dimension k).[In the literature, it is common to abuse terminology and refer to an equivalence class of k-dimensional cycles simply as a k-dimensional cycle.] The zeroth Betti number β_0 gives the number of connected components in a network, the first Betti number β_1 gives the number of inequivalent 1-dimensional cycles that enclose a void (i.e., it indicates loops), the second Betti number β_2 gives the number of inequivalent 2-dimensional cycles that enclose a void (i.e., it indicates cavities), and so on. 

Another useful way to examine the topological features that are determined by equivalence classes of k-dimensional cycles (i.e., components, loops, cavities, and so on) is to compute persistent homology (PH) of a network. For example, to compute PH for a weighted graph, one can first decompose it into a sequence of binary graphs. One way to do this is to begin with the empty graph and add one edge at a time in order of decreasing edge weights (see Fig. <ref>e). More formally and following <cit.>, this process can translate information about edge weights into a sequence of binary graphs as an example of what is called a filtration <cit.>. The sequence G_0 ⊂ G_1 ⊂…⊂ G_|ℰ| of unweighted graphs begins with the empty graph G_0, and one adds one edge at a time (or multiple edges, if some edges have the same weight) in order from largest edge weight to smallest edge weight. (One can also construct filtrations in other ways). Constructing a sequence of unweighted graphs in turn yields a sequence of clique complexes <cit.>, allowing one to examine equivalence classes of cycles as a function of the edge weight θ (or another filtration parameter). Important values of θ include the weight θ_birth associated with the first graph in which an equivalence class (i.e., a topological feature) occurs (i.e., its birth coordinate) and the edge weight θ_death associated with the first graph in which the feature disappears (i.e., its death coordinate), such as by being filled in with higher-dimensional simplices or by merging with an older feature. One potential marker of the relative importance of a particular feature (a component, a loop, and so on) in the clique complex is how long it persists, as quantified by its lifetime θ_birth - θ_death (although short-lived features can also be meaningful <cit.>). A large lifetime indicates robust features that persist over many values of a filtration parameter. Persistence diagrams (PDs) are one useful way to visualize the evolution of k-dimensional cycles with respect to a filtration parameter. PDs encode birth and death coordinates of features as a collection of persistence points (θ_birth,θ_death) in a planar region. One can construct a PD for each Betti number: a β_0 PD (denoted by PD_0) encodes the birth and death of components in a network, a β_1 PD (denoted by PD_1) encodes the birth and death of loops, and so on.



To demonstrate some key aspects of a filtration, the birth and death of topological features, and PDs, we borrow and adapt an example from <cit.>. Consider the small granular force network in Fig. <ref>a; the nodes represent particles in a 2D granular packing, and the colored edges represent the magnitude of the inter-particle forces (of which there are four distinct values) between contacting particles. In a 2D system like this one, the only relevant Betti numbers are β_0 and β_1, as all others are 0. In Fig. <ref>b, we show the flag complex (which is essentially the same as a clique complex <cit.>) of the granular network, where the color of a triangle indicates the value corresponding to the minimum force along any of its edges. Computing PH on a flag complex (which has been done in several studies of PH in granular force networks <cit.>) only counts loops that include 4 or more particles. That is, it does not count 3-particle loops (which are sometimes called “triangular loops”). Loops with 4 or more particles are associated with defects, because they would not exist in a collection of monosized disks that are packed perfectly and as densely as possible into a “crystalline" structure (which has only triangular loops) <cit.>. 

In Fig. <ref>c–f, we show the sequence of complexes that correspond to the filtration over the flag complex. One descends the four threshold levels (i.e., edge weights), beginning with the largest (θ_4) and ending with the smallest (θ_1). In Fig. <ref>g,h, we show the corresponding PDs for β_0 and β_1. It is helpful to discuss a few features of these diagrams. In PD_0, we observe four points that are born at θ_4; these points correspond to the four connected components that emerge at the first level of the filtration in Fig. <ref>c. Two of the components merge into one component at θ_3 (see Fig. <ref>d); this corresponds to the point at (θ_4,θ_3). A new component forms at θ_3 and dies at θ_2; this is represented by the point at (θ_3,θ_2) (see Fig. <ref>d). Additionally, two components born at θ_4 die at theta_2, corresponding to the two points at (θ_4,θ_2). One can continue this process until the end of the filtration, where there is just a single connected component (see Fig. <ref>f). This component is born at θ_4; it persists for all thresholds, and we use <cit.>'s convention to give it a death coordinate of -1; this yields the persistence point at (θ_4,-1). In PD_1, we observe that a loop emerges at θ_3 (see Fig. <ref>d), and it is then filled by triangles at θ_2 (see Fig. <ref>e), leading to the point at (θ_3,θ_2). Three more loops are born at θ_2 and never die (see Fig. <ref>e); using the convention in <cit.>, we assign these features a death coordinate of 0, so there are three persistence points at (θ_2,0). Finally, one more loop appears at θ_1 and does not die (see Fig. <ref>e); this is represented by a point at (θ_1,0).

<cit.> gave an in-depth exposition of how to apply PH to granular networks, and we refer interested readers to this paper for more information. Because studying PH is a general mathematical approach, it can be applied to different variations of force networks and can also be used on networks constructed from different types of experimental data (e.g., digital image data, particle-position data, or particle-interaction data). <cit.> also discussed a set of measures that can be used to compare and contrast the homology of force networks both within a single system (e.g., at two different packing fractions) and across different systems (e.g., if one uses particles of different sizes or shapes), and they explored the robustness of PH computations to noise and numerical errors. In Sec. <ref>, we further discuss applications of methods from algebraic and computational topology to granular materials.




 §.§ Some considerations when using network-based methods



Because there are many methods that one can use to analyze granular networks and many quantities that one can compute to measure properties of these networks, it is useful to discuss some relationships, similarities, and distinctions between them. Naturally, the meaning of any given network feature depends on how the network itself is defined, so we focus the present discussion on the most common representation of a granular system as a network. (See Sec. <ref> for discussions of other representations.) In this representation (see Fig. <ref>), nodes correspond to particles and edges correspond to contacts between particles. Edge weights can represent quantities such as normal or tangential forces between particles. In this type of granular network, it is important to be aware of which network quantities explicitly take into account spatial information or physical constraints in a system, which consider only network topology, and which consider only network geometry (i.e., both topology and edge weights, but not other information). Granular materials have a physical nature and are embedded in real space, so such considerations are extremely important. For discussions of how such issues manifest in spatial networks more generally, see <cit.>.

One way to explicitly include spatial or physical information into network analysis is to calculate quantities that are defined from some kind of distance (e.g., a Euclidean distance between nodes), whether directly or through a latent metric space, rather than a hop distance.
For example, as discussed in Sec. <ref>, one can define the edge length between two adjacent nodes from the physical distance between them, which allows quantities such as mean shortest path length, efficiency, and some centrality measures to directly incorporate spatial information. However, traditional network features such as degree and clustering coefficient depend only on network connectivity, although their values are influenced by spatial effects. In Sec. <ref>, we also saw that one can incorporate physical constraints from granular networks into community-detection methods by using a geographical null model, rather than the traditional NG null model, in modularity maximization.

Different computations in network analysis can also probe different spatial, topological, or geometrical scales. For example, measures such as degree, strength, and clustering coefficients are local measures that quantify information about the immediate neighborhood of a node. However, measures such as the mean shortest path length and global efficiency are global in nature, as they probe large-scale network organization. In between these extremes are mesoscale structures. A network-based framework can be very helpful for probing various types of intermediate-scale structures, ranging from very small ones (e.g., motifs, such as small cycles) to larger ones (e.g., communities), and tools such as PH were designed to reveal robust structural features across multiple scales. Crucially, although there are some clear qualitative similarities and differences between various network-analysis tools (and there are some known quantitative relationships between some of them <cit.>), it is a major open issue to achieve a precise understanding of the relationships between different network computations. Moreover, in spatially-embedded systems (as in any situation where there are additional constraints) one can also expect some ordinarily distinct quantities to become more closely related to each other <cit.>. Furthermore, the fact that a granular particle occupies a volume in space (volume exclusion) gives constraints beyond what arises from embeddedness in a low-dimensional space.



§ GRANULAR MATERIALS AS NETWORKS



We now review network-based models and approaches for studying granular materials. Over the past decade, network analysis has provided a novel view of the structure and dynamics of granular systems, insightfully complementing and extending traditional perspectives. See <cit.> for reviews of non-network approaches.

Perhaps the greatest advantages of using network representations and associated tools are their natural ability to (1) capture and quantify the complex and intrinsic heterogeneity that manifests in granular materials (e.g., in the form of force chains), and to (2) systematically and quantitatively investigate how the structure and organization of a granular system changes when subjected to external loads or perturbations (such as compression, shear, or tapping). In particular, network science and related subjects provide a set of tools that help quantify structure (and changes in structure) over a range of scales — including local, direct interactions between neighboring particles; larger, mesoscale collections of particles that can interact and reconfigure via more complicated patterns; and system-wide measurements of material (re)organization. It is thought that local, intermediate, and system-wide scales are all important for regulating emergent, bulk properties of granular systems. Because structure at each of these scales can play a role in processes such as acoustic transmission and heat transfer, it can be difficult to obtain a holistic, multiscale understanding of granular materials. For example, microscale particle-level approaches may not take into account collective organization that occurs on slightly larger scales, and continuum models and approaches that rely on averaging techniques may be insensitive to interesting and important material inhomogeneities <cit.>. 

Network representations also provide a flexible medium for modeling different types of granular materials (and other particulate matter). For example, network analysis is useful for both simulation and experimental data of granular materials, and methods from complex systems and network science can help improve understanding of both dense, quasistatically-deforming materials as well as granular flows. In any of these cases, one often seeks to understand how a system evolves over the course of an experiment or simulation. To study such dynamics, one can examine a network representation of a system as a function of a relevant physical quantity that parameterizes the system evolution. For example, for a granular system in which the packing fraction increases in small steps as the material is compressed quasistatically, one can extract a network representation of the system at each packing fraction during the compression process and then study how various features of that network representation change as the packing fraction increases. Even a particular type of granular system is amenable to multiple types of network representations, which can probe different aspects of the material and how it evolves under externally applied loads. For instance, one can build networks based only on knowledge of the locations of particles (which, in some cases, may be the only information available) or by considering the presence or absence of physical contacts between particles. 

If one knows additional information about the elements in a system or about their interactions, one can construct more complicated network representations of it. For example, it has long been known that granular materials exhibit highly heterogeneous patterns of force transmission, with a small subset of the particles carrying a majority of the force along force chains <cit.>. Recall from Sec. <ref> that, broadly speaking, a force chain (which is also sometimes called a force network) is a set of contacts that carry a load that is larger than the mean load <cit.>, and the mean orientation of a force chain often encodes the direction of the applied stress <cit.>. We illustrated an example of force chain structure in Fig. <ref>, and we further discuss force-chain organization in Sec. <ref>. Because of the nature of the distribution of force values and the interesting way in which forces are spatially distributed in a material, it is often very useful to consider network representations of granular materials that take into account information about inter-particle forces (see Sec. <ref>) and to use network-based methods that allow one to quantitatively investigate how the structure of a force network changes when one includes only contacts that carry at least some threshold force. See Sec. <ref> and Sec. <ref>.

In our ensuing discussion, we describe several network constructions that have been used to study granular materials, discuss how they have been investigated using many of the concepts and diagnostics introduced in Sec. <ref>, and review how these studies have improved scientific understanding of the underlying, complex physics of granular systems.



 §.§ Contact networks



A contact network is perhaps the simplest way to represent a granular system. Such networks (as well as the term “contact network”) were used to describe granular packings long before explicitly network science-based approaches were employed to study granular materials; see, for example, <cit.>. The structure of a contact network encodes important information about a material's mechanical properties. As its name suggests, a contact network embodies the physical connectivity and contact structure of the particles in a packing (see Fig. <ref>). In graph-theoretic terms, each particle in the packing is represented as a node, and an edge exists between any two particles that are in physical contact with one another. Note that it may not always be possible to experimentally determine which particles are in physical contact, and one may need to approximate contacts between particles using information about particle positions, radii, and inter-particle distances. (See Sec. <ref> for details.) By definition (and however it is constructed), a contact network is unweighted and undirected, and it can thus be described with an unweighted and undirected adjacency matrix (see Sec. <ref>):

    A_ij = {[ 1 ,  if particle i and j are in contact ,;                          0 ,  otherwise . ].


Because the organization of a contact network depends on and is constrained by the radii of the particles and their locations in Euclidean space, a contact network is a spatially-embedded graph <cit.>. In Sec. <ref>, we will see that this embedding in physical space has important consequences for the extraction of force-chain structures via community-detection techniques (see Sec. <ref>). In Fig. <ref>, we show an example of a contact network generated from a discrete-element-method (DEM) simulation (see Sec. <ref>) of biaxial compression <cit.>. The granular system in this figure is polydisperse, as it has more than two types of particles. (In this case, the particles have different sizes.) If all particles are identical in a granular system, it is called monodisperse; if there are two types of particles in a system, it is called bidisperse. In practice, although the presence or absence of a contact is definitive only in computer simulations, one can set reasonable thresholds and perform similar measurements in experiments <cit.> (see Sec. <ref>). 
It is also important to note that packing geometry and the resulting contact network do not completely define a granular system on their own. In particular, one can associate a given geometrical arrangement of particles with several configurations of inter-particle forces that satisfy force and torque balance constraints and the boundary conditions of a system <cit.>. This is a crucial concept to keep in mind when conducting investigations based only on contact networks, and it also motivates the inclusion of contact forces to construct more complete network representations of granular systems (see Sec. <ref>).

In the remainder of this subsection, we review some of the network-based approaches for characterizing contact networks of granular materials and how these approaches have been used to help understand the physical behavior of granular matter. We primarily label the following subsubsections according to the type of employed methodology. However, we also include some subsubsections about specific applications to certain systems.





  §.§.§ Coordination number and node degree.
 


One can study a contact network in several ways to investigate different features of a granular system. We begin our discussion by associating the mean node degree of a contact network with the familiar and well-studied coordination number (i.e., contact number) Z. Although early investigations of granular materials did not consciously make this connection, the mean degree and coordination number are synonymous quantities. The contact degree k_i of particle i is the number of particles with which i is directly in contact, and one can calculate it easily from an adjacency matrix [see Eq. (<ref>)]. A contact network is undirected, so its adjacency matrix A is symmetric, and its row sum and column sum each yield a vector of node degrees. The mean degree ⟨ k ⟩ of a contact network [ Eq. (<ref>)] is then what is usually known as the mean coordination number (i.e., contact number) Z, and it gives the mean number of contacts per particle. As we noted previously (see Sec. <ref>), Z is an important quantity in granular systems because of its connection with mechanical stability and rigidity — which, loosely speaking, is the ability of a system to withstand deformations — in these systems and its characterization of the jamming transition <cit.> and other mechanical properties. In particular, the condition for mechanical stability — i.e., the condition to have all translational and rotational degrees of freedom constrained so that there is force and torque balance — in a packing of frictionless spheres in d dimensions <cit.> is 

    Z ≥ 2d ≡ Z_iso .

The isostatic number Z_iso indicates the condition for isostaticity, which is defined as the minimum contact number that is needed for mechanical stability. One can use the coordination number (which is often tuned by changing the packing fraction ϕ) as an order parameter to describe the jamming transition for frictionless spheres in two and three dimensions <cit.>. Specifically, there is a critical packing fraction ϕ_c such that below ϕ_c, the contact number for these systems is Z = 0 (i.e. there are no load-bearing contacts), and at the critical packing fraction ϕ_c, the contact number jumps to the critical value Z_c = Z_iso = 2d. One can also generalize the use of the coordination number in order to examine mechanical stability and jamming in granular systems of frictional spheres. In these systems, the condition for stability is

    Z    ≥ Z^m_iso , 
    
    	Z^m_iso   ≡ (d+1) + 2N_m/d ,

where N_m is the mean number of contacts that have tangential forces f_t equal to the so-called Coulomb threshold — i.e., N_m is the mean number of contacts with f_t = μ f_n, where μ is the coefficient of friction and f_n is the normal force <cit.> — and Z^m_iso again designates the condition for isostaticity. Results from experimental systems have demonstrated that contact number also characterizes the jamming transition in frictional, photoelastic disks <cit.>. 

The coordination number has been studied for several years in the context of granular materials and jamming, and it is fruitful to connect it directly with ideas from network science. Several recent studies have formalized the notion of a contact network, and they deliberately modeled granular systems as such networks to take advantage of tools like those described in Sec. <ref>. Such investigations of contact networks allow one to go beyond the coordination number and further probe the rich behavior and properties of granular materials — including stability and the jamming transition <cit.>, force chains <cit.>, and acoustic propagation <cit.>.

Perhaps the simplest expansion of investigations into the role of the coordination number is the study of the degree distribution P(k) of the contact network of a packing. Calculating degree distributions can provide potential insights into possible generative mechanisms of a graph <cit.>, although one has to be very careful to avoid over-interpreting the results of such calculations <cit.>. In granular physics, it has been observed that the degree distribution of a contact network can track changes in network topology past the jamming transition in isotropically compressed simulations of a 2D granular system <cit.>. Specifically, the peak of P(k) shifts from a lower value of k to a higher value of k near the transition. Moreover, changes in the mean degree ⟨ k ⟩ and its standard deviation can anticipate the onset of different stages of deformation in DEM simulations (i.e., molecular-dynamics simulations) of granular systems under various biaxial compression tests <cit.>. 



  §.§.§ Investigating rigidity of a granular system using a contact network.
 


An important area of research in granular materials revolves around attempts to (1) understand how different types of systems respond when perturbed, and (2) determine what features of a system improve its integrity under such perturbations. As we noted in Sec. <ref>, it is well-known that coordination number (and hence node degree) is a key quantity for determining mechanical stability and understanding jamming in granular materials. However, contact networks obviously have many other structural features, and examining them can be very helpful for providing a more complete picture of these systems.

To the best of our knowledge, the stability of granular materials was first studied from a graph-theoretic standpoint in the context of structural rigidity <cit.>, and it has since been applied to amorphous solids more generally <cit.>. In structural rigidity theory, thought to have been studied originally by James Clerk Maxwell <cit.>, rods of fixed length are connected to one another by hinges, and one considers the conditions under which the associated structural graphs are able to resist deformations and support applied loads (see Fig. <ref>). A network is said to be minimally rigid (or isostatic) when it has exactly the number of bars needed for rigidity. This occurs when the number of constraints is equal to the number of degrees of freedom in the system (i.e., when Laman's minimal-rigidity criterion is satisfied). The network is flexible if there are too few rods, and it is overconstrained (i.e., self-stressed) if there are more rods than needed for minimal rigidity. Triangles are the smallest isostatic structures in two dimensions <cit.>; there are no allowed motions of the structure that preserve the lengths and connectivity of the bars, so triangles (i.e., 3-cycles) do not continuously deform due to an applied force. In comparison, a 4-cycle is structurally flexible and can continuously deform from one configuration to another while preserving the lengths and connectivity of the rods (see Fig. <ref>). 



Extending a traditional network of rods and hinges, concepts from structural rigidity yield interesting insights into contact networks of particulate matter. See <cit.> for a discussion of some of the earliest applications of such ideas to disordered systems and granular materials. <cit.> used structural rigidity theory to derive conditions for the isostaticity of a granular packing; and he tied the fact that random packings of particles tend to be isostatic to the origin of instabilities in granular piles. Later, similar concepts were used to show that in granular networks, cycles with an even number of edges allow contacting grains to roll without slipping when subject to shear; however, these relative rotations are “frustrated" in cycles with an odd number of edges, so such cycles can act as stabilizing structures in a network <cit.>. Several later studies (such as <cit.>) have confirmed that contact loops are often stabilizing mesoscale features in a contact network of a granular material. We specifically consider the role of cycles in granular contact networks in Sec. <ref>. 

Another type of network approach for understanding rigidity in granular systems is rigidity percolation <cit.> (see Sec. <ref>). <cit.> conducted an early investigation of an idealized version of bond percolation in a granular context. It is now known that hallmarks of this bond-percolation transition occur below isostaticity: <cit.> identified that a percolating (i.e., system-spanning) cluster of non-load-bearing contacts forms at a packing density below the jamming point. In modern contexts, the rigidity-percolation approach can be used to determine if a network is both percolating and rigid (see Sec. <ref>). Note that a rigid granular network is also percolating, but a percolating network need not be rigid. Rigidity percolation relies on tabulating local constraints via a pebble game <cit.>, which reveals connected, rigid regions (sometimes called “clusters”) in a network. In a series of papers <cit.> on simulated packings, Schwarz and coworkers went beyond Laman's minimal-rigidity criterion to investigate local versus global rigidity in a network, the size distribution of rigid clusters, the important role of spatial correlations, and the necessity of force balance. Building on the above work, <cit.> recently utilized a rigidity-percolation approach to identify floppy versus rigid regions in slowly-sheared granular materials and to characterize the nature of the phase transition from an underconstrained system to a rigid network. See also the recent paper <cit.>.



  §.§.§ Exploring the role of cycles.



We now consider the role of circuits (i.e., the conventional network notion of cycles, which we discussed in Sec. <ref>) in granular contact networks. Cycles in a contact network can play crucial stabilizing roles in several situations. Specifically, as we will discuss in detail in this section, simulations (and some experiments) suggest that (1) odd cycles (especially 3-cycles) can provide stability to granular materials by frustrating rotation among grains and by providing lateral support to surrounding particles, and that (2) a contact network loses these stabilizing structures as the corresponding granular system approaches failure. 

Noting that 3-cycles are the smallest arrangement of particles that can support (via force balance) a variety of 2D perturbations to a compressive load without deforming the contact structure, 
<cit.> studied the effects of friction and tilting on the evolution of contact-loop organization in a granular bed. In their simulations, they implemented tilting by incrementally increasing the angle of a gravity vector with respect to the vertical direction, while preserving the orientation of the granular bed and maintaining quasistatic conditions. In untilted granular packings, they observed that lowering inter-particle friction yields networks with a higher density of 3-cycles and 4-cycles, where they defined the “density” of an l-cycle to be the number of l-cycles divided by the total number of particles. By examining the contact network as a function of tilting angle, <cit.> also observed that the density of 4-cycles increases prior to failure — likely due to the fracture of stabilizing 3-cycles — and that this trend was distinguishable from changes in coordination number alone. 

Cycles have also been studied in the context of DEM simulations of dense, 2D granular assemblies subject to quasistatic, biaxial compression tests <cit.>. In many of these studies, the setup consists of a collection of disks in 2D that are compressed slowly at a constant strain rate in the vertical direction, while allowed to expand under constant confining pressure in the horizontal direction <cit.>. In another variation of boundary-driven biaxial compression, a sample can be compressed with a constant volume but a varying confining pressure <cit.>). Before describing specifics of the network analysis for these systems, it is important to note that for the previously described conditions, the axial strain increases in small increments (i.e., “steps”) as compression proceeds, and one can extract the inter-particle contacts and forces at each strain value during loading to examine the evolution of a system as a function of strain. Additionally, these systems undergo a change in behavior from a solid-like state to a liquid-like state and are characterized by different regimes of deformation as a function of increasing axial strain <cit.>. In particular, the granular material first undergoes a period of strain hardening, followed by strain softening, after which it enters a critical state. In the strain-hardening regime, the system is stable and the shear stress increases monotonically with axial strain up to a peak value. After the peak shear stress, strain softening sets in; this state is marked by a series of steep drops in the shear stress that indicate reduced load-carrying capacity. Finally, in the critical state, a persistent shear band has fully formed, and the shear stress fluctuates around a steady-state value. The shear band is a region of localized deformation and gives one signature of material failure <cit.>. Inside the shear band, force chains can both form and buckle <cit.>. One can also associate increases in the energy dissipation rate of the system with particle rearrangements (such as those that occur during force-chain buckling) and loss of stability in the material.



Examining the temporal evolution of cycles in an evolving granular contact network can reveal important information about changes that occur in a material during deformation. Using DEM simulations (see the previous paragraph), <cit.> computed the total number of cycles of different lengths in a minimal cycle basis (see Sec. <ref>) of a contact network at each strain state during loading, and they observed that there are many more 3-cycles and 4-cycles than longer cycles in the initial, solid-like state of the system. However, as axial strain increases and one approaches the maximum shear stress, the total number of 3-cycles falls off steeply. (The same is true for 4-cycles, though it is less dramatic.) Additionally, during axial-strain steps (i.e., axial-strain “intervals”) corresponding to drops in shear stress, <cit.> observed large increases in the number of 3-cycles and 4-cycles that open up to become longer cycles. In Fig. <ref>a, we show an example of the evolution of cycle organization with increasing axial strain for a subset of particles from a DEM simulation of a granular material under biaxial compression, carried out by <cit.>. The authors observed that in this system, both the global clustering coefficient C [Eq. (<ref>)] and the mean subgraph centrality Y decrease with increasing axial strain, drop sharply at peak shear stress, and then level out (see Fig. <ref>b,c). Recalling that C is a measure of triangle density in a graph and that subgraph centrality measures participation of nodes in closed walks (with more weight given to shorter paths), these results also imply that the loss of small cycles co-occurs with the deformation and failure of a system due to increasing load. <cit.> also computed the network bipartivity R <cit.> of the contact network to quantify the contribution to mean subgraph centrality Y from closed walks of even length [see Eq. <ref>]. They observed that R increases with increasing axial strain, revealing that closed walks of even length become more prevalent during loading (see Fig. <ref>c). The authors suggested that this trend may be due to a decrease in the prevalence of 3-cycles (which are stabilizing, as discussed in Sec. <ref> and elsewhere). <cit.> also examined the stability of cycles of various lengths in both DEM simulations and experimental data, and they observed that, during loading, 3-cycles tend to be more stable (as quantified by a measure of stability based on a structural-mechanics framework <cit.>) than cycles of other lengths in a minimal cycle basis of the network. 

Minimal cycle bases and the easier-to-compute subgraph centrality have also been used to examine fluctuations in kinetic energy in simulations of deforming sand. <cit.> computed a minimal cycle basis and then constructed cycle-participation vectors (see Sec. <ref>) from a contact network after each strain step (i.e., at each strain state) during loading. They observed that temporal changes in the cycle-participation vectors of the particles between consecutive strain steps are correlated positively with temporal changes in kinetic energy over those steps.
They also observed that large values in the temporal changes of particle cycle-participation vectors and particle subgraph centrality occur in the shear-band region. <cit.> also studied a minimal cycle basis and corresponding cycle-participation vectors to examine structural transitions in a 3D experimental granular system of hydrogel spheres under uniaxial compression. As pointed out in <cit.>, developing quantitative predictors that are based on topological information alone is extremely important for furthering understanding of how failure and rearrangements occur in systems in which energy or force measurements are not possible. 

Examining cycles in contact networks can also shed light on the behavior of force chains. The stability, load-bearing capacity, and buckling of force chains depend on neighboring particles (so-called spectator grains) to provide lateral support <cit.>. Because 3-cycles appear to be stabilizing features, it is interesting to consider the co-evolution of force chains and 3-cycles in a contact network. Such an investigation requires a precise definition of what constitutes a force chain, so that it is possible to (1) extract these structures from a given packing of particles and (2) characterize and quantify force-chain properties. Several definitions of force chains have been proposed; see, e.g., <cit.>. The studies that we describe in the next three paragraphs used a notion of “force chains” from <cit.>, in which force chain particles are identified based on their particle-load vectors (where each particle is associated with a single particle-load vector that indicates the main direction of force transmission). More specifically, a single chain is a set of three or more particles for which the magnitude of each of their particle-load vectors is larger than the mean particle-load vector magnitude over all particles, and for which the directions of the particle load vectors are, within some tolerance, aligned with one another (i.e., they are “quasilinear”.) We note that an important point for future work is to conduct network-based studies of force-chain structure for different definitions of force chains, and to investigate if there are qualitative differences in their associated network properties.

Using DEM simulations of a densely packed system of polydisperse disks under biaxial loading — i.e., compressed quasistatically at a constant strain rate in the vertical direction, while allowed to expand under constant confining pressure in the horizontal direction — <cit.> quantified the co-evolution of force chains and 3-cycles in several ways. For example, they computed a minimal cycle basis (see Sec. <ref>) of a contact network and then examined (1) the ratio of 3-cycles to the total number of cycles in which particles from a force chain participate and (2) the force chain's 3-cycle concentration, which is defined as the ratio of 3-cycles involving force-chain particles to the total number of particles in the force chain. When averaged over all force chains, the above two measures decrease rapidly with increased loading. Additionally, <cit.> observed that force chains that do not fail by buckling (see <cit.> for how “buckling" was defined) have a larger ratio of 3-cycle participation to total cycle participation than force chains that do buckle. <cit.> observed, in both DEM simulations of biaxial loading (see above) and 2D photoelastic disk experiments under pure shear, that a particular measure (developed by <cit.>) of structural stability of force chains is correlated positively with the mean of the local clustering coefficient [Eq. (<ref>)] over force-chain particles. Their results also suggest that 3-cycles are more stable structures than cycles of longer length during loading and that force chains with larger 3-cycle participation tend to be more structurally stable. These observations suggest that cycles — and especially 3-cycles — in contact networks are stabilizing structures that can provide lateral support to force chains. It would be interesting to study these ideas further, and to relate them to structural rigidity theory (see Fig. <ref> and Sec. <ref>), especially in light of the difference between 3-cycles (which are rigid) and deformable cycles (e.g., 4-cycles). 
	
DEM simulations of 3D, ellipsoidal-particle systems subject to triaxial compression also suggest that 3-cycles are important features in granular contact networks <cit.>. 
Similar to the aforementioned results from simulations of 2D systems with disk-shaped particles, the number of 3-cycles in a minimal cycle basis of the contact networks (and the global clustering coefficient [Eq. (<ref>)]) initially decrease and then saturate with increasing load, and particles in force chains have a larger number of 3-cycles per particle than particles that are not in force chains. <cit.> also observed that the set of 3-cycles that survive throughout loading tend to lie outside the strain-localization region (where force chains buckle). The dearth of 3-cycles in certain regions in a material may thus be a signature of strain-localization zones. 
Another paper to note is <cit.>, which examined and compared the temporal evolution of cycles (and several other contact-network diagnostics) in a set of DEM simulations using a variety of different material properties and boundary conditions.

In another interesting study, <cit.> examined the phenomenon of aging <cit.> — a process in which the shear strength and stiffness of a granular material increase with time — in collections of photoelastic disks subject to multiple cycles of pure shear under constant volume. Because aging is a slow process, it can be difficult both to uncover meaningful temporal changes in dynamics and to characterize important features in packing structure that accompany aging. To overcome these challenges, <cit.> first analyzed the time series of the stress ratio (using techniques from dynamical-systems theory) to uncover distinct temporal changes in the dynamics of the system. (See <cit.> for details.) After each small, quasistatic strain step, they also extracted the contact network of the packing at that time to relate aging to changes in topological features of the network structure. As one approaches the shear-jammed regime during prolonged cyclic shear, they observed on average that force chains are associated with more 3-cycles and 4-cycles from the minimal cycle basis. 



We have just discussed many papers that concern transitions in granular matter from a solid-like regime to a liquid-like regime. One can also use changes in the loop structure of a contact network to describe the opposite transition, in which a granular material goes from an underconstrained, flowing state to a solid-like state during the process known as jamming (see Sec. <ref>). Studying 2D frictional simulations of isotropically compressed granular packings, <cit.> examined a granular contact network as a function of packing fraction. They observed that the number of cycles (which were called polygons in <cit.>) in the contact network grows suddenly when the packing fraction approaches the critical value ϕ_c that marks the transition to a rigid state (see Fig. <ref>). They also observed that 3-cycles appear to be special: they continue to grow in number above the jamming point, whereas longer cycles slowly decrease in number after jamming. Although they observed a nonlinear relationship near the jamming point between Z (the contact number, which is the usual order parameter for the jamming transition) and the number of 3-cycles <cit.>, these quantities appear to depend linearly on each other after the transition. These results suggest that one can use the evolution of contact loops to understand the transition to a rigid state and to characterize subsequent changes in the system.



 

Application to tapped granular materials. 


Properties of contact networks have also been used to study tapped granular materials, in which a packing of grains is subject to external pulses of excitation against it. In most studies of tapped granular materials, the packing and pulses are both vertical. The intensity Γ of these mechanical perturbations (so-called “taps") is usually quantified as a dimensionless ratio of accelerations, such as the ratio of the peak acceleration of the excitation to the acceleration of gravity <cit.>. Tapped granular materials are interesting, because the packing fraction ϕ is not a monotonic function of the tapping intensity Γ <cit.>. It reaches a minimum value ϕ_min at an intensity of Γ_min, and it then increases as the tap intensity increases (see Fig. <ref>a). Consequently, one can achieve steady states with the same packing fraction by using different tap intensities (i.e., both a “low" tap intensity, which is smaller than Γ_min, and a “high" tap intensity, which is larger than Γ_min). These steady states are not equivalent to each other, as they have different force-moment tensors <cit.>, for example. An interesting question is thus the following: What features of a granular packing distinguish between states at the same packing fraction that are reached by using different tap intensities? 

Recent work has suggested that properties of contact networks — especially cycles (which, in this case, are particle contact loops) — can distinguish between steady-state configurations at the same packing fraction but that are generated from different tap intensities in simulated 2D granular packings subjected to tapping <cit.> (see Fig. <ref>b). For example, as Γ is increased in the regime Γ < Γ_min, the number of 3-cycles (i.e., triangles) and the number of 4-cycles (i.e., squares) both decrease. As Γ is increased in the regime Γ > Γ_min, the opposite trend occurs, so the numbers of 3-cycles and 4-cycles increase. This makes it possible to differentiate configurations at the same ϕ obtained from low and high tap intensities. (See Fig. <ref>e,f for a plot of the number of triangles versus Γ and ϕ.) However, geometrical measures like the pair-correlation function, distributions of Voronoi tessellation areas, or bond orientational order parameters do not seem to be as sensitive to differences in these two different states of the system (see Fig. <ref>c,d), perhaps because they quantify only local proximity rather than directly examining contacts. (See <cit.> and references therein for details of these descriptors.) These results suggest that topological features (e.g., mesoscale features) of a contact network can capture valuable information about the organization of granular packings.



  §.§.§ Other subgraphs in contact networks.



When studying contact networks, it can also be helpful to explore network motifs other than cycles. Recall from Sec. <ref> that motifs are subgraphs that appear often in a network (e.g., in comparison to a null model) and are thus often construed as being related to system function <cit.>. Network motifs, which traditionally refer to small subgraphs, are a type of mesoscale feature, and it can be insightful to examine how their prevalences change in a granular material as it deforms.

One system in which motifs and their dynamics have been studied is frictional, bidisperse, photoelastic disks subject to quasistatic cyclic shear <cit.>. After each small strain increment (i.e., strain step) in a shear cycle, the authors considered the contact network of the granular packing. For each particle i in the contact network, they extracted the subgraph of particles (nodes) and contacts (edges) formed by the central particle i and particle i's contacting neighbors. This process results in a set of N subgraphs[In some degenerate cases (e.g., when there is a network component that consists of a clique), this includes duplications of subgraphs, even when the nodes of a network are labeled.] (which, borrowing terminology from <cit.>, we call conformation subgraphs), where N is the number of particles in the network.

To examine packing rearrangements as a system is sheared, <cit.> represented each distinct conformation subgraph present at any time during loading as one “state" in a Markov transition matrix, and they studied transitions between the conformation subgraphs as a discrete-time Markov process. More specifically, each element in the n_c× n_c transition matrix (where n_c is the total number of unique conformation subgraphs and hence the number of unique states) captured the fraction of times that each conformation subgraph transformed to each other conformation subgraph (including itself) after four quasistatic steps of the shearing process. <cit.> reported that force-chain particles typically occur in network regions with high mean degree, high mean local clustering coefficients, and many 3-cycles. (Note that this study, as well as the others that we describe in the present and the following paragraph, define force-chains as in <cit.>.) Furthermore, when considering the conformation subgraphs of particles in force chains that fail by buckling (see <cit.> for details on the definition of “buckling"), the most likely transformations to occur tend either to maintain the topology of those conformation subgraphs or to involve transitions from conformation subgraphs in which the central particle has a larger degree or is part of more 3-cycles to conformation subgraphs in which the degree of the central particle is smaller or in which it participates in fewer 3-cycles. <cit.> also used force information to compute a measure of structural stability (based on a structural-mechanics framework <cit.> and summarized in a single number) for each conformation subgraph. They then split the full range of the stability values into several smaller “stability intervals” (i.e., small ranges of contiguous structural-stability values) and modeled transitions between stability intervals as a Markov chain. They examined the number of conformation subgraphs that occupy each stability interval and observed pronounced peaks in some intervals that persist during loading. They also reported that conformation subgraphs whose central particles belong to force chains tend to be more stable and that conformation subgraphs whose central particles are part of buckling force chains have a higher probability of transitioning from high-stability states to low-stability states than vice versa. (For details, see Fig. 7 of <cit.> and the corresponding discussions.) 



<cit.> used similar methods to study self-assembly in an almost-frictionless, 3D system of hydrogel spheres under quasistatic, cyclic uniaxial compression (see Fig. <ref>a–c). After every compression step, they constructed the contact network for the system and examined two types of subgraphs for each particle: (1) conformation subgraphs, which (as discussed earlier) consist of a single central particle i and that particle's contacts; and (2) the cycle-participation vector of each particle (see Sec. <ref>). <cit.> determined the set of all unique conformation subgraphs that exist during the above compression process. They then used each of those conformation subgraphs as one state in a transition matrix, the elements of which give the fraction of times (across the whole experiment) that a particle in one state transitions to any other state in consecutive compression steps. 
To focus on the presence or absence of a particle in an l-cycle (using cycle lengths up to l = 10), they binarized each element of the cycle-participation vectors. (The new vectors thus indicate, for each particle, whether a particle is part of at least one l-cycle.) They then constructed a transition matrix in which each state is a unique binarized cycle-participation vector that occurs during the experiment. The two transition matrices capture useful information about the most likely transformations that occur between different conformation subgraphs and cycle-participation vectors as one compresses or decompresses the granular system. For both types of mesoscale structures, <cit.> used their transition matrices to extract almost-invariant sets, which indicate sets of conformation subgraphs or cycle-participation vectors (i.e., states) that tend to transition among themselves more than to states in another almost-invariant set. (See <cit.> for details.) In Fig. <ref>d, we show the most common conformation subgraphs in each almost-invariant set of the 
subgraphs. The conformation subgraphs formed by force-chain particles belong mostly to Set 3 (see Fig. <ref>d), which consists of densely-connected conformation subgraphs in which there are many contacts between particles. To characterize structural changes that occur in a packing as it moves towards or away from a jammed configuration, <cit.> tracked the number of conformation subgraphs (and cycle-participation vectors) in each almost-invariant set across time. In Fig. <ref>e, we show the temporal evolution of the numbers of elements in the almost-invariant sets of the conformation subgraphs. <cit.> also proposed transition pathways (see Fig. <ref>f) that may be useful for thermo-micro-mechanical constitutive-modeling efforts <cit.>. (A transition pathway consists of a sequence of conformation subgraphs in different almost-invariant sets, and transitions between them.) 

Another way to study various types of subgraphs in granular materials is through the classification of superfamilies <cit.> (see Sec. <ref>). A recent investigation by <cit.> considered superfamilies that result from examining 4-particle subgraphs (see Fig. <ref>a) in a variety of different granular systems, including experimental packings of sand and photoelastic disks, and DEM simulations for different types of loading and in different dimensions. In their study, the authors defined a superfamily as a set of networks in which the prevalence of the different 4-node subgraphs have the same rank-ordering. (They did not consider whether the subgraph was a motif in the sense of occurring more frequently than in a random-graph null model.) Despite the diversity of system types, they observed several trends in the transitions between superfamilies that occur as a system transitions from a pre-failure regime to a failure regime. The most important change in the superfamilies appears to be a switch in relative prevalence of 4-edge motifs with 3 edges arranged as a triangle to acyclic 3-edge motifs (see Fig. <ref>). This observation highlights the important role that small mesoscale structures can play as building blocks in granular systems. It also suggests that examining the prevalence and temporal evolution of such motifs can (1) help characterize the macroscopic states of a granular system and (2) help quantify what structural changes occur as a system transitions between different states. 

Notably, although calculating the prevalence of cycles and small motifs can be useful for gaining insights into contact-network structure, it is also important to employ other types of network analysis that examine structure on larger scales. For example, in simulations of 2D packings of disks under isotropic compression, <cit.> observed that the mean shortest-path length (<ref>) of a contact network reflects changes in the organization of a packing as one approaches the jamming point and changes that occur after the jamming transition takes place. The path length appears to reach a maximum value at a packing fraction below ϕ_c. With further increases in ϕ below ϕ_c, the path length then decreases rapidly, likely due to the formation of new contacts that shorten the distance between particles as the system nears the jamming point. After the jamming transition, the path length decreases further.

Before moving on, we note that because it can be difficult to measure force information accurately in many experimental granular systems, continuing to develop relevant measures for studying contact topology (i.e., without incorporating weights) remains an important area of investigation.





 §.§ Force-weighted networks



Although studying only a contact network can be rather informative (see Sec. <ref>), it is important to incorporate more information into network representations to more fully capture the rich behavior of granular materials. Many of the approaches for quantifying unweighted networks can be generalized to weighted networks (see Sec. <ref>), although significant complications often arise (e.g., because typically there are numerous choices for precisely how one should do the generalizing). From both physics and network-science perspectives, it is sensible to construct weighted networks that incorporate information about the forces between particles. This can shed considerable light on phenomena that have eluded understanding from studying only unweighted contact networks.



One important physical motivation for incorporating information about inter-particle forces is that photoelastic-disk experiments and numerical simulations have both highlighted that, particularly just above isostaticity (see the bottom of Sec. <ref>), loads placed on granular systems are not shared evenly across all particles. Instead, forces are carried primarily by a backbone of force chains. It has often been claimed that the statistical distribution of the forces are approximately exponential <cit.>, but in fact, it depends on whether forces are recorded at a granular system's boundary or in its bulk <cit.>, as well as on the loading history <cit.>. Illuminating how force-chain structures arise provides crucial information for understanding how one can control the elastic modulus and mechanical stability <cit.> and acoustic transmission <cit.> in granular materials.

However, despite the ability of humans to see force chains easily in photoelastic images, it is difficult to characterize quantitatively what is or is not a force chain, and it can also be difficult to quantify how force chains evolve under compression or shear. Part of the challenge lies in the fact that force chains are spatially anisotropic, can exhibit long-range spatial correlations <cit.>, and can have complex temporal fluctuations <cit.>. Consequently, understanding emergent phenomena in granular systems is typically difficult using continuum theories or approaches based only on local structure. On the other hand, a network-theoretic perspective provides a fruitful way to explore interesting material properties and organization that arise from inter-particle contact forces in granular materials. Importantly, in addition to data from simulations, multiple techniques now exist for measuring inter-particle forces between grains in experiments; these include photoelasticity <cit.>, x-ray diffraction measurements of microscopic <cit.> or macroscopic <cit.> deformations, and fluorescence with light sheets <cit.>. As we will see, incorporating information about inter-particle forces into network-based investigations has yielded fascinating insights into the organization and collective structure in granular packings (and other particulate materials) for both numerically-simulated and experimental systems.

The most common method for constructing a network that captures the structure of forces in a granular system is to let a node represent a particle and then weight the edge between two particles that are in contact according to the value of the force between them. One can describe such a network with a weighted adjacency matrix (see Sec. <ref>) W with elements

    W_ij = {[ f_ij ,  if particles i and j are in contact ,;                              0 ,  otherwise , ].

where f_ij is the inter-particle force between particles i and j. Such a force network also encodes the information in the associated contact network, and one can recover a contact network from a force-weighted network by setting all non-zero weights in the latter to 1. Although most work has determined edge weights using the normal component of the inter-particle force, one could alternatively weight the edges by inter-particle tangential forces. With the advent of high-performance computational capabilities, one can determine inter-particle forces from DEM simulations <cit.> of hundreds to millions of particles. In experiments, it is possible to determine inter-particle forces using photoelastic disks (in 2D) <cit.> or x-ray tomography (in 3D) <cit.>, although typically these techniques are limited to systems of hundreds to thousands of particles.

We now review network-based approaches for investigating force-weighted networks constructed from granular materials and discuss the resulting insights that these approaches have provided into the physical behavior of such systems. We label most of the following subsubsections according to the type of employed methodology, although we also include a subsubsection about some specific applications to different systems.



  §.§.§ Examining weighted cycles and other structural features.






In Sec. <ref>, we discussed why examining cycles can be useful for studying granular contact networks. It is also useful to examine cycles when investigating force networks, which are weighted. For example, <cit.> studied the evolution of weighted contact loops in a simulation of a quasistatically tilted granular packing. They used topological information (i.e., which particles are in contact) to define the presence of a cycle, and they defined a notion of loop stability, 

    ξ_l = 1/f^l∏_i=1^lf^edge_i ,
	
to quantify the range of compressive loads that a given loop can support. In Eq. (<ref>), l is the number of edges in the loop (i.e., its length), f^edge_i is the contact force of the i^th edge, and f is the mean edge weight (i.e., mean force) over all of the edges in the loop. See Fig. <ref> for a schematic of this stability measure for a 3-cycle. For l = 3, the quantity ξ_3≈ 1 corresponds to having approximately equal contact forces on all edges and is the most stable configuration (see Fig. <ref>a). The value of ξ_3 approaches 0 as the contact force on one edge becomes much smaller than those on the other two edges. As illustrated in Fig. <ref>b, this situation is rather unstable. Both the density of 3-cycles (specifically, the number of 3-cycles in the system divided by the total number of particles) and a normalized 3-cycle loop stability 

    ξ^*_3 = ⟨ξ_3(θ_g) ⟩/⟨ξ_3(θ_g = 0) ⟩
	
(where the brackets denote means over all 3-cycles in a network) tend to decrease with increasing tilting angle θ_g (see Fig. <ref>c). <cit.> also reported that the effect of tilting on loop stability is largely independent from the effect of tilting on mean coordination number (i.e., mean degree).

<cit.> examined what they called a force cycle, which is a cycle of physically-connected particles in which each contact carries a force above the global mean. Using DEM simulations of a biaxially compressed, dense granular system — the sample was compressed quasistatically at a constant strain rate in the vertical direction, while allowed to expand under constant confining pressure in the horizontal direction — they studied the evolution of 3-force cycles (i.e., force cycles with 3 particles) in a minimal cycle basis of the contact network with respect to axial strain. They observed that 3-force cycles initially decrease in number during strain hardening, before increasing in number at the onset of force-chain buckling <cit.>, and finally leveling out in number in the critical-state regime. (See the third paragraph of Sec. <ref> for a brief description of these different regimes of deformation.) 
In Fig. <ref>, we show a plot of the number of 3-force cycles and the shear stress versus axial strain. The 3-force cycles that arise at the onset of buckling are often part of force chains (using the definition from <cit.>). Additionally, these 3-force cycles tend to concentrate in the region of the shear band, where they may act as stabilizing structures both by frustrating relative rotations and by providing strong lateral support to force chains. However, with increased loading, the system eventually fails, and <cit.> suggested that the increase in the number of 3-force cycles may be an indicator of failure onset. Qualitatively similar results have been observed when examining the evolution of 3-force cycles in three DEM simulations (each with slightly different material properties and boundary conditions) <cit.> and in DEM simulations of 3D ellipsoidal-particle packings subject to triaxial compression <cit.>. 



Using similar DEM simulations for biaxial compression as those described in the previous paragraph, <cit.> examined the evolution of force-weighted networks with axial strain using several of the network concepts that we discussed in Sec. <ref>. Unsurprisingly, they found that including contact forces in their network analysis yields a more complete characterization of a granular system than ignoring them.

One measure that they felt was particularly useful is a weighted version of subgraph centrality (see Sec. <ref>). From a contact network, <cit.> first extracted all conformation subgraphs. As described in Sec. <ref>, a conformation subgraph is a subgraph that consists of a given particle i and that particle's immediate contacts. (Each particle in a network thus yields one conformation subgraph.) To incorporate inter-particle force information, <cit.> generated force-weighted conformation subgraphs by weighting each edge in the conformation subgraphs by the magnitude of the normal force component along that contact. They then computed a weighted subgraph centrality Y^w_i for each force-weighted conformation subgraph, and they computed changes |Ỹ^w| in magnitude between consecutive strain steps of a suitably averaged version (Ỹ^w) of this quantity (see <cit.> for details). They observed that the temporal evolution of |Ỹ^w| with strain step is effective at tracking large changes in energy dissipation rate that can occur due to rearrangement events (e.g., force-chain buckling) associated with the loss of inter-particle contacts. They also observed that the central particle in the conformation subgraphs that undergo the largest changes in weighted subgraph centrality seems to be associated with locations with much dissipation, such as in the shear band and in buckling force chains or the neighboring particles of those force chains.

See <cit.> for the employed definition of force chains and <cit.> for the employed specification of force-chain buckling. <cit.> highlighted that network analysis — and especially examination of mesoscale features — can be helpful for gaining insights into mechanisms that regulate deformation in granular materials. Such studies can perhaps also help guide efforts in thermo-mechanical constitutive modeling <cit.>.



  §.§.§ Extracting multiscale architectures from a force network using community detection.



A major benefit of studying a network representation of a granular system (and using associated computational tools) is that it provides a natural framework in which to probe structure and dynamics across several spatial scales. One can examine different spatial scales in multiple ways, including both physically (e.g., using distance in Euclidean space or in some other metric space) or topologically (e.g., using the hop distance along edges in a network). In these studies, one can use network diagnostics and approaches like the ones discussed in Sec. <ref>. The ability to successfully study mesoscale architecture, which manifests at intermediate scales between particle-level and system-wide organization, is an especially important contribution of network analysis. One of the most common ways to examine mesoscale structures is with community detection (see Sec. <ref>), which one can use to extract sets of densely-connected nodes in a network <cit.>. One can also tune community-detection methods to examine sets of nodes of a range of sizes, from very small sets (with a lower limit of one node per set) to large sets (with an upper limit of all nodes in a single set).

By applying multiscale community-detection methods to force-weighted contact networks of photoelastic disks, <cit.> identified chain-like structures that are visually reminiscent of force chains in granular packings. Notably, the algorithmic extraction of these “distributed” mesoscale structures <cit.>, in contrast to “compact”, densely-connected geographical domains in a material <cit.>, required the development of a geographical null model, which can be used in modularity maximization and which encodes the fact that a given particle can exert force only on other particles with which it is in direct contact <cit.> (see Sec. <ref>). The different type of mesoscale organization extracted by this geographical null model highlights the fact that using physically motivated network-based approaches, which incorporate spatial and/or physical constraints, may give different information about a granular system than what one obtains when doing network calculations based only on network structure (i.e., without considering known context about a network, so that “one size fits all”). In a modularity-maximization approach to community detection, one can also tune a resolution parameter of a modularity objective function to identify and characterize network communities of different sizes. This can also be combined with inference procedures to determine particularly important scales.



One interesting result from <cit.> is that properties of force chain-like communities can distinguish frictional, laboratory packings from frictionless, simulated ones, allowing a quantification of structural differences between these systems. In later work, <cit.> used similar techniques to examine friction-dependence and pressure-dependence of community structure in 3D simulations of compressed granular materials. To further quantify such mesoscale organization and examine how it changes with compression, <cit.> extracted communities using the geographical null model, and used ideas from algebraic topology to define a topological compactness factor that quantifies the amount of branching — versus compact, densely-interconnected regions — in communities of a force network from 2D granular systems. The approach from <cit.> was extended to multilayer networks (see Sec. <ref>) in <cit.>, providing a way to link particulate communities across compression steps (rather than extracting new ones at each step) when examining how such communities reconfigure. These studies helped lay groundwork to improve understanding of how the multiscale nature of force-chain architecture impacts bulk material properties. 
Various community-detection approaches have also been used for identifying other types of inhomogeneities in granular matter <cit.>.

Before moving on, it is important to note that — although related — the definition of force-chain structure using the community-detection approaches that we described above <cit.> differs from the definitions of force chains that have been used in some other studies (e.g., see <cit.>). In future work, it is important to examine how the properties of force chains differ when they are defined in different ways. 



  §.§.§ Some applications.


Comminution processes. A network-based approach can give fascinating insights into comminution, the fragmentation of a material into smaller pieces. <cit.> used DEM simulations to study comminution in a granular material under uniaxial compression and reported that the degree distribution of the system's contact network (which we recall is unweighted) evolves towards a power law during this process. This is consistent with the development of a power-law grain-size distribution, in which large particles are hubs that have many smaller, neighboring particles, which make up the majority of a packing. <cit.> also examined several other features (such as measures of network efficiency, node betweenness, and cycle populations) of both contact networks and networks weighted by the normal force between particles as a function of increasing strain to examine what changes occur in a granular system during comminution. 

 

Heat transfer. Another problem that has been examined using network-based methods is heat transfer in granular matter. Using a heat-transport model on simulations of a compressed granular material, <cit.> probed the effects of heterogeneity in the force distribution and the spatial arrangements of forces in a system on heat transfer through the material. Specifically, they compared measures of transport in the (normal) force-weighted network of a granular system to two null-model networks with the same contact topology but with either (1) homogenous, uniform edge weights that were equal to the mean force of the packing; or (2) the same set of heterogeneous edge weights from the actual granular network, but assigned uniformly at random to the edges. <cit.> estimated the thermal diffusivity and effective conductivity from simulations on each network, and they observed that the real granular system has significantly higher diffusivity and effective conductivity than the homogenous null model. Additionally, comparing the results from the real material to the null model with randomly reassigned edge weights demonstrated that the qualitative differences between the real granular network and the homogenous null model could not be explained by the heterogeneity in the force distribution alone, as the authors observed that this second null model (with randomly reassigned edge weights) was also not a good medium for heat transfer. 

To investigate what features of a granular network facilitate efficient heat transfer, <cit.> defined a weighted network distance (see Sec. <ref>) between particles i and j as d^w_ij = 1/H_ij, where H_ij is the local heat-transfer coefficient, such that the network distance between two particles in contact is proportional to that contact's resistance to heat transfer. Note that H_ij∝ f_ij^ν, where f_ij is the magnitude of the normal force between i and j, and ν≥ 0 is a constant. They then defined a network-based (but physically-motivated) measure of heat transport efficiency as the weighted efficiency E^w (see Sec. <ref>) computed using the distances d^w_ij. In a comparison between the real granular system and the two null models, E^w gave the same quantitative results as the effective conductivity. In particular, the calculations in <cit.> revealed that the real granular system has a larger efficiency than that of either null model, suggesting that the spatial distribution of force-chain structure in the granular network appears to facilitate heat transport. Finally, iterative edge removals in decreasing order of geodesic edge betweenness centrality [Eq. (<ref>)] yield a faster decrease in effective conductivity than either edge removals done uniformly at random or edge removals in decreasing order of the local heat-transport coefficient, further illustrating the utility of network-theoretic measures for examining transport phenomena in granular systems.

 

Acoustic transmission. One can also examine the effect of network structure on properties such as electrical conductivity in systems composed of metallic particles or on other types of transport (such as sound propagation) through a particulate material. The transmission of acoustic signals through granular materials is poorly understood <cit.>, and it is particularly challenging to explain using continuum or particulate models <cit.>. A few years ago, <cit.> represented compressed, 2D packings of bidisperse, photoelastic disks as force-weighted contact networks, and found that some network diagnostics are able to identify injection versus scattering phases of acoustic signals transmitted through a granular material. Among the diagnostics that they computed, the authors observed that network efficiency (see Eq. (<ref>) in Sec. <ref>) is correlated positively with acoustic transmission during the signal-injection phase, suggesting that high-amplitude and high-energy signals are transmitted preferentially along short paths (and perhaps even shortest paths) through a force-weighted contact network. In contrast, low-amplitude and low-energy signals that reverberate through a packing during the subsequent signal-scattering phase correlate positively with the intra-community strength z-score, which characterizes how strongly a node connects to other nodes in its own community. These results suggest that one can use network diagnostics that probe diverse spatial scales in a system to describe different bulk properties. Because <cit.> did not use community-detection approaches informed by a geographical null model (see Sec. <ref>), it did not address (and it is not yet fully understood) how acoustic transmission depends on the multiscale architecture of chain-like structures reminiscent of force chains. This remains an open issue, and network-based approaches — e.g., using geographical null models and other ideas that pay attention to the role of space — are likely to be important in future work on this topic.



  §.§.§ Thresholded force networks.



Before researchers started using network-based methods as a common perspective for studying granular materials, <cit.> reported that under stress, the force network in simulations of granular matter organizes into two subsets: one set with “strong" contacts and another set with “weak" contacts. The strong subnetwork of forces forms a backbone of chain-like structures that carry most of a system's load and which tend to align approximately with the direction of compression. Between these strong force chains, there is a weak subnetwork of contacts that carry forces that are less than the mean. This weak subnetwork tends to have an anisotropy that is orthogonal to the compression, and it may provide support to the backbone of strong forces. Such heterogeneity in a force network is an interesting feature of granular materials, and network-based approaches provide a direct way to examine how strong and weak contacts can play important roles in material properties and stability.




Theses ideas have been explored using force-thresholded networks <cit.>, in which one retains only contacts that carry a force of at least some threshold f_th. That is,

    A^th_ij = {[ 1 ,  if particles i and j are in contact and  f_ij≥ f_th ,;                                            0 ,  otherwise. ].


The threshold f_th should be lower-bounded by the smallest non-zero contact force f_min in the system being considered. (Note that when f_th = f_min, one includes all contacts in the force-thresholded network.) It is common to use a system's mean force ⟨ f ⟩ as a reference point and to systematically vary f_th to be different fractions of the mean force. Whether one uses the normal or the tangential component of the force, varying f_th results in a series of thresholded networks that one can subsequently characterize using traditional network or percolation-based analyses (which we discuss in the following three paragraphs of this subsubsection) or using methods from computational algebraic topology (which we discuss in Sec. <ref>). See Fig. <ref> for an example depicting the contacts and forces in a simulation of a 2D granular material at two different values of the force threshold. There is also a connection between the idea of force-thresholded networks and carrying out modularity maximization with the geographical null model [Eq. (<ref>)] and a resolution parameter γ (see Sec. <ref>). Similar to how increasing the threshold f_th selects the subset of particles in a network for which inter-particle forces exceed or equal the threshold, in modularity maximization, increasing the resolution-parameter value γ tends to yield communities of particles such that within a community, the inter-particle forces are at least as large as the threshold γ⟨ f ⟩.

<cit.> examined several network diagnostics — e.g., mean degree, shortest-path length, diameter, LCC size, and component-size distributions — as a function of the force threshold f_th (and at different values of the packing fraction) in DEM simulations of 2D granular packings under isotropic compression. The computations in <cit.> suggest that the way in which many of these measures change as a function of f_th depends on packing fraction, and many of the measures can thus potentially be used to help understand differences in the organization of force networks above versus below the jamming point. For example, for packing fractions above the jamming point, the LCC size and the shortest-path length undergo qualitative changes in behavior near a threshold f_th≈⟨ f ⟩, signifying interesting structural changes in the organization of the force networks at that threshold. The relationship between the number of 3-particle contact cycles (i.e., triangles) and the force threshold was examined in <cit.>. In the jammed state, they observed a steep decline in the number of 3-cycles in the networks as they increased the threshold and considered only progressively larger forces. (As we discuss in Sec. <ref>, one can also use methods from computational algebraic topology to examine the evolution of cycle organization in dense granular materials.) This observation suggests that triangles (or at least one of their contacts) belong primarily to the weak subnetwork of forces that help support strong, filamentary force-chain structures. See Fig. <ref> and Secs. <ref> and  <ref> for other discussions of the roles of cycles and their relationship to force chains.

Another way to study force-thresholded granular networks is using a percolation-like approach (see Sec. <ref>). For example, one can examine the sizes and number of connected components in a thresholded network as a function of f_th <cit.>. For dense packings, the intuition is that when f_th is very large, a force-thresholded granular network splits up into several disconnected components; and as one decreases f_th, these components begin to merge until eventually all contacts are in a single component. In this type of bond percolation, which edges are included in the network thus depends on the force that the edges carry. One can look for a critical threshold f^c_th such that for f_th > f_th^c, the network fragments into many small components, but as one lowers the threshold towards f^c_th, a large, percolating cluster forms in the system (see Sec. <ref>). Quantities that are often investigated when studying this type of force-percolation transition include f^c_th and “critical exponents” for the transition (see Sec. <ref>). Because (both experimental and computational) granular systems are finite in practice, such investigations often use finite-size scaling techniques.

Several simulation-based studies of granular systems have deployed percolation analyses based on force-thresholded networks to quantify the organization of granular force networks, how such organization changes with increasing compression, and other phenomena <cit.>. For example, <cit.> studied force percolation in the q-model <cit.> of anisotropic granular force networks, in which there is a preferred direction of force propagation. They concluded that the asymmetry in the model has a significant effect on the percolation transition, and they found that the critical exponents differ from those of isotropically compressed granular force networks. <cit.> investigated force percolation in a variety of simulations of slowly compressed, 2D granular systems. They examined the effects of polydispersity and friction, finding that these factors can qualitatively influence various features of the percolation transition. Very recently, <cit.> also investigated the force-percolation transition in simulations of jammed granular packings at fixed pressures.



  §.§.§ Methods from computational algebraic topology.
 


In addition to traditional approaches for network analysis, one can also study the architecture of granular networks using ideas from algebraic topology. Persistent homology (PH) <cit.> (see Sec. <ref>) seems especially appropriate, and over the past several years, it has provided a fruitful perspective on the structure of compressed <cit.>, and tapped <cit.> granular materials. Very recently, it has also been used to study responses of granular materials to impact <cit.>. One way to characterize the organization and evolution of granular force networks is to examine how Betti numbers (see Sec. <ref>)) change as a function of a force threshold (and also as a function of packing fraction) in compressed granular systems. This includes studying the birth and death of components (determined by β_0) and loops (determined by β_1) as a function of a force threshold by computing and analyzing persistence diagrams (see Fig. <ref> of Sec. <ref>). Examining when and how long different features persist in a network provides a detailed characterization of the structure of granular force networks, and one can quantify differences between two networks by defining measures of “distance" between their associated persistence diagrams. These capabilities allow a PH framework to provide a distinct set of contributions that complement more traditional network-based analyses of components, cycle structure, and other features.

<cit.> investigated how simulated 2D granular force networks evolve under slow compression as they cross the jamming point. They first demonstrated that one can identify the jamming transition by a significant change in behavior of β_0 (specifically, there is an increase in the number of components at a force threshold approximately equal to the mean force ⟨ f ⟩), and that structural properties of the network — such as the size of the connected components — continue to change above jamming. <cit.> also demonstrated that β_0 and β_1 can quantitatively describe the effects of friction and polydispersity on the organization of force networks (and can distinguish how friction and polydispersity alter the structure of a force network). This work was extended in <cit.>, who examined numerical simulations of 2D, slowly compressed, dense granular materials using PH. In addition to examining the values of the Betti numbers, they also computed β_0 and β_1 persistence diagrams (PD_0 and PD_1, respectively) as the system was compressed to quantify the appearance, disappearance, and lifetimes of components and loops. In <cit.>, they defined a filtration over the clique complex of the networks (see Sec. <ref>), so only loops with four or more particles were counted. To extract useful information from the PDs, they binned the persistence points in each diagram into different regions corresponding to features that (1) are born at any force threshold but have relatively short lifetimes compared to those that are born at either (2) strong, (3) medium, or (4) weak forces and that persist for a large range of thresholds. Their persistence analysis led to several insights into the structure of a normal-force network as a granular system is compressed, as well as insights into differences in the structure of a normal-force network for systems with different amounts of friction and different polydispersities. For example, <cit.> observed that, near the jamming point, frictionless packings appear to have more “extreme” features than frictional packings, in the sense that frictionless packings have many more β_0 persistence points that are born at either weak or strong forces and that are relatively long-lived. Observations on the effects of polydispersity and friction may be difficult to observe using traditional measures such as the probability density function of the normal forces. 



<cit.> also used PH to study force networks from simulations of slowly compressed, polydisperse packings of disks in 2D as they traverse the jamming transition (see Fig. <ref>a). As they compressed the system 
through a range of different packing fractions ρ∈ [0.63,0.9], they extracted force information at approximately fixed time intervals during the simulation, and they then computed PDs of components and loops (i.e., PD_0 and PD_1, respectively) for each force network sampled during compression (see Fig. <ref>b). As in other studies, <cit.> used the clique complex of the force networks to avoid counting 3-particle loops. They then used the bottleneck distance d_B and two variants, d_W1 and d_W2, of the Wasserstein distance <cit.> to quantify differences between two PDs and thereby help quantify differences in local and global features between two granular force networks. As described in <cit.>, the bottleneck distance captures only the largest difference between two PDs, whereas the Wasserstein distances include all differences between two PDs, with d_W1 more sensitive than d_W2 to small changes. (See <cit.> for more details.) Calculating the various types of distance between the two β_0 PDs and the two β_1 PDs for consecutive samples (i.e., consecutive states) of a force network allows one to characterize different kinds of variations in the geometry of force networks as a function of packing fraction (see Fig. <ref>c). Using the employed distance measures, <cit.> observed that in the unjammed state, there can be significant (but localized) reorganization in force geometry as a packing is compressed. They also concluded that the jamming transition is characterized by rapid and dramatic global rearrangements of a force network, followed by smoother and less dramatic reconfiguration in the system above jamming, where the distances between consecutive states of a packing are much smaller than in the unjammed state. <cit.> also found that tangential-force networks seem to exhibit similar behavior (in most respects) to that of normal-force networks. They also observed that friction can have a significant impact on how the geometry of the forces reconfigures with compression. For example, they showed that the rate of change of loop features (as measured by the distances between consecutive PD_1s) is larger for a frictional system than for a frictionless one below the jamming point, but just before jamming and thereafter (i.e., during the entire jammed state) differences in loop structure between consecutive packing fractions are larger for frictionless systems than for frictional ones. In very recent work, <cit.> used PH to examine the temporal scales on which granular force networks evolve during slow compression. They simulated dense 2D granular materials and studied the influence of the externally-imposed time scale — set by the rate of compression — on how frequently one should sample a system during compression to be able to evaluate its dynamics near the jamming transition. By varying the sampling rate and carrying out a persistence analysis to quantify the distance between consecutive sampled states of the system, their results indicate that close to jamming, a force network evolves on a time scale that is much faster than the one imposed by the external compression. See <cit.> for further details.

One can also use PH to study force networks in tapped granular systems. <cit.> examined DEM simulations of two different 2D systems exposed to tapping. One type of packing consisted of monosized disk-shaped particles, and the other type consisted of monosized pentagon-shaped particles. <cit.>'s investigation suggested that particle shape can play an important role in mechanical responses, which is consistent with observations from classical investigations of granular materials <cit.>. More specifically, <cit.> computed β_0 and β_1 as a function of force threshold in both normal-force networks and tangential-force networks. They observed for both types of force-weighted networks (but particularly for the tangential one) that the first two Betti numbers are able to clearly distinguish between disks and pentagons, where β_0 (respectively, β_1) is consistently larger (respectively, smaller) for pentagons across a wide range of force thresholds. However, using only β_0 and β_1, <cit.> were unable to clearly differentiate states with similar packing fractions but that result from different tap intensities. 

In a follow-up investigation, <cit.> simulated a series of several taps to granular packings and used PH to examine how normal and tangential force-weighted networks vary between individual tap realizations. Specifically, they computed distances between PD_0s and between PD_1s of force networks associated either with individual realizations of tapping to the same system or with individual realizations of tapping to two different systems. In one part of their study, they examined systems of disks exposed to a series of taps at two different tap intensities. (See Sec. <ref> for a rough delineation of “low" tapping intensity versus “high" tapping intensity.) They observed that in terms of loop structure, the set of networks generated from a series of taps at low intensity differ far more substantially from each other — as quantified by the distribution of distances between the PD_1s for different realizations of the tapping — than do the set of force networks from a series of taps at high intensity. They also observed that the distances between different realizations of low-intensity tapping are as large as the distances between low-intensity tapping and high-intensity tapping realizations. Therefore, although the high-intensity tapping and low-intensity tapping regimes yield networks with approximately the same packing fraction (see Sec. <ref>), one can use methods from PH to help explain some of the differences between the packing structure in the two regimes. In another part of their study, <cit.> carried out a persistence analysis of tapped packings of disks and tapped packings of pentagons, and they observed clear distinctions between the two systems based on calculations of β_1 PDs. For example, for each system, they computed the PD_1s for a set of networks associated with several individual realizations of the same tapping intensity, and they then computed a distance between each pair of PD_1s for realizations within the same packing and across the two types of packings. They observed that the distribution of distances between the PD_1s of individual tapping realizations to the packing of pentagons is narrower and centered at a smaller value than the distribution of distances between individual realizations of taps for the packing of disks. They also observed that the distances between the disk and pentagon systems are much larger than those between different realizations of the disk system. Thus, <cit.> were able to distinguish clearly between tapped disk packings and tapped pentagon packings using PH, especially when considering properties of loop structures. Past work using 2D experiments has also been able to distinguish between the dynamics of disk and pentagon packings using conventional approaches <cit.>.

One can also use methods from computational algebraic topology to study granular networks in which one uses edge weights from something other than a force. For example, <cit.> used only particle positions (in the form of point clouds) and computed Betti numbers to distinguish states at the same density but that are at different mechanical equilibria. Using both experimental and simulated 2D granular packings of monodisperse particles, they constructed networks by locating the center of each particle and then introducing a filtration parameter δ, such that any two particles separated by a Euclidean distance less than or equal to δ are adjacent to each other in a graph. They considered δ∈ [d,1.12d], where d is the particle diameter and the domain for δ resembles the choices that are used for determining if particles are in physical contact with each other. The authors computed, as a function of δ, the first Betti number β_1 on the whole network to count the total number of loops at a given δ. They also computed β_1 on the flag
complex (see Sec. <ref>), thus counting the number of loops with four or more nodes at a given value of δ. For values of δ that are slightly larger than d, <cit.> were able to separate states at the same packing fraction that are generated by tapping at different intensities. They observed for a fixed packing fraction that states that arise from lower-intensity tapping have a larger value of β_1 when computed on the whole network and a lower value of β_1 when computed on the flag complex. Their results were robust to both noise and errors in particle-position data.



 §.§ Other network representations and approaches





  §.§.§ Network-flow models of force transmission.



Another technique for gaining insight into the organization of forces in deforming granular systems, and how microscale aspects of a force network lead to macroscale phenomena such as shear bands and material failure, is to view force transmission from the perspective of maximum-flow–minimum-cut and maximum-flow–minimum-cost problems (see Sec. <ref>). To examine a granular system using such a perspective, one can consider the “flow" of force through a contact network (with some contacts able to transmit more force than others), which in turn yields a “cost” to the system in terms of energy dissipation at the transmitting contacts <cit.>. One can calculate flow and costs in routes through a network (and hence determine bottlenecks in force transmission) to gain understanding of how contact structure relates to and constrains a system's ability to transmit forces in a material. For example, <cit.> constructed flow networks from DEM simulations of a system of polydisperse particles compressed quasistatically under a constant strain rate in the vertical direction and allowed to expand under constant confining pressure in the horizontal direction. At a given axial-strain value (i.e., “strain state”), they assigned uniform capacities u_ij to each edge of the contact network to reflect the maximum flow that can be transmitted through each contact, and they assigned costs c_ij to each edge to model dissipation of energy at each contact. After each axial-strain increment during loading, they then solved the maximum-flow–minimum-cost problem for the network at the given strain state, finding that edges in the minimum cut (yielding bottlenecks in the force-transmission networks) localize in the material's shear band. By using costs c_ij that reflected the type of inter-particle contact (specifically, elastic contacts versus various types of plastic contacts) <cit.> were able to track different stages of deformation (i.e., strain-hardening, strain-softening, and the critical-state regime). They also computed a minimal cycle basis and observed that a large majority of force-chain particles and particles in 3-cycles are involved in the set of contacts that comprise the maximum-flow–minimum-cost routes.

One can use the above approach with various definitions of force capacity and cost functions. Using simulations of the same type of system as that in the previous paragraph, <cit.> constructed networks — one for each strain state as a system was loaded until it failed — that incorporated information about both the inter-particle contacts at a given strain state and the particle displacements that occur between consecutive strain steps. Specifically, if nodes i and j are in contact at a given strain state, the weight of the edge between them is the inverse of the absolute value of the magnitude of the relative displacement vector of particles i and j, where one computes the displacement for each particle from the previous strain state to the current one. The distance between nodes i and j is 0 if the associated particles are not in contact. The intuition behind this capacity function is that, when there is more relative motion between a pair of particles, one expects those particles to have less capacity to transmit force to each other. <cit.> used capacities that incorporate 3-cycle memberships of edges in a study of minimum cuts of a flow network in two samples of 3D sand under triaxial compression and in a 3D DEM simulation of simple shear. Grains in the bottlenecks localize early during loading, and are indicative of subsequent shear-band formation. Other work <cit.> studied DEM simulations of compressed, 3D bonded granular materials (where bonded signifies that the grains are connected via solid bonds of some strength) and a system of 2D photoelastic disks under shear stress with the goal of testing the hypothesis that an appropriate maximum-flow–minimum-cost approach can identify experimentally-determined load-bearing particles and force-chain particles without relying on knowledge of contact forces. <cit.> examined different combinations of force-transmission capacity and cost functions, and they examined the fraction of force-chain particles that are part of the associated maximum-flow–minimum-cost network for a given capacity and cost function. In both cases, costs based on 3-cycle membership of edges seem to yield large values of these fractions, and <cit.> were able to successfully forecast most of the particles that eventually become part of force chains without using information about contact forces.

 


  §.§.§ Broken-link networks.



In previous discussions (see Secs. <ref>, <ref>), we have seen that one way to investigate the evolution of a granular system under an applied load is (1) to compute contact networks or force-weighted networks for the system as a function of packing fraction, strain, or some other control parameter, and then (2) to study how different features and properties of the networks emerge and change as one varies that parameter. One can also use other network constructions to explore different mesoscale features and examine system dynamics. For example, <cit.> designed a broken-link network (see Fig. <ref>) to study the dynamics of 3D granular flows. They conducted an experiment on a collection of acrylic beads immersed in a box of liquid medium, shearing the system at a constant rate (Ω≈ 1.05 × 10^-3 rad/s) by a rotating circular disk at the bottom of the box. (See <cit.> for details about their experiments.) First, they constructed proximity networks (a variant of a contact network) as a function of time. (Specifically, they constructed one network every 3-degree increment of rotation.) In their proximity network, they assigned a “contact" (edge) to each pair of particles whose distance from one another in the given frame was within a specified distance threshold, which they chose to be a conservative upper bound for particle contact. They then defined a broken link (relative to some reference frame) as an existing edge between two particles in the reference frame that was subsequently absent in at least two later time frames (due to the particles having moved apart). A broken-link network for the frame in which a pair of particles moved apart had an edge between the two particles that were initially in contact, and broken links were not allowed to reform later. In Fig. <ref>a, we illustrate this procedure for constructing a broken-link network. 



Studying the temporal evolution of a broken-link network provides a quantitative approach for examining particle rearrangement events in granular matter. To examine the temporal evolution of a granular system, <cit.> examined the size of the LCC in a sequence of broken-link networks as a function of applied shear, drawing an analogy between the fraction χ_b of broken links and the occupation probability in traditional percolation problems (see Sec. <ref>).  
They observed that the fraction s_g of nodes in the LCC of the broken-link network grows with χ_b in a way that suggests that there is a continuous phase transition in s_g, and they approximated the value of χ_b at which this transition occurs. (However, as we noted in Sec. <ref>, because these networks have finite sizes, one needs to be cautious regarding statements about percolation and phase transitions.) 
From a physical standpoint, this transition region corresponds to a characteristic deformation at which broken links — which are due to particle rearrangements — start to become globally connected, relative to a reference proximity network. By examining χ_b as they applied shear to the system, <cit.> approximated a characteristic amount of strain associated with this transition region (by mapping the value of χ_b associated with the transition region to a corresponding value of strain), and they suggested that the determined strain scale may be useful for identifying the onset of global reorganization in the system. 

In a later study, <cit.> used a similar approach to examine progression from reversible to irreversible dynamics in granular suspensions under oscillatory shear strain. In experiments similar to the one described above <cit.>, the authors considered a series of 20 shear cycles of the following form: for each cycle in one experiment, a suspension was sheared with a rotating disk up to a strain amplitude of θ_r, after which the direction of rotation was reversed and the disk was rotated by -θ_r back to its original location. The authors performed experiments for θ_r values of 2^∘, 4^∘, 10^∘, 20^∘, and 40^∘. In a qualitative sense, measuring “reversibility" of the dynamics in this setup entails determining the extent to which particles return to their original positions after some number of shear cycles. To quantify this idea, one can compute the mean square displacement (MSD) of the particles after each cycle, where smaller MSD values correspond to more reversible dynamics. To study this system, <cit.> adjusted the idea of a broken-link network to include “healing", so that broken links that are repaired in later frames do not contribute edges in the broken-link network, and they studied the temporal evolution of this broken-link network as a function of the cyclic shear and for different amplitudes θ_r. In addition to the extent of spatial (ir)reversibility measured by calculating the MSD, <cit.> also proposed a notion of topological (ir)reversibility by examining the temporal evolution of the size of the LCC in their broken-link networks. For low values of θ_r (specifically, for θ_r ≤ 20^∘), the system appears to be almost reversible: proximity-based contacts break, but they reform after shear reversal, and the fraction of particles in the LCC of the broken-link network thus grows before subsequently shrinking to almost 0 after reversal. However, for a higher shearing amplitude (specifically, for θ_r = 40^∘), the system shows signatures of irreversibility. Many broken links do not reform after a shear cycle, and after reversal, the LCC of the broken-link network remains at a value that constitutes a substantial fraction of the total system size.



  §.§.§ Constructing networks from time series of node properties or from kinematic data.



Another way to examine the organization of deforming granular materials is to construct networks based on the temporal evolution of particle properties <cit.>, an idea that draws from earlier work in complex systems on constructing networks from different kinds of time-series data. (See, for example, <cit.>.) In one type of construction, which yields what are sometimes called functional networks <cit.>, one records the time series of some property of each particle and places an edge (which can potentially be weighted and/or directed) between two particles according to some relationship (e.g., some type of distance or other measure of similarity) between the particle-property time series. When using a different property, the nodes are the same, but the edges and especially edge weights will in general be different. Once networks have been constructed, one can examine them using techniques such as those in Sec. <ref>. Some authors have also used the generic term particle-property networks to describe networks that they constructed in some way from particle properties.

<cit.> used particle-property networks to study the evolution of DEM simulations of a quasistatically deforming, 2D granular material under biaxial compression. They constructed time series for two features (as well as their coevolution) for each particle — membership in force chains (determined as in <cit.>) and membership in 3-cycles of a minimal cycle basis — by recording a 1 if a particle has the given property and a 0 if it does not. They then quantified the similarity of particle-property evolution using the Hamming distance <cit.>, and they added (undirected and unweighted) edges between each particle and its k closest (i.e., most similar) particles until eventually obtaining a single network with one connected component. <cit.> then extracted sets of particles that exhibit similar dynamic behavior by detecting communities (see Sec. <ref>) in the particle-property networks. This uncovered distinct regions in the material — including the shear band and different subnetworks composed of primarily force chain or primarily non-force chain particles — as well as interlaced regions in the shear band that alternate between jammed and unjammed configurations. See <cit.> for additional studies that used membership in cycles of length up to l = 7 for the construction of particle-property networks. 

One can also construct particle-property networks from data that do not rely on knowledge of contact or force networks. For example, <cit.> studied deformation in sand subject to plane-strain compression using measurements of grain-scale displacements from digital image correlations (DICs) <cit.>. From these data, they constructed kinematic networks, a type of network that arises from some measurement of motion (such as displacement or a rotation) over a small time increment. <cit.> considered each observation grid point of digital images of a sample to be a node, and they placed edges between nodes with similar displacement vectors during a small axial-strain increment. This yields a collection of (undirected and unweighted) time-ordered kinematic networks. In another study, <cit.> calculated particle rotations and displacements for triaxial compression tests on sand using x-ray micro-tomography scanning <cit.>. They generated an ordered set of networks from these data for several strain steps by treating each grain as a node and linking nodes with similar kinematics during the specified interval. More specifically, they represented the displacements and rotations of each particle as points in a state space, and they connected particles that are nearby in that state space according to Euclidean distance. In the network that they constructed, each particle is adjacent to k nearest neighbors in state space, where k is as small as possible so that the (unweighted and undirected) network is connected. Various notions of what it means to be “similar”, and thus how to quantitatively define edges, are discussed in <cit.>. To probe the collective dynamics of interacting groups of particles for the network in each strain step, <cit.> detected communities corresponding to mesoscale regions in the material that exhibit similar dynamic behavior. Calculating the mean shortest-path length between pairs of particles in the same community yields a potentially important intermediate spatial scale (that they concluded is consistent with the shear-band diameter) of a granular system. For each strain step, they computed a variant of closeness centrality — it is similar to the one in Sec. <ref>, but it is not exactly the same — of each particle in the corresponding network, and observed that particles with large closeness centrality localize in the region of the shear band early in loading (and, in particular, before the shear band develops). This study highlights the potential of network analysis to provide early warning to detect regions of failure in particulate systems. 

Methods from nonlinear time-series analysis have also been used in network-based studies of stick-slip dynamics in a granular packing sheared by a slider <cit.>. Using so-called phase-space networks (see <cit.> for a description of them) to construct networks from measurements of a slider time series, <cit.> associated network communities with slip events. 




 §.§ Comparing and contrasting different network representations and approaches


Network representations of granular and other particulate systems, in combination with methods from network science and related disciplines, provide a plethora of ways to analyze granular materials. For these tools to be optimally useful in furthering understanding of the physics of these systems, it is also important to draw connections between the many existing approaches. The application of network-based frameworks to the study of granular matter is relatively young, and little work thus far has focused specifically on exploring such connections, but it is crucial to determine (1) how conclusions from different approaches relate to one another and (2) how similarities and differences in these conclusions depend on the system being studied. In this short subsection, we point out a few relationships and places to begin thinking about such questions.

First, it is important to consider the network representation itself. We have discussed several different representations in this review, and some are more related to each other than others. Broadly speaking, one class of granular networks tries to encode the physical structure of a material — in the sense that edges exist only when there is an inter-particle contact — at a given point in an experiment or simulation. Such networks include contact networks (see Sec. <ref>) and force-weighted networks (see Sec. <ref>). One can obtain a contact network from an associated force-weighted network by discarding the edge weights and keeping information only about connectivity. Force-weighted networks thus contain much more information, and they allow one to investigate phenomena — such as force-chain organization — that may not be possible to probe with the contact network alone. However, it is still important to develop tools for the analysis of contact networks and understand what phenomena arises from features of the connectivity alone (or what one can learn and explain from connectivity alone), as force information may not always be available or may not be necessary to understand certain behaviors of a granular system. Generalizing some quantities (e.g., betweenness centralities) from unweighted networks to weighted networks also involves choices, and it is often desirable to conduct investigations that have as few of these potentially confounding factors as possible. Other types of granular networks do not encode physical connectivity, but instead directly represent something about the dynamics or changes that occur in a system during an experiment or simulation. Examples of such networks include broken-link networks (see Sec. <ref>) and particle-property networks (see Sec. <ref>). These different classes of network representations offer distinct ways of studying granular materials, and utilizing each of them should improve understanding of these systems. 

For a given network representation, it is reasonable to expect that conclusions that arise from similar network quantities or methods of analysis are related to one another, whereas conclusions that result from tools designed to probe very different kinds of organization in a network provide rather different information about the underlying granular system <cit.>. For example, some studies have suggested that in deforming granular materials, results based on calculations of clustering coefficients are similar to those from studying 3-cycles. This is intuitively reasonable, given that calculating a local clustering coefficient yields one type of local 3-cycle density. We have also observed that conclusions drawn from examinations of small subgraphs in a deforming granular system may also be related to whether or not those subgraphs contain cycles of certain lengths. As we discussed in Sec. <ref>, another way to draw connections between different approaches is to consider the spatial, topological, or other scales that are probed by the different approaches. For instance, in a force-weighted network, node strength is a particle-scale property, and it encompasses only very local information about a granular material. However, granular systems exhibit collective organization and dynamics on several larger scales that may be difficult to understand by exclusively computing local measures and distributions of such measures. To obtain an understanding of larger-scale structures, it is necessary to also employ different methods, such as community detection, persistent homology, and the examination of conformation subgraphs that are composed of more than just a single particle. Utilizing such approaches has provided insights into force-chain structure, shear band formation, and reconfiguration in granular systems under load that one may not be able to obtain by considering only local network quantities. It is thus important to continue to use and develop methods to analyze granular networks across multiple scales, as doing so can provide important and new information about a system. Finally, we note that even a single approach, depending on how it is used, can provide multiple types of information about a granular network. A good illustration of this is community detection. In Sec. <ref>, for example, we saw that using different null models in modularity maximization allows one to probe rather different types of mesoscale architecture in granular force networks.

We look forward to forthcoming studies that directly compare the results and assumptions in different approaches (both network-based and traditional ones) and different network representations of granular and other particulate systems. Conducting principled investigations into how conclusions from various network-based approaches are related is indeed an important direction for future work.



 §.§ Limitations and practicalities of simulations and experiments
 


Network-based studies of granular materials have examined inter-particle contact and force data (and associated dynamics, such as in the presence of external loading) from both experiments and simulations. A recent summary of the many experimental techniques available for obtaining data about inter-particle contacts is available in the focus issue <cit.>. Such techniques include laser-sheet scanning <cit.>, photoelasticity <cit.>, x-ray tomography <cit.>, and nuclear magnetic resonance <cit.>. Using each of the first three approaches, it is possible to measure both particle positions and inter-particle forces. If one is careful, it is sometimes possible to measure the forces as vectors (i.e., including both the normal and tangential components), but some techniques or systems do not have sufficient resolution to allow more than scalar or coarse-grained values. Determining the forces also helps experimentalists to confidently construct contact (i.e., unweighted) networks of particulate materials. In deciding whether or not two particles are actually in contact, rather than merely being adjacent in space, it is necessary to perform a detailed study of the effects of thresholding the data <cit.>. Any experimental technique will imperfectly report the presence versus absence of the weakest contacts in a system. Additionally, because of the difficulty of accessing the interior of granular materials, much more data is available for 2D force networks than for 3D force networks <cit.>.

The most widely-used simulation techniques are discrete element methods (DEMs) <cit.>, in which the dynamics of individual particles (usually spheres) are determined by their pairwise interactions under Newton's laws of motion. The normal forces are typically determined from a Hertzian-like contact law (see, e.g., the sidebar in <cit.> for an introduction to Hertzian contacts) via an energy penalty for the overlap of two particles. The tangential (frictional) forces are most commonly modeled using the Cundall–Strack method <cit.> of a spring and a dashpot, but they have also been modeled via the surface roughness created by a connected set of smaller spheres <cit.>. For a given application, it is not known whether these simplified models capture all of the salient features of inter-grain contacts, and the situation likely differs for different applications. For example, experimental measurements of sound propagation in photoelastic disks <cit.> suggest that the amplitude of sounds waves may be largest along a force chain network, an effect not observed in DEM simulations <cit.>. This is likely a consequence of real particles physically deforming their shape to create an increased contact area through which sound can be transmitted; existing DEM simulations do not account for this effect. Another important use of particle simulations is to provide a means to investigate the robustness of network-based analyses to various amounts of experimental error <cit.>. Simulations provide an important check on experimental uncertainties in the determination of force-weighted networks and other network representations of granular materials. Conversely, network-based approaches provide a means to compare how faithfully simulations are able to reproduce experimental results.




§ OPEN PROBLEMS AND FUTURE DIRECTIONS



We now discuss a few open problems and research directions for which we anticipate important progress in the near future. We divide our comments into three main areas: the construction of different types of networks that encode various physical relationships or other properties (see Sec. <ref>), the application of network analysis to additional types of materials (see Sec. <ref>), and the application of network-based approaches to the design of materials (see Sec. <ref>). Network tools can provide valuable insights — both explanatory and predictive — into particulate materials and their dynamics, and a lot of fascinating research is on the horizon.




 §.§ Network representations and computations



To briefly explore the potential of different approaches for constructing granular (and other particulate) networks for providing insights into the physics of granular materials (and particulate matter more generally), we discuss choices of nodes, choices of edges, edge-to-node dual networks, multilayer networks, and annotated networks.[Other ideas that are worth considering include memory networks <cit.>, adaptive networks <cit.>, and various representations of temporal networks <cit.>.] It is also worth thinking about what calculations to do once one has constructed a network representation of a particulate system, so we also briefly consider the important issue of developing physically-informed methods and diagnostics for network analysis. 




  §.§.§ Definitions of nodes and edges.
 

There are many choices — both explicit and implicit — for constructing a network <cit.>, and these choices can impact the physics that one can probe in granular networks <cit.>. Perhaps the most obvious choices lie in how one defines nodes and edges. 

In the study of granular materials, a common definition is to treat individual particles as nodes and to treat contacts as edges (often with weights from the inter-particle forces). A natural set of open questions lies in how contact network architectures depend on different features of the grains in a system. For example, there have been several recent studies on systems composed of particles that are neither spheres nor disks — including ones with U-shaped particles <cit.>, Z-shaped particles <cit.>, squares and rods <cit.>, dimers and ellipses <cit.>, and others <cit.>. It would be interesting to build network representations of these systems, examine how different grain geometries affect network organization, and investigate how that organization relates to the mechanical properties of a system <cit.>. It seems particularly important to develop an understanding of which (quantitative and qualitative) aspects of network structure depend on features of grains (such as shape, polydispersity, friction, cohesiveness, and so on <cit.>) and which are more universal.

One can also consider defining particulate networks in a variety of other ways. For example, when determining edges and edge weights, one can examine the tangential (rather than, or in addition to, the usual normal) component of the force between two grains. Such extensions may facilitate increasingly detailed investigations into a packing's organization <cit.>. It may also be useful to retain information about both the magnitude and direction of forces when defining edges. One may even wish to construct signed networks, for which edges can take either positive or negative values, thereby conveying further information about the relationship between nodes. In such studies, one can perhaps take advantage of advancements in community-detection techniques, such as by using signed null models <cit.>. Additionally, as we discussed in Sec. <ref>, particle-property networks <cit.> and networks constructed from particle-displacement information <cit.> are other informative ways to build networks for particulate systems. One can also construct edges (and determine edge weights) by incorporating information about inter-grain relationships based on similarities in particle properties such as orientation <cit.> (see Fig. <ref>), coefficient of friction <cit.>, or size <cit.>. Constructing networks whose edges are determined or weighted by inter-particle similarities may be particularly useful for achieving a better understanding of mesoscale physics in polydisperse packings, which are thought to depend on the spatial distributions of particles of different types <cit.>. A perhaps nonintuitive choice is to use a bipartite representation of a granular network, such as the approach used in <cit.>.


	
The above choices for network construction give a grain-centric view of the physics of particulate materials. One can also consider edge-to-node “duals” of such networks to provide a contact-centric perspective. In a contact-centric approach, one treats contacts between physical objects as nodes and the objects themselves as edges between contacts. (Compare this idea to the notion of a line graph <cit.> of a network G.)
A contact-centric network approach was used recently in the study of nanorod dispersions <cit.>. Contacts between rods were treated as nodes, and the effective conductance of the rod was treated as a weighted edge. The treatment of a grain or other physical object as an edge rather than as a node is also a particularly appealing way to describe networks of fibers in both human-made and natural systems. Recent examples of such fiber networks that have benefited from network analysis include collagen networks <cit.> (see Fig. <ref>), fibrin networks <cit.>, and axonal-fiber networks <cit.>.

Another way to study granular systems (especially porous ones) is to consider a network constructed from the pore space <cit.>, in which pores (i.e., empty volumes between contacting grains) are nodes and throats (i.e., flow pathways that connect pores) are edges (which can be weighted in various ways). Conveniently, there are several methods to precisely determine pores and grains <cit.>. Studying pore networks is a common way to examine flow through porous materials or to understand the responses of granular materials to external stresses, but only more recently have such networks been studied explicitly from a network-science perspective. See <cit.> for some recent examples, and see <cit.> for a study of force chains in porous media. Given how pore networks are formulated, we expect that they can also be studied using PH (see Sec. <ref>).



  §.§.§ Multilayer networks.
 

When considering different ways to construct a network, it is also natural to ask whether it is beneficial to combine two or more approaches in an integrated way to more accurately (and elaborately) represent a particulate system as a network. Ideally, doing network analysis using a more complicated representation will also lead to improved physical understanding. One way to incorporate heterogeneous node types, multiple types of relationships between nodes, and time-dependence into network representations is with multilayer networks <cit.>. Recall that <cit.> (see Sec. <ref>) studied one type of multilayer community structure in a compressed granular material. Another type of multilayer network that may be useful for the study of particulate systems are multiplex networks, in which one can encode different types of relationships between nodes in different layers of a network. For example, one can construct a multiplex network in which two particles are linked both by normal force-weighted contacts and also according to another relationship (such as the tangential force between them, their contact angle, or by a measure of similarity of one or several particle-properties, as discussed in Sec. <ref>). One can also envision using multilayer networks to study particulate systems with multiple types of particles. For example, if a system consists of particles of different shapes or sizes, one possibility is that each particle of a given shape (or size) is in a different layer, intralayer edges represent interactions between particles of the same type, and interlayer edges represent interactions between particles of different types. Another possibility is to let each layer represent a time window, perhaps with intralayer edges representing a mean interaction strength during that window. The study of multilayer networks is one of the most popular areas of network science, and we expect that they will be very illuminating for studies of particulate matter.



  §.§.§ Annotated graphs.



The network-construction techniques that we have been discussing throughout this review represent data in the form of an adjacency matrix (for an ordinary graph) or an adjacency tensor (for a multilayer network) <cit.>. However, it may be desirable to encode not only relationships between grains, but also features of the grains and/or their interactions. One option is to use annotated graphs (one can also use multilayer networks) to encode inter-node relationships, individual-node properties, and individual-edge properties <cit.>. One can use annotated graphs — also sometimes called labeled graphs or tagged graphs <cit.> — to study properties of interaction patterns (such as force-weighted contacts) that may be driven by local particle properties (e.g., size, shape, spatial position in a material, membership in cycles or force chains, and so on). Available tools for annotated graphs include clustering techniques in the form of stochastic block models (SBMs) that combine information from both connectivity and annotations to infer densely-connected sets of nodes (e.g., including when there is sparse or missing data) <cit.>.




  §.§.§ Beyond pairwise interactions.



It is also desirable to explicitly examine interactions between three or more particles, rather than restricting one's study to representations of pairwise interactions. We discussed this idea briefly in Sec. <ref> in the context of simplicial complexes <cit.>, and we note that one can also encode edges among three or more nodes by using hypergraphs <cit.>.




  §.§.§ Physically-informed network calculations.


In addition to exploring different choices of how to construct particulate networks, it is also important to consider ways to generalize the tools of network science to incorporate physical constraints <cit.>. A natural way to begin is to build null models of spatially-embedded graphs that obey simple mechanical constraints like force balance <cit.>. One can also develop network diagnostics that incorporate physical geometry, spatial embedding, and latent spatial structure of networks <cit.>. See <cit.> for an example of a kinetic approach. It should also be useful to incorporate ideas of flow and other forms of dynamics into community detection <cit.>. Additionally, it is desirable to develop methods to more explicitly characterize the relationship between a network's structure and the geometry of a space in which it is embedded <cit.> (as well as latent geometry) and to use such techniques to better understand the plurality of packings that are consistent with a single force distribution via a force-network ensemble, which is the space of allowed force configurations for a fixed contact geometry <cit.>. We also point out that a spatial embedding can induce correlations between network diagnostics <cit.>, and it is therefore important to develop a better understanding of the extent to which this occurs in networks that arise from particulate systems.



 §.§ Beyond granular materials



Although we have focused primarily on network analysis of the canonical type of granular systems, which consist of discrete, macroscopic particles that interact via contact forces, one can potentially also use network-based approaches to characterize granular materials with more complex interactions as well as soft materials more broadly <cit.>. As reviewed recently in <cit.>, these materials include colloids, glasses, polymers, and gels. In these systems (and others) the particles (or other entities) can interact via various (and often complicated) means.[One can also examine particulate networks of hard materials that admit Hamiltonian descriptions <cit.>.] For example, system components can have attractive and/or repulsive long-range interactions <cit.>, can be cohesive <cit.>, and/or can interact with one another via chemical gradients <cit.>, electric charges <cit.>, or through other molecular or mechanical processes <cit.>. In each of these cases, interaction strengths can yield edge weights, either as the sole relationship studied between entities in a graph or as one of several inter-entity relationships in a multilayer network. 

One particularly interesting avenue for future work may be to study polymer and fiber networks <cit.>, which are important in biological systems. For example, in biopolymer assemblies, cross-linking can glue filaments together into large-scale web-like structures (e.g., as in the cytoskeleton of a cell) <cit.>. Such cross-linked actin filaments are critical for cellular function and integrity <cit.>, and it is thus important to understand the structural organization of these kinds of networks, their mechanical properties, and how force transmission is regulated in them <cit.>. Indeed, there has already been some work examining network representations of gels and polymers, and employing graph-theoretic analyses to quantify the structural properties of these systems (e.g.<cit.>). One can also use network analysis to study systems at smaller spatial scales. Because traditional network-based approaches are agnostic to the physical scale of a system, off-the-shelf calculations and algorithms are directly portable to microscale, nanoscale, and even smaller-scale systems <cit.>. However, despite the technical portability of common tools, the investigation of much smaller-scale systems should benefit from the extension of classical network-based tools in ways that incorporate additional underlying physics. For example, we have described extensions of network tools to assimilate ideas and constraints from classical physics (e.g., spatial embeddedness) <cit.>. For such investigations, we expect that ideas from the study of random geometric graphs (RGGs) and their extensions will be helpful <cit.>. One can also consider employing ideas that take into account principles from quantum physics <cit.> or other areas. The study of quantum networks is a new and exciting area of network science, and there are many opportunities for progress <cit.>.




 §.§ Implications for material design



As is the case with mathematics more generally, network analysis gives a flexible approach for studying many systems, because it is agnostic to many details of their physical (or other) nature <cit.>. Such versatility supports the application of network-science techniques to both living and non-living materials to examine the architectures and dynamics of these systems and to gain insights into relationships between structure and function. The tools and approaches of network science also have the potential to inform the design of new materials. For example, it should be possible to use network theory (e.g., via the tuning of a system's network architecture) to provide guidance for how to engineer a material to exhibit specific mechanical, electrical, or other properties. Material design has become increasingly popular with recent advances in the study and development of metamaterials <cit.>. Metamaterials can take advantage of precisely-defined component shapes, geometries, orientations, arrangements, and connectivity patterns (rather than specific material and physical characteristics of individual units) to produce tailored mechanical <cit.>, acoustic <cit.>, and electromagnetic <cit.> properties. The control of a single unit or component is relatively straightforward, but the question of how to link many components in ways that yield complex material properties is a very challenging one. Approaches that use ideas from network science have the potential to offer guidance for construction patterns of material units that support desired bulk properties. 

There are likely many ways to use network-based approaches to inform the design of new materials. One reasonable possibility is to employ evolutionary and genetic computer algorithms <cit.> and other tools from algorithmic game theory <cit.>. For example, the combination of multi-objective functions and Pareto optimality <cit.> can offer a targeted path through the space of possible network architectures (i.e., through a network morphospace <cit.>). If exact simulations are not computationally tractable, one can use machine-learning techniques to offer fast estimates of material properties and behavior <cit.>. One can perhaps begin with a single material network structure that is physically realizable and then rewire the initial network architecture with a cost function that drives the system towards an arrangement that enhances a desired property. One can perform such a rewiring explicitly along Pareto-optimal fronts using a set of rewiring rules that preserve physical constraints and laws. This approach, which builds on prior applications to other types of spatially-embedded networks <cit.>, selects for physically-feasible network designs that purposely maximize specific properties. A network analysis of these “evolved” (and evolving) networks may help elucidate relationships between structural features of system architecture (e.g., clustering-coefficient values) and material properties (e.g., stability under load), providing a link between structure and function. Such an evolutionary-design approach can complement recent efforts to identify optimal shapes with which to construct a packing <cit.>, to design rules with which to pack those shapes <cit.>, and to construct “allosteric materials” with specific functionalities via evolution according to fitness functions <cit.>.

One set of problems for which network-based tools may be useful are those related to rigidity and mechanical responses of disordered material networks, which have been studied previously using a rigidity-percolation framework <cit.>. In terms of material design, a particularly interesting line of future work may be to use network analysis in conjunction with methods such as tuning-by-pruning, in which individual bonds are selectively removed from a disordered material network to engineer a specific property into a system <cit.>. For example, beginning with simulated, disordered spring networks (derived from jammed particle packings), <cit.> used tuning-by-pruning to design networks with different ratios of the shear modulus to the bulk modulus. Motivated by allosteric responses in proteins, <cit.> developed an approach that allows careful control of local mechanical responses in disordered elastic networks. In particular, with the removal of very few bonds, one can tune the local strain in one region of a network via a response to an applied strain in a distant region of a system. <cit.> studied three different spring networks and continuously tuned their mechanical rigidity (measured by different parameters, such as the distance above isostaticity) to examine the effect of such tuning on material failure. They observed that for a fixed amount of disorder (which the authors measured using a scalar, following the approach in <cit.>), the width of the failure zone in a stressed material increases as the rigidity decreases. Recently, <cit.> studied a model of amorphous networks that incorporates angle-bending forces, and employed pruning-based methods to design auxetic metamaterials. In light of these findings, it is natural to ask what network quantities and methods may be related to the above types of global or local mechanical properties, and whether they can inform which bonds to remove (or rearrange <cit.>) to invoke particular types of functional responses.

In developing a network-based approach for designing and building new materials, it is desirable to capitalize on the ability of network analysis to quantify multiscale structures (with respect to space, time, and network architectures) in a wide variety of systems (regardless of the exact details of their composition). For example, network analysis has revealed mesoscale architectures that are often crucial for determining material properties in disordered media, but such heterogeneities also appear to be important for biological materials, including networks of collagen fibers <cit.>, tissues and tendons <cit.>, muscle fibers <cit.>, and axonal fibers <cit.>. Tools from network science should be useful (and, in principle, flexible) for designing nontrivial structural organizations that yield desired material functions. One can imagine using a network-theoretic framework to design localized, mesoscale, and/or system-level properties, to design and manipulate human-made and natural (including biological) materials, and to precisely control both static and dynamic material properties. 



§ CONCLUSIONS
 


Network science is an interdisciplinary subject — drawing on methods from physics, mathematics, statistics, computer science, social science, and many other disciplines — that has been used successfully to help understand the structure and function of many complex systems. Much recent work on networks has yielded fascinating insights into granular materials, which consist of collections of discrete, macroscopic particles whose contact interactions give rise to many interesting behaviors and intricate organization on multiple spatial and temporal scales. These insights have increased scientific understanding of the structure and dynamics of the heterogeneous material architecture that manifests in granular matter (as well as network-based approaches to quantify such architecture) and the response of granular systems to external perturbations such as compression, shear, tapping, and tilting. In this paper, we have reviewed the increasingly fertile intersection of network science and granular materials. Future efforts should help provide a better understanding of the physics of particulate matter more generally, elucidate the roles of mesoscale interaction patterns in mechanical failure, inform the design of new materials with desired properties, and further scientific understanding of numerous important problems and applications in granular physics, soft-matter physics more generally, and even biophysics. 



§ ACKNOWLEDGEMENTS


We thank Alejandro J. Martínez, Ann E. Sizemore, Konstantin Mischaikow, Jen Schwarz, and an anonymous referee for helpful comments. We also thank our collaborators, students, and mentors, who have shaped our views on the subjects of this review. KD is grateful for support from the National Science Foundation (DMR-0644743, DMR-1206808) and the James S. McDonnell Foundation. LP is grateful to the National Science Foundation for a Graduate Research Fellowship. DSB is grateful to the Alfred P. Sloan Foundation, the John D. and Catherine T. MacArthur Foundation, the Paul G. Allen Foundation, and to the National Science Foundation (PHY-1554488). The content is solely the responsibility of the authors and does not necessarily represent the official views of any of the funding agencies. 



389
	
	
	urlstyle[Jaeger et al.(1996)Jaeger, Nagel, and Behringer]Jaeger1996
	H M Jaeger, S R Nagel, and R P Behringer.
	Granular solids, liquids, and gases.
	Rev Mod Phys, 680 (4):0 1259–1273, 1996.
	
	[Duran(1999)]duran1999sands
	J Duran.
	Sands, Powders, and Grains: An Introduction to the Physics of
		Granular Materials.
	Springer-Verlag, 1999.
	
	[Mehta(2007)]Mehta2007
	A Mehta.
	Granular Physics.
	Cambridge University Press, 2007.
	
	[Franklin and Shattuck(2015)]Franklin2015
	S V Franklin and M D Shattuck, editors.
	Handbook of Granular Materials.
	CRC Press, 2015.
	
	[Andreotti et al.(2013)Andreotti, Forterre, and
	Pouliquen]Andreotti2013
	B Andreotti, Y Forterre, and O Pouliquen.
	Granular Media: Between Solid and Fluid.
	Cambridge University Press, 2013.
	
	[Nagel(2017)]Nagel:2017a
	S R Nagel.
	Experimental soft-matter science.
	Rev Mod Phys, 890 (2):0 025002, 2017.
	
	[Mort et al.(2015)Mort, Michaels, Behringer, Campbell, Kondic,
	Kheiripour Langroudi, Shattuck, Tang, Tardos, and Wassgren]Mort2015
	P Mort, J N Michaels, R P Behringer, C S Campbell, L Kondic,
	M Kheiripour Langroudi, M Shattuck, J Tang, G I Tardos, and C Wassgren.
	Dense granular flow — A collaborative study.
	Powder Technol, 284:0 571–584, 2015.
	
	[Liu et al.(1995)Liu, Nagel, Schecter, Coppersmith, Majumdar, Narayan,
	and Witten]Liu:1995aa
	C-H Liu, S R Nagel, D A Schecter, S N Coppersmith, S Majumdar, O Narayan, and
	T A Witten.
	Force fluctuations in bead packs.
	Science, 2690 (5223):0 513–515, 1995.
	
	[Mueth et al.(1998)Mueth, Jaeger, and Nagel]Mueth:1998a
	D M Mueth, H M Jaeger, and S R Nagel.
	Force distribution in a granular medium.
	Phys Rev E, 570 (3):0 3164–3169, 1998.
	
	[Coppersmith et al.(1996)Coppersmith, Liu, Majumdar, Narayan, and
	Witten]Coppersmith:1996a
	S N Coppersmith, C-h Liu, S Majumdar, O Narayan, and T A Witten.
	Model for force fluctuations in bead packs.
	Phys Rev E, 530 (5):0 4673–4685, 1996.
	
	[Claudin et al.(1998)Claudin, Bouchaud, Cates, and
	Wittmer]Claudin:1998a
	P Claudin, J-P Bouchaud, M E Cates, and J P Wittmer.
	Models of stress fluctuations in granular media.
	Phys Rev E, 570 (4):0 4441–4457, 1998.
	
	[Sexton et al.(1999)Sexton, Socolar, and Schaeffer]Sexton:1999a
	M G Sexton, J E S Socolar, and D G Schaeffer.
	Force distribution in a scalar model for noncohesive granular
	material.
	Phys Rev E, 600 (2):0 1999–2008, 1999.
	
	[Socolar et al.(2002)Socolar, Schaeffer, and Claudin]Socolar:2002a
	J E S Socolar, D G Schaeffer, and P Claudin.
	Directed force chain networks and stress response in static granular
	materials.
	Eur Phys J E, 70 (4):0 353–370, 2002.
	
	[Peters et al.(2005)Peters, Muthuswamy, Wibowo, and
	Tordesillas]Peters:2005aa
	J F Peters, M Muthuswamy, J Wibowo, and A Tordesillas.
	Characterization of force chains in granular material.
	Phys Rev E, 720 (4):0 041307, 2005.
	
	[Howell et al.(1999)Howell, Behringer, and Veje]Howell1999
	D Howell, R P Behringer, and C Veje.
	Stress fluctuations in a 2D granular Couette experiment: A
	continuous transition.
	Phys Rev Lett, 820 (26):0 5241–5244, 1999.
	
	[Majmudar and Behringer(2005)]Majmudar:2005aa
	T S Majmudar and R P Behringer.
	Contact force measurements and stress-induced anisotropy in granular
	materials.
	Nature, 4350 (1079):0 1079–1082, 2005.
	
	[Geng et al.(2001)Geng, Howell, Longhi, Behringer, Reydellet, Vanel,
	Clément, and Luding]Geng:2001a
	J Geng, D Howell, E Longhi, R P Behringer, G Reydellet, L Vanel, E Clément,
	and S Luding.
	Footprints in sand: The response of a granular material to local
	perturbations.
	Phys Rev Lett, 870 (3):0 035506, 2001.
	
	[Radjai et al.(1998)Radjai, Wolf, Jean, and Moreau]Radjai:1998aa
	F Radjai, D E Wolf, M Jean, and J-J Moreau.
	Bimodal character of stress transmission in granular packings.
	Phys Rev Lett, 800 (1):0 61–64, 1998.
	
	[Cates et al.(1999)Cates, Wittmer, Bouchaud, and
	Claudin]cates1999jamming
	M E Cates, J P Wittmer, J-P Bouchaud, and P Claudin.
	Jamming and static stress transmission in granular materials.
	Chaos, 90 (3):0 511–522, 1999.
	
	[Bassett et al.(2012)Bassett, Owens, Daniels, and
	Porter]bassett2012influence
	D S Bassett, E T Owens, K E Daniels, and M A Porter.
	Influence of network topology on sound propagation in granular
	materials.
	Phys Rev E, 860 (4):0 041306, 2012.
	
	[Richard et al.(2005)Richard, Nicodemi, Delannay, Ribière, and
	Bideau]richard2005slow
	P Richard, M Nicodemi, R Delannay, P Ribière, and D Bideau.
	Slow relaxation and compaction of granular systems.
	Nat Mater, 40 (2):0 121–128, 2005.
	
	[Owens and Daniels(2011)]Owens:2011
	E T Owens and K E Daniels.
	Sound propagation and force chains in granular materials.
	EPL (Europhysics Letters), 940 (5):0 54005,
	2011.
	
	[Smart et al.(2007)Smart, Umbanhowar, and Ottino]Smart:2007a
	A Smart, P Umbanhowar, and J Ottino.
	Effects of self-organization on transport in granular matter: A
	network-based approach.
	EPL (Europhysics Letters), 790 (2):0 24002,
	2007.
	
	[Gervois et al.(1989)Gervois, Ammi, Travers, Bideau, Messager, and
	Troadec]Gervois:1989a
	A Gervois, M Ammi, T Travers, D Bideau, J-C Messager, and J-P Troadec.
	Importance of disorder in the conductivity of packings under
	compression.
	Physica A, 1570 (1):0 565–569, 1989.
	
	[Combe et al.(2015)Combe, Richefeu, Stasiak, and
	Atman]combe2015experimental
	G Combe, V Richefeu, M Stasiak, and A P F Atman.
	Experimental validation of a nonextensive scaling law in confined
	granular media.
	Phys Rev Lett, 1150 (23):0 238301, 2015.
	
	[Behringer et al.(2014)Behringer, Bi, Chakraborty, Clark, Dijksman,
	Ren, and Zhang]Behringer2014:Statistical
	R P Behringer, D Bi, B Chakraborty, A Clark, J Dijksman, J Ren, and J Zhang.
	Statistical properties of granular materials near jamming.
	Journal of Statistical Mechanics: Theory and Experiment,
	20140 (6):0 P06004, 2014.
	
	[Bassett et al.(2015a)Bassett, Owens, Porter, Manning, and
	Daniels]bassett2015extraction
	D S Bassett, E T Owens, M A Porter, M L Manning, and K E Daniels.
	Extraction of force-chain network architecture in granular materials
	using community detection.
	Soft Matter, 110 (14):0 2731–2744,
	2015a.
	
	[Herrera et al.(2011)Herrera, McCarthy, Slotterback, Cephas, Losert,
	and Girvan]Herrera:2011aa
	M Herrera, S McCarthy, S Slotterback, E Cephas, W Losert, and M Girvan.
	Path to fracture in granular flows: Dynamics of contact networks.
	Phys Rev E, 830 (6):0 061303, 2011.
	
	[Newman(2010)]newman2010networks
	M E J Newman.
	Networks: An Introduction.
	Oxford University Press, 2010.
	
	[Bollobás(1998)]Bollobas:1998a
	B Bollobás.
	Modern Graph Theory.
	Springer-Verlag, 1998.
	
	[Fortunato and Hric(2016)]Fortunato2016
	S Fortunato and D Hric.
	Community detection in networks: A user guide.
	Phys Rep, 659:0 1–44, 2016.
	
	[Porter et al.(2009)Porter, Onnela, and Mucha]Porter2009
	M A Porter, J-P Onnela, and P J Mucha.
	Communities in networks.
	Not Amer Math Soc, 560 (9):0 1082–1097,
	1164–1166, 2009.
	
	[Fortunato(2010)]Fortunato2010
	S Fortunato.
	Community detection in graphs.
	Phys Rep, 4860 (3–5):0 75–174, 2010.
	
	[Csermely et al.(2013)Csermely, London, Wu, and Uzzi]csermely2013
	P Csermely, A London, L-Y Wu, and B Uzzi.
	Structure and dynamics of core–periphery networks.
	Journal of Complex Networks Complex Networks, 10
	(2):0 93–123, 2013.
	
	[Newman(2011)]newman2011complex
	M E J Newman.
	Complex systems: A survey.
	Am J Phys, 79:0 800–810, 2011.
	
	[Kivelä et al.(2014)Kivelä, Arenas, Barthelemy, Gleeson, Moreno,
	and Porter]kivela2014multilayer
	M Kivelä, A Arenas, M Barthelemy, J P Gleeson, Y Moreno, and M A Porter.
	Multilayer networks.
	Journal of Complex Networks, 20 (3):0
	203–271, 2014.
	
	[Barthélemy(2011)]barth2011
	M Barthélemy.
	Spatial networks.
	Phys Rep, 4990 (1):0 1–101, 2011.
	
	[Cruz Hidalgo et al.(2002)Cruz Hidalgo, Grosse, Kun, Reinhardt, and
	Herrmann]Hidalgo2002
	R Cruz Hidalgo, C U Grosse, F Kun, H W Reinhardt, and H J Herrmann.
	Evolution of Percolating Force Chains in Compressed Granular Media.
	Phys Rev Lett, 890 (20):0 205501, 2002.
	
	[Candelier et al.(2009)Candelier, Dauchot, and
	Biroli]candelier2009building
	R Candelier, O Dauchot, and G Biroli.
	Building blocks of dynamical heterogeneities in dense granular media.
	Phys Rev Lett, 1020 (8):0 088001, 2009.
	
	[Mehta et al.(2008)Mehta, Barker, and Luck]mehta2008heterogeneities
	A Mehta, G C Barker, and J M Luck.
	Heterogeneities in granular dynamics.
	Proc Natl Acad Sci, 1050 (24):0 8244–8249,
	2008.
	
	[Keys et al.(2007)Keys, Abate, Glotzer, and Durian]Keys:2007a
	A S Keys, A R Abate, S C Glotzer, and D J Durian.
	Measurement of growing dynamical length scales and prediction of the
	jamming transition in a granular material.
	Nat Phys, 30 (4):0 260–264, 2007.
	
	[Digby(1981)]Digby:1981
	PJ Digby.
	The effective elastic moduli of porous granular rocks.
	J Appl Mech, 480 (4):0 803–808, 1981.
	
	[Velický and Caroli(2002)]Velicky:2002
	B Velický and C Caroli.
	Pressure dependence of the sound velocity in a two-dimensional
	lattice of Hertz–Mindlin balls: Mean-field description.
	Phys Rev E, 650 (2):0 021307, 2002.
	
	[Goddard(1990)]Goddard:1990
	J D Goddard.
	Nonlinear elasticity and pressure-dependent wave speeds in granular
	media.
	Proc R Soc A, 4300 (1878):0 105–131, 1990.
	
	[Makse et al.(1999)Makse, Gland, Johnson, and Schwartz]Makse-1999-WEM
	H A Makse, N Gland, D L Johnson, and L M Schwartz.
	Why effective medium theory fails in granular materials.
	Phys Rev Lett, 830 (24):0 5070–5073, 1999.
	
	[Goldenberg and Goldhirsch(2005)]Goldenberg-2005-FEE
	C Goldenberg and I Goldhirsch.
	Friction Enhances Elasticity In Granular Solids.
	Nature, 435:0 188–191, 2005.
	
	[Smart and Ottino(2008a)]Smart:2008b
	A Smart and J M Ottino.
	Granular matter and networks: Three related examples.
	Soft Matter, 40 (11):0 2125–2131,
	2008a.
	
	[Holme and Saramäki(2012)]Holme2011
	P Holme and J Saramäki.
	Temporal networks.
	Phys Rep, 5190 (3):0 97–125, 2012.
	
	[Sayama et al.(2013)Sayama, Pestov, Schmidt, Bush, Wong, Yamanoi, and
	Gross]thilo-adaptive
	H Sayama, I Pestov, J Schmidt, B J Bush, C Wong, J Yamanoi, and T Gross.
	Modeling complex systems with adaptive networks.
	Computers & Mathematics with Applications, 650
	(10):0 1645–1664, 2013.
	
	[Newman et al.(2006)Newman, Barabási, and Watts]Newman:2006aa
	M E J Newman, A-L Barabási, and D J Watts.
	The Structure and Dynamics of Networks.
	Princeton University Press, 2006.
	
	[Giusti et al.(2016a)Giusti, Ghrist, and
	Bassett]giusti2016twos
	C Giusti, R Ghrist, and D S Bassett.
	Two's company, three (or more) is a simplex: Algebraic-topological
	tools for understanding higher-order structure in neural data.
	J Comput Neurosci, 410 (1):0 1–14,
	2016a.
	
	[Porter and Gleeson(2016)]porter2016dynamical
	M A Porter and J P Gleeson.
	Dynamical systems on networks: A tutorial.
	In Frontiers in Applied Dynamical Systems: Reviews and
		Tutorials, volume 4. Springer-Verlag, 2016.
	
	[Liben-Nowell and Kleinberg(2008)]liben2008tracing
	D Liben-Nowell and J Kleinberg.
	Tracing information flow on a global scale using internet
	chain-letter data.
	Proc Natl Acad Sci, 1050 (12):0 4633–4638,
	2008.
	
	[Scott(2012)]Scott:2012a
	J Scott.
	Social Network Analysis.
	Sage Publications, 2012.
	
	[Hurd et al.(2017)Hurd, Gleeson, and Melnik]hurd2017framework
	T R Hurd, J P Gleeson, and S Melnik.
	A framework for analyzing contagion in assortative banking networks.
	PLoS One, 120 (2):0 e0170579, 2017.
	
	[Sporns(2013)]Sporns:2013aa
	O Sporns.
	Structure and function of complex brain networks.
	Dialogues Clin Neurosci, 150 (3):0 247–262,
	2013.
	
	[Bassett and Sporns(2017)]Bassett:2017a
	D S Bassett and O Sporns.
	Network neuroscience.
	Nat Neurosci, 200 (3):0 353–364, 2017.
	
	[Walker and Tordesillas(2010)]Walker:2010aa
	D M Walker and A Tordesillas.
	Topological evolution in dense granular materials: A complex
	networks perspective.
	Int J Solids Struct, 470 (5):0 624–639, 2010.
	
	[Barrat et al.(2004)Barrat, Barthélemy, Pastor-Satorras, and
	Vespignani]Barrat:2004aa
	A Barrat, M Barthélemy, R Pastor-Satorras, and A Vespignani.
	The architecture of complex weighted networks.
	Proc Natl Acad Sci, 1010 (11):0 3747–3752,
	2004.
	
	[Newman(2004)]Newman:2004aa
	M E J Newman.
	Analysis of weighted networks.
	Phys Rev E, 700 (5):0 056131, 2004.
	
	[Alexander(1998)]Alexander:1998aa
	S Alexander.
	Amorphous solids: Their structure, lattice dynamics and elasticity.
	Phys Rep, 2960 (2–4):0 65–236, 1998.
	
	[Wyart(2005)]Wyart:2005a
	M Wyart.
	On the rigidity of amorphous solids.
	Annales De Physique, 300 (3):0 1, 2005.
	
	[Liu and Nagel(2010)]Liu:2010aa
	A J Liu and S R Nagel.
	The jamming transition and the marginally jammed solid.
	Annu Rev Condens Matter Phys, 10 (1):0
	347–369, 2010.
	
	[van Hecke(2010)]VanHecke2010
	M van Hecke.
	Jamming of soft particles: Geometry, mechanics, scaling and
	isostaticity.
	J Phys Condens Matter, 220 (3):0 033101, 2010.
	
	[Masuda et al.(2017)Masuda, Porter, and Lambiotte]masuda2016
	N Masuda, M A Porter, and R Lambiotte.
	Random walks and diffusion on networks.
	Phys Rep, 716-717:0 1–58, 2017.
	
	[Krioukov et al.(2010)Krioukov, Papadopoulos, Kitsak, Vahdat, and
	Boguñá]boguna2010
	D Krioukov, F Papadopoulos, M Kitsak, A Vahdat, and M Boguñá.
	Hyperbolic geometry of complex networks.
	Phys Rev E, 820 (3):0 036106, 2010.
	
	[Skiena(2008)]skiena2008algorithm
	S Skiena.
	The Algorithm Design Manual.
	Springer-Verlag, 2008.
	
	[Watts and Strogatz(1998)]Watts:1998aa
	D J Watts and S H Strogatz.
	Collective dynamics of `small-world' networks.
	Nature, 3930 (6684):0 440–442, 1998.
	
	[Latora and Marchiori(2001)]latora2001efficient
	V Latora and M Marchiori.
	Efficient behavior of small-world networks.
	Phys Rev Lett, 870 (19):0 198701, 2001.
	
	[Rubinov and Sporns(2009)]Rubinov2009
	M Rubinov and O Sporns.
	Complex network measures of brain connectivity: Uses and
	interpretations.
	NeuroImage, 520 (3):0 1059–1069, 2009.
	
	[Latora and Marchiori(2003)]Latora2003
	V Latora and M Marchiori.
	Economic small-world behavior in weighted networks.
	Eur Phys J B, 32:0 249–263, 2003.
	
	[Estrada and Hatano(2016)]estrada2016
	E Estrada and N Hatano.
	Communicability angle and the spatial efficiency of networks.
	SIAM Review, 580 (4):0 692–715, 2016.
	
	[Arévalo et al.(2010a)Arévalo, Zuriguel, and
	Maza]Arevalo:2010aa
	R Arévalo, I Zuriguel, and D Maza.
	Topology of the force network in the jamming transition of an
	isotropically compressed granular packing.
	Phys Rev E, 810 (4):0 041302,
	2010a.
	
	[Gross and Yellen(2005)]Gross:2005a
	J L Gross and J Yellen.
	Graph theory and its Applications.
	CRC Press, 2005.
	
	[Kavitha et al.(2009)Kavitha, Liebchen, Mehlhorn, Michail, Rizzi,
	Ueckerdt, and Zweig]Kavitha:2009a
	T Kavitha, C Liebchen, K Mehlhorn, D Michail, R Rizzi, T Ueckerdt, and K A
	Zweig.
	Cycle bases in graphs characterization, algorithms, complexity, and
	applications.
	Comp Sci Rev, 30 (4):0 199–243, 2009.
	
	[Griffin(2017)]Griffin:2017a
	C Griffin.
	Graph theory: Penn State Math 485 lecture notes.
	Webpage, 2017.
	URL <http://www.personal.psu.edu/cxg286/Math485.pdf>.
	(This manuscript includes contributions by S Shekhar.).
	
	[Horton(1987)]Horton:1987a
	J D Horton.
	A polynomial-time algorithm to find the shortest cycle basis of a
	graph.
	SIAM Journal on Computing, 160 (2):0 358–366,
	1987.
	
	[Mehlhorn and Michail(2007)]Mehlhorn:2006a
	K Mehlhorn and D Michail.
	Implementing minimum cycle basis algorithms.
	J Exp Algorithmics, 11, 2007.
	
	[Walker et al.(2014a)Walker, Tordesillas, and
	Froyland]Tordesillas:2014a
	D M Walker, A Tordesillas, and G Froyland.
	Mesoscale and macroscale kinetic energy fluxes from granular fabric
	evolution.
	Phys Rev E, 890 (3):0 032205,
	2014a.
	
	[Walker et al.(2015a)Walker, Tordesillas, Brodu, Dijksman,
	Behringer, and Froyland]Walker:2015b
	D M Walker, A Tordesillas, N Brodu, J A Dijksman, R P Behringer, and
	G Froyland.
	Self-assembly in a near-frictionless granular material:
	Conformational structures and transitions in uniaxial cyclic compression of
	hydrogel spheres.
	Soft Matter, 11:0 2157–2173, 2015a.
	
	[Smart and Ottino(2008b)]Smart:2008aa
	A G Smart and J M Ottino.
	Evolving loop structure in gradually tilted two-dimensional granular
	packings.
	Phys Rev E, 770 (4):0 041307,
	2008b.
	
	[Arévalo et al.(2009)Arévalo, Zuriguel, and Maza]Arevalo:2009aa
	R Arévalo, I Zuriguel, and D Maza.
	Topological properties of the contact network of granular materials.
	Int J Bifurc Chaos, 190 (02):0 695–702, 2009.
	
	[Arévalo et al.(2010b)Arévalo, Zuriguel, Trevijano,
	and Maza]Arevalo:2010ba
	R Arévalo, I Zuriguel, S A Trevijano, and D Maza.
	Third order loops of contacts in a granular force network.
	Int J Bifurc Chaos, 200 (03):0 897–903,
	2010b.
	
	[Tordesillas et al.(2010a)Tordesillas, Walker, and
	Lin]Tordesillas:2010aa
	A Tordesillas, D M Walker, and Q Lin.
	Force cycles and force chains.
	Phys Rev E, 810 (1):0 011302,
	2010a.
	
	[Rivier(2006)]rivier2006extended
	N Rivier.
	Extended constraints, arches and soft modes in granular materials.
	J Non-Cryst Solids, 3520 (42–49):0
	4505–4508, 2006.
	
	[Newman(2003)]Newman:2003aa
	M E J Newman.
	The structure and function of complex networks.
	SIAM Review, 450 (2):0 167–256, 2003.
	
	[Barrat and Weigt(2000)]Barrat:2000aa
	A Barrat and M Weigt.
	On the properties of small-world network models.
	Eur Phys J B, 130 (3):0 547–560, 2000.
	
	[Saramäki et al.(2007)Saramäki, Kivelä, Onnela, Kaski, and
	Kertész]jari2007
	J Saramäki, M Kivelä, J-P Onnela, K Kaski, and J Kertész.
	Generalizations of the clustering coefficient to weighted complex
	networks.
	Phys Rev E, 750 (2):0 027105, 2007.
	
	[Onnela et al.(2005)Onnela, Saramäki, Kertész, and
	Kaski]onnela2005
	J-P Onnela, J Saramäki, J Kertész, and K Kaski.
	Intensity and coherence of motifs in weighted complex networks.
	Phys Rev E, 710 (6):0 065103, 2005.
	
	[Zhang and Horvath(2005)]zhang2005general
	B Zhang and S Horvath.
	A general framework for weighted gene co-expression network analysis.
	Stat Appl Genet Mol Biol, 40 (1):0 17, 2005.
	
	[Freeman(1977)]Freeman:1977aa
	L C Freeman.
	A set of measures of centrality based on betweenness.
	Sociometry, 400 (1):0 35–41, 1977.
	
	[Girvan and Newman(2002)]Girvan2002
	M Girvan and M E J Newman.
	Community structure in social and biological networks.
	Proc Natl Acad Sci, 990 (12):0 7821–7826,
	2002.
	
	[Estrada and Rodríguez-Velázquez(2005)]Estrada:2005aa
	E Estrada and J A Rodríguez-Velázquez.
	Subgraph centrality in complex networks.
	Phys Rev E, 710 (5):0 056103, 2005.
	
	[Estrada et al.(2012)Estrada, Hatano, and Benzi]estrada2012
	E Estrada, N Hatano, and M Benzi.
	The physics of communicability in complex networks.
	Phys Rep, 5140 (3):0 89–119, 2012.
	
	[Estrada and Rodríguez-Velázquez(2005)]Estrada:2005ab
	E Estrada and J A Rodríguez-Velázquez.
	Spectral measures of bipartivity in complex networks.
	Phys Rev E, 720 (4):0 046105, 2005.
	
	[Milo et al.(2002)Milo, Shen-Orr, Itzkovitz, Kashtan, Chklovskii, and
	Alon]Milo:2002a
	R Milo, S S Shen-Orr, S Itzkovitz, N Kashtan, D Chklovskii, and U Alon.
	Network motifs: Simple building blocks of complex networks.
	Science, 2980 (5594):0 824–827, 2002.
	
	[Shen-Orr et al.(2002)Shen-Orr, Milo, Mangan, and
	Alon]shenorr2002network
	S S Shen-Orr, R Milo, S Mangan, and U Alon.
	Network motifs in the transcriptional regulation network of
	Escherichia coli.
	Nat Genet, 310 (1):0 64–68, 2002.
	
	[Milo et al.(2004)Milo, Itzkovitz, Kashtan, Levitt, Shen-Orr,
	Ayzenshtat, Sheffer, and Alon]milo2004superfamilies
	R Milo, S Itzkovitz, N Kashtan, R Levitt, S Shen-Orr, I Ayzenshtat, M Sheffer,
	and U Alon.
	Superfamilies of evolved and designed networks.
	Science, 3030 (5663):0 1538–1542, 2004.
	
	[Alon(2007)]alon2007network
	U Alon.
	Network motifs: Theory and experimental approaches.
	Nat Rev Genet, 80 (6):0 450–461, 2007.
	
	[Schreiber and Schwobbermeyer(2005)]schreiber2005frequency
	F Schreiber and H Schwobbermeyer.
	Frequency concepts and pattern detection for the analysis of motifs
	in networks.
	Transactions on Computational Systems Biology, III:0
	89–104, 2005.
	
	[Wernicke(2006)]wernicke2006efficient
	S Wernicke.
	Efficient detection of network motifs.
	IEEE/ACM Trans Comput Biol Bioinform, 30 (4):0
	347–359, 2006.
	
	[Grochow and Kellis(2007)]grochow2007network
	J A Grochow and M Kellis.
	Network motif discovery using sub-graph enumeration and
	symmetry-breaking.
	RECOMB, pages 92–106, 2007.
	
	[Omidi et al.(2009)Omidi, Schreiber, and Masoudi-Nejad]omidi2009moda
	S Omidi, F Schreiber, and A Masoudi-Nejad.
	MODA: An efficient algorithm for network motif discovery in
	biological networks.
	Genes Genet Syst, 840 (5):0 385–395, 2009.
	
	[Kashani et al.(2009)Kashani, Ahrabian, Elahi, Nowzari-Dalini, Ansari,
	Asadi, Mohammadi, Schreiber, and Masoudi-Nejad]kashani2009kavosh
	Z R Kashani, H Ahrabian, E Elahi, A Nowzari-Dalini, E S Ansari, S Asadi,
	S Mohammadi, F Schreiber, and A Masoudi-Nejad.
	Kavosh: A new algorithm for finding network motifs.
	BMC Bioinformatics, 100 (318):0 318, 2009.
	
	[Paulau et al.(2015)Paulau, Feenders, and Blasius]paulau2015motif
	P V Paulau, C Feenders, and B Blasius.
	Motif analysis in directed ordered networks and applications to food
	webs.
	Sci Rep, 5:0 11926, 2015.
	
	[Sporns and Kotter(2004)]sporns2004motifs
	O Sporns and R Kotter.
	Motifs in brain networks.
	PLoS Biol, 20 (11):0 e369, 2004.
	
	[Xu et al.(2008)Xu, Zhang, and Small]Xu:2008aa
	X Xu, J Zhang, and M Small.
	Superfamily phenomena and motifs of networks induced from time
	series.
	Proc Natl Acad Sci, 1050 (50):0 19601–19605,
	2008.
	
	[Walker et al.(2014b)Walker, Tordesillas, Small,
	Behringer, and Tse]Walker:2014aa
	D M Walker, A Tordesillas, M Small, R P Behringer, and C K Tse.
	A complex systems analysis of stick-slip dynamics of a laboratory
	fault.
	Chaos, 240 (1):0 013132, 2014b.
	
	[Walker et al.(2015b)Walker, Tordesillas, Zhang,
	Behringer, Andò, Viggiani, Druckrey, and Alshibli]Walker:2015c
	D M Walker, A Tordesillas, J Zhang, R P Behringer, E Andò, G Viggiani,
	A Druckrey, and K Alshibli.
	Structural templates of disordered granular media.
	Int J Solids Struct, 54:0 20–30, 2015b.
	
	[Tordesillas et al.(2012a)Tordesillas, Walker, Froyland,
	Zhang, and Behringer]Tordesillas:2012aa
	A Tordesillas, D M Walker, G Froyland, J Zhang, and R P Behringer.
	Transition dynamics and magic-number-like behavior of frictional
	granular clusters.
	Phys Rev E, 860 (1):0 011306,
	2012a.
	
	[Peixoto(2017)]peixoto2017
	T P Peixoto.
	Bayesian stochastic blockmodeling.
	arXiv, arXiv:1705.10225 [stat.ML], 2017.
	
	[Giusti et al.(2016b)Giusti, Papadopoulos, Owens, Daniels,
	and Bassett]Giusti:2016a
	C Giusti, L Papadopoulos, E T Owens, K E Daniels, and D S Bassett.
	Topological and geometric measurements of force-chain structure.
	Phys Rev E, 940 (3):0 032909,
	2016b.
	
	[Papadopoulos et al.(2016)Papadopoulos, Puckett, Daniels, and
	Bassett]papadopoulos2016evolution
	L Papadopoulos, J G Puckett, K E Daniels, and D S Bassett.
	Evolution of network architecture in a granular material under
	compression.
	Phys Rev E, 940 (3):0 032908, 2016.
	
	[Walker and Tordesillas(2012)]Walker:2012aa
	D M Walker and A Tordesillas.
	Taxonomy of granular rheology from grain property networks.
	Phys Rev E, 850 (1):0 011304, 2012.
	
	[Tordesillas et al.(2013a)Tordesillas, Walker, Andò,
	and Viggiani]Tordesillas:2013a
	A Tordesillas, D M Walker, E Andò, and G Viggiani.
	Revisiting localized deformation in sand with complex systems.
	Proc Math Phys Eng Sci, 4690 (2152),
	2013a.
	
	[Walker and Tordesillas(2014)]Walker:2014ba
	D M Walker and A Tordesillas.
	Examining overlapping community structures within grain property
	networks.
	In 2014 IEEE International Symposium on Circuits and Systems
		(ISCAS), pages 1275–1278, 2014.
	
	[Walker et al.(2012)Walker, Tordesillas, Pucilowski, Lin, Rechenmacher,
	and Abedi]Walker:2012a
	D M Walker, A Tordesillas, S Pucilowski, Q Lin, A L Rechenmacher, and S Abedi.
	Analysis of grain-scale measurements of sand using kinematical
	complex networks.
	Int J Bifurc Chaos, 220 (12):0 1230042, 2012.
	
	[Newman and Girvan(2004)]NG2004
	M E J Newman and M Girvan.
	Finding and evaluating community structure in networks.
	Phys Rev E, 690 (2):0 026113, 2004.
	
	[Newman(2006a)]Newman2006b
	M E J Newman.
	Finding community structure in networks using the eigenvectors of
	matrices.
	Phys Rev E, 740 (3):0 036104,
	2006a.
	
	[Rosvall and Bergstrom(2008)]rosvall2008maps
	M Rosvall and C T Bergstrom.
	Maps of random walks on complex networks reveal community structure.
	Proc Natl Acad Sci, 1050 (4):0 1118–1123,
	2008.
	
	[Clauset(2005)]clauset2005finding
	A Clauset.
	Finding local community structure in networks.
	Phys Rev E, 720 (2):0 026132, 2005.
	
	[Jeub et al.(2015)Jeub, Balachandran, Porter, Mucha, and
	Mahoney]jeub2015
	L G S Jeub, P Balachandran, M A Porter, P J Mucha, and W M Mahoney.
	Think locally, act locally: The detection of small, medium-sized,
	and large communities in large networks.
	Phys Rev E, 910 (1):0 012821, 2015.
	
	[Ahn et al.(2010)Ahn, Bagrow, and Lehmann]ahn2010link
	Y Y Ahn, J P Bagrow, and S Lehmann.
	Link communities reveal multiscale complexity in networks.
	Nature, 4660 (7307):0 761–764, 2010.
	
	[Good et al.(2010)Good, de Montjoye, and Clauset]Good2010
	B H Good, Y A de Montjoye, and A Clauset.
	Performance of modularity maximization in practical contexts.
	Phys Rev E, 810 (4):0 046106, 2010.
	
	[Fortunato and Barthélemy(2007)]Fortunato2007
	S. Fortunato and M. Barthélemy.
	Resolution limit in community detection.
	Proc Natl Acad Sci, 1040 (1):0 36–41, 2007.
	
	[Bassett et al.(2013)Bassett, Porter, Wymbs, Grafton, Carlson, and
	Mucha]bassett2013robust
	D S Bassett, M A Porter, N F Wymbs, S T Grafton, J M Carlson, and P J Mucha.
	Robust detection of dynamic community structure in networks.
	Chaos, 230 (1):0 013142, 2013.
	
	[Newman(2006b)]Newman2006
	M E J Newman.
	Modularity and community structure in networks.
	Proc Natl Acad Sci, 1030 (23):0 8577–8582,
	2006b.
	
	[Bazzi et al.(2016)Bazzi, Porter, Williams, McDonald, Fenn, and
	Howison]Bazzi2016
	M Bazzi, M A Porter, S Williams, M McDonald, D J Fenn, and S D Howison.
	Community detection in temporal multilayer networks, with an
	application to correlation networks.
	Multiscale Model Simul, 140 (1):0 1–41, 2016.
	
	[Sarzynska et al.(2016)Sarzynska, Leicht, Chowell, and
	Porter]Sarzynska:2015aa
	M Sarzynska, E A Leicht, G Chowell, and M A Porter.
	Null models for community detection in spatially embedded, temporal
	networks.
	Journal of Complex Networks, 4:0 363–406, 2016.
	
	[Brandes et al.(2008)Brandes, Delling, Gaertler, Görke, Hoefer,
	Nikoloski, and Wagner]Brandes2008
	U Brandes, D Delling, M Gaertler, R Görke, M Hoefer, Z Nikoloski, and
	D Wagner.
	On modularity clustering.
	IEEE Trans on Knowl Data Eng, 200 (2):0
	172–188, 2008.
	
	[Blondel et al.(2008)Blondel, Guillaume, Lambiotte, and
	Lefebvre]Blondel2008
	V D Blondel, J L Guillaume, R Lambiotte, and E Lefebvre.
	Fast unfolding of community hierarchies in large networks.
	Journal of Statistical Mechanics: Theory and Experiment,
	20080 (10):0 P10008, 2008.
	
	[Jeub et al.(2011–2016)Jeub, Bazzi, Jutla, and Mucha]genlouvain2016
	L G S Jeub, M Bazzi, I S Jutla, and P J Mucha.
	A generalized Louvain method for community detection implemented in
	matlab, 2011–2016.
	URL <https://github.com/GenLouvain/GenLouvain>.
	
	[Lancichinetti and Fortunato(2012)]Lancichinetti2012
	A Lancichinetti and S Fortunato.
	Consensus clustering in complex networks.
	Sci Rep, 2:0 336, 2012.
	
	[Jeub et al.(2017)Jeub, Sporns, and Fortunato]jeub-santo2017
	L G S Jeub, O Sporns, and S Fortunato.
	Multiresolution consensus clustering in networks.
	arXiv, arXiv:1710.02249 [cs.SI], 2017.
	
	[Boccaletti et al.(2014)Boccaletti, Bianconi, Criado, Del Genio,
	Gómez-Gardeñes, Romance, Sendiña-Nadal, Wang, and
	Zanin]Boccaletti2014
	S Boccaletti, G Bianconi, R Criado, C I Del Genio, J Gómez-Gardeñes,
	M Romance, I Sendiña-Nadal, Z Wang, and M Zanin.
	The structure and dynamics of multilayer networks.
	Phys Rep, 5440 (1):0 1–122, 2014.
	
	[Mucha et al.(2010)Mucha, Richardson, Macon, Porter, and
	Onnela]Mucha2010
	P J Mucha, T Richardson, K Macon, M A Porter, and J-P Onnela.
	Community structure in time-dependent, multiscale, and multiplex
	networks.
	Science, 3280 (5980):0 876–878, 2010.
	
	[De Domenico et al.(2013)De Domenico, Solé-Ribalta, Cozzo, Kivelä,
	Moreno, Porter, Gomez, and Arenas]domenico2013mathematical
	M De Domenico, A Solé-Ribalta, E Cozzo, M Kivelä, Y Moreno, M A Porter,
	S Gomez, and A Arenas.
	Mathematical formulation of multi-layer networks.
	Phys Rev X, 30 (4):0 041022, 2013.
	
	[Kolda and Bader(2009)]Kolda2009Tensor
	T G Kolda and B W Bader.
	Tensor decompositions and applications.
	SIAM Rev, 510 (3):0 455–500, 2009.
	
	[Cranmer et al.(2015)Cranmer, Menninga, and Mucha]Cranmer:2014ut
	S J Cranmer, E J Menninga, and P J Mucha.
	Kantian fractionalization predicts the conflict propensity of the
	international system.
	Proc Nat Acad Sci, 1120 (38):0 11812–11816,
	2015.
	
	[Danchev and Porter(2017)]danchev2016
	V Danchev and M A Porter.
	Neither global nor local: Heterogeneous connectivity in spatial
	network structures of world migration.
	Social Networks, in press;
	<https://doi.org/10.1016/j.socnet.2017.06.003>, 2017.
	
	[Bassett et al.(2011)Bassett, Wymbs, Porter, Mucha, Carlson, and
	Grafton]Bassett2011b
	D S Bassett, N F Wymbs, M A Porter, P J Mucha, J M Carlson, and S T Grafton.
	Dynamic reconfiguration of human brain networks during learning.
	Proc Natl Acad Sci, 1080 (18):0 7641–7646,
	2011.
	
	[Bassett et al.(2015b)Bassett, Yang, Wymbs, and
	Grafton]bassett2015learning
	D S Bassett, M Yang, N F Wymbs, and S T Grafton.
	Learning-induced autonomy of sensorimotor systems.
	Nat Neurosci, 180 (5):0 744–751,
	2015b.
	
	[Braun et al.(2015)Braun, Schafer, Walter, Erk, Romanczuk-Seiferth,
	Haddad, Schweiger, Grimm, Heinz, Tost, Meyer-Lindenberg, and
	Bassett]braun2015dynamic
	U Braun, A Schafer, H Walter, S Erk, N Romanczuk-Seiferth, L Haddad, J I
	Schweiger, O Grimm, A Heinz, H Tost, A Meyer-Lindenberg, and D S Bassett.
	Dynamic reconfiguration of frontal brain networks during executive
	cognition in humans.
	Proc Natl Acad Sci, 1120 (37):0 11678–11683,
	2015.
	
	[Blinder et al.(2013)Blinder, Tsai, Kaufhold, Knutsen, Suhl, and
	Kleinfeld]Blinder:2013aa
	P Blinder, P S Tsai, J P Kaufhold, P M Knutsen, H Suhl, and D Kleinfeld.
	The cortical angiome: An interconnected vascular network with
	noncolumnar patterns of blood flow.
	Nat Neurosci, 160 (7):0 889–897, 2013.
	
	[Katifori et al.(2010)Katifori, Szöll ősi, and Magnasco]Katifori:2010a
	E Katifori, G J Szöll ősi, and M O
	Magnasco.
	Damage and fluctuations induce loops in optimal transport networks.
	Phys Rev Lett, 1040 (4):0 048704, 2010.
	
	[Lee et al.(2017)Lee, Fricker, and Porter]shl2017
	S H Lee, M D Fricker, and M A Porter.
	Mesoscale analyses of fungal networks as an approach for quantifying
	phenotypic traits.
	Journal of Complex Networks, 50 (1):0
	145–159, 2017.
	
	[Bebber et al.(2007)Bebber, Hynes, Darrah, Boddy, and
	Fricker]Bebber:2007a
	D P Bebber, J Hynes, P R Darrah, L Boddy, and M D Fricker.
	Biological solutions to transport network design.
	Proc R Soc Lond B Biol Sci, 2740 (1623):0
	2307–2315, 2007.
	
	[Banavar et al.(2000)Banavar, Colaiori, Flammini, Maritan, and
	Rinaldo]Banavar:2000a
	J R Banavar, F Colaiori, A Flammini, A Maritan, and A Rinaldo.
	Topology of the fittest transportation network.
	Phys Rev Lett, 840 (20):0 4745–4748, 2000.
	
	[Gastner and Newman(2006)]Newman:2006a
	M T Gastner and M E J Newman.
	Optimal design of spatial distribution networks.
	Phys Rev E, 740 (1):0 016117, 2006.
	
	[Kurant and Thiran(2006)]Kurant:2006a
	M Kurant and P Thiran.
	Extraction and analysis of traffic and topologies of transportation
	networks.
	Phys Rev E, 740 (3):0 036114, 2006.
	
	[Bertsekas(1998)]Bertsekas:1998
	D P Bertsekas.
	Network Optimization: Continuous and Discrete Models.
	Athena Scientific, 1998.
	
	[Ahuja et al.(1993)Ahuja, Magnanti, and Orlin]Ahuja:flow
	R K Ahuja, T L Magnanti, and J B Orlin.
	Network Flows: Theory, Algorithms, and Applications.
	Prentice-Hall, Inc., 1993.
	
	[Kaczynski et al.(2004)Kaczynski, Mischaikow, and
	Mrozek]Kaczynski:2004aa
	T Kaczynski, K Mischaikow, and M Mrozek.
	Computational Homology.
	Springer-Verlag, 2004.
	
	[Kesten(2006)]kesten-whatis
	H Kesten.
	What is ... percolation?
	Not Am Math Soc, 530 (5):0 572–573, 2006.
	
	[Stauffer and Aharony(1994)]Stauffer:1994a
	D Stauffer and A Aharony.
	Introduction to Percolation Theory.
	CRC Press, 1994.
	
	[Saberi(2015)]saberi2015
	A A Saberi.
	Recent advances in percolation theory and its applications.
	Phys Rep, 578:0 1–32, 2015.
	
	[Broadbent and Hammersley(1957)]broadbent1957percolation
	S Broadbent and J Hammersley.
	Percolation processes I. Crystals and mazes.
	Proc Camb Philos Soc, 53:0 629–641, 1957.
	
	[Albert and Barabási(2002)]Albert2002
	R Albert and A-L Barabási.
	Statistical mechanics of complex networks.
	Rev Mod Phys, 740 (1):0 47–98, 2002.
	
	[Erdős and Reńyi(1959)]Erdos:1959aa
	P Erdős and A Reńyi.
	On random graphs I.
	Publicationes Mathematicae Debrecen, 6:0 290–297,
	1959.
	
	[Erdős and Reńyi(1960)]Erdos:1960aa
	P Erdős and A Reńyi.
	On the evolution of random graphs.
	Publications of the Mathematical Institute of the Hungarian
		Academy of Sciences, 5:0 17–61, 1960.
	
	[Stauffer(1979)]Stauffer1979a
	D. Stauffer.
	Scaling theory of percolation clusters.
	Phys Rep, 540 (1):0 1–74, 1979.
	
	[Slotterback et al.(2012)Slotterback, Mailman, Ronaszegi, van Hecke,
	Girvan, and Losert]Slotterback:2012aa
	S Slotterback, M Mailman, K Ronaszegi, M van Hecke, M Girvan, and W Losert.
	Onset of irreversibility in cyclic shear of granular packings.
	Phys Rev E, 850 (2):0 021309, 2012.
	
	[Kondic et al.(2012)Kondic, Goullet, O'Hern, Kramár, Mischaikow, and
	Behringer]Kondic:2012aa
	L Kondic, A Goullet, C S O'Hern, M Kramár, K Mischaikow, and R P Behringer.
	Topology of force networks in compressed granular media.
	EPL (Europhysics Letters), 970 (5):0 54001,
	2012.
	
	[Kramár et al.(2013)Kramár, Goullet, Kondic, and
	Mischaikow]Kramar:2013aa
	M Kramár, A Goullet, L Kondic, and K Mischaikow.
	Persistence of force networks in compressed granular media.
	Phys Rev E, 870 (4):0 042207, 2013.
	
	[Kramár et al.(2014a)Kramár, Goullet, Kondic, and
	Mischaikow]Kramar:2014b
	M Kramár, A Goullet, L Kondic, and K Mischaikow.
	Quantifying force networks in particulate systems.
	Physica D, 283:0 37–55, 2014a.
	
	[Kramár et al.(2014b)Kramár, Goullet, Kondic, and
	Mischaikow]Kramar:2014aa
	M Kramár, A Goullet, L Kondic, and K Mischaikow.
	Evolution of force networks in dense particulate media.
	Phys Rev E, 900 (5):0 052203,
	2014b.
	
	[Ardanza-Trevijano et al.(2014)Ardanza-Trevijano, Zuriguel, Arévalo,
	and Maza]Ardanza-Trevijano:2014aa
	S Ardanza-Trevijano, I Zuriguel, R Arévalo, and D Maza.
	Topological analysis of tapped granular media using persistent
	homology.
	Phys Rev E, 890 (5):0 052212, 2014.
	
	[Kondic et al.(2016)Kondic, Kramár, Pugnaloni, Carlevaro, and
	Mischaikow]Kondic:2016a
	L Kondic, M Kramár, L A Pugnaloni, C M Carlevaro, and K Mischaikow.
	Structure of force networks in tapped particulate systems of disks
	and pentagons. II. Persistence analysis.
	Phys Rev E, 930 (6):0 062903, 2016.
	
	[Pugnaloni et al.(2016)Pugnaloni, Carlevaro, Kramár, Mischaikow, and
	Kondic]Pugnaloni:2016a
	L A Pugnaloni, C M Carlevaro, M Kramár, K Mischaikow, and L Kondic.
	Structure of force networks in tapped particulate systems of disks
	and pentagons. I. Clusters and loops.
	Phys Rev E, 930 (6):0 062902, 2016.
	
	[Feng(1985)]Feng1985
	S Feng.
	Percolation properties of granular elastic networks in two
	dimensions.
	Phys Rev B, 320 (1):0 510–513, 1985.
	
	[Moukarzel and Duxbury(1995)]Moukarzel:1995a
	C Moukarzel and P M Duxbury.
	Stressed backbone and elasticity of random central-force systems.
	Phys Rev Lett, 750 (22):0 4055–4058, 1995.
	
	[Jacobs and Thorpe(1995)]Jacobs1995
	D Jacobs and M F Thorpe.
	Generic rigidity percolation: The pebble game.
	Phys Rev Lett, 750 (22):0 4051–4054, 1995.
	
	[Aharonov and Sparks(1999)]Aharonov:1999a
	E Aharonov and D Sparks.
	Rigidity phase transition in granular packings.
	Phys Rev E, 600 (6):0 6890–6896, 1999.
	
	[Lois et al.(2008)Lois, Blawzdziewicz, and O'Hern]Lois:2008a
	G Lois, J Blawzdziewicz, and C S O'Hern.
	Jamming transition and new percolation universality classes in
	particulate systems with attraction.
	Phys Rev Lett, 1000 (2):0 028001, 2008.
	
	[Shen et al.(2012)Shen, O'Hern, and Shattuck]Shen:2012a
	T Shen, C S O'Hern, and M D Shattuck.
	Contact percolation transition in athermal particulate systems.
	Phys Rev E, 850 (1):0 011308, 2012.
	
	[Kovalcinova et al.(2015)Kovalcinova, Goullet, and
	Kondic]Kovalcinova:2015a
	L Kovalcinova, A Goullet, and L Kondic.
	Percolation and jamming transitions in particulate systems with and
	without cohesion.
	Phys Rev E, 920 (3):0 032204, 2015.
	
	[Henkes et al.(2016)Henkes, Quint, Fily, and Schwarz]Henkes2016
	S Henkes, D A Quint, Y Fily, and J M Schwarz.
	Rigid Cluster Decomposition Reveals Criticality in Frictional
		Jamming.
	Phys Rev Lett, 1160 (2):0 028301, 2016.
	
	[Thorpe(1985)]Thorpe1985a
	M F Thorpe.
	Rigidity Percolation.
	Springer-Verlag, 1985.
	
	[Thorpe and Duxbury(1999)]Thorpe:1999aa
	M F Thorpe and P M Duxbury, editors.
	Rigidity Theory and Applications.
	Springer-Verlag, 1999.
	
	[Kovalcinova et al.(2016)Kovalcinova, Goullet, and
	Kondic]Kovalcinova:2016a
	L Kovalcinova, A Goullet, and L Kondic.
	Scaling properties of force networks for compressed particulate
	systems.
	Phys Rev E, 930 (4):0 042903, 2016.
	
	[Pastor-Satorras and Miguel(2012)]PastorSatorras2012
	R Pastor-Satorras and M-C Miguel.
	Percolation analysis of force networks in anisotropic granular
	matter.
	Journal of Statistical Mechanics: Theory and Experiment,
	20120 (02):0 P02008, 2012.
	
	[Pathak et al.(2017)Pathak, Esposito, Coniglio, and
	Ciamarra]pathak2017force
	S N Pathak, V Esposito, A Coniglio, and M P Ciamarra.
	Force percolation transition of jammed granular systems.
	Phys Rev E, 960 (4):0 042901, 2017.
	
	[Edelsbrunner(2010)]Edelsbrunner:2010a
	H Edelsbrunner.
	Computational Topology: An Introduction.
	American Mathematical Society, 2010.
	
	[Carlsson(2009)]Carlsson:2009a
	G Carlsson.
	Topology and data.
	Bull Am Math Soc, 460 (2):0 255–308, 2009.
	
	[Ghrist(2014)]ghrist2014elementary
	R Ghrist.
	Elementary Applied Topology.
	Createspace, 2014.
	Available at <https://www.math.upenn.edu/ ghrist/notes.html>.
	
	[Dlotko et al.(2012)Dlotko, Juda, Mrozek, and
	Ghrist]dlotko2012distributed
	P Dlotko, M Juda, M Mrozek, and R Ghrist.
	Distributed computation of coverage in sensor networks by homological
	methods.
	Appl Algebr Eng Comm, 23:0 29–58, 2012.
	
	[Taylor et al.(2015)Taylor, Klimm, Harrington, Kramár, Mischaikow,
	Porter, and Mucha]taylor2015topological
	D Taylor, F Klimm, H A Harrington, M Kramár, K Mischaikow, M A Porter, and
	P J Mucha.
	Topological data analysis of contagion maps for examining spreading
	processes on networks.
	Nat Commun, 6:0 7723, 2015.
	
	[Sizemore et al.(2017a)Sizemore, Giusti, and
	Bassett]sizemore2016classification
	A Sizemore, C Giusti, and D S Bassett.
	Classification of weighted networks through mesoscale homological
	features.
	Journal of Complex Networks, 50 (2):0 245,
	2017a.
	
	[Sizemore et al.(2017b)Sizemore, Giusti, Kahn, Betzel, and
	Bassett]sizemore2016cliques
	A Sizemore, C Giusti, A E Kahn, R F Betzel, and D S Bassett.
	Cliques and cavities in the human connectome.
	Journal of Computational Neuroscience:
		https://doi.org/10.1007/s10827-017-0672-6, 2017b.
	
	[Otter et al.(2017)Otter, Porter, Tullmann, Grindrod, and
	Harrington]otter2015
	N Otter, M A Porter, U Tullmann, P Grindrod, and H A Harrington.
	A roadmap for the computation of persistent homology.
	EPJ Data Science, 60 (1):0 17, 2017.
	
	[Patania et al.(2017)Patania, Vaccarino, and Petri]Patania2017
	A Patania, F Vaccarino, and G Petri.
	Topological analysis of data.
	EPJ Data Science, 60 (1):0 7, 2017.
	
	[Stolz et al.(2017)Stolz, Harrington, and Porter]Stolz2016
	B J Stolz, H A Harrington, and M A Porter.
	Persistent homology of time-dependent functional networks constructed
	from coupled time series.
	Chaos, 270 (4):0 047410, 2017.
	
	[Kozlov(2007)]kozlov2007combinatorial
	D Kozlov.
	Combinatorial Algebraic Topology, volume 21.
	Springer-Verlag, 2007.
	
	[Nanda and Sazdanović(2014)]nanda2014simplicial
	V Nanda and R Sazdanović.
	Simplicial Models and Topological Inference in Biological
		Systems, pages 109–141.
	Springer-Verlag, 2014.
	
	[Petri et al.(2013)Petri, Scolamiero, Donato, and
	Vaccarino]Petri2013a
	G Petri, M Scolamiero, I Donato, and F Vaccarino.
	Topological strata of weighted complex networks.
	PLoS One, 80 (6):0 1–8, 2013.
	
	[Dantu(1957)]Dantu1957
	P Dantu.
	Contribution l'étude méchanique et
		géométrique des milieux pulvérulents.
	In Proceedings of the fourth International Conference on Soil
		Mechanics and Foundation Engineering, London, pages 144–148, 1957.
	
	[Drescher and de Josselin de Jong(1972)]Drescher1972
	A Drescher and G de Josselin de Jong.
	Photoelastic Verification Of A Mechanical Model For Flow Of A
		Granular Material.
	Journal Of The Mechanics And Physics Of Solids, 200
	(5):0 337–340, 1972.
	
	[Luding(1997)]luding1997stress
	S Luding.
	Stress distribution in static two-dimensional granular model media in
	the absence of friction.
	Phys Rev E, 550 (4):0 4720, 1997.
	
	[Silbert et al.(2002)Silbert, Grest, and Landry]silbert2002statistics
	L E Silbert, G S Grest, and J W Landry.
	Statistics of the contact network in frictional and frictionless
	granular packings.
	Phys Rev E, 660 (6):0 061303, 2002.
	
	[Tordesillas(2007)]Tordesillas:2007aa
	A Tordesillas.
	Force chain buckling, unjamming transitions and shear banding in
	dense granular assemblies.
	Philos Mag, 870 (32):0 4987–5016, 2007.
	
	[Majmudar et al.(2007)Majmudar, Sperl, Luding, and
	Behringer]Majmudar:2007aa
	T S Majmudar, M Sperl, S Luding, and R P Behringer.
	Jamming transition in granular systems.
	Phys Rev Lett, 980 (5):0 058001, 2007.
	
	[Snoeijer et al.(2004a)Snoeijer, Vlugt, van Hecke, and van
	Saarloos]Snoeijer-2004-FNE
	J H Snoeijer, T J H Vlugt, M van Hecke, and W van Saarloos.
	Force network ensemble: A new approach to static granular matter.
	Phys Rev Lett, 920 (5):0 54302,
	2004a.
	
	[Snoeijer et al.(2004b)Snoeijer, Vlugt, Ellenbroek,
	Van Hecke, and Van Leeuwen]Snoeijer-2004-ETF
	J H Snoeijer, T J H Vlugt, W G Ellenbroek, M Van Hecke, and J M J Van Leeuwen.
	Ensemble theory for force networks in hyperstatic granular matter.
	Phys Rev E, 700 (6):0 61306,
	2004b.
	
	[Tighe et al.(2010)Tighe, Snoeijer, Vlugt, and van Hecke]Tighe2010a
	B P Tighe, J H Snoeijer, T J H Vlugt, and M van Hecke.
	The force network ensemble for granular packings.
	Soft Matter, 60 (13):0 2908–2917, 2010.
	
	[Kollmer and Daniels(2017)]Kollmer2017
	J E Kollmer and K E Daniels.
	An experimental investigation of the force network ensemble.
	In Powders and Grains 2017, volume 140, page 02024, 2017.
	
	[Liu et al.(2011)Liu, Nagel, Van Saarloos, and Wyart]Liu:2011aa
	A J Liu, S R Nagel, W Van Saarloos, and M Wyart.
	The jamming scenario — An introduction and outlook.
	Oxford University Press, 2011.
	
	[Henkes et al.(2010)Henkes, van Hecke, and van
		Saarloos]Henkes:2010aa
	S Henkes, M van Hecke, and W van Saarloos.
	Critical jamming of frictional grains in the generalized isostaticity
	picture.
	EPL (Europhysics Letters), 900 (1):0 14003,
	2010.
	
	[Shundyak et al.(2007)Shundyak, van Hecke, and van
	Saarloos]Shundyak:2007aa
	K Shundyak, M van Hecke, and W van Saarloos.
	Force mobilization and generalized isostaticity in jammed packings of
	frictional grains.
	Phys Rev E, 750 (1):0 010301, 2007.
	
	[Stumpf and Porter(2012)]stumpf2012critical
	M P Stumpf and M A Porter.
	Mathematics. critical truths about power laws.
	Science, 3350 (6069):0 665–666, 2012.
	
	[Tordesillas et al.(2010b)Tordesillas, O'Sullivan, Walker,
	and Paramitha]Tordesillas:2010a
	A Tordesillas, P O'Sullivan, D M Walker, and Paramitha.
	Evolution of functional connectivity in contact and force chain
	networks: Feature vectors, k-cores and minimal cycles.
	Comptes Rendus Mécanique, 3380 (10):0
	556–569, 2010b.
	
	[Duxbury et al.(1999)Duxbury, Jacobs, Thorpe, and
	Moukarzel]Duxbury1999
	P Duxbury, D Jacobs, M Thorpe, and C Moukarzel.
	Floppy modes and the free energy: Rigidity and connectivity
		percolation on Bethe lattices.
	Phys Rev E, 590 (2):0 2084–2092, 1999.
	
	[Maxwell(1864)]Maxwell:1864aa
	J C Maxwell.
	On the calculation of the equilibrium and stiffness of frames.
	Philosophical Magazine Series 4, 270 (182):0
	294–299, 1864.
	
	[Laman(1970)]Laman:1970aa
	G Laman.
	On graphs and rigidity of plane skeletal structures.
	J Eng Math, 40 (4):0 331–340, 1970.
	
	[Asimow and Roth(1978)]Asimow:1978aa
	L Asimow and B Roth.
	The rigidity of graphs.
	Trans Am Math Soc, 245:0 279–289, 1978.
	
	[Crapo(1979)]Crapo:1979aa
	H Crapo.
	Structural rigidity.
	Structural Topology, 1:0 26–45, 1979.
	
	[Guyon et al.(1990)Guyon, Roux, Hansen, Bideau, Troadec, and
	Crapo]Guyon:1990aa
	E Guyon, S Roux, A Hansen, D Bideau, J-P Troadec, and H Crapo.
	Non-local and non-linear problems in the mechanics of disordered
	systems: Application to granular media and rigidity problems.
	Rep Prog Phys, 530 (4):0 373, 1990.
	
	[Moukarzel(1998)]Moukarzel:1998aa
	C F Moukarzel.
	Isostatic phase transition and instability in stiff granular
	materials.
	Phys Rev Lett, 810 (8):0 1634–1637, 1998.
	
	[Tordesillas et al.(2011)Tordesillas, Lin, Zhang, Behringer, and
	Shi]Tordesillas:2011aa
	A Tordesillas, Q Lin, J Zhang, R P Behringer, and J Shi.
	Structural stability and jamming of self-organized cluster
	conformations in dense granular materials.
	J Mech Phys Solids, 590 (2):0 265–296, 2011.
	
	[Tordesillas et al.(2012b)Tordesillas, Pucilowski, Walker,
	Peters, and Hopkins]Tordesillas:2010b
	A Tordesillas, S Pucilowski, D M Walker, J Peters, and M Hopkins.
	A complex network analysis of granular fabric evolution in
	three-dimensions.
	Dynam Cont Dis Ser B, 190 (4–5):0 417–495,
	2012b.
	
	[Walker et al.(2014c)Walker, Tordesillas, Ren, Dijksman,
	and Behringer]Walker:2014c
	D M Walker, A Tordesillas, J Ren, J A Dijksman, and R P Behringer.
	Uncovering temporal transitions and self-organization during slow
	aging of dense granular media in the absence of shear bands.
	EPL (Europhysics Letters), 1070 (1):0 18005,
	2014c.
	
	[Jeng and Schwarz(2008)]Jeng2008
	M Jeng and J M Schwarz.
	On the Study of Jamming Percolation.
	J Stat Phys, 1310 (4):0 575–595, 2008.
	
	[Jeng and Schwarz(2010)]Jeng2010
	M Jeng and J M Schwarz.
	Force-balance percolation.
	Phys Rev E, 810 (1):0 011134, 2010.
	
	[Cao and Schwarz(2012)]Cao2012
	L Cao and J M Schwarz.
	Correlated percolation and tricriticality.
	Phys Rev E, 860 (6):0 061131, 2012.
	
	[Lopez et al.(2013)Lopez, Cao, and Schwarz]Lopez:2013aa
	J H Lopez, L Cao, and J M Schwarz.
	Jamming graphs: A local approach to global mechanical rigidity.
	Phys Rev E, 880 (6):0 062130, 2013.
	
	[Heroy et al.(2017)Heroy, Taylor, Shi, Forest, and Mucha]heroy2017
	S Heroy, D Taylor, F B Shi, M G Forest, and P J Mucha.
	Rigid graph compression: Motif-based rigidity analysis for
	disordered fiber networks.
	arXiv, arXiv:1711.05790 [cond-mat.dis-nn], 2017.
	
	[Oda and Kazama(1998)]Oda:1998a
	M Oda and H Kazama.
	Microstructure of shear bands and its relation to the mechanisms of
	dilatancy and failure of dense granular soils.
	Géotechnique, 480 (4):0 465–481, 1998.
	
	[Tordesillas et al.(2009)Tordesillas, Zhang, and
	Behringer]Tordesillas:2009b
	A Tordesillas, J Zhang, and R Behringer.
	Buckling force chains in dense granular assemblies: Physical and
	numerical experiments.
	Geomech Geoeng, 40 (1):0 3–16, 2009.
	
	[Bagi(2007)]Bagi:2007a
	K Bagi.
	On the concept of jammed configurations from a structural mechanics
	perspective.
	Granular Matter, 90 (1):0 109–134, 2007.
	
	[Tordesillas and Muthuswamy(2009)]Tordesillas:2009aa
	A Tordesillas and M Muthuswamy.
	On the modeling of confined buckling of force chains.
	J Mech Phys Solids, 570 (4):0 706–727, 2009.
	
	[Cates et al.(1998)Cates, Wittmer, Bouchaud, and Claudin]Cates:1998aa
	M E Cates, J P Wittmer, J-P Bouchaud, and P Claudin.
	Jamming, force chains, and fragile matter.
	Phys Rev Lett, 810 (9):0 1841–1844, 1998.
	
	[Muthuswamy and Tordesillas(2006)]Muthuswamy2006a
	M Muthuswamy and A Tordesillas.
	How do interparticle contact friction, packing density and degree of
	polydispersity affect force propagation in particulate assemblies?
	Journal of Statistical Mechanics: Theory and Experiment,
	20060 (09):0 P09003, 2006.
	
	[Kob and Barrat(1997)]kob1997aging
	W Kob and J L Barrat.
	Aging effects in a Lennard-Jones glass.
	Phys Rev Lett, 780 (24):0 4581–4584, 1997.
	
	[Kabla and Debregeas(2004)]kabla2004contact
	A Kabla and G Debregeas.
	Contact dynamics in a gently vibrated granular pile.
	Phys Rev Lett, 920 (3):0 35501, 2004.
	
	[Steinhardt et al.(1983)Steinhardt, Nelson, and
	Ronchetti]Steinhardt:1983a
	P J Steinhardt, D R Nelson, and M Ronchetti.
	Bond-orientational order in liquids and glasses.
	Phys Rev B, 280 (2):0 784–805, 1983.
	
	[Arévalo et al.(2013)Arévalo, Pugnaloni, Zuriguel, and
	Maza]Arevalo:2013aa
	R Arévalo, L A Pugnaloni, I Zuriguel, and D Maza.
	Contact network topology in tapped granular media.
	Phys Rev E, 870 (2):0 022203, 2013.
	
	[Nowak et al.(1998)Nowak, Knight, Ben-Naim, Jaeger, and
	Nagel]Nowak:1998a
	E R Nowak, J B Knight, E Ben-Naim, H M Jaeger, and S R Nagel.
	Density fluctuations in vibrated granular materials.
	Phys Rev E, 570 (2):0 1971–1982, 1998.
	
	[Pugnaloni et al.(2010)Pugnaloni, Sánchez, Gago, Damas, Zuriguel, and
	Maza]Pugnaloni:2010aa
	L A Pugnaloni, I Sánchez, P A Gago, J Damas, I Zuriguel, and D Maza.
	Towards a relevant set of state variables to describe static granular
	packings.
	Phys Rev E, 820 (5):0 050301, 2010.
	
	[Pugnaloni et al.(2008)Pugnaloni, Mizrahi, Carlevaro, and
	Vericat]Pugnaloni:2008aa
	L A Pugnaloni, M Mizrahi, C M Carlevaro, and F Vericat.
	Nonmonotonic reversible branch in four model granular beds subjected
	to vertical vibration.
	Phys Rev E, 780 (5):0 051305, 2008.
	
	[Gago et al.(2009)Gago, Bueno, and Pugnaloni]Gago:2009aa
	P A Gago, N E Bueno, and L A Pugnaloni.
	High intensity tapping regime in a frustrated lattice gas model of
	granular compaction.
	Granular Matter, 110 (6):0 365–369, 2009.
	
	[Carlevaro and Pugnaloni(2011)]Carlevaro:2011a
	C M Carlevaro and L A Pugnaloni.
	Steady state of tapped granular polygons.
	Journal of Statistical Mechanics: Theory and Experiment,
	20110 (01):0 P01007, 2011.
	
	[Arévalo et al.(2013)Arévalo, Pugnaloni, Maza, and
	Zuriguel]Arevalo:2013ba
	R Arévalo, L A Pugnaloni, D Maza, and I Zuriguel.
	Tapped granular packings described as complex networks.
	Philosophical Magazine, 930 (31-33):0
	4078–4089, 2013.
	
	[Itzkovitz and Alon(2005)]itzkovitz2005subgraphs
	S Itzkovitz and U Alon.
	Subgraphs and network motifs in geometric networks.
	Phys Rev E, 710 (2):0 026117, 2005.
	
	[Shoval and Alon(2010)]shoval2010snap
	O Shoval and U Alon.
	SnapShot: Network motifs.
	Cell, 1430 (2):0 326–326.e1, 2010.
	
	[Brodu et al.(2015)Brodu, Dijksman, and Behringer]Brodu2015
	N Brodu, J A Dijksman, and R P Behringer.
	Spanning the scales of granular materials through microscopic force
	imaging.
	Nat Commun, 6:0 6361, 2015.
	
	[Dijksman et al.(2017)Dijksman, Brodu, and Behringer]Dijksman2017
	J A Dijksman, N Brodu, and R P Behringer.
	Refractive index matched scanning and detection of soft particles.
	Rev Sci Instrum, 880 (5):0 051807, 2017.
	
	[Sepiani and Ghazavi(2009)]sepiani2009thermo
	H A Sepiani and A Ghazavi.
	A thermo-micro-mechanical modeling for smart shape memory alloy woven
	composite under in-plane biaxial deformation.
	International Journal of Mechanics and Materials in Design,
	50 (2):0 111, 2009.
	
	[Tighe and Vlugt(2011)]Tighe:2011a
	B P Tighe and T J H Vlugt.
	Stress fluctuations in granular force networks.
	Journal of Statistical Mechanics: Theory and Experiment,
	20110 (04):0 P04002, 2011.
	
	[Daniels et al.(2017)Daniels, Kollmer, and Puckett]Daniels2017
	K E Daniels, J E Kollmer, and J G Puckett.
	Photoelastic force measurements in granular materials.
	Rev Sci Instrum, 880 (5):0 051808, 2017.
	
	[Hurley et al.(2016)Hurley, Hall, Andrade, and Wright]Hurley:2016a
	R C Hurley, S A Hall, J E Andrade, and J Wright.
	Quantifying interparticle forces and heterogeneity in 3d granular
	materials.
	Phys Rev Lett, 1170 (9):0 098005, 2016.
	
	[Mukhopadhyay and Peixinho(2011)]Mukhopadhyay2011
	S Mukhopadhyay and J Peixinho.
	Packings of deformable spheres.
	Phys Rev E, 840 (1):0 011302, 2011.
	
	[Saadatfar et al.(2012)Saadatfar, Sheppard, Senden, and
	Kabla]Saadatfar2012
	M Saadatfar, A P Sheppard, T J Senden, and A J Kabla.
	Mapping forces in a 3D elastic assembly of grains.
	J Mech Phys Solids, 600 (1):0 55–66, 2012.
	
	[Pöschel and Schwager(2005)]Poschel2005
	T Pöschel and T Schwager.
	Computational Granular Dynamics: Models and Algorithms.
	Springer-Verlag, 2005.
	
	[Weis and Schröter(2017)]weis:17
	S Weis and M Schröter.
	Analyzing X-ray tomographies of granular packings.
	Rev Sci Instrum, 880 (5):0 051809, 2017.
	
	[Tordesillas and Muthuswamy(2008)]Tordesillas:2008a
	A Tordesillas and M Muthuswamy.
	A thermomicromechanical approach to multiscale continuum modeling of
	dense granular materials.
	Acta Geotechnica, 30 (3):0 225–240, 2008.
	
	[Huang and Daniels(2016)]Huang:2016a
	Y Huang and K E Daniels.
	Friction and pressure-dependence of force chain communities in
	granular materials.
	Granular Matter, 180 (4):0 85, 2016.
	
	[Navakas et al.(2014)Navakas, Džiugys, and Peters]Navakas:2014aa
	R Navakas, A Džiugys, and B Peters.
	A community-detection based approach to identification of
	inhomogeneities in granular matter.
	Physica A, 407:0 312–331, 2014.
	
	[Walker et al.(2011)Walker, Tordesillas, Einav, and
	Small]Walker:2011a
	D M Walker, A Tordesillas, I Einav, and M Small.
	Complex networks in confined comminution.
	Phys Rev E, 840 (2):0 021301, 2011.
	
	[Radjai et al.(1999)Radjai, Roux, and Moreau]Radjai:1999aa
	F Radjai, S Roux, and J J Moreau.
	Contact forces in a granular packing.
	Chaos, 90 (4):0 544–550, 1999.
	
	[Peña et al.(2009)Peña, Herrmann, and Lind]Pena:2009aa
	A A Peña, H J Herrmann, and P G Lind.
	Force chains in sheared granular media of irregular particles.
	AIP Conference Proceedings, 11450 (1):0
	321–324, 2009.
	
	[Ostojic et al.(2007)Ostojic, Vlugt, and Nienhuis]Ostojic:2007a
	S Ostojic, T J H Vlugt, and B Nienhuis.
	Universal anisotropy in force networks under shear.
	Phys Rev E, 750 (3):0 030301, 2007.
	
	[Kondic et al.(2017)Kondic, Kramár, Kovalčinová, and
	Mischaikow]Kondic:2017a
	L Kondic, M Kramár, L Kovalčinová, and K Mischaikow.
	Evolution of force networks in dense granular matter close to
	jamming.
	EPJ Web Conf, 140:0 15014, 2017.
	
	[Tadanaga et al.(2017)Tadanaga, Clark, Majmudar, and
	Kondic]kondic-impact2017
	T Tadanaga, A H Clark, T Majmudar, and L Kondic.
	Granular response to impact: Topology of the force networks.
	arXiv, arXiv:1709.06957 [cond-mat.soft], 2017.
	
	[Lim and Behringer(2017)]lim-impact2017
	M X Lim and R P Behringer.
	Topology of force networks in granular media under impact.
	arXiv, arXiv:1709.01884 [cond-mat.soft], 2017.
	
	[Tordesillas et al.(2015a)Tordesillas, Tobin, Cil,
	Alshibli, and Behringer]Tordesillas:2015a
	A Tordesillas, S T Tobin, M Cil, K Alshibli, and R P Behringer.
	Network flow model of force transmission in unbonded and bonded
	granular media.
	Phys Rev E, 910 (6):0 062204,
	2015a.
	
	[Tordesillas et al.(2015b)Tordesillas, Pucilowski, Tobin,
	Kuhn, Andò, Viggiani, Druckrey, and Alshibli]Tordesillas:2015b
	A Tordesillas, S Pucilowski, S Tobin, M R Kuhn, E Andò, G Viggiani,
	A Druckrey, and K Alshibli.
	Shear bands as bottlenecks in force transmission.
	EPL (Europhysics Letters), 1100 (5):0 58005,
	2015b.
	
	[Tordesillas et al.(2013b)Tordesillas, Cramer, and
	Walker]Tordesillas:2013b
	A Tordesillas, A Cramer, and D M Walker.
	Minimum cut and shear bands.
	AIP Conference Proceedings, 15420 (1):0
	507–510, 2013b.
	
	[Lin and Tordesillas(2013)]Lin:2013a
	Q Lin and A Tordesillas.
	Constrained optimisation in granular network flows: Games with a
	loaded dice.
	AIP Conference Proceedings, 15420 (1):0
	547–550, 2013.
	
	[Lin and Tordesillas(2014)]Lin:2014a
	Q Lin and A Tordesillas.
	Towards an optimization theory for deforming dense granular
	materials: Minimum cost maximum flow solutions.
	J Ind Manag Optim, 100 (1):0 337–362, 2014.
	
	[Walker and Tordesillas(2013)]Walker:2013b
	D M Walker and A Tordesillas.
	Understanding multi-scale structural evolution in granular systems
	through gmems.
	AIP Conference Proceedings, 15420 (1):0
	145–148, 2013.
	
	[Zhang and Small(2006)]Zhang:2006aa
	J Zhang and M Small.
	Complex network from pseudoperiodic time series: Topology versus
	dynamics.
	Phys Rev Lett, 960 (23):0 238701, 2006.
	
	[Yang and Yang(2008)]Yang:2008a
	Y Yang and H Yang.
	Complex network-based time series analysis.
	Physica A, 3870 (5–6):0 1381–1386, 2008.
	
	[Lacasa et al.(2008)Lacasa, Luque, Ballesteros, Luque, and
	Nuño]Lacasa:2008a
	L Lacasa, B Luque, F Ballesteros, J Luque, and J C Nuño.
	From time series to complex networks: The visibility graph.
	Proc Natl Acad Sci, 1050 (13):0 4972–4975,
	2008.
	
	[Gao and Jin(2009)]Gao:2009a
	Z Gao and N Jin.
	Complex network from time series based on phase space reconstruction.
	Chaos, 190 (3):0 033137, 2009.
	
	[Marwan et al.(2009)Marwan, Donges, Zou, Donner, and
	Kurths]Marwan:2009a
	N Marwan, J F Donges, Y Zou, R V Donner, and J Kurths.
	Complex network approach for recurrence analysis of time series.
	Phys Lett A, 3730 (46):0 4246–4254, 2009.
	
	[MacKay(2003)]MacKay:2003
	D J C MacKay.
	Information Theory, Inference and Learning Algorithms.
	Cambridge University Press, 2003.
	
	[Rechenmacher(2006)]Rechenmacher:2006a
	A L Rechenmacher.
	Grain-scale processes governing shear band initiation and evolution
	in sands.
	J Mech Phys Solids, 540 (1):0 22–45, 2006.
	
	[Rechenmacher et al.(2010)Rechenmacher, Abedi, and
	Chupin]Rechenmacher:2010a
	A L Rechenmacher, S Abedi, and O Chupin.
	Evolution of force chains in shear bands in sands.
	Géotechnique, 600 (5):0 343–351, 2010.
	
	[Andò et al.(2012)Andò, Hall, Viggiani, Desrues, and
	Bésuelle]Ando:2012a
	E Andò, S A Hall, G Viggiani, J Desrues, and P Bésuelle.
	Grain-scale experimental investigation of localised deformation in
	sand: A discrete particle tracking approach.
	Acta Geotechnica, 70 (1):0 1–13, 2012.
	
	[Amon et al.(2017)Amon, Born, Daniels, Dijksman, Huang, Parker,
	Schröter, Stannarius, and Wierschem]amon_focus_2017
	A Amon, P Born, K E Daniels, J A Dijksman, K Huang, D Parker, M Schröter,
	R Stannarius, and A Wierschem.
	Preface: Focus on imaging methods in granular physics.
	Rev Sci Instrum, 880 (5):0 051701, 2017.
	
	[Stannarius(2017)]stannarius:17
	R Stannarius.
	Magnetic resonance imaging of granular materials.
	Rev Sci Instrum, 880 (5):0 051806, 2017.
	
	[Porter et al.(2015)Porter, Kevrekidis, and Daraio]PT2015
	M A Porter, P G Kevrekidis, and C Daraio.
	Granular crystals: Nonlinear dynamics meets materials engineering.
	Physics Today, 680 (11):0 44, 2015.
	
	[Cundall and Strack(1979)]Cundall1979
	P A Cundall and O D L Strack.
	Discrete Numerical-model For Granular Assemblies.
	Geotechnique, 290 (1):0 47–65, 1979.
	
	[Papanikolaou et al.(2013)Papanikolaou, O'Hern, and
	Shattuck]Papanikolaou2013
	S Papanikolaou, C S O'Hern, and M D Shattuck.
	Isostaticity at frictional jamming.
	Phys Rev Lett, 1100 (19):0 198002, 2013.
	
	[Somfai et al.(2005)Somfai, Roux, Snoeijer, van Hecke, and van
	Saarloos]Somfai2005
	E Somfai, J N Roux, J H Snoeijer, M van Hecke, and W van Saarloos.
	Elastic wave propagation in confined granular systems.
	Phys Rev E, 720 (2):0 21301, 2005.
	
	[Rosvall et al.(2014)Rosvall, Esquivel, Lancichinetti, West, and
	Lambiotte]Rosvall2014NatComm
	M Rosvall, A V Esquivel, A Lancichinetti, J D West, and R Lambiotte.
	Memory in network flows and its effects on spreading dynamics and
		community detection.
	Nat Commun, 5:0 4630, 2014.
	
	[Holme(2015)]Holme2015EurPhysJB
	P Holme.
	Modern temporal network theory: A colloquium.
	Eur Phys J B, 880 (9):0 234, 2015.
	
	[Butts(2009)]butts2009
	C T Butts.
	Revisiting the foundations of network analysis.
	Science, 3250 (5939):0 414–6, 2009.
	
	[Gravish et al.(2012)Gravish, Franklin, Hu, and
	Goldman]gravish2012entangled
	N Gravish, S V Franklin, D L Hu, and D I Goldman.
	Entangled granular media.
	Phys Rev Lett, 1080 (20):0 208001, 2012.
	
	[Murphy et al.(2016)Murphy, Reiser, Choksy, Singer, and
	Jaeger]murphy2015freestanding
	K A Murphy, N Reiser, D Choksy, C E Singer, and H M Jaeger.
	Freestanding loadbearing structures with Z-shaped particles.
	Granular Matter, 180 (2):0 26, 2016.
	
	[Hidalgo et al.(2009)Hidalgo, Zuriguel, Maza, and
	Pagonabarraga]Hidalgo:2009aa
	R Cruz Hidalgo, I Zuriguel, D Maza, and I Pagonabarraga.
	Role of particle shape on the stress propagation in granular
	packings.
	Phys. Rev. Lett., 1030 (11):0 118001, 2009.
	
	[Trepanier and Franklin(2010)]trepanier2010column
	M Trepanier and S V Franklin.
	Column collapse of granular rods.
	Phys Rev E, 820 (1):0 011308, 2010.
	
	[Schreck et al.(2010)Schreck, Xu, and O'Hern]Schreck:2010a
	C F Schreck, N Xu, and C S O'Hern.
	A comparison of jamming behavior in systems composed of dimer- and
	ellipse-shaped particles.
	Soft Matter, 60 (13):0 2960–2969, 2010.
	
	[Athanassiadis et al.(2014)Athanassiadis, Miskin, Kaplan, Rodenberg,
	Lee, Merritt, Brown, Amend, Lipson, and Jaeger]Athanassiadis:2014a
	A G. Athanassiadis, M Z. Miskin, P Kaplan, N Rodenberg, S H Lee, J Merritt,
	E Brown, J Amend, H Lipson, and H M Jaeger.
	Particle shape effects on the stress response of granular packings.
	Soft Matter, 100 (1):0 48–59, 2014.
	
	[Harrington and Durian(2017)]durian2017
	M Harrington and D J Durian.
	Anisotropic particles strengthen granular pillars under compression.
	arXiv, arXiv:1709.09511 [cond-mat.soft], 2017.
	
	[Azéma et al.(2013)Azéma, Radjai, and Dubois]Azema:2013aa
	E Azéma, F Radjai, and F Dubois.
	Packings of irregular polyhedral particles: Strength, structure,
	and effects of angularity.
	Phys Rev E, 870 (6):0 062203, 2013.
	
	[Gómez et al.(2009)Gómez, Jensen, and Arenas]Gomez:2009a
	S Gómez, P Jensen, and A Arenas.
	Analysis of community structure in networks of correlated data.
	Phys Rev E, 800 (1):0 016114, 2009.
	
	[Traag and Bruggeman(2009)]Traag:2009a
	V A Traag and J Bruggeman.
	Community detection in networks with positive and negative links.
	Phys Rev E, 800 (3):0 036115, 2009.
	
	[Zhang et al.(2016)Zhang, Bassett, and Winkelstein]Zhang:2015
	S Zhang, D S Bassett, and B A Winkelstein.
	Stretch-induced network reconfiguration of collagen fibers in the
	human facet capsular ligament.
	J R Soc Interface, 130 (114):0 20150883, 2016.
	
	[Puckett and Daniels(2013)]Puckett2013
	J G Puckett and K E Daniels.
	Equilibrating temperaturelike variables in jammed granular
	subsystems.
	Phys Rev Lett, 1100 (5):0 058001, 2013.
	
	[Shaebani et al.(2012)Shaebani, Madadi, Luding, and
	Wolf]shaebani2012influence
	M R Shaebani, M Madadi, S Luding, and D E Wolf.
	Influence of polydispersity on micromechanics of granular materials.
	Phys Rev E, 850 (1):0 011301, 2012.
	
	[Kumar et al.(2016)Kumar, Magnanimo, Ramaioli, and Luding]Kumar:2016a
	N Kumar, V Magnanimo, M Ramaioli, and S Luding.
	Tuning the bulk properties of bidisperse granular mixtures by small
	amount of fines.
	Powder Technol, 293:0 94–112, 2016.
	
	[Slanina(2017)]slanina2017
	F Slanina.
	Localization in random bipartite graphs: Numerical and empirical
	study.
	Phys Rev E, 950 (5):0 052149, 2017.
	
	[Harary(1972)]harary1972
	F Harary.
	Graph Theory.
	Addison-Wesley, 1972.
	
	[Shi et al.(2013)Shi, Wang, Forest, and Mucha]Shi:2013aa
	F Shi, S Wang, M G Forest, and P J Mucha.
	Percolation-induced exponential scaling in the large current tails of
	random resistor networks.
	Multiscale Model Sim, 110 (4):0 1298–1310,
	2013.
	
	[Shi et al.(2014)Shi, Wang, Forest, Mucha, and Zhou]Shi:2014aa
	F Shi, S Wang, M G Forest, P J Mucha, and R Zhou.
	Network-based assessments of percolation-induced current
	distributions in sheared rod macromolecular dispersions.
	Multiscale Model Sim, 120 (1):0 249–264,
	2014.
	
	[Abhilash et al.(2014)Abhilash, Baker, Trappmann, Chen, and
	Shenoy]abhilash2014remodeling
	A S Abhilash, B M Baker, B Trappmann, C S Chen, and V B Shenoy.
	Remodeling of fibrous extracellular matrices by contractile cells:
	Predictions from discrete fiber network simulations.
	Biophys J, 1070 (8):0 1829–1840, 2014.
	
	[Purohit et al.(2011)Purohit, Litvinov, Brown, Discher, and
	Weisel]purohit2011protein
	P K Purohit, R I Litvinov, A E Brown, D E Discher, and J W Weisel.
	Protein unfolding accounts for the unusual mechanical behavior of
	fibrin networks.
	Acta Biomater, 70 (6):0 2374–2783, 2011.
	
	[Bullmore and Sporns(2009)]Bullmore2009
	E Bullmore and O Sporns.
	Complex brain networks: Graph theoretical analysis of structural
	and functional systems.
	Nat Rev Neurosci, 100 (3):0 186–198, 2009.
	
	[Blunt(2001)]Blunt:2001a
	M J Blunt.
	Flow in porous media — pore-network models and multiphase flow.
	Curr Opin Colloid Interface Sci, 60 (3):0
	197–207, 2001.
	
	[Al-Raoush et al.(2003)Al-Raoush, Thompson, and
	Willson]Al-Raoush:2003a
	R Al-Raoush, K Thompson, and C S Willson.
	Comparison of network generation techniques for unconsolidated porous
	media.
	Soil Sci Soc Am J, 670 (6):0 1687–1700, 2003.
	
	[Vo et al.(2013)Vo, Walker, and Tordesillas]Vo:2013
	K Vo, D M Walker, and A Tordesillas.
	Transport pathways within percolating pore space networks of granular
	materials.
	AIP Conf Proc, 15420 (1):0 551–554, 2013.
	
	[Walker et al.(2013)Walker, Vo, and Tordesillas]Walker:2013a
	D M Walker, K Vo, and A Tordesillas.
	On reynolds' dilatancy and shear band evolution: A new perspective.
	Int J Bifurc Chaos, 230 (09):0 1330034, 2013.
	
	[van der Linden et al.(2016)van der Linden, Narsilio, and
	Tordesillas]Linden:2016a
	J H van der Linden, G A Narsilio, and A Tordesillas.
	Machine learning framework for analysis of transport through complex
	networks in porous, granular media: A focus on permeability.
	Phys Rev E, 940 (2):0 022904, 2016.
	
	[Russell et al.(2016)Russell, Walker, and Tordesillas]Russell:2016a
	S Russell, D M Walker, and A Tordesillas.
	A characterization of the coupled evolution of grain fabric and pore
	space using complex networks: Pore connectivity and optimized flows in the
	presence of shear bands.
	J Mech Phys Solids, 88:0 227–251, 2016.
	
	[Jimenez-Martinez and Negre(2017)]Martinez:2017a
	J Jimenez-Martinez and C F A Negre.
	Eigenvector centrality for geometric and topological characterization
	of porous media.
	Phys Rev E, 960 (1):0 013310, 2017.
	
	[Laubie et al.(2017)Laubie, Radjai, Pellenq, and Ulm]porous2017
	H Laubie, F Radjai, R Pellenq, and F J Ulm.
	Stress transmission and failure in disordered porous media.
	Phys Rev Lett, 1190 (7):0 075501, 2017.
	
	[Newman and Clauset(2016)]newman2016structure
	M E J Newman and A Clauset.
	Structure and inference in annotated networks.
	Nat Commun, 7:0 11863, 2016.
	
	[Hric et al.(2016)Hric, Peixoto, and Fortunato]tiago2016
	D Hric, T P Peixoto, and S Fortunato.
	Network structure, metadata, and the prediction of missing nodes and
	annotations.
	Phys Rev X, 60 (3):0 031038, 2016.
	
	[Palla et al.(2008)Palla, Farkas, Pollner, Derenyi, and
	Vicsek]palla2008fundamental
	G Palla, I J Farkas, P Pollner, I Derenyi, and T Vicsek.
	Fundamental statistical features and self-similar properties of
	tagged networks.
	New J Phys, 100 (12):0 123026, 2008.
	
	[Edelsbrunner and Harer(2010)]edels2010
	H Edelsbrunner and J Harer.
	Computational Topology: An Introduction.
	American Mathematical Society, 2010.
	
	[Ramola and Chakraborty(2017)]Ramola:2017a
	K Ramola and B Chakraborty.
	Stress response of granular systems.
	J Stat Phys, 1690 (1):0 1–17, 2017.
	
	[Taylor-King et al.(2017)Taylor-King, Basanta, Chapman, and
	Porter]jkt2017
	J P Taylor-King, D Basanta, S J Chapman, and M A Porter.
	Mean-field approach to evolving spatial networks, with an application
	to osteocyte network formation.
	Phys Rev E, 960 (1):0 012301, 2017.
	
	[Beguerisse-Diaz et al.(2014)Beguerisse-Diaz, Garduno-Hernandez,
	Vangelov, Yaliraki, and Barahona]beguerisse2014interest
	M Beguerisse-Diaz, G Garduno-Hernandez, B Vangelov, S N Yaliraki, and
	M Barahona.
	Interest communities and flow roles in directed networks: The
	Twitter network of the UK riots.
	J R Soc Interface, 110 (101):0 20140940, 2014.
	
	[Bassett et al.(2010)Bassett, Greenfield, Meyer-Lindenberg, Weinberger,
	Moore, and Bullmore]Bassett2010
	D S Bassett, D L Greenfield, A Meyer-Lindenberg, D R Weinberger, S Moore, and
	E Bullmore.
	Efficient physical embedding of topologically complex information
	processing networks in brains and computer circuits.
	PLoS Comput Biol, 60 (4):0 e1000748, 2010.
	
	[Modes et al.(2016)Modes, Magnasco, and Katifori]Modes:2016a
	C D Modes, M O Magnasco, and E Katifori.
	Extracting hidden hierarchies in 3D distribution networks.
	Phys Rev X, 60 (3):0 031009, 2016.
	
	[Tighe et al.(2008)Tighe, van Eerd, and Vlugt]Tighe2008
	B P Tighe, A R T van Eerd, and T J H Vlugt.
	Entropy maximization in the force network ensemble for granular
	solids.
	Phys Rev Lett, 1000 (23):0 238001, 2008.
	
	[Ronhovde et al.(2012)Ronhovde, Chakrabarty, Sahu, Sahu, Kelton, Mauro,
	and Nussinov]ronhovde2011detection
	P Ronhovde, S Chakrabarty, M Sahu, K K Sahu, K F Kelton, N Mauro, and
	Z Nussinov.
	Detection of hidden structures for arbitrary scales in complex
	physical systems.
	Sci Rep, 2:0 329, 2012.
	
	[Agarwala and Shenoy(2017)]amorphous2017
	A Agarwala and V B Shenoy.
	Topological insulators in amorphous systems.
	Phys Rev Lett, 1180 (23):0 236402, 2017.
	
	[Muller and Luding(2009)]muller2009homogeneous
	M K Muller and S Luding.
	Homogeneous cooling with repulsive and attractive long‚Äêrange
	interactions.
	AIP Conf Proc, 11450 (1):0 697–700, 2009.
	
	[Mitarai and Nori(2006)]Mitarai:2006a
	N Mitarai and F Nori.
	Wet granular materials.
	Adv Phys, 550 (1-2):0 1–45, 2006.
	
	[Wrobel and Sundararaghavan(2014)]wrobel2014directed
	M R Wrobel and H G Sundararaghavan.
	Directed migration in neural tissue engineering.
	Tissue Eng Part B Rev, 200 (2):0 93–105,
	2014.
	
	[Huttenlocher and Poznansky(2008)]huttenlocher2008reverse
	A Huttenlocher and M C Poznansky.
	Reverse leukocyte migration can be attractive or repulsive.
	Trends in Cell Biol, 180 (6):0 298–306, 2008.
	
	[Hartveit and Veruki(2012)]hartveit2012electrical
	E Hartveit and M L Veruki.
	Electrical synapses between aii amacrine cells in the retina:
	Function and modulation.
	Brain Res, 1487:0 160–172, 2012.
	
	[Nualart-Marti et al.(2013)Nualart-Marti, Solsona, and
	Fields]nualart2013biochim
	A Nualart-Marti, C Solsona, and R D Fields.
	Gap junction communication in myelinating glia.
	Biochim Biophys Acta, 18280 (1):0 69–78,
	2013.
	
	[Pahtz et al.(2010)Pahtz, Herrmann, and Shinbrot]Pahtz:2010a
	T Pahtz, H J Herrmann, and T Shinbrot.
	Why do particle clouds generate electric charges?
	Nat Phys, 60 (5):0 364–368, 2010.
	
	[Ladoux et al.(2015)Ladoux, Nelson, Yan, and
	Mege]ladoux2015mechanotransduction
	B Ladoux, W J Nelson, J Yan, and R M Mege.
	The mechanotransduction machinery at work at adherens junctions.
	Integr Biol (Camb), 70 (10):0 1109–1119,
	2015.
	
	[Bausch and Kroy(2006)]Bausch:2006aa
	A R Bausch and K Kroy.
	A bottom-up approach to cell mechanics.
	Nat Phys, 20 (4):0 231–238, 2006.
	
	[Broedersz and MacKintosh(2014)]Broedersz:2014aa
	C P Broedersz and F C MacKintosh.
	Modeling semiflexible polymer networks.
	Rev Mod Phys, 860 (3):0 995–1036, 2014.
	
	[Lieleg et al.(2009)Lieleg, Schmoller, Claessens, and
	Bausch]Lieleg:2009a
	O Lieleg, K M Schmoller, M M Claessens, and A R Bausch.
	Cytoskeletal polymer networks: Viscoelastic properties are
	determined by the microscopic interaction potential of cross-links.
	Biophys J, 960 (11):0 4725–4732, 2009.
	
	[Fletcher and Mullins(2010)]Fletcher:2010aa
	D A Fletcher and R D Mullins.
	Cell mechanics and the cytoskeleton.
	Nature, 4630 (7280):0 485–492, 2010.
	
	[Mizuno et al.(2007)Mizuno, Tardin, Schmidt, and
	MacKintosh]Mizuno2007
	D Mizuno, C Tardin, C F Schmidt, and F C MacKintosh.
	Nonequilibrium mechanics of active cytoskeletal networks.
	Science, 3150 (5810):0 370–373, 2007.
	
	[Gardel et al.(2008)Gardel, Kasza, Brangwynne, Liu, and
	Weitz]Gardel:2008aa
	M L Gardel, K E Kasza, C P Brangwynne, J Liu, and D A Weitz.
	Mechanical response of cytoskeletal networks.
	Methods in Cell Biol, 89:0 487–519, 2008.
	
	[Majumdar et al.(2017)Majumdar, Foucard, Levine, and Gardel]actin2017
	S Majumdar, L C Foucard, A J Levine, and M L Gardel.
	Encoding mechano-memories in actin networks.
	arXiv, arXiv:1706.05336 [cond-mat.soft], 2017.
	
	[Billen et al.(2009)Billen, Wilson, Rabinovitch, and
	Baljon]Billen2009
	J Billen, M Wilson, A Rabinovitch, and A R C Baljon.
	Topological changes at the gel transition of a reversible polymeric
	network.
	EPL (Europhysics Letters), 870 (6):0 68003,
	2009.
	
	[Kim et al.(2014)Kim, Litvinov, Weisel, and Alber]Kim2014
	O V Kim, R I Litvinov, J W Weisel, and M S Alber.
	Structural basis for the nonlinear mechanics of fibrin networks under
	compression.
	Biomaterials, 350 (25):0 6739–6749, 2014.
	
	[Gavrilov et al.(2015)Gavrilov, Komarov, and Khalatur]Gavrilov2015
	A A Gavrilov, P V Komarov, and P G Khalatur.
	Thermal properties and topology of epoxy networks: A multiscale
	simulation methodology.
	Macromolecules, 480 (1):0 206–212, 2015.
	
	[Liang et al.(2016)Liang, Jones, Chen, Sun, and Jiao]Liang:2016
	L Liang, C Jones, S Chen, B Sun, and Y Jiao.
	Heterogeneous force network in 3D cellularized collagen networks.
	Phys Biol, 130 (6):0 066001, 2016.
	
	[Venkatesan et al.(2017)Venkatesan, Vivek-Ananth, Sreejith,
	Mangalapandi, Hassanali, and Samal]samal2017
	S Venkatesan, R P Vivek-Ananth, R P Sreejith, P Mangalapandi, A A Hassanali,
	and A Samal.
	Network approach towards understanding the crazing in glassy
	amorphous polymers.
	arXiv, arXiv:1710.01996 [cond-mat.soft], 2017.
	
	[Bouzid and Del Gado(2017)]Bouzid2017
	M Bouzid and E Del Gado.
	Network topology in soft gels: Hardening and softening materials
	network topology in soft gels: Hardening and softening materials network
	topology in soft gels: Hardening and softening materials.
	Langmuir, 2017.
	URL <http://dx.doi.org/10.1021/acs.langmuir.7b02944>.
	
	[Ahnert et al.(2017)Ahnert, Grant, and Pickard]Ahnert:2017aa
	S E Ahnert, W P Grant, and C J Pickard.
	Revealing and exploiting hierarchical material structure through
	complex atomic networks.
	npj Computational Materials, 30 (1):0 35,
	2017.
	
	[Setford(2014)]setford2014
	J Setford.
	Models of granular networks in two and three dimensions.
	Undergraduate Thesis, Department of Physics, University of Oxford
	(available at
	<http://www.math.ucla.edu/ mason/research/setford-final.pdf>), 2014.
	
	[Mülken and Blumen(2011)]Mulken:2011aa
	O Mülken and A Blumen.
	Continuous-time quantum walks: Models for coherent transport on
	complex networks.
	Phys Rep, 5020 (2–3):0 37–87, 2011.
	
	[Bianconi(2015)]Bianconi:2015a
	G Bianconi.
	Interdisciplinary and physics challenges of network theory.
	EPL (Europhysics Letters), 1110 (5):0 56001,
	2015.
	
	[Biamonte et al.(2017)Biamonte, Faccin, and
	De Domenico]Biamonte:2017a
	J Biamonte, M Faccin, and M De Domenico.
	Complex networks: From classical to quantum.
	arXiv, arXiv:1702.08459 [quant-ph], 2017.
	
	[Boccaletti et al.(2006)Boccaletti, Latora, Moreno, Chavez, and
	Hwang]bocca2006
	S Boccaletti, V Latora, Y Moreno, M Chavez, and D-U Hwang.
	Complex networks: Structure and dynamics.
	Phys Rep, 4240 (4):0 175–308, 2006.
	
	[Liu and Zhang(2011)]Liu:2011
	Y Liu and X Zhang.
	Metamaterials: A new frontier of science and technology.
	Chem Soc Rev, 40:0 2494–2507, 2011.
	
	[Turpin et al.(2014)Turpin, Bossard, Morgan, Werner, and
	Werner]Turpin:2014
	J P Turpin, J A Bossard, K L Morgan, D H Werner, and P L Werner.
	Reconfigurable and tunable metamaterials: A review of the theory
	and applications.
	International Journal of Antennas and Propagation,
	20140 (429837), 2014.
	
	[Lee et al.(2012)Lee, Singer, and Thomas]Lee2012
	J H Lee, J P Singer, and E L Thomas.
	Micro-/nanostructured mechanical metamaterials.
	Adv Mater, 240 (36):0 4782–4810, 2012.
	
	[Greaves et al.(2011)Greaves, Greer, Lakes, and Rouxel]Greaves:2011a
	G N Greaves, A L Greer, R S Lakes, and T Rouxel.
	Poisson's ratio and modern materials.
	Nat Mater, 100 (11):0 823–837, 2011.
	
	[Rocklin et al.(2017)Rocklin, Zhou, Sun, and Mao]Rocklin:2015
	D Z Rocklin, S Zhou, K Sun, and X Mao.
	Transformable topological mechanical metamaterials.
	Nat Commun, 8:0 14201, 2017.
	
	[Fang et al.(2006)Fang, Xi, Xu, Ambati, Srituravanich, Sun, and
	Zhang]Fang:2006
	N Fang, D Xi, J Xu, M Ambati, W Srituravanich, C Sun, and X Zhang.
	Ultrasonic metamaterials with negative modulus.
	Nat Mater, 50 (6):0 452–456, 2006.
	
	[Nicolaou and Motter(2012)]Nicolaou:2012aa
	Z G Nicolaou and A E Motter.
	Mechanical metamaterials with negative compressibility transitions.
	Nat Mater, 110 (7):0 608–613, 2012.
	
	[Simovski et al.(2012)Simovski, Belov, Atrashchenko, and
	Kivshar]Simovski2012
	C R Simovski, P A Belov, A V Atrashchenko, and Y S Kivshar.
	Wire metamaterials: Physics and applications.
	Adv Mater, 240 (31):0 4229–4248, 2012.
	
	[Smith et al.(2004)Smith, Pendry, and Wiltshire]Smith:2004aa
	D. R. Smith, J. B. Pendry, and M. C. K. Wiltshire.
	Metamaterials and negative refractive index.
	Science, 3050 (5685):0 788–792, 2004.
	
	[Eiben and Smith(2015)]eiben2015from
	A E Eiben and J Smith.
	From evolutionary computation to the evolution of things.
	Nature, 5210 (7553):0 476–482, 2015.
	
	[Díaz-Manríquez et al.(2016)Díaz-Manríquez, Toscano,
	Barron-Zambrano, and Tello-Leal]diaz2016review
	A Díaz-Manríquez, G Toscano, J H Barron-Zambrano, and E Tello-Leal.
	A review of surrogate assisted multiobjective evolutionary
	algorithms.
	Comput Intell Neurosci, 20160 (9420460), 2016.
	
	[Papadimitriou(2014)]papadimitriou2014algorithms
	C Papadimitriou.
	Algorithms, complexity, and the sciences.
	Proc Natl Acad Sci, 1110 (45):0 15881–15887,
	2014.
	
	[Goldberg(1989)]Goldberg:1989
	D E Goldberg.
	Genetic Algorithms in Search, Optimization and Machine
		Learning.
	Addison-Wesley Longman Publishing Co., Inc., 1989.
	
	[McGhee(1999)]McGhee:1997
	G R McGhee.
	Theoretical Morphology: The Concept and its Applications.
	Columbia University Press, 1999.
	
	[Valera et al.(2017)Valera, Guo, Kelly, Matz, Cantu, Percus, Hyman,
	Srinivasan, and Viswanathan]valera2017machine
	M Valera, Z Guo, P Kelly, S Matz, A Cantu, A G Percus, J D Hyman, G Srinivasan,
	and H S Viswanathan.
	Machine learning for graph-based representations of three-dimensional
	discrete fracture networks.
	arXiv, arXiv:1705.09866 [physics.geo-ph], 2017.
	
	[Avena-Koenigsberger et al.(2014a)Avena-Koenigsberger,
	Goñi, Solé, and Sporns]Avena-Koenigsberger:2014
	A Avena-Koenigsberger, J Goñi, R Solé, and O Sporns.
	Network morphospace.
	J R Soc Interface, 120 (103), 2014a.
	
	[Avena-Koenigsberger et al.(2014b)Avena-Koenigsberger,
	Goñi, Betzel, van den Heuvel, Griffa, Hagmann, Thiran, and
	Sporns]Avena-Koenigsberger:2013
	A Avena-Koenigsberger, J Goñi, R F Betzel, M P van den Heuvel, A Griffa,
	P Hagmann, J-P Thiran, and O Sporns.
	Using pareto optimality to explore the topology and dynamics of the
	human connectome.
	Philos Trans R Soc Lond B Biol Sci, 3690 (1653),
	2014b.
	
	[Goñi et al.(2013)Goñi, Avena-Koenigsberger, de Mendizabal,
	van den Heuvel, Betzel, and Sporns]Goni:2013
	J Goñi, A Avena-Koenigsberger, NV de Mendizabal, M P van den Heuvel,
	R F Betzel, and O Sporns.
	Exploring the morphospace of communication efficiency in complex
	networks.
	PLoS One, 80 (3):0 e58070, 2013.
	
	[Jaeger and de Pablo(2016)]jaeger2016evolutionary
	H M Jaeger and J J de Pablo.
	Perspective: Evolutionary design of granular media and block
	copolymer patterns.
	APL Materials, 40 (5):0 053209, 2016.
	
	[Miskin and Jaeger(2013)]miskin2013adapting
	M Z Miskin and H M Jaeger.
	Adapting granular materials through artificial evolution.
	Nat Mater, 12:0 326–331, 2013.
	
	[Miskin and Jaeger(2014)]miskin2014evolving
	M Z Miskin and H M Jaeger.
	Evolving design rules for the inverse granular packing problem.
	Soft Matter, 10:0 3708–3715, 2014.
	
	[Roth and Jaeger(2016)]roth2016optimizing
	L K Roth and H M Jaeger.
	Optimizing packing fraction in granular media composed of overlapping
	spheres.
	Soft Matter, 12:0 1107–1115, 2016.
	
	[Yan et al.(2017)Yan, Ravasio, Brito, and Wyart]Yan:2017a
	L Yan, R Ravasio, C Brito, and M Wyart.
	Architecture and coevolution of allosteric materials.
	Proc Natl Acad Sci, 1140 (10):0 2526–2531,
	2017.
	
	[Ellenbroek et al.(2015)Ellenbroek, Hagh, Kumar, Thorpe, and van
	Hecke]Ellenbroek:2015aa
	W G Ellenbroek, V F Hagh, A Kumar, M F Thorpe, and M van Hecke.
	Rigidity loss in disordered systems: Three scenarios.
	Phys Rev Lett, 1140 (13):0 135501, 2015.
	
	[Goodrich et al.(2015)Goodrich, Liu, and Nagel]Goodrich:2015a
	C P Goodrich, A J Liu, and S R Nagel.
	The principle of independent bond-level response: Tuning by pruning
	to exploit disorder for global behavior.
	Phys Rev Lett, 1140 (22):0 225501, 2015.
	
	[Rocks et al.(2017)Rocks, Pashine, Bischofberger, Goodrich, Liu, and
	Nagel]Rocks:2017a
	J W Rocks, N Pashine, I Bischofberger, C P Goodrich, A J Liu, and S R Nagel.
	Designing allostery-inspired response in mechanical networks.
	Proc Natl Acad Sci, 1140 (10):0 2520–2525,
	2017.
	
	[Driscoll et al.(2016)Driscoll, Chen, Beuman, Ulrich, Nagel, and
	Vitelli]Driscoll2016
	M M Driscoll, B G-G Chen, T H Beuman, S Ulrich, S R Nagel, and V Vitelli.
	The role of rigidity in controlling material failure.
	Proc Natl Acad Sci U S A, 1130 (39):0
	10813–10817, 2016.
	
	[Shekhawat et al.(2013)Shekhawat, Zapperi, and
	Sethna]shekhawat2013damage
	A Shekhawat, S Zapperi, and J P Sethna.
	From damage percolation to crack nucleation through finite size
	criticality.
	Phys Rev Lett, 1100 (18):0 185505, 2013.
	
	[Reid et al.(2017)Reid, Pashine, Jaeger, Liu, Nagel, and
	de Pablo]ReidAuxeticMetamaterials:2017
	D R Reid, N Pashine, H M Jaeger, A J Liu, S R Nagel, and J J de Pablo.
	Auxetic metamaterials from disordered networks.
	arXiv, arXiv:1710.02493 [cond-mat.soft], 2017.
	
	[Quinn and Winkelstein(2011)]Quinn:2011
	K P Quinn and B A Winkelstein.
	Preconditioning is correlated with altered collagen fiber alignment
	in ligament.
	J Biomech Eng, 1330 (6):0 064506–064506,
	2011.
	
	[Zhao et al.(2014)Zhao, Chen, and Reich]Zhao:2014
	R Zhao, C S Chen, and D H Reich.
	Force-driven evolution of mesoscale structure in engineered 3D
	microtissues and the modulation of tissue stiffening.
	Biomaterials, 350 (19):0 5056–5064, 2014.
	
	[Han et al.(2013)Han, Heo, Driscoll, Smith, Mauck, and
	Elliott]Han:2013
	W M Han, S-J Heo, T P Driscoll, L J Smith, R L Mauck, and D M Elliott.
	Macro- to microscale strain transfer in fibrous tissues is
	heterogeneous and tissue-specific.
	Biophys J, 1050 (3):0 807–817, 2013.
	
	[Pong et al.(2011)Pong, Adams, Bray, Feinberg, Sheehy, Werdich, and
	Parker]Pong:2011
	T Pong, W J Adams, M-A Bray, A W Feinberg, S P Sheehy, A A Werdich, and K K
	Parker.
	Hierarchical architecture influences calcium dynamics in engineered
	cardiac muscle.
	Exp Biol Med, 2360 (3):0 366–373, 2011.
	
	[Sporns(2014)]Sporns:2014aa
	O Sporns.
	Towards network substrates of brain disorders.
	Brain, 1370 (8):0 2117–2118, 2014.
	



