Coverage for /pythoncovmergedfiles/medio/medio/usr/local/lib/python3.9/dist-packages/networkx/algorithms/centrality/flow_matrix.py: 26%

1# Helpers for current-flow betweenness and current-flow closeness 

2# Lazy computations for inverse Laplacian and flow-matrix rows. 

3import networkx as nx 

4 

5 

6@nx._dispatch(edge_attrs="weight") 

7def flow_matrix_row(G, weight=None, dtype=float, solver="lu"): 

8 # Generate a row of the current-flow matrix 

9 import numpy as np 

10 

11 solvername = { 

12 "full": FullInverseLaplacian, 

13 "lu": SuperLUInverseLaplacian, 

14 "cg": CGInverseLaplacian, 

15 } 

16 n = G.number_of_nodes() 

17 L = nx.laplacian_matrix(G, nodelist=range(n), weight=weight).asformat("csc") 

18 L = L.astype(dtype) 

19 C = solvername[solver](L, dtype=dtype) # initialize solver 

20 w = C.w # w is the Laplacian matrix width 

21 # row-by-row flow matrix 

22 for u, v in sorted(sorted((u, v)) for u, v in G.edges()): 

23 B = np.zeros(w, dtype=dtype) 

24 c = G[u][v].get(weight, 1.0) 

25 B[u % w] = c 

26 B[v % w] = -c 

27 # get only the rows needed in the inverse laplacian 

28 # and multiply to get the flow matrix row 

29 row = B @ C.get_rows(u, v) 

30 yield row, (u, v) 

31 

32 

33# Class to compute the inverse laplacian only for specified rows 

34# Allows computation of the current-flow matrix without storing entire 

35# inverse laplacian matrix 

36class InverseLaplacian: 

37 def __init__(self, L, width=None, dtype=None): 

38 global np 

39 import numpy as np 

40 

41 (n, n) = L.shape 

42 self.dtype = dtype 

43 self.n = n 

44 if width is None: 

45 self.w = self.width(L) 

46 else: 

47 self.w = width 

48 self.C = np.zeros((self.w, n), dtype=dtype) 

49 self.L1 = L[1:, 1:] 

50 self.init_solver(L) 

51 

52 def init_solver(self, L): 

53 pass 

54 

55 def solve(self, r): 

56 raise nx.NetworkXError("Implement solver") 

57 

58 def solve_inverse(self, r): 

59 raise nx.NetworkXError("Implement solver") 

60 

61 def get_rows(self, r1, r2): 

62 for r in range(r1, r2 + 1): 

63 self.C[r % self.w, 1:] = self.solve_inverse(r) 

64 return self.C 

65 

66 def get_row(self, r): 

67 self.C[r % self.w, 1:] = self.solve_inverse(r) 

68 return self.C[r % self.w] 

69 

70 def width(self, L): 

71 m = 0 

72 for i, row in enumerate(L): 

73 w = 0 

74 x, y = np.nonzero(row) 

75 if len(y) > 0: 

76 v = y - i 

77 w = v.max() - v.min() + 1 

78 m = max(w, m) 

79 return m 

80 

81 

82class FullInverseLaplacian(InverseLaplacian): 

83 def init_solver(self, L): 

84 self.IL = np.zeros(L.shape, dtype=self.dtype) 

85 self.IL[1:, 1:] = np.linalg.inv(self.L1.todense()) 

86 

87 def solve(self, rhs): 

88 s = np.zeros(rhs.shape, dtype=self.dtype) 

89 s = self.IL @ rhs 

90 return s 

91 

92 def solve_inverse(self, r): 

93 return self.IL[r, 1:] 

94 

95 

96class SuperLUInverseLaplacian(InverseLaplacian): 

97 def init_solver(self, L): 

98 import scipy as sp 

99 

100 self.lusolve = sp.sparse.linalg.factorized(self.L1.tocsc()) 

101 

102 def solve_inverse(self, r): 

103 rhs = np.zeros(self.n, dtype=self.dtype) 

104 rhs[r] = 1 

105 return self.lusolve(rhs[1:]) 

106 

107 def solve(self, rhs): 

108 s = np.zeros(rhs.shape, dtype=self.dtype) 

109 s[1:] = self.lusolve(rhs[1:]) 

110 return s 

111 

112 

113class CGInverseLaplacian(InverseLaplacian): 

114 def init_solver(self, L): 

115 global sp 

116 import scipy as sp 

117 

118 ilu = sp.sparse.linalg.spilu(self.L1.tocsc()) 

119 n = self.n - 1 

120 self.M = sp.sparse.linalg.LinearOperator(shape=(n, n), matvec=ilu.solve) 

121 

122 def solve(self, rhs): 

123 s = np.zeros(rhs.shape, dtype=self.dtype) 

124 s[1:] = sp.sparse.linalg.cg(self.L1, rhs[1:], M=self.M, atol=0)[0] 

125 return s 

126 

127 def solve_inverse(self, r): 

128 rhs = np.zeros(self.n, self.dtype) 

129 rhs[r] = 1 

130 return sp.sparse.linalg.cg(self.L1, rhs[1:], M=self.M, atol=0)[0]