Changeset 7505

Timestamp:

Sep 10, 2009, 2:14:40 PM (15 years ago)

Author:

steve

Message:

Added a new Time_space_boundary. Just provide the domain and a function of

t,x,y

Location:

anuga_core

Files:

• source/anuga/abstract_2d_finite_volumes/generic_boundary_conditions.py (modified) (3 diffs)
• source/anuga/abstract_2d_finite_volumes/test_generic_boundary_conditions.py (modified) (1 diff)
• source/anuga_parallel/test_parallel_distribute_domain.py (modified) (1 diff)
• source/anuga_parallel/test_parallel_sw_flow.py (modified) (1 diff)
• test_all.py (modified) (2 diffs)

anuga_core/source/anuga/abstract_2d_finite_volumes/generic_boundary_conditions.py

@@ -110 +110 @@
         self.verbose = verbose
 
-        
+        if domain is None:
+            raise Exception, 'You must specify a domain to Time_boundary'
+
         # FIXME: Temporary code to deal with both f and function
         if function is not None and f is not None:
@@ -142 +144 @@
 
         d = len(domain.conserved_quantities)
-        msg = 'Return value for f must be a list or an array of length %d' %d
+        msg = 'Return value for function must be a list or an array of length %d' %d
         assert len(q) == d, msg
 
@@ -186 +188 @@
         return res
 
-
+class Time_space_boundary(Boundary):
+    """Time and spatially dependent boundary returns values for the
+    conserved quantities as a function of time and space.
+    Must specify domain to get access to model time and a function of t,x,y
+    which must return conserved quantities at specified time and location.
+    
+    Example:
+      B = Time_space_boundary(domain, 
+                              function=lambda t,x,y: [(60<t<3660)*2, 0, 0])
+      
+      This will produce a boundary condition with is a 2m high square wave
+      starting 60 seconds into the simulation and lasting one hour.
+      Momentum applied will be 0 at all times.
+                        
+    """
+
+
+    def __init__(self, domain=None,
+                 function=None,
+                 default_boundary=None,
+                 verbose=False):
+        Boundary.__init__(self)
+
+        self.default_boundary = default_boundary
+        self.default_boundary_invoked = False    # Flag
+        self.domain = domain
+        self.verbose = verbose
+
+        
+        if function is None:
+            raise Exception, 'You must specify a function to Time_space_boundary'
+            
+        try:
+            q = function(0.0, 0.0, 0.0)
+        except Exception, e:
+            msg = 'Function for time_space_boundary could not be executed:\n%s' %e
+            raise msg
+
+
+        try:
+            q = num.array(q, num.float)
+        except:
+            msg = 'Return value from time_space_boundary function could '
+            msg += 'not be converted into a numeric array of floats.\n'
+            msg += 'Specified function should return either list or array.\n'
+            msg += 'I got %s' %str(q)
+            raise msg
+
+        msg = 'ERROR: Time_space_boundary function must return a 1d list or array '
+        assert len(q.shape) == 1, msg
+
+        d = len(domain.conserved_quantities)
+        msg = 'Return value for function must be a list or an array of length %d' %d
+        assert len(q) == d, msg
+
+        self.function = function
+        self.domain = domain
+
+    def __repr__(self):
+        return 'Time space boundary'
+
+    def evaluate(self, vol_id=None, edge_id=None):
+        # FIXME (Ole): I think this should be get_time(), see ticket:306
+
+        [x,y] = self.domain.get_edge_midpoint_coordinate(vol_id, edge_id)
+
+        try:
+            res = self.function(self.domain.get_time(), x, y)
+        except Modeltime_too_early, e:
+            raise Modeltime_too_early, e
+        except Modeltime_too_late, e:
+            if self.default_boundary is None:
+                raise Exception, e # Reraise exception
+            else:
+                # Pass control to default boundary
+                res = self.default_boundary.evaluate(vol_id, edge_id)
+                
+                # Ensure that result cannot be manipulated
+                # This is a real danger in case the 
+                # default_boundary is a Dirichlet type 
+                # for instance. 
+                res = res.copy() 
+                
+                if self.default_boundary_invoked is False:
+                    if self.verbose:                
+                        # Issue warning the first time
+                        msg = '%s' %str(e)
+                        msg += 'Instead I will use the default boundary: %s\n'\
+                            %str(self.default_boundary) 
+                        msg += 'Note: Further warnings will be supressed'
+                        log.critical(msg)
+               
+                    # FIXME (Ole): Replace this crude flag with
+                    # Python's ability to print warnings only once.
+                    # See http://docs.python.org/lib/warning-filter.html
+                    self.default_boundary_invoked = True
+
+        return res
 
 class File_boundary(Boundary):

anuga_core/source/anuga/abstract_2d_finite_volumes/test_generic_boundary_conditions.py

@@ -47 +47 @@
         q = Bd.evaluate()
         assert num.allclose(q, x)
+
+
+    def test_time(self):
+        from quantity import Quantity
+
+        a = [0.0, 0.0]
+        b = [0.0, 2.0]
+        c = [2.0,0.0]
+        d = [0.0, 4.0]
+        e = [2.0, 2.0]
+        f = [4.0,0.0]
+
+        points = [a, b, c, d, e, f]
+
+        #bac, bce, ecf, dbe
+        elements = [ [1,0,2], [1,2,4], [4,2,5], [3,1,4] ]
+
+        domain = Domain(points, elements)
+        domain.check_integrity()
+
+        domain.conserved_quantities = ['stage', 'ymomentum']
+        domain.quantities['stage'] =\
+                                   Quantity(domain, [[1,2,3], [5,5,5],
+                                                     [0,0,9], [-6, 3, 3]])
+
+        domain.quantities['ymomentum'] =\
+                                   Quantity(domain, [[2,3,4], [5,5,5],
+                                                     [0,0,9], [-6, 3, 3]])
+
+
+        domain.check_integrity()
+
+        #Test time bdry, you need to provide a domain and function
+        try:
+            T = Time_boundary(domain)
+        except:
+            pass
+        else:
+            raise 'Should have raised exception'
+
+        #Test time bdry, you need to provide a function
+        try:
+            T = Time_boundary()
+        except:
+            pass
+        else:
+            raise 'Should have raised exception'
+
+
+        def function(t):
+            return [1.0, 0.0]
+        
+        T = Time_boundary(domain, function)
+
+        from anuga.config import default_boundary_tag
+        domain.set_boundary( {default_boundary_tag: T} )
+
+
+        #FIXME: should not necessarily be true always.
+        #E.g. with None as a boundary object.
+        assert len(domain.boundary) == len(domain.boundary_objects)
+
+        q = T.evaluate(0, 2)  #Vol=0, edge=2
+
+        assert num.allclose(q, [1.0, 0.0])
+
+
+    def test_time_space_boundary(self):
+        from quantity import Quantity
+
+        a = [0.0, 0.0]
+        b = [0.0, 2.0]
+        c = [2.0,0.0]
+        d = [0.0, 4.0]
+        e = [2.0, 2.0]
+        f = [4.0,0.0]
+
+        points = [a, b, c, d, e, f]
+
+        #bac, bce, ecf, dbe
+        elements = [ [1,0,2], [1,2,4], [4,2,5], [3,1,4] ]
+
+        domain = Domain(points, elements)
+        domain.check_integrity()
+
+        domain.conserved_quantities = ['stage', 'ymomentum']
+        domain.quantities['stage'] =\
+                                   Quantity(domain, [[1,2,3], [5,5,5],
+                                                     [0,0,9], [-6, 3, 3]])
+
+        domain.quantities['ymomentum'] =\
+                                   Quantity(domain, [[2,3,4], [5,5,5],
+                                                     [0,0,9], [-6, 3, 3]])
+
+
+        domain.check_integrity()
+
+        #Test time space bdry, you need to provide a domain and function
+        try:
+            T = Time_space_boundary(domain)
+        except:
+            pass
+        else:
+            raise 'Should have raised exception'
+
+        #Test time bdry, you need to provide a function
+        try:
+            T = Time_space_boundary()
+        except:
+            pass
+        else:
+            raise 'Should have raised exception'
+
+
+        def function(t,x,y):
+            return [x,y]
+        
+        T = Time_space_boundary(domain, function)
+
+        from anuga.config import default_boundary_tag
+        domain.set_boundary( {default_boundary_tag: T} )
+
+
+        #FIXME: should not necessarily be true always.
+        #E.g. with None as a boundary object.
+        assert len(domain.boundary) == len(domain.boundary_objects)
+
+        q = T.evaluate(0, 2)  #Vol=0, edge=2
+        assert num.allclose(q, domain.get_edge_midpoint_coordinate(0,2))
+
+
+        q = T.evaluate(1, 1)  #Vol=0, edge=2
+        assert num.allclose(q, domain.get_edge_midpoint_coordinate(1,1))        
+
+
+
 
 

anuga_core/source/anuga_parallel/test_parallel_distribute_domain.py

@@ -48 +48 @@
 quantity = 'stage'
 nprocs = 4
-verbose = True
+verbose = False
 
 #--------------------------------------------------------------------------

anuga_core/source/anuga_parallel/test_parallel_sw_flow.py

@@ -37 +37 @@
 N = 25
 M = 25
-verbose =  True
+verbose = False
 
 #---------------------------------

anuga_core/test_all.py

@@ -8 +8 @@
 os.chdir('source')
 os.chdir('anuga')
-    
+print
+print '======================= anuga tests ================================='    
 print 'Changing to', os.getcwd() # This is now different from buildroot   
 execfile('test_all.py')
@@ -16 +17 @@
 os.chdir('pymetis')
 print
-print '*********************************************************************'
+print '===================== pymetis tests ================================='
 print 'Changing to', os.getcwd()
 execfile('test_all.py')
 
-#os.chdir('..')
-#os.chdir('anuga_parallel')
-#print
-#print '*********************************************************************'
-#print 'Changing to', os.getcwd()
-#execfile('test_all.py')
+os.chdir('..')
+os.chdir('anuga_parallel')
+print
+print '===================== anuga_parallel tests =========================='
+print 'Changing to', os.getcwd()
+execfile('test_all.py')
 
 print