[Scipy-svn] r3261 - trunk/scipy/sandbox/multigrid

scipy-svn@scip... scipy-svn@scip...
Fri Aug 24 10:55:07 CDT 2007


Author: wnbell
Date: 2007-08-24 10:55:04 -0500 (Fri, 24 Aug 2007)
New Revision: 3261

Modified:
   trunk/scipy/sandbox/multigrid/multilevel.py
   trunk/scipy/sandbox/multigrid/utils.py
Log:
switched inline code to sparse.extract_diagonal()



Modified: trunk/scipy/sandbox/multigrid/multilevel.py
===================================================================
--- trunk/scipy/sandbox/multigrid/multilevel.py	2007-08-23 20:52:36 UTC (rev 3260)
+++ trunk/scipy/sandbox/multigrid/multilevel.py	2007-08-24 15:55:04 UTC (rev 3261)
@@ -154,7 +154,7 @@
 
 if __name__ == '__main__':
     from scipy import *
-    A = poisson_problem2D(200).T
+    A = poisson_problem2D(200)
     asa = smoothed_aggregation_solver(A)
     #asa = ruge_stuben_solver(A)
     x = rand(A.shape[0])

Modified: trunk/scipy/sandbox/multigrid/utils.py
===================================================================
--- trunk/scipy/sandbox/multigrid/utils.py	2007-08-23 20:52:36 UTC (rev 3260)
+++ trunk/scipy/sandbox/multigrid/utils.py	2007-08-24 15:55:04 UTC (rev 3261)
@@ -3,8 +3,9 @@
 import numpy,scipy,scipy.sparse,scipy.weave
 from numpy import ravel,arange
 from scipy.sparse import isspmatrix,isspmatrix_csr,isspmatrix_csc, \
-                        csr_matrix,csc_matrix
+                        csr_matrix,csc_matrix,extract_diagonal
 
+
 def inf_norm(A):
     """
     Infinity norm of a sparse matrix (maximum absolute row sum).  This serves 
@@ -26,31 +27,8 @@
        - return a csr_matrix with A on the diagonal
     """
     
-    if isspmatrix_csr(A) or isspmatrix_csc(A):
-        n_row = len(A.indptr) - 1
-        data,indices,indptr = A.data,A.indices,A.indptr
-
-        diag = numpy.zeros(n_row,dtype=A.dtype)
-        
-        code = """
-        #line 33 "sparse.py"
-
-        for(int i = 0; i < n_row; i++){
-          for(int jj = indptr(i); jj < indptr(i+1); jj++){
-            if(indices(jj) == i){
-              diag(i) = data(jj);
-            }
-          }
-        }
-        """
-
-        err = scipy.weave.inline(code,
-                                 ['data', 'indices', 'indptr', 'n_row', 'diag'],
-                                 type_converters = scipy.weave.converters.blitz,
-                                 compiler = 'gcc')
-        return diag
-    elif isspmatrix(A): 
-        return ravel(array([float(A[i,i]) for i in range(min(A.shape))]))
+    if isspmatrix(A):
+        return extract_diagonal(A)
     else:
         return csr_matrix((A,arange(len(A)),arange(len(A)+1)),(len(A),len(A)))
 



More information about the Scipy-svn mailing list