Statistics
| Branch: | Revision:

iof-tools / networkxMiCe / networkx-master / networkx / algorithms / centrality / flow_matrix.py @ 5cef0f13

History | View | Annotate | Download (4.27 KB)

1
# Helpers for current-flow betweenness and current-flow closness
2
# Lazy computations for inverse Laplacian and flow-matrix rows.
3
import networkx as nx
4

    
5

    
6
def flow_matrix_row(G, weight=None, dtype=float, solver='lu'):
7
    # Generate a row of the current-flow matrix
8
    import numpy as np
9
    from scipy import sparse
10
    from scipy.sparse import linalg
11
    solvername = {"full": FullInverseLaplacian,
12
                  "lu": SuperLUInverseLaplacian,
13
                  "cg": CGInverseLaplacian}
14
    n = G.number_of_nodes()
15
    L = laplacian_sparse_matrix(G, nodelist=range(n), weight=weight,
16
                                dtype=dtype, format='csc')
17
    C = solvername[solver](L, dtype=dtype)  # initialize solver
18
    w = C.w  # w is the Laplacian matrix width
19
    # row-by-row flow matrix
20
    for u, v in sorted(sorted((u, v)) for u, v in G.edges()):
21
        B = np.zeros(w, dtype=dtype)
22
        c = G[u][v].get(weight, 1.0)
23
        B[u % w] = c
24
        B[v % w] = -c
25
        # get only the rows needed in the inverse laplacian
26
        # and multiply to get the flow matrix row
27
        row = np.dot(B, C.get_rows(u, v))
28
        yield row, (u, v)
29

    
30

    
31
# Class to compute the inverse laplacian only for specified rows
32
# Allows computation of the current-flow matrix without storing entire
33
# inverse laplacian matrix
34
class InverseLaplacian(object):
35
    def __init__(self, L, width=None, dtype=None):
36
        global np
37
        import numpy as np
38
        (n, n) = L.shape
39
        self.dtype = dtype
40
        self.n = n
41
        if width is None:
42
            self.w = self.width(L)
43
        else:
44
            self.w = width
45
        self.C = np.zeros((self.w, n), dtype=dtype)
46
        self.L1 = L[1:, 1:]
47
        self.init_solver(L)
48

    
49
    def init_solver(self, L):
50
        pass
51

    
52
    def solve(self, r):
53
        raise nx.NetworkXError("Implement solver")
54

    
55
    def solve_inverse(self, r):
56
        raise nx.NetworkXError("Implement solver")
57

    
58
    def get_rows(self, r1, r2):
59
        for r in range(r1, r2 + 1):
60
            self.C[r % self.w, 1:] = self.solve_inverse(r)
61
        return self.C
62

    
63
    def get_row(self, r):
64
        self.C[r % self.w, 1:] = self.solve_inverse(r)
65
        return self.C[r % self.w]
66

    
67
    def width(self, L):
68
        m = 0
69
        for i, row in enumerate(L):
70
            w = 0
71
            x, y = np.nonzero(row)
72
            if len(y) > 0:
73
                v = y - i
74
                w = v.max() - v.min() + 1
75
                m = max(w, m)
76
        return m
77

    
78

    
79
class FullInverseLaplacian(InverseLaplacian):
80
    def init_solver(self, L):
81
        self.IL = np.zeros(L.shape, dtype=self.dtype)
82
        self.IL[1:, 1:] = np.linalg.inv(self.L1.todense())
83

    
84
    def solve(self, rhs):
85
        s = np.zeros(rhs.shape, dtype=self.dtype)
86
        s = np.dot(self.IL, rhs)
87
        return s
88

    
89
    def solve_inverse(self, r):
90
        return self.IL[r, 1:]
91

    
92

    
93
class SuperLUInverseLaplacian(InverseLaplacian):
94
    def init_solver(self, L):
95
        from scipy.sparse import linalg
96
        self.lusolve = linalg.factorized(self.L1.tocsc())
97

    
98
    def solve_inverse(self, r):
99
        rhs = np.zeros(self.n, dtype=self.dtype)
100
        rhs[r] = 1
101
        return self.lusolve(rhs[1:])
102

    
103
    def solve(self, rhs):
104
        s = np.zeros(rhs.shape, dtype=self.dtype)
105
        s[1:] = self.lusolve(rhs[1:])
106
        return s
107

    
108

    
109
class CGInverseLaplacian(InverseLaplacian):
110
    def init_solver(self, L):
111
        global linalg
112
        from scipy.sparse import linalg
113
        ilu = linalg.spilu(self.L1.tocsc())
114
        n = self.n - 1
115
        self.M = linalg.LinearOperator(shape=(n, n), matvec=ilu.solve)
116

    
117
    def solve(self, rhs):
118
        s = np.zeros(rhs.shape, dtype=self.dtype)
119
        s[1:] = linalg.cg(self.L1, rhs[1:], M=self.M)[0]
120
        return s
121

    
122
    def solve_inverse(self, r):
123
        rhs = np.zeros(self.n, self.dtype)
124
        rhs[r] = 1
125
        return linalg.cg(self.L1, rhs[1:], M=self.M)[0]
126

    
127

    
128
# graph laplacian, sparse version, will move to linalg/laplacianmatrix.py
129
def laplacian_sparse_matrix(G, nodelist=None, weight=None, dtype=None,
130
                            format='csr'):
131
    import numpy as np
132
    import scipy.sparse
133
    A = nx.to_scipy_sparse_matrix(G, nodelist=nodelist, weight=weight,
134
                                  dtype=dtype, format=format)
135
    (n, n) = A.shape
136
    data = np.asarray(A.sum(axis=1).T)
137
    D = scipy.sparse.spdiags(data, 0, n, n, format=format)
138
    return D - A