Skip to content
This repository
Browse code

fix laplacian LIL-matrix bug

  • Loading branch information...
commit 98d2758620930dd84c67159969207ba4840413b5 1 parent a89be18
Jake Vanderplas authored July 10, 2012 rgommers committed July 17, 2012
20  scipy/sparse/csgraph/_laplacian.py
@@ -8,7 +8,8 @@
8 8
 # License: BSD
9 9
 
10 10
 import numpy as np
11  
-from scipy.sparse import isspmatrix
  11
+from scipy.sparse import isspmatrix, coo_matrix
  12
+
12 13
 
13 14
 ###############################################################################
14 15
 # Graph laplacian
@@ -31,7 +32,7 @@ def laplacian(csgraph, normed=False, return_diag=False):
31 32
     -------
32 33
     lap: ndarray
33 34
         The N x N laplacian matrix of graph.
34  
-    
  35
+
35 36
     diag: ndarray
36 37
         The length-N diagonal of the laplacian matrix.
37 38
         diag is returned only if return_diag is True.
@@ -70,7 +71,7 @@ def laplacian(csgraph, normed=False, return_diag=False):
70 71
                    or np.issubdtype(csgraph.dtype, np.uint)):
71 72
         csgraph = csgraph.astype(np.float)
72 73
 
73  
-    if isspmatrix(csgraph):            
  74
+    if isspmatrix(csgraph):
74 75
         return _laplacian_sparse(csgraph, normed=normed,
75 76
                                  return_diag=return_diag)
76 77
     else:
@@ -89,14 +90,14 @@ def _laplacian_sparse(graph, normed=False, return_diag=False):
89 90
         # The sparsity pattern of the matrix has holes on the diagonal,
90 91
         # we need to fix that
91 92
         diag_idx = lap.row[diag_mask]
92  
-
93  
-        lap = lap.tolil()
94  
-
95 93
         diagonal_holes = list(set(range(n_nodes)).difference(
96 94
                                 diag_idx))
97  
-        lap[diagonal_holes, diagonal_holes] = 1
98  
-        lap = lap.tocoo()
  95
+        new_data = np.concatenate([lap.data, np.ones(len(diagonal_holes))])
  96
+        new_row = np.concatenate([lap.row, diagonal_holes])
  97
+        new_col = np.concatenate([lap.col, diagonal_holes])
  98
+        lap = coo_matrix((new_data, (new_row, new_col)), shape=lap.shape)
99 99
         diag_mask = (lap.row == lap.col)
  100
+
100 101
     lap.data[diag_mask] = 0
101 102
     w = -np.asarray(lap.sum(axis=1)).squeeze()
102 103
     if normed:
@@ -115,7 +116,7 @@ def _laplacian_sparse(graph, normed=False, return_diag=False):
115 116
 
116 117
 def _laplacian_dense(graph, normed=False, return_diag=False):
117 118
     n_nodes = graph.shape[0]
118  
-    lap = -graph.copy()
  119
+    lap = -np.asarray(graph)  # minus sign leads to a copy
119 120
 
120 121
     # set diagonal to zero
121 122
     lap.flat[::n_nodes + 1] = 0
@@ -132,4 +133,3 @@ def _laplacian_dense(graph, normed=False, return_diag=False):
132 133
     if return_diag:
133 134
         return lap, w
134 135
     return lap
135  
-
2  scipy/sparse/csgraph/tests/test_graph_laplacian.py
... ...
@@ -1,4 +1,5 @@
1 1
 # Author: Gael Varoquaux <gael.varoquaux@normalesup.org>
  2
+#         Jake Vanderplas <vanderplas@astro.washington.edu>
2 3
 # License: BSD
3 4
 
4 5
 import numpy as np
@@ -12,6 +13,7 @@ def test_graph_laplacian():
12 13
                 np.ones((7, 7)),
13 14
                 np.eye(19),
14 15
                 np.vander(np.arange(4)) + np.vander(np.arange(4)).T,
  16
+                sparse.diags([1, 1], [-1, 1], shape=(4, 4)).toarray(),
15 17
                ):
16 18
         sp_mat = sparse.csr_matrix(mat)
17 19
         for normed in (True, False):

0 notes on commit 98d2758

Please sign in to comment.
Something went wrong with that request. Please try again.