Add adjacency tensor function (#664)

* feat: added adj tensor

* test: added + fix support for empty hypergraph

* fix: test

* docs: added new function

* style: isort + black

* review changes

* review: added normalized option

* fix: dtype

* fix: coldict

docs/source/api/linalg/xgi.linalg.hypergraph_matrix.rst

@@ -8,6 +8,7 @@
    .. rubric:: Functions
 
    .. autofunction:: adjacency_matrix
+   .. autofunction:: adjacency_tensor
    .. autofunction:: clique_motif_matrix
    .. autofunction:: degree_matrix
    .. autofunction:: incidence_matrix

tests/linalg/test_matrix.py

@@ -756,3 +756,48 @@ def test_empty():
     assert xgi.boundary_matrix(H).shape == (0, 0)
 
     assert xgi.hodge_laplacian(H).shape == (0, 0)
+
+
+def test_adjacency_tensor(edgelist1):
+    el1 = edgelist1
+    H1 = xgi.Hypergraph(el1)
+
+    A11, node_dict = xgi.adjacency_tensor(H1, order=1, index=True)
+    node_dict1 = {k: v for v, k in node_dict.items()}
+
+    assert type(A11) == np.ndarray
+    assert A11.shape == (8, 8)
+    assert A11.sum() == 2
+    assert list(np.unique(A11)) == [0, 1]
+    assert np.all(A11 == xgi.adjacency_matrix(H1, order=1, sparse=False))
+
+    A12 = xgi.adjacency_tensor(H1, order=2)
+    assert type(A12) == np.ndarray
+    assert A12.shape == (8, 8, 8)
+    assert A12.sum() == 12
+    assert list(np.unique(A12)) == [0, 1]
+
+    assert A12[node_dict1[1], node_dict1[2], node_dict1[3]] == 1
+    assert A12[node_dict1[1], node_dict1[3], node_dict1[2]] == 1
+    assert A12[node_dict1[2], node_dict1[3], node_dict1[1]] == 1
+    assert A12[node_dict1[2], node_dict1[1], node_dict1[3]] == 1
+    assert A12[node_dict1[3], node_dict1[1], node_dict1[2]] == 1
+    assert A12[node_dict1[3], node_dict1[2], node_dict1[1]] == 1
+
+    assert A12[node_dict1[6], node_dict1[7], node_dict1[8]] == 1
+    assert A12[node_dict1[6], node_dict1[8], node_dict1[7]] == 1
+    assert A12[node_dict1[7], node_dict1[8], node_dict1[6]] == 1
+    assert A12[node_dict1[7], node_dict1[6], node_dict1[8]] == 1
+    assert A12[node_dict1[8], node_dict1[6], node_dict1[7]] == 1
+    assert A12[node_dict1[8], node_dict1[7], node_dict1[6]] == 1
+
+    # normalization
+    A11_norm = xgi.adjacency_tensor(H1, order=1, normalized=True)
+    assert np.allclose(A11_norm, A11)
+
+    A12_norm = xgi.adjacency_tensor(H1, order=2, normalized=True)
+    assert np.allclose(A12_norm, A12 / 2)
+
+    A13 = xgi.adjacency_tensor(H1, order=3)
+    A13_norm = xgi.adjacency_tensor(H1, order=3, normalized=True)
+    assert np.allclose(A13_norm, A13 / 6)

xgi/linalg/hypergraph_matrix.py

@@ -44,6 +44,8 @@
 
 """
 
+from itertools import permutations
+from math import factorial
 from warnings import catch_warnings, warn
 
 import numpy as np
@@ -55,6 +57,7 @@
     "intersection_profile",
     "degree_matrix",
     "clique_motif_matrix",
+    "adjacency_tensor",
 ]
 
 
@@ -294,3 +297,63 @@ def clique_motif_matrix(H, sparse=True, index=False):
     W, rowdict = adjacency_matrix(H, sparse=sparse, weighted=True, index=True)
 
     return (W, rowdict) if index else W
+
+
+def adjacency_tensor(H, order, normalized=False, index=False):
+    """
+    Compute the order-d adjacency tensor of a hypergraph from its incidence matrix.
+
+    The order-d adjacency tensor B is defined as B[i1, i2, ..., id] = 1 if there exists
+    a hyperedge containing nodes (i1, i2, ..., id), and 0 otherwise.
+
+    Parameters
+    ----------
+    H : xgi.Hypergraph
+        The input hypergraph.
+    order : int
+        The order of interactions to consider.
+    index: bool, default: False
+        Specifies whether to output disctionaries mapping the node IDs to indices
+
+    Returns
+    -------
+    B : np.ndarray
+        Adjacency tensor
+
+    References
+    ----------
+    Galuppi, F., Mulas, R., & Venturello, L. (2023).   
+    "Spectral theory of weighted hypergraphs via tensors"  
+    Linear and Multilinear Algebra, 71(3), 317-347.  
+    """
+
+    I, rowdict, coldict = incidence_matrix(H, order=order, sparse=False, index=True)
+
+    if I.shape == (0, 0):
+        if not rowdict:
+            B = np.empty((0,) * (order + 1))
+        if not coldict:
+            shape = (H.num_nodes,) * (order + 1)
+            B = np.zeros(shape, dtype=int)
+        return (B, {}) if index else B
+
+    # find nodes participating in each hyperedge
+    hyperedges = H.edges.filterby("order", order).members()
+
+    nodedict = {v: k for (k, v) in rowdict.items()}
+    N = H.num_nodes
+
+    # initialize adjacency tensor
+    dtype = float if normalized else int
+    B = np.zeros((N,) * (order + 1), dtype=dtype)
+
+    # populate adjacency tensor
+    for edge in hyperedges:
+        edge_node_ids = [nodedict[node] for node in edge]
+        for node_idx in permutations(edge_node_ids, order + 1):
+            B[node_idx] = 1
+
+    if normalized:
+        B /= factorial(order)
+
+    return (B, rowdict) if index else B