Changeset 4437

Timestamp:

May 15, 2007, 4:31:03 PM (18 years ago)

Author:

ole

Message:

Enabled test using new tight slope limiters. It failed before.

Location:

anuga_core/source/anuga

Files:

• abstract_2d_finite_volumes/domain.py (modified) (3 diffs)
• shallow_water/test_data_manager.py (modified) (17 diffs)
• shallow_water/test_shallow_water_domain.py (modified) (1 diff)

anuga_core/source/anuga/abstract_2d_finite_volumes/domain.py

@@ -827 +827 @@
 
         if finaltime is not None and duration is not None:
-            print 'F', finaltime, duration
+            #print 'F', finaltime, duration
             msg = 'Only one of finaltime and duration may be specified'
             raise msg
@@ -864 +864 @@
 
         while True:
-
             #Compute fluxes across each element edge
             self.compute_fluxes()
@@ -983 +982 @@
                     timestep = min_timestep  #Try enforcing min_step
 
-                    #raise msg
+
+                    print self.timestepping_statistics(track_speeds=True)
+
+
+                    raise Exception
                 else:
                     #Try to overcome situation by switching to 1 order

anuga_core/source/anuga/shallow_water/test_data_manager.py

@@ -43 +43 @@
         domain = Domain(points, vertices, boundary)
         domain.default_order = 2
-
 
         #Set some field values
@@ -342 +341 @@
         sww.store_connectivity()
         sww.store_timestep('stage')
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -396 +396 @@
         sww.store_connectivity()
         sww.store_timestep('stage')
-
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -449 +449 @@
         self.domain.store = True
         self.domain.beta_h = 0
-
+        #self.domain.limit2007 = 1
 
         #Evolution
         for t in self.domain.evolve(yieldstep = 1.0, finaltime = 4.0):
+            
             #########self.domain.write_time(track_speeds=True)
             stage = self.domain.quantities['stage'].vertex_values
@@ -459 +460 @@
             fid = NetCDFFile(self.domain.writer.filename, 'r')
             stage_file = fid.variables['stage']
-
+            
             if t == 0.0:
                 assert allclose(stage, self.initial_stage)
@@ -491 +492 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1118 +1120 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
+
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1291 +1295 @@
         sww.store_connectivity()
         sww.store_timestep('stage')
-
+        
+        domain.limit2007 = 1
         domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1476 +1481 @@
         sww.store_timestep('stage')
 
+        #domain.limit2007 = 1
         domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1624 +1630 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1735 +1742 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -1993 +2001 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -2092 +2101 @@
         sww.store_timestep('stage')
 
+        #self.domain.limit2007 = 1
         self.domain.evolve_to_end(finaltime = 0.01)
         sww.store_timestep('stage')
@@ -2753 +2763 @@
         self.domain.reduction = mean
         self.domain.set_datadir('.')
+        #self.domain.limit2007 = 1        
 
 
@@ -2775 +2786 @@
         assert allclose(ymin, 0.0)
         assert allclose(ymax, 1.0)
-        assert allclose(stagemin, -0.85)
+        assert allclose(stagemin, -0.85), 'stagemin=%.4f' %stagemin
         assert allclose(stagemax, 0.15)
 
@@ -2824 +2835 @@
         domain.check_integrity()
         #Evolution
+        #domain.limit2007 = 1
         for t in domain.evolve(yieldstep = yiel, finaltime = 0.05):
             #domain.write_time()
@@ -6022 +6034 @@
         pass
     suite = unittest.makeSuite(Test_Data_Manager,'test')
-    runner = unittest.TextTestRunner() #verbosity=2)
+    runner = unittest.TextTestRunner()#(verbosity=2)
     runner.run(suite)
 

anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

@@ -4198 +4198 @@
 
 
+
+    def test_tight_slope_limiters(self):
+        """Test that new slope limiters (Feb 2007) don't induce extremely
+        small timesteps. This test actually reveals the problem as it
+        was in March-April 2007 
+        """
+
+        import time, os
+        from Numeric import array, zeros, allclose, Float, concatenate
+        from Scientific.IO.NetCDF import NetCDFFile
+        from data_manager import get_dataobject, extent_sww
+        from mesh_factory import rectangular
+
+        
+        #Create basic mesh
+        points, vertices, boundary = rectangular(2, 2)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.default_order = 2
+        domain.limit2007 = 1
+        domain.H0 = 0.01
+        
+
+
+        #Set some field values
+        domain.set_quantity('elevation', lambda x,y: -x)
+        domain.set_quantity('friction', 0.03)
+
+
+        ######################
+        # Boundary conditions
+        B = Transmissive_boundary(domain)
+        domain.set_boundary( {'left': B, 'right': B, 'top': B, 'bottom': B})
+
+
+        ######################
+        #Initial condition - with jumps
+
+
+        bed = domain.quantities['elevation'].vertex_values
+        stage = zeros(bed.shape, Float)
+
+        h = 0.3
+        for i in range(stage.shape[0]):
+            if i % 2 == 0:
+                stage[i,:] = bed[i,:] + h
+            else:
+                stage[i,:] = bed[i,:]
+
+        domain.set_quantity('stage', stage)
+
+
+        domain.distribute_to_vertices_and_edges()               
+
+        
+
+        domain.set_name('tight_limiters')
+        domain.format = 'sww'
+        domain.smooth = True
+        domain.reduction = mean
+        domain.set_datadir('.')
+        domain.smooth = False
+        domain.store = True
+        domain.beta_h = 0
+        
+
+        #Evolution
+        for t in domain.evolve(yieldstep = 0.1, finaltime = 0.3):
+            
+            #domain.write_time(track_speeds=True)
+            stage = domain.quantities['stage'].vertex_values
+
+            #Get NetCDF
+            fid = NetCDFFile(domain.writer.filename, 'r')
+            stage_file = fid.variables['stage']
+            
+            fid.close()
+
+        os.remove(domain.writer.filename)
+
+
     def test_pmesh2Domain(self):
          import os