Index: /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/boundary_regions.py
===================================================================
--- /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/boundary_regions.py	(revision 8997)
+++ /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/boundary_regions.py	(revision 8997)
@@ -0,0 +1,167 @@
+"""region.py - Classes for implementing region conditions
+
+NOTE: I've only set it up for testing X values of constants and functions.
+Although it should work for vectors/arrays of values, that implies
+knowing info about the actual triangles, and that is not how the user should
+operate.
+"""
+# FIXME (DSG-DSG) add better comments
+
+import numpy as num
+
+
+class Region:
+    """Base class for modifying quantities based on a region.
+    """   
+
+    def __init__(self, location='vertices'):
+        self.location = location
+        
+    def __call__(self, tag, elements, domain):
+        msg = 'Generic class Boundary must be subclassed'
+        raise Exception(msg)
+
+
+    def build_indices(self, elements, domain):
+        """
+        Return a list of triangle_id or vertex_id, depending on the location
+        """
+        if self.location == 'unique vertices':
+            return domain.get_unique_vertices(elements)
+        else:
+            return elements
+                
+class Set_region(Region):
+    
+    def __init__(self, tag, quantity, X, location='vertices'):
+        """ 
+        name: Name of quantity
+        X: const or function 
+        location: Where values are to be stored.
+        Permissible options are: vertices, centroid
+        """
+        
+        Region.__init__(self)
+        self.tag = tag 
+        self.quantity = quantity
+        self.location = location
+        self.X = X
+
+    def __repr__(self):
+        pass
+    
+    def __call__(self, tag, elements, domain):
+        """
+        """    
+        if tag == self.tag:
+            domain.set_quantity(self.quantity,
+                                self.X,
+                                location=self.location,
+                                indices=self.build_indices(elements, domain))
+
+        
+class Add_value_to_region(Region):
+    """
+    Will add a value to the current quantity value.
+    
+    quantity = initial_quantity + X.
+    
+    This method does not work with functions
+    
+    Use average to add X to a mean average of the quantity values.
+    eg if applying this to elevation this will give a flat roof.
+    """
+    
+    def __init__(self, tag, quantity, X, location='vertices', initial_quantity=None, average=False):
+        #I have to get this going!
+        #Region.__init__(self)
+        self.tag = tag 
+        self.quantity_answer = quantity
+        self.location = location
+        self.X = X
+        self.average = average
+        if initial_quantity is None:
+            self.quantity_initial_value = quantity
+        else:
+            self.quantity_initial_value = initial_quantity
+        if callable(X):
+            raise Exception('This class does not work with functions')
+
+    def __repr__(self):
+        pass
+    
+    def __call__(self, tag, elements, domain):
+        """
+        """    
+        if tag == self.tag:
+            #new_values = domain.get_quantity(self.quantity_initial_value,
+            #              indices=self.build_indices(elements, domain),
+            #              location=self.location) + self.X
+            Q = domain.get_quantity(self.quantity_initial_value)
+            if self.average is True:
+                # Average the points, and then the verts in the averaged point.  
+                values = Q.get_values(indices=self.build_indices(elements, domain),
+                                      location=self.location)
+                av = num.average(values)
+                if self.location == "vertices":
+                    av = num.average(av)
+                new_values = av + self.X    
+            else:
+                new_values = Q.get_values(indices=self.build_indices(elements, domain),
+                                      location=self.location) + self.X    
+            domain.set_quantity(self.quantity_answer, new_values,
+                                indices=self.build_indices(elements, domain),
+                                location=self.location)
+
+class Add_quantities(Region):
+    """
+    Will add a quantity to the current quantity value.
+    """
+    
+    def __init__(self, tag, quantity_answer, adding_quantity, location='vertices'):
+        #I have to get this going!
+        #Region.__init__(self)
+        self.tag = tag 
+        self.quantity_answer = quantity_answer
+        self.adding_quantity = adding_quantity
+        self.location = location
+
+    def __repr__(self):
+        pass
+    
+    def __call__(self, tag, elements, domain):
+        """
+        """    
+        if tag == self.tag:
+
+            #new_values = domain.get_quantity(self.quantity_answer,
+            #              indices=self.build_indices(elements, domain),
+            #              location=self.location) \
+            #              + domain.get_quantity(self.adding_quantity,
+            #              indices=self.build_indices(elements, domain),
+            #              location=self.location)
+
+
+
+            indices = self.build_indices(elements, domain)
+            location = self.location
+            Q1 = domain.get_quantity(self.quantity_answer)
+            Q2 = domain.get_quantity(self.adding_quantity)            
+            new_values = Q1.get_values(indices=indices, location=location) +\
+                         Q2.get_values(indices=indices, location=location)
+
+            
+            domain.set_quantity(self.quantity_answer, new_values,
+                                indices=self.build_indices(elements, domain),
+                                location=self.location)
+
+
+class Stage_no_less_than_elevation(Region):
+    """
+    Will set the stage to not be less than the elevation.
+    This would be good, but it's not region dependent.
+    Wait for it to become a default for pyvolution.
+    """
+    
+    def __init__(self):
+        pass
Index: /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_boundary_regions.py
===================================================================
--- /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_boundary_regions.py	(revision 8997)
+++ /trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_boundary_regions.py	(revision 8997)
@@ -0,0 +1,258 @@
+#!/usr/bin/env python
+
+import unittest
+from math import sqrt
+
+from mesh_factory import rectangular
+from anuga.shallow_water.shallow_water_domain import Domain
+from generic_domain import Generic_Domain
+from boundary_region import *
+#from anuga.config import epsilon
+
+import numpy as num
+
+
+"""
+This is what the mesh in these tests look like;
+
+3---7
+|5 /|
+| /6|
+2---6
+|3 /|
+| /2|
+1---5
+|1 /|
+| /0|
+0---4
+"""
+
+def add_x_y(x, y):
+    return x+y
+
+def give_me_23(x, y):
+    return 23.0
+
+class Test_Region(unittest.TestCase):
+    def setUp(self):
+        pass
+
+
+    def tearDown(self):
+        pass
+
+
+    def test_region_tags(self):
+        """get values based on triangle lists."""
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(1, 3)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.build_tagged_elements_dictionary({'bottom': [0,1],
+                                                 'top': [4,5],
+                                                 'all': [0,1,2,3,4,5]})
+
+        #Set friction
+        manning = 0.07
+        domain.set_quantity('friction', manning)
+
+        a = Set_region('bottom', 'friction', 0.09)
+        b = Set_region('top', 'friction', 1.0)
+        domain.set_region([a, b])
+
+        expected = [[ 0.09,  0.09,  0.09],
+                    [ 0.09,  0.09,  0.09],
+                    [ 0.07,  0.07,  0.07],
+                    [ 0.07,  0.07,  0.07],
+                    [ 1.0,   1.0,   1.0],
+                    [ 1.0,   1.0,   1.0]]
+        msg = ("\ndomain.quantities['friction']=%s\nexpected value=%s"
+               % (str(domain.quantities['friction'].get_values()),
+                  str(expected)))
+        assert num.allclose(domain.quantities['friction'].get_values(),
+                            expected), msg
+
+        #c = Add_Value_To_region('all', 'friction', 10.0)
+        domain.set_region(Add_value_to_region('all', 'friction', 10.0))
+        #print domain.quantities['friction'].get_values()
+        assert num.allclose(domain.quantities['friction'].get_values(),
+                            [[ 10.09, 10.09, 10.09],
+                             [ 10.09, 10.09, 10.09],
+                             [ 10.07, 10.07, 10.07],
+                             [ 10.07, 10.07, 10.07],
+                             [ 11.0,  11.0,  11.0],
+                             [ 11.0,  11.0,  11.0]])
+
+        # trying a function
+        domain.set_region(Set_region('top', 'friction', add_x_y))
+        #print domain.quantities['friction'].get_values()
+        assert num.allclose(domain.quantities['friction'].get_values(),
+                            [[ 10.09, 10.09, 10.09],
+                             [ 10.09, 10.09, 10.09],
+                             [ 10.07, 10.07, 10.07],
+                             [ 10.07, 10.07, 10.07],
+                             [ 5./3,  2.0,  2./3],
+                             [ 1.0,  2./3,  2.0]])
+
+        domain.set_quantity('elevation', 10.0)
+        domain.set_quantity('stage', 10.0)
+        domain.set_region(Add_value_to_region('top', 'stage', 1.0,initial_quantity='elevation'))
+        #print domain.quantities['stage'].get_values()
+        assert num.allclose(domain.quantities['stage'].get_values(),
+                            [[ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 11.0,  11.0,  11.0],
+                             [ 11.0,  11.0,  11.0]])
+
+        
+        domain.set_quantity('elevation', 10.0)
+        domain.set_quantity('stage', give_me_23)
+        #this works as well, (is cleaner, but doesn't work for regions)
+        #domain.set_quantity('stage',
+        #                    domain.quantities['stage'].vertex_values+ \
+        #                    domain.quantities['elevation'].vertex_values)
+        domain.set_region(Add_quantities('top', 'elevation','stage'))
+        #print domain.quantities['stage'].get_values()
+        assert num.allclose(domain.quantities['elevation'].get_values(),
+                            [[ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 10., 10., 10.],
+                             [ 33.,  33.0,  33.],
+                             [ 33.0,  33.,  33.]])
+        
+    def test_unique_vertices(self):
+        """get values based on triangle lists."""
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(1, 3)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.build_tagged_elements_dictionary({'bottom':[0,1],
+                                                 'top':[4,5],
+                                                 'all':[0,1,2,3,4,5]})
+
+        #Set friction
+        manning = 0.07
+        domain.set_quantity('friction', manning)
+
+        a = Set_region('bottom', 'friction', 0.09, location = 'unique vertices')
+        domain.set_region(a)
+        assert num.allclose(domain.quantities['friction'].get_values(),
+                            [[ 0.09,  0.09,  0.09],
+                             [ 0.09,  0.09,  0.09],
+                             [ 0.09,  0.07,  0.09],
+                             [ 0.07,  0.09,  0.07],
+                             [ 0.07,  0.07,  0.07],
+                             [ 0.07,  0.07,  0.07]])
+
+
+    def test_unique_verticesII(self):
+        """
+        get values based on triangle lists.
+        """
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(1, 3)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.build_tagged_elements_dictionary({'bottom':[0,1],
+                                                 'top':[4,5],
+                                                 'all':[0,1,2,3,4,5]})
+
+        #Set friction
+        manning = 0.07
+        domain.set_quantity('friction', manning)
+
+        domain.set_region(Add_value_to_region('bottom', 'friction', 1.0,initial_quantity='friction', location = 'unique vertices'))
+
+        #print domain.quantities['friction'].get_values()
+        assert num.allclose(domain.quantities['friction'].get_values(),\
+                            [[ 1.07,  1.07,  1.07],
+                             [ 1.07,  1.07,  1.07],
+                             [ 1.07,  0.07,  1.07],
+                             [ 0.07,  1.07,  0.07],
+                             [ 0.07,  0.07,  0.07],
+                             [ 0.07,  0.07,  0.07]])
+                         
+    def test_unique_vertices_average_loc_vert(self):
+        """Get values based on triangle lists."""
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(1, 3)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.build_tagged_elements_dictionary({'bottom': [0, 1],
+                                                 'top': [4, 5],
+                                                 'not_bottom': [2, 3, 4, 5]})
+
+        #Set friction
+        domain.set_quantity('friction', add_x_y)
+        av_bottom = 2.0 / 3.0
+        add = 60.0
+        calc_frict = av_bottom + add
+        domain.set_region(Add_value_to_region('bottom', 'friction', add,
+                                              initial_quantity='friction',
+                                              location='vertices',
+                                              average=True))
+
+        frict_points = domain.quantities['friction'].get_values()
+
+        expected = [calc_frict, calc_frict, calc_frict]
+        msg = ('frict_points[0]=%s\nexpected=%s'
+               % (str(frict_points[0]), str(expected)))
+        assert num.allclose(frict_points[0], expected), msg
+
+        msg = ('frict_points[1]=%s\nexpected=%s'
+               % (str(frict_points[1]), str(expected)))
+        assert num.allclose(frict_points[1], expected), msg
+  
+    def test_unique_vertices_average_loc_unique_vert(self):
+        """
+        get values based on triangle lists.
+        """
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(1, 3)
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.build_tagged_elements_dictionary({'bottom':[0,1],
+                                                 'top':[4,5],
+                                                 'not_bottom':[2,3,4,5]})
+
+        #Set friction
+        domain.set_quantity('friction', add_x_y)
+        av_bottom = 2.0/3.0
+        add = 60.0
+        calc_frict = av_bottom + add
+        domain.set_region(Add_value_to_region('bottom', 'friction', add,
+                          initial_quantity='friction',
+                           location='unique vertices',
+                           average=True
+                          ))
+
+        #print domain.quantities['friction'].get_values()
+        frict_points = domain.quantities['friction'].get_values()
+        assert num.allclose(frict_points[0],\
+                            [ calc_frict, calc_frict, calc_frict])
+        assert num.allclose(frict_points[1],\
+                            [ calc_frict, calc_frict, calc_frict])
+        assert num.allclose(frict_points[2],\
+                            [ calc_frict, 1.0 + 2.0/3.0, calc_frict])
+        assert num.allclose(frict_points[3],\
+                            [ 2.0/3.0,calc_frict, 1.0 + 2.0/3.0])
+                                                
+                         
+#-------------------------------------------------------------
+
+if __name__ == "__main__":
+    suite = unittest.makeSuite(Test_Region, 'test')    
+    runner = unittest.TextTestRunner()
+    runner.run(suite)
