1from functools import wraps
2
3import numpy as np
4import scipy.sparse as sp
5from sklearn.utils.extmath import row_norms
6
7from orangecontrib.network import Network
8
9__all__ = ["path", "cycle", "complete", "complete_bipartite", "barbell",
10           "ladder", "circular_ladder", "grid", "hypercube",
11           "lollipop", "star", "geometric"]
12
13def from_csr(f):
14    @wraps(f)
15    def wrapped(*args):
16        edges = f(*args)
17        return Network(
18            range(edges.shape[0]), edges,
19            name=f"{f.__name__}{args}".replace(",)", ")"))
20    return wrapped
21
22
23def from_row_col(f):
24    @wraps(f)
25    def wrapped(*args):
26        row, col, *n = f(*args)
27        n = n[0] if n else max(np.max(row), np.max(col)) + 1
28        edges = sp.csr_matrix((np.ones(len(row)), (row, col)), shape=(n, n))
29        return Network(
30            range(n), edges,
31            name=f"{f.__name__}{args}".replace(",)", ")"))
32    return wrapped
33
34
35def from_dense(f):
36    @wraps(f)
37    def wrapped(*args):
38        m = f(*args)
39        return Network(
40            range(len(m)), sp.csr_matrix(m),
41            name=f"{f.__name__}{args}".replace(",)", ")"))
42    return wrapped
43
44
45@from_row_col
46def path(n):
47    return np.arange(n - 1), np.arange(n - 1) + 1, n
48
49
50@from_row_col
51def cycle(n):
52    return np.arange(n), np.arange(1, n + 1) % n
53
54
55@from_csr
56def complete(n):
57    if n < 5_000:
58        return sp.triu(np.ones((n, n)), k=1, format="csr")
59
60    # manually construct sparse matrix elements to avoid storing a dense (n x n) matrix in memory
61    indptr = [0]
62    for n_above_diag in range(n - 1, 0, -1):
63        indptr.append(indptr[-1] + n_above_diag)
64    indptr.append(indptr[-1])
65    indptr = np.array(indptr)
66    indices = np.concatenate([np.arange(row, n) for row in range(1, n)])
67    data = np.ones(n * (n - 1) // 2)
68
69    return sp.csr_matrix((data, indices, indptr), shape=(n, n))
70
71
72@from_row_col
73def complete_bipartite(n, m):
74    return np.repeat(np.arange(n), m), n + np.arange(n * m) % m
75
76
77@from_dense
78def barbell(n, m):
79    e = np.zeros((n + m, n + m))
80    e[:n, :n] = 1  # complete
81    e[-m:, -m:] = 1  # complete
82    e[n, -m - 1] = e[n - 1, -m] = 1  # bridge
83    e -= np.eye(n + m, n + m)  # no loops
84    return np.triu(e)
85
86
87def ladder(n):
88    return grid(n, 2)
89
90
91def circular_ladder(n):
92    return grid(n, 2, True)
93
94
95@from_csr
96def grid(n, m=4, circular=False):
97    t = n * m
98    row = np.arange(2 * t - n) % t
99    col = np.zeros(2 * t - n, dtype=int)
100    data = np.ones(2 * t - n)
101
102    col[:t] = np.arange(1, t + 1)  # horizontal edges
103    col[n - 1:t:n] = np.arange(m)  # circle to prevent index out of range
104    if not circular:
105        data[n - 1:t:n] = 0  # eliminate circle
106
107    col[t:] = np.arange(n, t)  # vertical edges
108    edges = sp.csr_matrix((data, (row, col)), shape=(t, t))
109    edges.eliminate_zeros()
110    return edges
111
112
113# https://math.stackexchange.com/questions/2328139/adjacency-matrix-for-n-dimensional-hypercube-graph
114def _hypercube(ndim):
115    """Recursively construct the edge-connectivity of a hypercube
116
117    Parameters
118    ----------
119    ndim : int
120        Dimension of the hypercube
121
122    Returns
123    -------
124    ndarray, [2**ndim, 2**ndim], bool
125        connectivity pattern of the hypercube
126    """
127    if ndim == 0:
128        return np.array([[0]])
129    else:
130        d = _hypercube(ndim-1)
131        i = np.eye(len(d))
132        return np.block([
133            [d, i],
134            [i, d],
135        ])
136
137
138@from_dense
139def hypercube(ndim):
140    return _hypercube(ndim)
141
142
143@from_row_col
144def star(n):
145    return np.zeros(n), np.arange(1, n + 1)
146
147
148@from_dense
149def lollipop(n, m):
150    e = np.zeros((n + m, n + m))
151    e[:n, :n] = 1 - np.eye(n)
152    e[n, -m - 1] = e[n - 1, -m] = 1  # bridge
153    for i in range(m, n + m):
154        e[i - 1, i] = 1
155    return np.triu(e)
156
157
158def geometric(n_nodes, n_edges):
159    n_pairs = n_nodes * (n_nodes - 1) // 2
160    if n_edges > n_pairs:
161        raise ValueError(
162            f"There are only {n_pairs} (< {n_edges}) possible edges between "
163            f"{n_nodes} points")
164    xy = np.random.random((n_nodes, 2))
165    xx = row_norms(xy, squared=True)[:, np.newaxis]
166    distances = np.dot(xy, xy.T)
167    distances *= -2
168    distances += xx
169    distances += xx.T
170    ur = np.triu_indices(n_nodes, k=1)
171    # skip zeros and repetitions
172    dist_threshold = np.partition(distances[ur], n_edges)[n_edges]
173    mask = distances <= dist_threshold
174    mask[np.tril_indices(n_nodes)] = False
175    row, col = mask.nonzero()
176    edges = sp.csr_matrix((np.ones(len(row)), (row, col)),
177                          shape=(n_nodes, n_nodes))
178    return Network(
179        range(n_nodes), edges,
180        name=f"geometric({n_nodes},{n_edges})",
181        coordinates=xy
182    )
183