Changeset 5729 for anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

Timestamp:

Sep 4, 2008, 2:45:27 PM (16 years ago)

Author:

ole

Message:

Refactored get_flow_through_cross_section and added a runtime version to
shallow_water_domain. It still needs some work and testing.

File:

• anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py (modified) (3 diffs)

anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

@@ -1406 +1406 @@
             
         
+        
+    def test_get_flow_through_cross_section_with_geo(self):
+        """test_get_flow_through_cross_section(self):
+
+        Test that the total flow through a cross section can be
+        correctly obtained at run-time from the ANUGA domain.
+        
+        This test creates a flat bed with a known flow through it and tests
+        that the function correctly returns the expected flow.
+
+        The specifics are
+        e = -1 m
+        u = 2 m/s
+        h = 2 m
+        w = 3 m (width of channel)
+
+        q = u*h*w = 12 m^3/s
+
+
+        This run tries it with georeferencing and with elevation = -1
+        
+        """
+
+        import time, os
+        from Numeric import array, zeros, allclose, Float, concatenate
+        from Scientific.IO.NetCDF import NetCDFFile
+
+        # Setup
+        from mesh_factory import rectangular
+
+        # Create basic mesh (20m x 3m)
+        width = 3
+        length = 20
+        t_end = 1
+        points, vertices, boundary = rectangular(length, width,
+                                                 length, width)
+
+        # Create shallow water domain
+        domain = Domain(points, vertices, boundary,
+                        geo_reference=Geo_reference(56,308500,6189000))
+
+        domain.default_order = 2
+        domain.set_quantities_to_be_stored(None)
+
+
+        e = -1.0
+        w = 1.0
+        h = w-e
+        u = 2.0
+        uh = u*h
+
+        Br = Reflective_boundary(domain)     # Side walls
+        Bd = Dirichlet_boundary([w, uh, 0])  # 2 m/s across the 3 m inlet: 
+
+
+        # Initial conditions
+        domain.set_quantity('elevation', e)
+        domain.set_quantity('stage', w)
+        domain.set_quantity('xmomentum', uh)
+        domain.set_boundary( {'left': Bd, 'right': Bd, 'top': Br, 'bottom': Br})
+        
+        
+        # Interpolation points down the middle
+        I = [[0, width/2.],
+             [length/2., width/2.],
+             [length, width/2.]]
+        interpolation_points = domain.geo_reference.get_absolute(I)        
+        
+        # Shortcuts to quantites
+        stage = domain.get_quantity('stage')        
+        xmomentum = domain.get_quantity('xmomentum')       
+        ymomentum = domain.get_quantity('ymomentum')               
+
+        for t in domain.evolve(yieldstep=0.1, finaltime = t_end):
+            # Check that quantities are they should be in the interior
+
+            w_t = stage.get_values(interpolation_points)            
+            uh_t = xmomentum.get_values(interpolation_points)
+            vh_t = ymomentum.get_values(interpolation_points)            
+            
+            assert allclose(w_t, w)
+            assert allclose(uh_t, uh)            
+            assert allclose(vh_t, 0.0)                        
+            
+            
+            # Check flows through the middle
+            for i in range(5):
+                x = length/2. + i*0.23674563 # Arbitrary
+                cross_section = [[x, 0], [x, width]]
+
+                cross_section = domain.geo_reference.get_absolute(cross_section)            
+                Q = domain.get_flow_through_cross_section(cross_section,
+                                                          verbose=False)
+
+                assert allclose(Q, uh*width)
+
+
+            
+        
+        
 
     def test_another_runup_example(self):
@@ -1431 +1531 @@
         domain = Domain(points, vertices, boundary) # Create domain
         domain.set_quantities_to_be_stored(None)
-        domain.set_maximum_allowed_speed(100) #
+        domain.set_maximum_allowed_speed(100) #FIXME (Ole): try to remove this
         
         # FIXME (Ole): Need tests where this is commented out
-        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                 
+        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)
         domain.H0 = 0 # Backwards compatibility (6/2/7)
         domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
@@ -5805 +5905 @@
 if __name__ == "__main__":
 
-    suite = unittest.makeSuite(Test_Shallow_Water,'test')
-
-    #suite = unittest.makeSuite(Test_Shallow_Water,'test_bedslope_problem_first_order_moresteps')    
+    #suite = unittest.makeSuite(Test_Shallow_Water,'test')
+
+    suite = unittest.makeSuite(Test_Shallow_Water,'test_get_flow_through_cross_section_with_geo')    
     #suite = unittest.makeSuite(Test_Shallow_Water,'test_fitting_using_shallow_water_domain')    
     #suite = unittest.makeSuite(Test_Shallow_Water,'test_tight_slope_limiters')