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

scipy-svn@scip... scipy-svn@scip...
Mon Jan 14 06:58:20 CST 2008


Author: wnbell
Date: 2008-01-14 06:58:10 -0600 (Mon, 14 Jan 2008)
New Revision: 3829

Added:
   trunk/scipy/sandbox/multigrid/gallery/
   trunk/scipy/sandbox/multigrid/gallery/__init__.py
   trunk/scipy/sandbox/multigrid/gallery/elasticity.py
   trunk/scipy/sandbox/multigrid/gallery/info.py
   trunk/scipy/sandbox/multigrid/gallery/poisson.py
   trunk/scipy/sandbox/multigrid/gallery/setup.py
   trunk/scipy/sandbox/multigrid/gallery/tests/
   trunk/scipy/sandbox/multigrid/gallery/tests/test_elasticity.py
   trunk/scipy/sandbox/multigrid/gallery/tests/test_poisson.py
Modified:
   trunk/scipy/sandbox/multigrid/__init__.py
   trunk/scipy/sandbox/multigrid/multilevel.py
   trunk/scipy/sandbox/multigrid/setup.py
   trunk/scipy/sandbox/multigrid/tests/test_sa.py
Log:
added matrix example gallery
added linear elasticity examples
moved poisson examples


Modified: trunk/scipy/sandbox/multigrid/__init__.py
===================================================================
--- trunk/scipy/sandbox/multigrid/__init__.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/__init__.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -3,6 +3,7 @@
 from info import __doc__
 
 from multilevel import *
+from gallery import *
 
 __all__ = filter(lambda s:not s.startswith('_'),dir())
 from scipy.testing.pkgtester import Tester

Added: trunk/scipy/sandbox/multigrid/gallery/__init__.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/__init__.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/__init__.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,10 @@
+"Matrix Gallery for Multigrid Solvers"
+
+from info import __doc__
+
+from poisson import *
+from elasticity import *
+
+__all__ = filter(lambda s:not s.startswith('_'),dir())
+from scipy.testing.pkgtester import Tester
+test = Tester().test

Added: trunk/scipy/sandbox/multigrid/gallery/elasticity.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/elasticity.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/elasticity.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,144 @@
+"""Linear Elasticity Examples"""
+
+__all__ = [ 'linear_elasticity' ]
+
+from scipy import array, matrix, ones, zeros, arange, empty, \
+        hstack, vstack, tile, ravel, mgrid, concatenate, \
+        cumsum
+from scipy.linalg import inv, det
+from scipy.sparse import coo_matrix, bsr_matrix
+   
+
+def linear_elasticity( grid, spacing=None, E=1e5, nu=0.3, format=None):
+    if len(grid) == 2:
+        return q12d( grid, spacing=spacing, E=E, nu=nu, format=format )
+    else:
+        raise NotImplemented,'no support for grid=%s' % str(grid)
+
+def q12d( grid, spacing=None, E = 1e5, nu = 0.3, dirichlet_boundary = True, format=None):
+    """ Q1 elements in 2 dimensions """
+    X,Y = grid
+    
+    if X < 1 or Y < 1:
+        raise ValueError,'invalid grid shape'
+
+    if dirichlet_boundary:
+        X += 1
+        Y += 1
+
+    pts = mgrid[0:X+1,0:Y+1]
+    pts = hstack((pts[0].T.reshape(-1,1) - X/2.0, pts[1].T.reshape(-1,1) - Y/2.0)) 
+    
+    if spacing is None:
+        DX,DY = 1,1
+    else:
+        DX,DY = spacing
+        pts = [DX,DY]
+
+    #compute local stiffness matrix
+    lame = E * nu / ((1 + nu) * (1 - 2*nu)) # Lame's first parameter
+    mu   = E / (2 + 2*nu)                   # shear modulus
+    vertices = array([[0,0],[DX,0],[DX,DY],[0,DY]])
+    K = stima4(vertices, lame, mu ) 
+    
+    nodes = arange( (X+1)*(Y+1) ).reshape(X+1,Y+1)
+    LL = nodes[:-1,:-1]
+    I  = (2*LL).repeat( K.size ).reshape(-1,8,8)
+    J  = I.copy()
+    I += tile( [0,1,2,3, 2*X + 4, 2*X + 5, 2*X + 2, 2*X + 3], (8,1) )
+    J += tile( [0,1,2,3, 2*X + 4, 2*X + 5, 2*X + 2, 2*X + 3], (8,1) ).T
+    V  = tile( K, (X*Y,1) )
+     
+    I = ravel(I)
+    J = ravel(J)
+    V = ravel(V)
+
+    A = coo_matrix( (V,(I,J)), shape=(pts.size,pts.size) ).tocsr() #sum duplicates
+    A = A.tobsr(blocksize=(2,2))
+
+    del I,J,V,LL,nodes
+
+    B = zeros(( 2 * (X+1)*(Y+1), 3))
+    B[0::2,0] = 1
+    B[1::2,1] = 1
+    B[0::2,2] = -pts[:,1]
+    B[1::2,2] =  pts[:,0]
+
+    if dirichlet_boundary:
+        mask = zeros((X+1, Y+1),dtype='bool')
+        mask[1:-1,1:-1] = True
+        mask = ravel(mask)
+        data = zeros( ((X-1)*(Y-1),2,2) )
+        data[:,0,0] = 1
+        data[:,1,1] = 1
+        indices = arange( (X-1)*(Y-1) )
+        indptr = concatenate((array([0]),cumsum(mask)))
+        P = bsr_matrix((data,indices,indptr),shape=(2*(X+1)*(Y+1),2*(X-1)*(Y-1)))
+        Pt = P.T
+        A = P.T * A * P
+
+        B = Pt * B
+
+    return A.asformat(format),B 
+
+def stima4(vertices, lame, mu):
+    """local stiffness matrix for two dimensional elasticity on a square element
+
+    Material Parameters:
+        - lame  : Lame's first parameter
+        - mu : shear modulus
+
+    Note:
+        Vertices should be listed in counter-clockwise order:
+            [3]----[2]
+             |      |
+             |      |
+            [0]----[1]
+        
+        Degrees of freedom are enumerated as follows:
+            [x=6,y=7]----[x=4,y=5]
+                |            |
+                |            |
+            [x=0,y=1]----[x=2,y=3]
+    """
+
+    M    = lame + 2*mu # P-wave modulus
+
+    R_11 = matrix([[  2, -2, -1,  1],
+                   [ -2,  2,  1, -1],
+                   [ -1,  1,  2, -2],
+                   [  1, -1, -2,  2]]) / 6.0
+    
+    R_12 = matrix([[  1,  1, -1, -1],
+                   [ -1, -1,  1,  1],
+                   [ -1, -1,  1,  1],
+                   [  1,  1, -1, -1]]) / 4.0
+    
+    R_22 = matrix([[  2,  1, -1, -2],
+                   [  1,  2, -2, -1],
+                   [ -1, -2,  2,  1],
+                   [ -2, -1,  1,  2]]) / 6.0
+     
+    F = inv( vstack( (vertices[1] - vertices[0], vertices[3] - vertices[0]) ) )
+    
+    K = zeros((8,8)) # stiffness matrix
+    
+    E = F.T * matrix([[M, 0],[0, mu]]) * F
+    K[0::2,0::2] = E[0,0] * R_11 + E[0,1] * R_12 + E[1,0] * R_12.T + E[1,1] * R_22
+    
+    E = F.T * matrix([[mu, 0],[0, M]]) * F
+    K[1::2,1::2] = E[0,0] * R_11 + E[0,1] * R_12 + E[1,0] * R_12.T + E[1,1] * R_22
+    
+    E = F.T * matrix([[0, mu],[lame, 0]]) * F;
+    K[1::2,0::2] = E[0,0] * R_11 + E[0,1] * R_12 + E[1,0] * R_12.T + E[1,1] * R_22
+    
+    K[0::2,1::2] = K[1::2,0::2].T
+    
+    K /= det(F)
+
+    return K
+
+
+   
+   
+   

Added: trunk/scipy/sandbox/multigrid/gallery/info.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/info.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/info.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,12 @@
+"""Matrix Gallery for Multigrid Solvers
+
+Functions
+=========
+
+    - poisson() : Poisson problem using Finite Differences
+    - linear_elasticity() : Linear Elasticity using Finite Elements
+
+"""
+
+postpone_import = 1
+

Added: trunk/scipy/sandbox/multigrid/gallery/poisson.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/poisson.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/poisson.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,55 @@
+__all__ = ['poisson']
+
+from scipy import array, empty 
+from scipy.sparse import dia_matrix
+
+def poisson(N, stencil='5pt', dtype=float, format=None):
+    """Finite Difference approximations to the Poisson problem
+
+    TheDirichlet boundary conditions are
+   
+    
+    Parameters
+    ==========
+        - N : integer 
+            - grid size
+        - stencil : one of the following strings
+            - '3pt' : 3-point Finite Difference stencil in 1 dimension
+            - '5pt' : 5-point Finite Difference stencil in 2 dimensions
+            - '8pt' : NotImplemented
+            - '27pt' : NotImplemented
+
+    """
+    if N < 1:
+        raise ValueError,'invalid grid size %s' % N
+
+    if stencil == '3pt':
+        if N == 1:
+            diags   = array( [[2]], dtype=dtype)
+            return dia_matrix((diags,[0]), shape=(1,1)).asformat(format)
+        else:
+            data = empty((3,N),dtype=dtype)
+            data[0,:] = 2 #main diagonal
+            data[1,:] = -1 
+            data[2,:] = -1
+    
+            return dia_matrix((data,[0,-1,1]),shape=(N,N)).asformat(format)
+    elif stencil == '5pt':
+        if N == 1:
+            data = array( [[4]], dtype=dtype)
+            return dia_matrix((diags,[0]), shape=(1,1)).asformat(format)
+        else:
+            diags = array([0,-N,N,-1,1])
+
+            data = empty((5,N**2),dtype=dtype)
+
+            data[0]  =  4 #main diagonal
+            data[1::,:] = -1
+            data[3,N-1::N] = 0 
+            data[4,N::N]   = 0    
+            
+            return dia_matrix((data,diags),shape=(N**2,N**2)).asformat(format)
+    else:
+        raise NotImplementedError,'unsupported stencil=%s' % stencil
+
+

Added: trunk/scipy/sandbox/multigrid/gallery/setup.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/setup.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/setup.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,17 @@
+#!/usr/bin/env python
+
+import sys
+
+def configuration(parent_package='',top_path=None):
+    import numpy
+    from numpy.distutils.misc_util import Configuration
+
+    config = Configuration('gallery',parent_package,top_path)
+
+    config.add_data_dir('tests')
+
+    return config
+
+if __name__ == '__main__':
+    from numpy.distutils.core import setup
+    setup(**configuration(top_path='').todict())


Property changes on: trunk/scipy/sandbox/multigrid/gallery/setup.py
___________________________________________________________________
Name: svn:executable
   + *

Added: trunk/scipy/sandbox/multigrid/gallery/tests/test_elasticity.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/tests/test_elasticity.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/tests/test_elasticity.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,115 @@
+from scipy.testing import *
+
+from scipy import matrix, array
+from scipy.sparse import coo_matrix
+
+from scipy.sandbox.multigrid.gallery.elasticity import linear_elasticity, \
+        stima4
+
+class TestLinearElasticity(TestCase):
+    def test_1x1(self):
+        A_expected = matrix ([[ 230769.23076923,       0.        ],
+                              [      0.        ,  230769.23076923]])
+        B_expected = array([[1, 0, 0],
+                            [0, 1, 0]])
+
+        A,B = linear_elasticity( (1,1), E=1e5, nu=0.3 )
+        
+        assert_almost_equal(A.todense(),A_expected)
+        assert_almost_equal(B,B_expected)
+    
+    def test_1x1(self):
+
+        data = array([ 230769.23076923,  -76923.07692308,   19230.76923077,
+                       -28846.15384615,  -24038.46153846,  230769.23076923,
+                        19230.76923077,  -76923.07692308,  -24038.46153846,
+                       -28846.15384615,  -76923.07692308,  230769.23076923,
+                       -28846.15384615,   24038.46153846,   19230.76923077,
+                        19230.76923077,  230769.23076923,   24038.46153846,
+                       -28846.15384615,  -76923.07692308,   19230.76923077,
+                       -28846.15384615,   24038.46153846,  230769.23076923,
+                       -76923.07692308,  -76923.07692308,   24038.46153846,
+                       -28846.15384615,  230769.23076923,   19230.76923077,
+                       -28846.15384615,  -24038.46153846,   19230.76923077,
+                       -76923.07692308,  230769.23076923,  -24038.46153846,
+                       -28846.15384615,  -76923.07692308,   19230.76923077,
+                       230769.23076923])
+        row = array([0, 2, 4, 6, 7, 1, 3, 5, 6, 7, 0, 2, 4, 5, 6, 1, 3, 4, 5, 7, 0, 2, 3, 4, 6, 1, 2, 3, 5, 7, 0, 1, 2, 4, 6, 0, 1, 3, 5, 7])
+        col = array([0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7])
+
+        A_expected = coo_matrix((data,(row,col)), shape=(8,8)).todense()
+        B_expected = array([[ 1. ,  0. ,  0.5],
+                            [ 0. ,  1. , -0.5],
+                            [ 1. ,  0. ,  0.5],
+                            [ 0. ,  1. ,  0.5],
+                            [ 1. ,  0. , -0.5],
+                            [ 0. ,  1. , -0.5],
+                            [ 1. ,  0. , -0.5],
+                            [ 0. ,  1. ,  0.5]])
+        
+        A,B = linear_elasticity( (2,2), E=1e5, nu=0.3 )
+
+        assert_almost_equal(A.todense(),A_expected)
+        assert_almost_equal(B,B_expected)
+
+class TestLocalStiffnessMatrix(TestCase):
+    def test_stima4(self):
+        L = matrix([[  4,  3, -4,  3, -2, -3,  2, -3],
+                    [  3,  4, -3,  2, -3, -2,  3, -4],
+                    [ -4, -3,  4, -3,  2,  3, -2,  3],
+                    [  3,  2, -3,  4, -3, -4,  3, -2],
+                    [ -2, -3,  2, -3,  4,  3, -4,  3],
+                    [ -3, -2,  3, -4,  3,  4, -3,  2],
+                    [  2,  3, -2,  3, -4, -3,  4, -3],
+                    [ -3, -4,  3, -2,  3,  2, -3,  4]]) / 12.0
+        
+        M = matrix([[  4,  1, -2, -1, -2, -1,  0,  1],
+                    [  1,  4,  1,  0, -1, -2, -1, -2],
+                    [ -2,  1,  4, -1,  0, -1, -2,  1],
+                    [ -1,  0, -1,  4,  1, -2,  1, -2],
+                    [ -2, -1,  0,  1,  4,  1, -2, -1],
+                    [ -1, -2, -1, -2,  1,  4,  1,  0],
+                    [  0, -1, -2,  1, -2,  1,  4, -1],
+                    [  1, -2,  1, -2, -1,  0, -1,  4]]) / 4.0
+        
+        vertices = matrix([[ 0, 0],
+                           [ 1, 0],
+                           [ 1, 1],
+                           [ 0, 1]])
+        
+        assert_almost_equal( stima4(vertices, 1, 0) , L)
+        assert_almost_equal( stima4(vertices, 0, 1) , M)
+        assert_almost_equal( stima4(vertices, 1, 1) , L + M)
+           
+        
+           
+        L = matrix([[ 2,  3, -2,  3, -1, -3,  1, -3],
+                    [ 3,  8, -3,  4, -3, -4,  3, -8],
+                    [-2, -3,  2, -3,  1,  3, -1,  3],
+                    [ 3,  4, -3,  8, -3, -8,  3, -4],
+                    [-1, -3,  1, -3,  2,  3, -2,  3],
+                    [-3, -4,  3, -8,  3,  8, -3,  4],
+                    [ 1,  3, -1,  3, -2, -3,  2, -3],
+                    [-3, -8,  3, -4,  3,  4, -3,  8]]) / 12.0
+        
+        M = matrix([[ 4,  1,  0, -1, -2, -1, -2,  1],
+                    [ 1,  6,  1,  2, -1, -3, -1, -5],
+                    [ 0,  1,  4, -1, -2, -1, -2,  1],
+                    [-1,  2, -1,  6,  1, -5,  1, -3],
+                    [-2, -1, -2,  1,  4,  1,  0, -1],
+                    [-1, -3, -1, -5,  1,  6,  1,  2],
+                    [-2, -1, -2,  1,  0,  1,  4, -1],
+                    [ 1, -5,  1, -3, -1,  2, -1,  6]]) / 4.0
+        
+        vertices = matrix([[ 0, 0],
+                           [ 2, 0],
+                           [ 2, 1],
+                           [ 0, 1]])
+        
+        assert_almost_equal( stima4(vertices, 1, 0) , L)
+        assert_almost_equal( stima4(vertices, 0, 1) , M)
+        assert_almost_equal( stima4(vertices, 1, 1) , L + M)
+
+
+if __name__ == '__main__':
+    unittest.main()

Added: trunk/scipy/sandbox/multigrid/gallery/tests/test_poisson.py
===================================================================
--- trunk/scipy/sandbox/multigrid/gallery/tests/test_poisson.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/gallery/tests/test_poisson.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -0,0 +1,26 @@
+from scipy.testing import *
+
+from scipy import matrix
+
+from scipy.sandbox.multigrid.gallery.poisson import *
+
+class TestPoisson(TestCase):
+    def test_3pt(self):
+        cases = []
+
+        cases.append( (1,matrix([[2]])) )
+        cases.append( (2,matrix([[ 2,-1],
+                                 [-1, 2]])) )
+        cases.append( (4,matrix([[ 2,-1, 0, 0],
+                                 [-1, 2,-1, 0],
+                                 [ 0,-1, 2,-1],
+                                 [ 0, 0,-1, 2]])) )
+
+
+        for N,expected in cases:
+            assert_equal( poisson(N,stencil='3pt').todense(), expected )
+
+
+
+if __name__ == '__main__':
+    unittest.main()

Modified: trunk/scipy/sandbox/multigrid/multilevel.py
===================================================================
--- trunk/scipy/sandbox/multigrid/multilevel.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/multilevel.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -1,5 +1,4 @@
-__all__ = ['poisson_problem1D','poisson_problem2D',
-           'ruge_stuben_solver','smoothed_aggregation_solver',
+__all__ = ['ruge_stuben_solver','smoothed_aggregation_solver',
            'multilevel_solver']
 
 import scipy
@@ -15,43 +14,7 @@
 from utils import symmetric_rescaling, diag_sparse
 
 
-def poisson_problem1D(N):
-    """
-    Return a sparse CSR matrix for the 1d poisson problem
-    with standard 3-point finite difference stencil on a
-    grid with N points.
-    """
-    D = 2*ones(N)
-    O =  -ones(N)
-    return scipy.sparse.spdiags([D,O,O],[0,-1,1],N,N).tocoo().tocsr() #eliminate explicit zeros
 
-def poisson_problem2D(N, epsilon=1.0, dtype='d', format=None):
-    """
-    Return a sparse matrix for the 2d poisson problem
-    with standard 5-point finite difference stencil on a
-    square N-by-N grid.
-    """
-    if N == 1:
-        diags   = asarray( [[4]],dtype=dtype)
-        return dia_matrix((diags,[0]), shape=(1,1)).asformat(format)
-
-    offsets = array([0,-N,N,-1,1])
-
-    diags = empty((5,N**2),dtype=dtype)
-
-    diags[0]  =  (2 + 2*epsilon) #main diagonal
-    diags[1,:] = -1
-    diags[2,:] = -1
-    
-    diags[3,:] = -epsilon #first lower diagonal 
-    diags[4,:] = -epsilon #first upper diagonal 
-    diags[3,N-1::N] = 0  
-    diags[4,N::N]   = 0    
-    
-    return dia_matrix((diags,offsets),shape=(N**2,N**2)).tocoo().tocsr() #eliminate explicit zeros
-    #return dia_matrix((diags,offsets),shape=(N**2,N**2)).asformat(format)
-
-
 def ruge_stuben_solver(A,max_levels=10,max_coarse=500):
     """
     Create a multilevel solver using Ruge-Stuben coarsening (Classical AMG)

Modified: trunk/scipy/sandbox/multigrid/setup.py
===================================================================
--- trunk/scipy/sandbox/multigrid/setup.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/setup.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -9,6 +9,8 @@
 
     config = Configuration('multigrid',parent_package,top_path)
 
+    config.add_subpackage('gallery')
+
     config.add_data_dir('tests')
     config.add_data_dir(join('tests','sample_data'))
 

Modified: trunk/scipy/sandbox/multigrid/tests/test_sa.py
===================================================================
--- trunk/scipy/sandbox/multigrid/tests/test_sa.py	2008-01-14 08:33:50 UTC (rev 3828)
+++ trunk/scipy/sandbox/multigrid/tests/test_sa.py	2008-01-14 12:58:10 UTC (rev 3829)
@@ -19,12 +19,11 @@
 from scipy.sandbox.multigrid.sa import sa_strong_connections, sa_constant_interpolation, \
                                         sa_interpolation, sa_fit_candidates, \
                                         sa_smoothed_prolongator
-from scipy.sandbox.multigrid.multilevel import poisson_problem1D,poisson_problem2D, \
-                                        smoothed_aggregation_solver
+from scipy.sandbox.multigrid.multilevel import smoothed_aggregation_solver
 from scipy.sandbox.multigrid.utils import diag_sparse
 
+from scipy.sandbox.multigrid.gallery import poisson
 
-
 #def sparsity(A):
 #    A = A.copy()
 #
@@ -49,9 +48,9 @@
 
         # poisson problems in 1D and 2D
         for N in [2,3,5,7,10,11,19]:
-            self.cases.append( poisson_problem1D(N) )
+            self.cases.append( poisson(N,stencil='3pt',format='csr') )
         for N in [2,3,5,7,10,11]:
-            self.cases.append( poisson_problem2D(N) )
+            self.cases.append( poisson(N,stencil='5pt',format='csr') )
 
 
     def test_sa_strong_connections(self):
@@ -75,7 +74,7 @@
                 #assert_array_equal(S_result.todense(),S_expected.todense())
 
         # two aggregates in 1D
-        A = poisson_problem1D(6)
+        A = poisson(6,stencil='3pt')
         AggOp = csr_matrix((ones(6),array([0,0,0,1,1,1]),arange(7)),shape=(6,2))
         candidates = ones((6,1))
 
@@ -107,13 +106,13 @@
         user_cases = []
 
         #simple 1d example w/ two aggregates
-        A = poisson_problem1D(6)
+        A = poisson(6, stencil='3pt', format='csr')
         AggOp = csr_matrix((ones(6),array([0,0,0,1,1,1]),arange(7)),shape=(6,2))
         candidates = ones((6,1))
         user_cases.append((A,AggOp,candidates))
 
         #simple 1d example w/ two aggregates (not all nodes are aggregated)
-        A = poisson_problem1D(6)
+        A = poisson(6, stencil='3pt', format='csr')
         AggOp = csr_matrix((ones(4),array([0,0,1,1]),array([0,1,1,2,3,3,4])),shape=(6,2))
         candidates = ones((6,1))
         user_cases.append((A,AggOp,candidates))
@@ -184,8 +183,8 @@
     def setUp(self):
         self.cases = []
 
-        self.cases.append((poisson_problem1D(100),None))
-        self.cases.append((poisson_problem2D(50),None))
+        self.cases.append(( poisson(100, stencil='3pt', format='csr'), None))
+        self.cases.append(( poisson(100, stencil='5pt', format='csr'), None))
         # TODO add unstructured tests
 
 



More information about the Scipy-svn mailing list