Index: inundation/ga/storm_surge/pyvolution-1d/domain.py
===================================================================
--- inundation/ga/storm_surge/pyvolution-1d/domain.py	(revision 278)
+++ inundation/ga/storm_surge/pyvolution-1d/domain.py	(revision 279)
@@ -26,15 +26,15 @@
 
         #Allocate space for geometric quantities
-        self.vertex_coordinates   = zeros((N, 2), Float)
-        self.centroid_coordinates = zeros(N, Float)
-        self.areas                = zeros(N, Float)
+        self.vertices  = zeros((N, 2), Float)
+        self.centroids = zeros(N, Float)
+        self.areas     = zeros(N, Float)
         for i in range(N):
             xl = self.coordinates[i]
             xr = self.coordinates[i+1]
-            self.vertex_coordinates[i,0] = xl
-            self.vertex_coordinates[i,1] = xr
+            self.vertices[i,0] = xl
+            self.vertices[i,1] = xr
             
             centroid = (xl+xr)/2
-            self.centroid_coordinates[i] = centroid
+            self.centroids[i] = centroid
 
             msg = 'Coordinates should be ordered, smallest to largest'
@@ -43,17 +43,22 @@
             self.areas[i] = (xr-xl)
 
+##         print 'N', N
+##         print 'Centroid', self.centroids
+##         print 'Areas', self.areas
+##         print 'Vertex_Coordinates', self.vertices
+        
+    def get_centroids(self):
+        """Return all coordinates of centroids
+        Return x coordinate of centroid for each element as a N array
+        """
+        
+        return self.centroids
 
-        print 'N', N
-        print 'Coordinates', self.coordinates
-        print 'Centroid', self.centroid_coordinates
-        print 'Areas', self.areas
-        print 'Vertex_Coordinates', self.vertex_coordinates
-        
-    def get_centroid_coordinates(self):
+    def get_vertices(self):
         """Return all coordinates of centroids
         Return x coordinate of centroid for each element as a N array
         """
 
-        return self.centroid_coordinates
+        return self.vertices 
 
     def get_coordinate(self, elem_id, vertex=None):
@@ -64,8 +69,14 @@
 
         if vertex is None:
-            return self.centroid_coordinates[elem_id]
+            return self.centroids[elem_id]
         else:
-            return self.vertex_coordinates[elem_id,vertex]
-        
+            return self.vertices[elem_id,vertex]
+
+    def get_area(self, elem_id):
+        """Return area of element id
+        """
+
+        return self.areas[elem_id]
+    
 
 if __name__ == "__main__":
Index: inundation/ga/storm_surge/pyvolution-1d/quantity.py
===================================================================
--- inundation/ga/storm_surge/pyvolution-1d/quantity.py	(revision 278)
+++ inundation/ga/storm_surge/pyvolution-1d/quantity.py	(revision 279)
@@ -50,9 +50,7 @@
         self.centroid_values = zeros(N, Float)
 
-        #Intialise centroid and edge_values
+        #Intialise centroid values
         self.interpolate()
-
-
-
+        
     
     def interpolate(self):
@@ -68,21 +66,4 @@
             self.centroid_values[i] = (v0 + v1)/2
 
-if __name__ == "__main__":
-    from domain import Domain
-
-    points1 = [0.0, 1.0, 2.0, 3.0]
-    vertex_values = [[1.0,2.0],[4.0,5.0],[-1.0,2.0]]
-    
-    D1 = Domain(points1)
-
-    Q1 = Quantity(D1, vertex_values)
-
-    print Q1.vertex_values
-    print Q1.centroid_values
-    
-
-
-            
-            
     def set_values(self, X, location='vertices'):
         """Set values for quantity
@@ -95,5 +76,5 @@
         In case of location == 'centroid' the dimension values must
         be a list of a Numerical array of length N, N being the number
-        of elements in the mesh. Otherwise it must be of dimension Nx2
+        of elements in the mesh. Otherwise it must be of dimension Nx3
 
         The values will be stored in elements following their
@@ -102,5 +83,5 @@
         If values are described a function, it will be evaluated at specified points
 
-        If selected location is vertices, values for centroid
+        If selected location is vertices, values for centroid and edges
         will be assigned interpolated values.
         In any other case, only values for the specified locations
@@ -109,5 +90,5 @@
 
         if location not in ['vertices', 'centroids']:
-            msg = 'Invalid location: %s' %location
+            msg = 'Invalid location: %s, (possible choices vertices, centroids)' %location
             raise msg
 
@@ -126,4 +107,5 @@
             else:
                 self.vertex_values[:] = X                
+
         else:
             #Use array specific method
@@ -134,4 +116,5 @@
             self.interpolate()
             
+
 
 
@@ -150,55 +133,86 @@
         else:
             #Vertices
-            P = self.domain.get_vertex_coordinates()
+            P = self.domain.get_vertices()
             for i in range(2):
-                self.vertex_values[:,i] = f(P[:,2*i])
-
-            
-##     def set_array_values(self, values, location='vertices'):
-##         """Set values for quantity
-
-##         values: Numeric array
-##         location: Where values are to be stored.
-##                   Permissible options are: vertices, edges, centroid
-##                   Default is "vertices"
-
-##         In case of location == 'centroid' the dimension values must
-##         be a list of a Numerical array of length N, N being the number
-##         of elements in the mesh. Otherwise it must be of dimension Nx3
-
-##         The values will be stored in elements following their
-##         internal ordering.
-
-##         If selected location is vertices, values for centroid and edges
-##         will be assigned interpolated values.
-##         In any other case, only values for the specified locations
-##         will be assigned and the others will be left undefined.
-##         """
-
-##         from Numeric import array, Float
-
-##         values = array(values).astype(Float)
-
-##         N = self.centroid_values.shape[0]
-        
-##         msg = 'Number of values must match number of elements'
-##         assert values.shape[0] == N, msg
-        
-##         if location == 'centroids':
-##             assert len(values.shape) == 1, 'Values array must be 1d'
-##             self.centroid_values = values                        
-##         elif location == 'edges':
-##             assert len(values.shape) == 2, 'Values array must be 2d'
-##             msg = 'Array must be N x 3'            
-##             assert values.shape[1] == 3, msg
-            
-##             self.edge_values = values
-##         else:
-##             assert len(values.shape) == 2, 'Values array must be 2d'
-##             msg = 'Array must be N x 3'            
-##             assert values.shape[1] == 3, msg
-            
-##             self.vertex_values = values
-
+                self.vertex_values[:,i] = f(P[:,i])
+
+
+    def set_array_values(self, values, location='vertices'):
+        """Set values for quantity
+
+        values: Numeric array
+        location: Where values are to be stored.
+                  Permissible options are: vertices, centroid
+                  Default is "vertices"
+
+        In case of location == 'centroid' the dimension values must
+        be a list of a Numerical array of length N, N being the number
+        of elements in the mesh. Otherwise it must be of dimension Nx2
+
+        The values will be stored in elements following their
+        internal ordering.
+
+        If selected location is vertices, values for centroid
+        will be assigned interpolated values.
+        In any other case, only values for the specified locations
+        will be assigned and the others will be left undefined.
+        """
+
+        from Numeric import array, Float
+
+        values = array(values).astype(Float)
+
+        N = self.centroid_values.shape[0]
+        
+        msg = 'Number of values must match number of elements'
+        assert values.shape[0] == N, msg
+        
+        if location == 'centroids':
+            assert len(values.shape) == 1, 'Values array must be 1d'
+            self.centroid_values = values                        
+        else:
+            assert len(values.shape) == 2, 'Values array must be 2d'
+            msg = 'Array must be N x 2'            
+            assert values.shape[1] == 2, msg
+            
+            self.vertex_values = values
+
+
+
+if __name__ == "__main__":
+    from domain import Domain
+
+    points1 = [0.0, 1.0, 2.0, 3.0]
+    vertex_values = [[1.0,2.0],[4.0,5.0],[-1.0,2.0]]
+    
+    D1 = Domain(points1)
+
+    Q1 = Quantity(D1, vertex_values)
+
+    print Q1.vertex_values
+    print Q1.centroid_values
+
+    new_vertex_values = [[2.0,1.0],[3.0,4.0],[-2.0,4.0]]
+
+    Q1.set_values(new_vertex_values)
+    
+   
+    print Q1.vertex_values
+    print Q1.centroid_values
+
+    new_centroid_values = [20,30,40]
+    Q1.set_values(new_centroid_values,'centroids')
+            
+    print Q1.vertex_values
+    print Q1.centroid_values
+
+    def fun(x):
+
+        return x**2
+
+    Q1.set_values(fun,'vertices')    
+            
+    print Q1.vertex_values
+    print Q1.centroid_values
 
 
Index: inundation/ga/storm_surge/pyvolution-1d/test_domain.py
===================================================================
--- inundation/ga/storm_surge/pyvolution-1d/test_domain.py	(revision 278)
+++ inundation/ga/storm_surge/pyvolution-1d/test_domain.py	(revision 279)
@@ -10,5 +10,13 @@
 class TestCase(unittest.TestCase):
     def setUp(self):
-        
+        pass
+
+
+    def tearDown(self):
+        pass
+
+
+    def test_construct(self):
+
         a = 0.0
         b = 2.0
@@ -20,18 +28,26 @@
         points = [a, b, c, d, e, f]
 
-        self.points = [a, b, c, d, e, f]
-        self.domain = Domain(self.points)
-
-        
-    def tearDown(self):
-        pass
+        domain = Domain(points) 
+        assert allclose(domain.get_centroids(),[(a+b)*0.5, (b+c)*0.5, (c+d)*0.5, (d+e)*0.5, (e+f)*0.5])
 
 
-    def test_construct(self):
+    def test_area(self):
 
-        domain = self.domain
-        centroids = domain.get_centroids()
-        assert allclose(centroids, [1.0, 2.25, 2.75, 6.5, 11.0])
-            
+        a = 0.0
+        b = 2.0
+        c = 2.5
+        d = 3.0
+        e = 10.0
+        f = 12.0
+
+        points = [a, b, c, d, e, f]
+
+        domain = Domain(points) 
+        assert allclose(domain.get_area(0),b-a)
+        assert allclose(domain.get_area(1),c-b)
+        assert allclose(domain.get_area(2),d-c)
+        assert allclose(domain.get_area(3),e-d)
+        assert allclose(domain.get_area(4),f-e)
+        
     def test_in_order(self):
         a = 0.0
@@ -62,10 +78,10 @@
         domain = Domain(points)
 
-        assert allclose(domain.get_point(0), 0.5*(a+b))
-        assert allclose(domain.get_point(1), 0.5*(b+c))
-        assert allclose(domain.get_point(2), 0.5*(c+d))
+        assert allclose(domain.get_coordinate(0), 0.5*(a+b))
+        assert allclose(domain.get_coordinate(1), 0.5*(b+c))
+        assert allclose(domain.get_coordinate(2), 0.5*(c+d))
         
         try:
-            allclose(domain.get_point(3), 0.5*(a+b))
+            allclose(domain.get_coordinate(3), 0.5*(a+b))
         except:
             pass
@@ -74,11 +90,11 @@
             raise msg
         
-        assert allclose(domain.get_point(0,0), a)
-        assert allclose(domain.get_point(1,0), b)
-        assert allclose(domain.get_point(2,0), c)
+        assert allclose(domain.get_coordinate(0,0), a)
+        assert allclose(domain.get_coordinate(1,0), b)
+        assert allclose(domain.get_coordinate(2,0), c)
 
-        assert allclose(domain.get_point(0,1), b)
-        assert allclose(domain.get_point(1,1), c)
-        assert allclose(domain.get_point(2,1), d)
+        assert allclose(domain.get_coordinate(0,1), b)
+        assert allclose(domain.get_coordinate(1,1), c)
+        assert allclose(domain.get_coordinate(2,1), d)
 
 
Index: inundation/ga/storm_surge/pyvolution-1d/test_quantity.py
===================================================================
--- inundation/ga/storm_surge/pyvolution-1d/test_quantity.py	(revision 278)
+++ inundation/ga/storm_surge/pyvolution-1d/test_quantity.py	(revision 279)
@@ -85,97 +85,69 @@
 
 
-##     def test_set_values(self):
-##         quantity = Quantity(self.mesh4)
-
-
-##         quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]],
-##                             location = 'vertices')
-##         assert allclose(quantity.vertex_values,
-##                         [[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]])
-##         assert allclose(quantity.centroid_values, [2., 5., 3., 0.]) #Centroid
-##         assert allclose(quantity.edge_values, [[2.5, 2.0, 1.5],
-##                                                [5., 5., 5.],
-##                                                [4.5, 4.5, 0.],
-##                                                [3.0, -1.5, -1.5]])
-
-
-##         #Test default
-##         quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]])
-##         assert allclose(quantity.vertex_values,
-##                         [[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]])        
-##         assert allclose(quantity.centroid_values, [2., 5., 3., 0.]) #Centroid
-##         assert allclose(quantity.edge_values, [[2.5, 2.0, 1.5],
-##                                                [5., 5., 5.],
-##                                                [4.5, 4.5, 0.],
-##                                                [3.0, -1.5, -1.5]])
-
-##         #Test centroids
-##         quantity.set_values([1,2,3,4], location = 'centroids')
-##         assert allclose(quantity.centroid_values, [1., 2., 3., 4.]) #Centroid
-
-##         #Test edges
-##         quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]],
-##                             location = 'edges')
-##         assert allclose(quantity.edge_values,
-##                         [[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]])        
-
-##         #Test exceptions
-##         try:
-##             quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]],
-##                                 location = 'bas kamel tuba')
-##         except:
-##             pass
-
-
-##         try:
-##             quantity.set_values([[1,2,3], [0,0,9]])
-##         except AssertionError:
-##             pass
-##         except:
-##             raise 'should have raised Assertionerror'
-
-
-
-##     def test_set_values_const(self):
-##         quantity = Quantity(self.mesh4)
-
-##         quantity.set_values(1.0, location = 'vertices')
-##         assert allclose(quantity.vertex_values,
-##                         [[1,1,1], [1,1,1], [1,1,1], [1, 1, 1]])
-##         assert allclose(quantity.centroid_values, [1, 1, 1, 1]) #Centroid
-##         assert allclose(quantity.edge_values, [[1, 1, 1],
-##                                                [1, 1, 1],
-##                                                [1, 1, 1],
-##                                                [1, 1, 1]])
-
-
-##         quantity.set_values(2.0, location = 'centroids')
-##         assert allclose(quantity.centroid_values, [2, 2, 2, 2])
-
-##         quantity.set_values(3.0, location = 'edges')
-##         assert allclose(quantity.edge_values, [[3, 3, 3],
-##                                                [3, 3, 3],
-##                                                [3, 3, 3],
-##                                                [3, 3, 3]])        
-
-
-##     def test_set_values_func(self):
-##         quantity = Quantity(self.mesh4)
-
-##         def f(x, y):
-##             return x+y
-
-##         quantity.set_values(f, location = 'vertices')
-##         assert allclose(quantity.vertex_values,
-##                         [[2,0,2], [2,2,4], [4,2,4], [4,2,4]])        
-##         assert allclose(quantity.centroid_values,
-##                         [4.0/3, 8.0/3, 10.0/3, 10.0/3])
-##         assert allclose(quantity.edge_values,
-##                         [[1,2,1], [3,3,2], [3,4,3], [3,4,3]])                
-
-        
-##         quantity.set_values(f, location = 'centroids')
-##         assert allclose(quantity.centroid_values,
-##                         [4.0/3, 8.0/3, 10.0/3, 10.0/3])        
+    def test_set_values(self):
+        quantity = Quantity(self.domain5)
+
+
+        quantity.set_values([[1,2], [5,5], [0,0], [-6, 3], [-2,4]],
+                            location = 'vertices')
+        assert allclose(quantity.vertex_values,
+                        [[1,2], [5,5], [0,0], [-6, 3], [-2,4]])
+        assert allclose(quantity.centroid_values, [1.5, 5., 0., -1.5, 1.0]) #Centroid
+
+        #Test default
+        quantity.set_values([[1,2], [5,5], [0,0], [-6, 3], [-2,4]])
+        assert allclose(quantity.vertex_values,
+                        [[1,2], [5,5], [0,0], [-6, 3], [-2,4]])        
+        assert allclose(quantity.centroid_values, [1.5, 5., 0., -1.5, 1.0]) #Centroid
+
+        #Test centroids
+        quantity.set_values([1,2,3,4,5], location = 'centroids')
+        assert allclose(quantity.centroid_values, [1., 2., 3., 4., 5.]) #Centroid
+
+        #Test exceptions
+        try:
+            quantity.set_values([[1,2], [5,5], [0,0], [-6, 3], [-2,4]],
+                                location = 'bas kamel tuba')
+        except:
+            pass
+
+
+        try:
+            quantity.set_values([[1,2], [0,0]])
+        except AssertionError:
+            pass
+        except:
+            raise 'should have raised AssertionError'
+
+
+
+    def test_set_values_const(self):
+        quantity = Quantity(self.domain5)
+
+        quantity.set_values(1.0, location = 'vertices')
+        assert allclose(quantity.vertex_values,
+                        [[1,1], [1,1], [1,1], [1,1], [1,1]])
+        assert allclose(quantity.centroid_values, [1, 1, 1, 1, 1]) #Centroid
+
+
+        quantity.set_values(2.0, location = 'centroids')
+        assert allclose(quantity.centroid_values, [2, 2, 2, 2, 2])
+
+
+    def test_set_values_func(self):
+        quantity = Quantity(self.domain5)
+
+        def f(x):
+            return x**2
+
+        quantity.set_values(f, location = 'vertices')
+        assert allclose(quantity.vertex_values,
+                        [[0,1], [1,4], [4,6.25], [6.25,9.61], [9.61,16]])        
+        assert allclose(quantity.centroid_values,
+                        [0.5, 2.5, 5.125, 7.93, 12.805 ])
+        
+        quantity.set_values(f, location = 'centroids')
+        assert allclose(quantity.centroid_values,
+                        [0.25, 1.5**2, 2.25**2, 2.8**2, 3.55**2])        
 
 
@@ -223,5 +195,5 @@
 ##         assert b[0] == 0.0
 ##         #The others are OK
-##         for i in range(1,4):
+##      163.968025   for i in range(1,4):
 ##             assert allclose(a[i], 3.0)
 ##             assert allclose(b[i], 1.0)
@@ -303,5 +275,5 @@
 ## #         assert v2.conserved_quantities_vertex2 <= qmax
 ## #         assert v2.conserved_quantities_vertex2 >= qmin                    
-
+163.968025
 
 ## #         #Check that volume has been preserved    
