Index: /inundation/ga/storm_surge/pyvolution/sparse.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/sparse.py	(revision 473)
+++ /inundation/ga/storm_surge/pyvolution/sparse.py	(revision 474)
@@ -2,21 +2,59 @@
 """
 
-#from scipy_base import * #Hardly worth importing scipy just to get isscalar.
-from cg_solve import conjugate_gradient, VectorShapeError
 
 class Sparse:
 
-    def __init__(self, M, N):
-        """Set dimensions
-        """
-
-        self.M = M
-        self.N = N
-        self.shape = (M,N)
+    def __init__(self, *args):
+        """Create sparse matrix.
+        There are two construction forms
+        Usage:
+
+        Sparse(A)     #Creates sparse matrix from dense matrix A 
+        Sparse(M, N)  #Creates empty MxN sparse matrix
+
+
+        
+        """
+
         self.A = {}
+            
+        if len(args) == 1:
+            from Numeric import array
+            try:
+                A = array(args[0])
+            except:
+                raise 'Input must be convertable to a Numeric array'
+
+            assert len(A.shape) == 2, 'Input must be a 2d matrix'
+            
+            self.M, self.N = A.shape
+            for i in range(self.M):
+                for j in range(self.N):
+                    if A[i, j] != 0.0:
+                        self.A[i, j] = A[i, j]
+                
+            
+        elif len(args) == 2:
+            self.M = args[0]
+            self.N = args[1]
+        else:
+            raise 'Invalid construction'
+            
+        self.shape = (self.M, self.N) 
+
 
     def __repr__(self):
-        return '%d X %d sparse matrix:\n' %(self.M, self.N) + `self.A` 
-
+        return '%d X %d sparse matrix:\n' %(self.M, self.N) + `self.A`
+
+    def __len__(self):
+        """Return number of nonzeros of A
+        """
+        return len(self.A)
+
+    def nonzeros(self):
+        """Return number of nonzeros of A
+        """        
+        return len(self)
+    
     def __setitem__(self, key, x):
 
@@ -43,5 +81,5 @@
 
     def copy(self):
-
+        #FIXME: Use the copy module instead
         new = Sparse(self.M,self.N)
 
@@ -81,11 +119,10 @@
         #Assume numeric types from now on
         R = zeros((self.M,), Float) #Result
-        
-        if len(B.shape) == 1:
+
+        if len(B.shape) == 0:
+            #Scalar - use __rmul__ method
+            R = B*self
+        elif len(B.shape) == 1:
             #Vector
-
-##             print 'B.shape      ',B.shape
-##             print 'self.shape ',self.shape
-            
             assert B.shape[0] == self.N, 'Mismatching dimensions'
 
@@ -95,9 +132,11 @@
 
                 R[i] += self.A[key]*B[j]
-
-        else:
-            raise ValueError, 'Numeric matrix not yet implemented'
+        elif len(B.shape) == 2:
+            raise ValueError, 'Numeric matrix not yet implemented'            
+        else:
+            raise ValueError, 'Dimension too high: d=%d' %len(B.shape)
 
         return R
+    
 
     def __add__(self, other):
@@ -108,12 +147,8 @@
         
         new = other.copy()
-
- #       print 'self.shape',self.shape
- #       print 'other.shape',other.shape
-        
         for key in self.A.keys():
             i, j = key
 
-            new[i,j] = new[i,j] + self.A[key]
+            new[i,j] += self.A[key]
 
         return new
@@ -129,5 +164,6 @@
             other = float(other)
         except:
-            raise 'only right multiple with scalar implemented'
+            msg = 'Sparse matrix can only "right-multiply" onto a scalar'
+            raise TypeError, msg
         else:
             new = self.copy()
@@ -142,5 +178,5 @@
 
     def trans_mult(self, other):
-        """Multiply the transpose of  matrix with 'other' which can be
+        """Multiply the transpose of matrix with 'other' which can be
         a Numeric vector.
         """
@@ -155,6 +191,4 @@
 
         #Assume numeric types from now on
-
-        
         if len(B.shape) == 1:
             #Vector
Index: /inundation/ga/storm_surge/pyvolution/test_all.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/test_all.py	(revision 473)
+++ /inundation/ga/storm_surge/pyvolution/test_all.py	(revision 474)
@@ -17,5 +17,5 @@
 #exclude = ['test_least_squares.py', 'test_cg_solve.py',
           # 'test_interpolate_sww.py']
-exclude = ['test_cg_solve.py']
+#exclude = ['test_cg_solve.py']
 
 
@@ -34,7 +34,8 @@
 
 
-    for file in exclude:
-        files.remove(file)
-        print 'WARNING: File '+ file + ' excluded from testing'
+    if globals().has_key('exclude'):    
+        for file in exclude:
+            files.remove(file)
+            print 'WARNING: File '+ file + ' excluded from testing'
         
 
Index: /inundation/ga/storm_surge/pyvolution/test_sparse.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/test_sparse.py	(revision 474)
+++ /inundation/ga/storm_surge/pyvolution/test_sparse.py	(revision 474)
@@ -0,0 +1,164 @@
+#!/usr/bin/env python
+
+import unittest
+from math import sqrt
+
+
+from sparse import *
+from Numeric import allclose, array, transpose, Float 
+        
+class TestCase(unittest.TestCase):
+
+    def setUp(self):
+        pass
+        
+    def tearDown(self):
+        pass
+
+
+    def test_init1(Self):
+        """Test initialisation from dimensions
+        """        
+        A = Sparse(3,3)
+        A[1,1] = 4
+
+        for i in range(3):
+            for j in range(3):
+                if i==1 and j==1:
+                    assert A[i,j] == 4.0
+                else:
+                    assert A[i,j] == 0.0                    
+                    
+    def test_init2(Self):
+        """Test initialisation from dense matrix
+        """
+        
+        A = Sparse(3,3)
+        A[1,1] = 4
+        A[2,0] = 6
+        A[0,1] = 13
+        A[0,2] = -6
+        A[2,2] = 1                        
+
+        B = A.todense()
+
+        C = Sparse(B)
+        
+        assert allclose(C.todense(), B)
+
+                    
+    def test_dense(self):
+        A = Sparse(3,3)
+        A[1,1] = 4        
+                     
+        assert allclose(A.todense(), [[0,0,0], [0,4,0], [0,0,0]])
+
+
+    def test_reset_to_zero_possible(self):
+        """Test that nonzero element can be reset to zero
+        (This was not possible when we tried sparse from scipy).
+        """
+
+        A = Sparse(3,3)
+        A[1,1] = 4        
+        A[1,1] = 0
+
+        assert len(A) == 0
+        assert allclose(A.todense(), [[0,0,0], [0,0,0], [0,0,0]])
+
+        #Set an existing zero element to zero
+        A[1,2] = 0
+        assert len(A) == 0
+        assert allclose(A.todense(), [[0,0,0], [0,0,0], [0,0,0]])        
+
+    def test_sparse_multiplication(self):
+        A = Sparse(3,3)
+        
+        A[0,0] = 3
+        A[1,1] = 2
+        A[1,2] = 2
+        A[2,2] = 1
+
+        #Right hand side vector
+        v = [2,3,4]
+
+        
+        u = A*v
+        assert allclose(u, [6,14,4])
+
+
+
+        #Right hand side column
+        v = array([[2,4],[3,4],[4,4]])
+
+        u = A*v[:,0]
+        assert allclose(u, [6,14,4])
+
+        u = A*v[:,1]
+        assert allclose(u, [12,16,4])        
+
+
+    def test_sparse_transpose_multiplication(self):        
+        A = Sparse(3,3)
+        
+        A[0,0] = 3
+        A[1,1] = 2
+        A[1,2] = 2
+        A[2,2] = 1
+
+        #Right hand side vector
+        v = [2,3,4]
+
+        u = A.trans_mult(v)
+        assert allclose(u, [6,6,10])
+
+
+    def test_scalar_multiplication(self):
+        """Test method __rmul__
+        """
+        
+        A = Sparse(3,3)
+
+        A[0,0] = 3
+        A[1,1] = 2
+        A[1,2] = 2
+        A[2,2] = 1
+
+        B = 3*A
+        assert allclose(B.todense(), 3*A.todense())
+
+        B = A*3
+        assert allclose(B.todense(), 3*A.todense())        
+
+        try:
+            B = 'a'*A
+        except TypeError:
+            pass
+        else:
+            raise 'Should have failed'
+
+    def test_sparse_addition(self):
+        
+        A = Sparse(3,3)
+
+        A[0,0] = 3
+        A[1,1] = 2
+        A[1,2] = 2
+        A[2,2] = 1
+
+
+        B = 3*A
+        B[1,0] = 2
+
+        C = A+B
+
+        assert allclose(C.todense(), [[12,0,0], [2,8,8], [0,0,4]])        
+        
+
+
+#-------------------------------------------------------------
+if __name__ == "__main__":
+    suite = unittest.makeSuite(TestCase, 'test')
+    runner = unittest.TextTestRunner()
+    runner.run(suite)
+
