[Scipy-svn] r3875 - in trunk/scipy/sandbox/multigrid: . tests

scipy-svn@scip... scipy-svn@scip...
Tue Jan 29 19:17:50 CST 2008


Author: wnbell
Date: 2008-01-29 19:17:40 -0600 (Tue, 29 Jan 2008)
New Revision: 3875

Added:
   trunk/scipy/sandbox/multigrid/tests/test_rs.py
Modified:
   trunk/scipy/sandbox/multigrid/__init__.py
   trunk/scipy/sandbox/multigrid/multilevel.py
   trunk/scipy/sandbox/multigrid/rs.py
Log:
added RS unittest
moved RS to rs.py


Modified: trunk/scipy/sandbox/multigrid/__init__.py
===================================================================
--- trunk/scipy/sandbox/multigrid/__init__.py	2008-01-29 22:24:05 UTC (rev 3874)
+++ trunk/scipy/sandbox/multigrid/__init__.py	2008-01-30 01:17:40 UTC (rev 3875)
@@ -2,7 +2,9 @@
 
 from info import __doc__
 
-from multilevel import *
+from multilevel import multilevel_solver
+from rs import ruge_stuben_solver
+from sa import smoothed_aggregation_solver
 from gallery import *
 
 __all__ = filter(lambda s:not s.startswith('_'),dir())

Modified: trunk/scipy/sandbox/multigrid/multilevel.py
===================================================================
--- trunk/scipy/sandbox/multigrid/multilevel.py	2008-01-29 22:24:05 UTC (rev 3874)
+++ trunk/scipy/sandbox/multigrid/multilevel.py	2008-01-30 01:17:40 UTC (rev 3875)
@@ -1,51 +1,20 @@
-__all__ = ['ruge_stuben_solver','multilevel_solver']
+__all__ = ['multilevel_solver']
 
 import scipy
 import numpy
-from numpy import ones,zeros,zeros_like,array,asarray,empty
+from numpy import ones, zeros, zeros_like, array, asarray, empty
 from numpy.linalg import norm
 from scipy.splinalg import spsolve
-from scipy.sparse import dia_matrix
 
-from rs import rs_interpolation
 from relaxation import gauss_seidel,jacobi,sor
 from utils import symmetric_rescaling, diag_sparse
 
 
-
-def ruge_stuben_solver(A,max_levels=10,max_coarse=500):
-    """
-    Create a multilevel solver using Ruge-Stuben coarsening (Classical AMG)
-
-        References:
-            "Multigrid"
-                Trottenberg, U., C. W. Oosterlee, and Anton Schuller.
-                San Diego: Academic Press, 2001.
-                Appendix A
-
-    """
-    As = [A]
-    Ps = []
-
-    while len(As) < max_levels  and A.shape[0] > max_coarse:
-        P = rs_interpolation(A)
-
-        A = (P.T.tocsr() * A) * P     #galerkin operator
-
-        As.append(A)
-        Ps.append(P)
-
-    return multilevel_solver(As,Ps)
-
-
-
-
-
 class multilevel_solver:
-    def __init__(self,As,Ps,Rs=None,preprocess=None,postprocess=None):
+    def __init__(self, As, Ps, Rs=None, preprocess=None, postprocess=None):
         self.As = As
         self.Ps = Ps
-        self.preprocess = preprocess
+        self.preprocess  = preprocess
         self.postprocess = postprocess
 
         if Rs is None:

Modified: trunk/scipy/sandbox/multigrid/rs.py
===================================================================
--- trunk/scipy/sandbox/multigrid/rs.py	2008-01-29 22:24:05 UTC (rev 3874)
+++ trunk/scipy/sandbox/multigrid/rs.py	2008-01-30 01:17:40 UTC (rev 3875)
@@ -1,9 +1,40 @@
 from scipy.sparse import csr_matrix,isspmatrix_csr
 
+from multilevel import multilevel_solver
 import multigridtools
 
-__all__ = ['rs_strong_connections','rs_interpolation']
+__all__ = ['ruge_stuben_solver','rs_strong_connections','rs_interpolation']
 
+
+def ruge_stuben_solver(A, max_levels=10, max_coarse=500):
+    """
+    Create a multilevel solver using Ruge-Stuben coarsening (Classical AMG)
+
+        References:
+            "Multigrid"
+                Trottenberg, U., C. W. Oosterlee, and Anton Schuller.
+                San Diego: Academic Press, 2001.
+                Appendix A
+
+    """
+    As = [A]
+    Ps = []
+    Rs = []
+
+    while len(As) < max_levels  and A.shape[0] > max_coarse:
+        P = rs_interpolation(A)
+        R = P.T.tocsr()
+
+        A = R * A * P     #galerkin operator
+
+        As.append(A)
+        Ps.append(P)
+        Rs.append(R)
+
+    return multilevel_solver(As,Ps,Rs=Rs)
+
+
+
 def rs_strong_connections(A,theta):
     """Return a strength of connection matrix using the method of Ruge and Stuben
 

Added: trunk/scipy/sandbox/multigrid/tests/test_rs.py
===================================================================
--- trunk/scipy/sandbox/multigrid/tests/test_rs.py	2008-01-29 22:24:05 UTC (rev 3874)
+++ trunk/scipy/sandbox/multigrid/tests/test_rs.py	2008-01-30 01:17:40 UTC (rev 3875)
@@ -0,0 +1,38 @@
+from scipy.testing import *
+from scipy import rand
+from scipy.sparse import spdiags,csr_matrix,lil_matrix, \
+                         isspmatrix_csr,isspmatrix_csc,isspmatrix_coo, \
+                         isspmatrix_lil
+import numpy
+
+from scipy.sandbox.multigrid.rs import ruge_stuben_solver
+from scipy.sandbox.multigrid.gallery import poisson
+
+class TestRugeStubenSolver(TestCase):
+    def test_poisson(self):
+        cases = []
+        
+        cases.append( (1000,) )
+        cases.append( (500,500) )
+        cases.append( (50,50,50) )
+
+        for case in cases:
+            A = poisson( case, format='csr' )
+
+            numpy.random.seed(0) #make tests repeatable
+
+            x = rand(A.shape[0])
+            b = A*rand(A.shape[0]) #zeros_like(x)
+
+            ml = ruge_stuben_solver(A)
+
+            x_sol,residuals = ml.solve(b,x0=x,maxiter=20,tol=1e-12,return_residuals=True)
+
+            avg_convergence_ratio = (residuals[-1]/residuals[0])**(1.0/len(residuals))
+            
+            print "avg",avg_convergence_ratio
+
+            assert(avg_convergence_ratio < 0.15)
+
+if __name__ == '__main__':
+    nose.run(argv=['', __file__])



More information about the Scipy-svn mailing list