[Feat] Random expanders utilities (#6761)

Adds functionality to both generate and assess whether graphs are regular expanders.

* add maybe_regular_expander, is_regular_expander, random_regular_expander utilities

---------

Co-authored-by: Dan Schult <dschult@colgate.edu>
Co-authored-by: Ross Barnowski <rossbar@berkeley.edu>

doc/reference/generators.rst

@@ -54,6 +54,9 @@ Expanders
    margulis_gabber_galil_graph
    chordal_cycle_graph
    paley_graph
+   maybe_regular_expander
+   is_regular_expander
+   random_regular_expander_graph
 
 Lattice
 -------

networkx/conftest.py

@@ -219,6 +219,7 @@ def add_nx(doctest_namespace):
     "algorithms/node_classification.py",
     "algorithms/non_randomness.py",
     "algorithms/shortest_paths/dense.py",
+    "generators/expanders.py",
     "linalg/bethehessianmatrix.py",
     "linalg/laplacianmatrix.py",
     "utils/misc.py",
@@ -245,6 +246,7 @@ def add_nx(doctest_namespace):
     "convert_matrix.py",
     "drawing/layout.py",
     "generators/spectral_graph_forge.py",
+    "generators/expanders.py",
     "linalg/algebraicconnectivity.py",
     "linalg/attrmatrix.py",
     "linalg/bethehessianmatrix.py",

networkx/generators/expanders.py

@@ -5,7 +5,14 @@
 
 import networkx as nx
 
-__all__ = ["margulis_gabber_galil_graph", "chordal_cycle_graph", "paley_graph"]
+__all__ = [
+    "margulis_gabber_galil_graph",
+    "chordal_cycle_graph",
+    "paley_graph",
+    "maybe_regular_expander",
+    "is_regular_expander",
+    "random_regular_expander_graph",
+]
 
 
 # Other discrete torus expanders can be constructed by using the following edge
@@ -204,3 +211,253 @@ def paley_graph(p, create_using=None):
             G.add_edge(x, (x + x2) % p)
     G.graph["name"] = f"paley({p})"
     return G
+
+
+@nx.utils.decorators.np_random_state("seed")
+def maybe_regular_expander(n, d, *, create_using=None, max_tries=100, seed=None):
+    r"""Utility for creating a random regular expander.
+
+    Returns a random $d$-regular graph on $n$ nodes which is an expander
+    graph with very good probability.
+
+    Parameters
+    ----------
+    n : int
+      The number of nodes.
+    d : int
+      The degree of each node.
+    create_using : Graph Instance or Constructor
+      Indicator of type of graph to return.
+      If a Graph-type instance, then clear and use it.
+      If a constructor, call it to create an empty graph.
+      Use the Graph constructor by default.
+    max_tries : int. (default: 100)
+      The number of allowed loops when generating each independent cycle
+    seed : (default: None)
+      Seed used to set random number generation state. See :ref`Randomness<randomness>`.
+
+    Notes
+    -----
+    The nodes are numbered from $0$ to $n - 1$.
+
+    The graph is generated by taking $d / 2$ random independent cycles.
+
+    Joel Friedman proved that in this model the resulting
+    graph is an expander with probability
+    $1 - O(n^{-\tau})$ where $\tau = \lceil (\sqrt{d - 1}) / 2 \rceil - 1$. [1]_
+
+    Examples
+    --------
+    >>> G = nx.maybe_regular_expander(n=200, d=6)
+
+    Returns
+    -------
+    G : graph
+        The constructed undirected graph.
+
+    Raises
+    ------
+    NetworkXError
+        If $d % 2 != 0$ as the degree must be even.
+        If $n - 1$ is less than $ 2d $ as the graph is complete at most.
+        If max_tries is reached
+
+    See Also
+    --------
+    is_regular_expander
+    random_regular_expander_graph
+
+    References
+    ----------
+    .. [1] Joel Friedman,
+       A Proof of Alon’s Second Eigenvalue Conjecture and Related Problems, 2004
+       https://arxiv.org/abs/cs/0405020
+
+    """
+
+    import numpy as np
+
+    if n < 1:
+        raise nx.NetworkXError("n must be a positive integer")
+
+    if not (d >= 2):
+        raise nx.NetworkXError("d must be greater than or equal to 2")
+
+    if not (d % 2 == 0):
+        raise nx.NetworkXError("d must be even")
+
+    if not (n - 1 >= d):
+        raise nx.NetworkXError(
+            f"Need n-1>= d to have room for {d//2} independent cycles with {n} nodes"
+        )
+
+    G = nx.empty_graph(n, create_using)
+
+    if n < 2:
+        return G
+
+    cycles = []
+    edges = set()
+
+    # Create d / 2 cycles
+    for i in range(d // 2):
+        iterations = max_tries
+        # Make sure the cycles are independent to have a regular graph
+        while len(edges) != (i + 1) * n:
+            iterations -= 1
+            # Faster than random.permutation(n) since there are only
+            # (n-1)! distinct cycles against n! permutations of size n
+            cycle = np.concatenate((seed.permutation(n - 1), [n - 1]))
+
+            new_edges = {
+                (u, v)
+                for u, v in nx.utils.pairwise(cycle, cyclic=True)
+                if (u, v) not in edges and (v, u) not in edges
+            }
+            # If the new cycle has no edges in common with previous cycles
+            # then add it to the list otherwise try again
+            if len(new_edges) == n:
+                cycles.append(cycle)
+                edges.update(new_edges)
+
+            if iterations == 0:
+                raise nx.NetworkXError("Too many iterations in maybe_regular_expander")
+
+    G.add_edges_from(edges)
+
+    return G
+
+
+@nx.utils.not_implemented_for("directed")
+@nx.utils.not_implemented_for("multigraph")
+def is_regular_expander(G, *, epsilon=0):
+    r"""Determines whether the graph G is a regular expander. [1]_
+
+    An expander graph is a sparse graph with strong connectivity properties.
+
+    More precisely, this helper checks whether the graph is a
+    regular $(n, d, \lambda)$-expander with $\lambda$ close to
+    the Alon-Boppana bound and given by
+    $\lambda = 2 \sqrt{d - 1} + \epsilon$. [2]_
+
+    In the case where $\epsilon = 0 $ then if the graph successfully passes the test
+    it is a Ramanujan graph. [3]_
+
+    A Ramanujan graph has spectral gap almost as large as possible, which makes them
+    excellent expanders.
+
+    Parameters
+    ----------
+    G : NetworkX graph
+    epsilon : int, float, default=0
+
+    Returns
+    -------
+    bool
+        Whether the given graph is a regular $(n, d, \lambda)$-expander
+        where $\lambda = 2 \sqrt{d - 1} + \epsilon$.
+
+    Examples
+    --------
+    >>> G = nx.random_regular_expander_graph(20, 4)
+    >>> nx.is_regular_expander(G)
+    True
+
+    See Also
+    --------
+    maybe_regular_expander
+    random_regular_expander_graph
+
+    References
+    ----------
+    .. [1] Expander graph, https://en.wikipedia.org/wiki/Expander_graph
+    .. [2] Alon-Boppana bound, https://en.wikipedia.org/wiki/Alon%E2%80%93Boppana_bound
+    .. [3] Ramanujan graphs, https://en.wikipedia.org/wiki/Ramanujan_graph
+
+    """
+
+    import numpy as np
+    from scipy.sparse.linalg import eigsh
+
+    if epsilon < 0:
+        raise nx.NetworkXError("epsilon must be non negative")
+
+    if not nx.is_regular(G):
+        return False
+
+    _, d = nx.utils.arbitrary_element(G.degree)
+
+    A = nx.adjacency_matrix(G)
+    lams = eigsh(A.asfptype(), which="LM", k=2, return_eigenvectors=False)
+
+    # lambda2 is the second biggest eigenvalue
+    lambda2 = min(lams)
+
+    return abs(lambda2) < 2 ** np.sqrt(d - 1) + epsilon
+
+
+def random_regular_expander_graph(n, d, *, epsilon=0, create_using=None, max_tries=100):
+    r"""Returns a random regular expander graph on $n$ nodes with degree $d$.
+
+    An expander graph is a sparse graph with strong connectivity properties. [1]_
+
+    More precisely the returned graph is a $(n, d, \lambda)$-expander with
+    $\lambda = 2 \sqrt{d - 1} + \epsilon$, close to the Alon-Boppana bound. [2]_
+
+    In the case where $\epsilon = 0 $ it returns a Ramanujan graph.
+    A Ramanujan graph has spectral gap almost as large as possible,
+    which makes them excellent expanders. [3]_
+
+    Parameters
+    ----------
+    n : int
+      The number of nodes.
+    d : int
+      The degree of each node.
+    epsilon : int, float, default=0
+    max_tries : int, (default: 100)
+      The number of allowed loops, also used in the maybe_regular_expander utility
+
+    Raises
+    ------
+    NetworkXError
+        If max_tries is reached
+
+    Examples
+    --------
+    >>> G = nx.random_regular_expander_graph(20, 4)
+    >>> nx.is_regular_expander(G)
+    True
+
+    Notes
+    -----
+    This loops over `maybe_regular_expander` and can be slow when
+    $n$ is too big or $\epsilon$ too small.
+
+    See Also
+    --------
+    maybe_regular_expander
+    is_regular_expander
+
+    References
+    ----------
+    .. [1] Expander graph, https://en.wikipedia.org/wiki/Expander_graph
+    .. [2] Alon-Boppana bound, https://en.wikipedia.org/wiki/Alon%E2%80%93Boppana_bound
+    .. [3] Ramanujan graphs, https://en.wikipedia.org/wiki/Ramanujan_graph
+
+    """
+    G = maybe_regular_expander(n, d, create_using=create_using, max_tries=max_tries)
+    iterations = max_tries
+
+    while not is_regular_expander(G, epsilon=epsilon):
+        iterations -= 1
+        G = maybe_regular_expander(
+            n=n, d=d, create_using=create_using, max_tries=max_tries
+        )
+
+        if iterations == 0:
+            raise nx.NetworkXError(
+                "Too many iterations in random_regular_expander_graph"
+            )
+
+    return G

networkx/generators/tests/test_expanders.py

@@ -63,6 +63,45 @@ def test_paley_graph(p):
             assert (v, u) in G.edges
 
 
+@pytest.mark.parametrize("d, n", [(2, 7), (4, 10), (4, 16)])
+def test_maybe_regular_expander(d, n):
+    pytest.importorskip("numpy")
+    G = nx.maybe_regular_expander(n, d)
+
+    assert len(G) == n, "Should have n nodes"
+    assert len(G.edges) == n * d / 2, "Should have n*d/2 edges"
+    assert nx.is_k_regular(G, d), "Should be d-regular"
+
+
+@pytest.mark.parametrize("n", (3, 5, 6, 10))
+def test_is_regular_expander(n):
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+    G = nx.complete_graph(n)
+
+    assert nx.is_regular_expander(G) == True, "Should be a regular expander"
+
+
+@pytest.mark.parametrize("d, n", [(2, 7), (4, 10), (4, 16)])
+def test_random_regular_expander(d, n):
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+    G = nx.random_regular_expander_graph(n, d)
+
+    assert len(G) == n, "Should have n nodes"
+    assert len(G.edges) == n * d / 2, "Should have n*d/2 edges"
+    assert nx.is_k_regular(G, d), "Should be d-regular"
+    assert nx.is_regular_expander(G) == True, "Should be a regular expander"
+
+
+def test_random_regular_expander_explicit_construction():
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+    G = nx.random_regular_expander_graph(d=4, n=5)
+
+    assert len(G) == 5 and len(G.edges) == 10, "Should be a complete graph"
+
+
 @pytest.mark.parametrize("graph_type", (nx.Graph, nx.DiGraph, nx.MultiDiGraph))
 def test_margulis_gabber_galil_graph_badinput(graph_type):
     with pytest.raises(
@@ -84,3 +123,42 @@ def test_paley_graph_badinput():
         nx.NetworkXError, match="`create_using` cannot be a multigraph."
     ):
         nx.paley_graph(3, create_using=nx.MultiGraph)
+
+
+def test_maybe_regular_expander_badinput():
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+
+    with pytest.raises(nx.NetworkXError, match="n must be a positive integer"):
+        nx.maybe_regular_expander(n=-1, d=2)
+
+    with pytest.raises(nx.NetworkXError, match="d must be greater than or equal to 2"):
+        nx.maybe_regular_expander(n=10, d=0)
+
+    with pytest.raises(nx.NetworkXError, match="Need n-1>= d to have room"):
+        nx.maybe_regular_expander(n=5, d=6)
+
+
+def test_is_regular_expander_badinput():
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+
+    with pytest.raises(nx.NetworkXError, match="epsilon must be non negative"):
+        nx.is_regular_expander(nx.Graph(), epsilon=-1)
+
+
+def test_random_regular_expander_badinput():
+    pytest.importorskip("numpy")
+    pytest.importorskip("scipy")
+
+    with pytest.raises(nx.NetworkXError, match="n must be a positive integer"):
+        nx.random_regular_expander_graph(n=-1, d=2)
+
+    with pytest.raises(nx.NetworkXError, match="d must be greater than or equal to 2"):
+        nx.random_regular_expander_graph(n=10, d=0)
+
+    with pytest.raises(nx.NetworkXError, match="Need n-1>= d to have room"):
+        nx.random_regular_expander_graph(n=5, d=6)
+
+    with pytest.raises(nx.NetworkXError, match="epsilon must be non negative"):
+        nx.random_regular_expander_graph(n=4, d=2, epsilon=-1)