Index: inundation/ga/storm_surge/pyvolution/sparse.py
===================================================================
--- inundation/ga/storm_surge/pyvolution/sparse.py	(revision 471)
+++ 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
