Changeset 1134

Timestamp:

Mar 23, 2005, 11:21:51 AM (20 years ago)

Author:

prow

Message:

testing sww2domain

File:

• inundation/ga/storm_surge/pyvolution/test_data_manager.py (modified) (4 diffs)

inundation/ga/storm_surge/pyvolution/test_data_manager.py

@@ -728 +728 @@
 
         from coordinate_transforms.redfearn import redfearn
-
+        import os
         fid1 = NetCDFFile('test_ha.nc','w')
         fid2 = NetCDFFile('test_ua.nc','w')
@@ -836 +836 @@
                    origin = (56, 0, 0))
 
+        os.remove('test_va.nc')
+        os.remove('test_ua.nc')
+        os.remove('test_ha.nc')
+        os.remove('test_e.nc')
 
         #Read output file 'test.sww'
@@ -859 +863 @@
 
         #Cleanup
-        import os
         os.remove('test.sww')
-        os.remove('test_va.nc')
-        os.remove('test_ua.nc')
-        os.remove('test_ha.nc')
 
 
@@ -939 +939 @@
         os.remove('small.sww')
 
+    def test_sww2domain(self):
+        ################################################
+        #Create a test domain, and evolve and save it.
+        ################################################
+        from mesh_factory import rectangular
+        from shallow_water import Domain, Reflective_boundary, Dirichlet_boundary,\
+             Constant_height, Time_boundary, Transmissive_boundary
+        from Numeric import array
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(2,2)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.smooth = False
+        domain.visualise = False
+        domain.store = True
+        domain.filename = 'bedslope'
+        domain.default_order=2
+        #Bed-slope and friction
+        domain.set_quantity('elevation', lambda x,y: -x/3)
+        domain.set_quantity('friction', 0.1)
+        # Boundary conditions
+        from math import sin, pi
+        Br = Reflective_boundary(domain)
+        Bt = Transmissive_boundary(domain)
+        Bd = Dirichlet_boundary([0.2,0.,0.])
+        Bw = Time_boundary(domain=domain,
+                           f=lambda t: [(0.1*sin(t*2*pi)), 0.0, 0.0])
+
+        domain.set_boundary({'left': Bd, 'right': Br, 'top': Br, 'bottom': Br})
+        domain.quantities_to_be_stored.extend(['xmomentum','ymomentum'])
+        #Initial condition
+        h = 0.05
+        elevation = domain.quantities['elevation'].vertex_values
+        domain.set_quantity('stage', elevation + h)
+        #elevation = domain.get_quantity('elevation',location='unique vertices')
+        #domain.set_quantity('stage', elevation + h,location='unique vertices')
+
+        domain.check_integrity()
+        dir(domain)
+        #Evolution
+        for t in domain.evolve(yieldstep = 1, finaltime = 2.0):
+        #    domain.write_time()
+            pass
+
+
+        ##########################################
+        #Import the example's file as a new domain
+        ##########################################
+        from data_manager import sww2domain
+        from Numeric import allclose
+
+        filename = domain.datadir+'\\'+domain.filename+'.sww'
+
+        domain2 = sww2domain(filename,fail_if_NaN=False,verbose = False)
+
+        ###################
+        ##NOW TEST IT!!!
+        ##################
+
+        bits = ['xllcorner','yllcorner','vertex_coordinates','time','starttime']
+
+        for quantity in ['elevation']+domain.quantities_to_be_stored:
+            bits.append('get_quantity("%s")'%quantity)
+
+        for bit in bits:
+        #    print 'testing that domain.'+bit+' has been restored'
+            assert allclose(eval('domain.'+bit),eval('domain2.'+bit))
+
+        #print 'passed'
+
+
+    def test_sww2domain2(self):
+        ##################################################################
+        #Same as previous test, but this checks how NaNs are handled.
+        ##################################################################
+
+
+        from mesh_factory import rectangular
+        from shallow_water import Domain, Reflective_boundary, Dirichlet_boundary,\
+             Constant_height, Time_boundary, Transmissive_boundary
+        from Numeric import array
+
+        #Create basic mesh
+        points, vertices, boundary = rectangular(2,2)
+
+        #Create shallow water domain
+        domain = Domain(points, vertices, boundary)
+        domain.smooth = False
+        domain.visualise = False
+        domain.store = True
+        domain.filename = 'bedslope'
+        domain.default_order=2
+        domain.quantities_to_be_stored=['stage']
+
+        domain.set_quantity('elevation', lambda x,y: -x/3)
+        domain.set_quantity('friction', 0.1)
+
+        from math import sin, pi
+        Br = Reflective_boundary(domain)
+        Bt = Transmissive_boundary(domain)
+        Bd = Dirichlet_boundary([0.2,0.,0.])
+        Bw = Time_boundary(domain=domain,
+                           f=lambda t: [(0.1*sin(t*2*pi)), 0.0, 0.0])
+
+        domain.set_boundary({'left': Bd, 'right': Br, 'top': Br, 'bottom': Br})
+
+        h = 0.05
+        elevation = domain.quantities['elevation'].vertex_values
+        domain.set_quantity('stage', elevation + h)
+        #elevation = domain.get_quantity('elevation',location='unique vertices')
+        #domain.set_quantity('stage', elevation + h,location='unique vertices')
+
+        domain.check_integrity()
+        
+        for t in domain.evolve(yieldstep = 1, finaltime = 2.0):
+            pass
+            #domain.write_time()
+
+        ##################################
+        #Import the file as a new domain
+        ##################################
+        from data_manager import sww2domain
+        from Numeric import allclose
+
+        filename = domain.datadir+'\\'+domain.filename+'.sww'
+
+        #Fail because NaNs are present
+        try:
+            domain2 = sww2domain(filename,fail_if_NaN=True,verbose=False)
+            assert True == False
+        except:
+            #Now import it, filling NaNs to be 0
+            filler = 0
+            domain2 = sww2domain(filename,fail_if_NaN=False,NaN_filler = filler,verbose=False)
+
+        bits = ['xllcorner','yllcorner','vertex_coordinates','time','starttime']
+
+        for quantity in ['elevation']+domain.quantities_to_be_stored:
+            bits.append('get_quantity("%s")'%quantity)
+
+        for bit in bits:
+        #    print 'testing that domain.'+bit+' has been restored'
+            assert allclose(eval('domain.'+bit),eval('domain2.'+bit))
+
+        #print max(max(domain2.get_quantity('xmomentum')))
+        #print min(min(domain2.get_quantity('xmomentum')))
+        #print max(max(domain2.get_quantity('ymomentum')))
+        #print min(min(domain2.get_quantity('ymomentum')))
+
+        assert max(max(domain2.get_quantity('xmomentum')))==filler
+        assert min(min(domain2.get_quantity('xmomentum')))==filler
+        assert max(max(domain2.get_quantity('ymomentum')))==filler
+        assert min(min(domain2.get_quantity('ymomentum')))==filler
+
+        #print 'passed'
+
+        #cleanup
+        #import os
+        #os.remove(domain.datadir+'/'+domain.filename+'.sww')
 
 #-------------------------------------------------------------