Changeset 620 for inundation/ga/storm_surge

Timestamp:

Nov 24, 2004, 12:25:21 PM (20 years ago)

Author:

ole

Message:

Added functionality for storing momentum (x and y) as well as water levels in sww files.

Location:

inundation/ga/storm_surge

Files:

• analytical solutions/Analytical_solution_contracting_channel_import_mesh.py (modified) (1 diff)
• pyvolution/data_manager.py (modified) (3 diffs)
• pyvolution/shallow_water.py (modified) (2 diffs)

inundation/ga/storm_surge/analytical solutions/Analytical_solution_contracting_channel_import_mesh.py

@@ -45 +45 @@
 domain.filename = 'contracting_channel_first-order'
 
+#Decide which quantities are to be stored at each timestep
+domain.quantities_to_be_stored = ['level', 'xmomentum', 'ymomentum']
+
 #Reduction operation for get_vertex_values              
 from util import mean

inundation/ga/storm_surge/pyvolution/data_manager.py

@@ -215 +215 @@
                                 'number_of_points'))
 
+            fid.createVariable('uh', self.precision,
+                               ('number_of_timesteps',
+                                'number_of_points'))
+
+            fid.createVariable('vh', self.precision,
+                               ('number_of_timesteps',
+                                'number_of_points'))
 
         #Close
@@ -240 +247 @@
         y = fid.variables['y']
         z = fid.variables['z']
+        
+        volumes = fid.variables['volumes']
+
+        # Get X, Y and bed elevation Z
+        Q = domain.quantities['elevation']
+        X,Y,Z,V = Q.get_vertex_values(xy=True,
+                                      precision = self.precision)
+                                      
+
+
+        x[:] = X.astype(self.precision)
+        y[:] = Y.astype(self.precision)
+        z[:] = Z.astype(self.precision)
+
+        volumes[:] = V
+        
+        #Close
+        fid.close()
+
+
+    def store_timestep(self, names):
+        """Store time and named quantities to file
+        """
+        from Scientific.IO.NetCDF import NetCDFFile
+        import types
+        from time import sleep
+
+                
+        #Get NetCDF
+        retries = 0
+        file_open = False
+        while not file_open and retries < 10:
+            try:        
+                fid = NetCDFFile(self.filename, 'a') #Open existing file
+            except IOError:
+                #This could happen if someone was reading the file.
+                #In that case, wait a while and try again
+                msg = 'Warning (store_timestep): File %s could not be opened'\
+                %self.filename
+                msg += ' - trying again'
+                print msg
+                retries += 1
+                sleep(1)
+            else:
+               file_open = True
+               
+        if not file_open:
+            msg = 'File %s could not be opened for append' %self.filename
+            raise msg
+                
+            
+        domain = self.domain
+        
+        # Get the variables
+        time = fid.variables['time']
+        stage = fid.variables['stage']
+        uh = fid.variables['uh']
+        vh = fid.variables['vh']                        
+        i = len(time)
+
+        #Store time
+        time[i] = self.domain.time
+
+
+        if type(names) not in [types.ListType, types.TupleType]:
+            names = [names]
+
+        for name in names:    
+            # Get quantity
+            Q = domain.quantities[name]
+            A,V = Q.get_vertex_values(xy=False,
+                                      precision = self.precision)
+
+            #FIXME: Make this more general and rethink naming
+            if name == 'level':
+                stage[i,:] = A.astype(self.precision)
+            elif name == 'xmomentum':
+                uh[i,:] = A.astype(self.precision)
+            elif name == 'ymomentum':
+                vh[i,:] = A.astype(self.precision)                        
+            
+
+        #Flush and close
+        fid.sync()
+        fid.close()
+
+
+#Class for handling checkpoints data
+class Data_format_cpt(Data_format):
+    """Interface to native NetCDF format (.cpt)
+    """
+
+        
+    def __init__(self, domain, mode = 'w'):
+        from Scientific.IO.NetCDF import NetCDFFile
+        from Numeric import Int, Float, Float 
+
+        self.precision = Float #Use full precision
+        
+        Data_format.__init__(self, domain, 'sww', mode)
+
+
+        # NetCDF file definition
+        fid = NetCDFFile(self.filename, mode)
+        
+        if mode == 'w':
+            #Create new file
+            fid.institution = 'Geoscience Australia'
+            fid.description = 'Checkpoint data'
+            #fid.smooth = domain.smooth                
+            fid.order = domain.default_order    
+
+            # dimension definitions
+            fid.createDimension('number_of_volumes', self.number_of_volumes)  
+            fid.createDimension('number_of_vertices', 3)
+
+            #Store info at all vertices (no smoothing)
+            fid.createDimension('number_of_points', 3*self.number_of_volumes)
+            fid.createDimension('number_of_timesteps', None) #extensible
+
+            # variable definitions
+
+            #Mesh
+            fid.createVariable('x', self.precision, ('number_of_points',))
+            fid.createVariable('y', self.precision, ('number_of_points',))
+            #fid.createVariable('z', self.precision, ('number_of_points',))
+        
+            fid.createVariable('volumes', Int, ('number_of_volumes',
+                                                'number_of_vertices'))
+
+            fid.createVariable('time', self.precision,
+                               ('number_of_timesteps',))
+
+            #Allocate space for all quantities
+            for name in domain.quantities.keys():
+                fid.createVariable(name, self.precision,
+                                   ('number_of_timesteps',
+                                    'number_of_points'))
+
+        #Close
+        fid.close()                
+
+
+    def store_checkpoint(self):
+        """
+        Write x,y coordinates of triangles.
+        Write connectivity (
+        constituting
+        the bed elevation.
+        """
+
+        from Scientific.IO.NetCDF import NetCDFFile
+
+        from Numeric import concatenate
+
+        domain = self.domain
+
+        #Get NetCDF
+        fid = NetCDFFile(self.filename, 'a')  #Open existing file for append
+        
+        # Get the variables
+        x = fid.variables['x']
+        y = fid.variables['y']
         
         volumes = fid.variables['volumes']
@@ -311 +481 @@
 
 
-#Class for handling checkpoints data
-class Data_format_cpt(Data_format):
-    """Interface to native NetCDF format (.cpt)
-    """
-
-        
-    def __init__(self, domain, mode = 'w'):
-        from Scientific.IO.NetCDF import NetCDFFile
-        from Numeric import Int, Float, Float 
-
-        self.precision = Float #Use full precision
-        
-        Data_format.__init__(self, domain, 'sww', mode)
-
-
-        # NetCDF file definition
-        fid = NetCDFFile(self.filename, mode)
-        
-        if mode == 'w':
-            #Create new file
-            fid.institution = 'Geoscience Australia'
-            fid.description = 'Checkpoint data'
-            #fid.smooth = domain.smooth                
-            fid.order = domain.default_order    
-
-            # dimension definitions
-            fid.createDimension('number_of_volumes', self.number_of_volumes)  
-            fid.createDimension('number_of_vertices', 3)
-
-            #Store info at all vertices (no smoothing)
-            fid.createDimension('number_of_points', 3*self.number_of_volumes)
-            fid.createDimension('number_of_timesteps', None) #extensible
-
-            # variable definitions
-
-            #Mesh
-            fid.createVariable('x', self.precision, ('number_of_points',))
-            fid.createVariable('y', self.precision, ('number_of_points',))
-            #fid.createVariable('z', self.precision, ('number_of_points',))
-        
-            fid.createVariable('volumes', Int, ('number_of_volumes',
-                                                'number_of_vertices'))
-
-            fid.createVariable('time', self.precision,
-                               ('number_of_timesteps',))
-
-            #Allocate space for all quantities
-            for name in domain.quantities.keys():
-                fid.createVariable(name, self.precision,
-                                   ('number_of_timesteps',
-                                    'number_of_points'))
-
-        #Close
-        fid.close()                
-
-
-    def store_checkpoint(self):
-        """
-        Write x,y coordinates of triangles.
-        Write connectivity (
-        constituting
-        the bed elevation.
-        """
-
-        from Scientific.IO.NetCDF import NetCDFFile
-
-        from Numeric import concatenate
-
-        domain = self.domain
-
-        #Get NetCDF
-        fid = NetCDFFile(self.filename, 'a')  #Open existing file for append
-        
-        # Get the variables
-        x = fid.variables['x']
-        y = fid.variables['y']
-        
-        volumes = fid.variables['volumes']
-
-        # Get X, Y and bed elevation Z
-        Q = domain.quantities['elevation']
-        X,Y,Z,V = Q.get_vertex_values(xy=True,
-                                      precision = self.precision)
-                                      
-
-
-        x[:] = X.astype(self.precision)
-        y[:] = Y.astype(self.precision)
-        z[:] = Z.astype(self.precision)
-
-        volumes[:] = V
-        
-        #Close
-        fid.close()
-
-
-    def store_timestep(self, name):
-        """Store time and named quantity to file
-        """
-        from Scientific.IO.NetCDF import NetCDFFile
-        from time import sleep
-                
-        #Get NetCDF
-        retries = 0
-        file_open = False
-        while not file_open and retries < 10:
-            try:        
-                fid = NetCDFFile(self.filename, 'a') #Open existing file
-            except IOError:
-                #This could happen if someone was reading the file.
-                #In that case, wait a while and try again
-                msg = 'Warning (store_timestep): File %s could not be opened'\
-                %self.filename
-                msg += ' - trying again'
-                print msg
-                retries += 1
-                sleep(1)
-            else:
-               file_open = True
-               
-        if not file_open:
-            msg = 'File %s could not be opened for append' %self.filename
-            raise msg
-                
-            
-        domain = self.domain
-        
-        # Get the variables
-        time = fid.variables['time']
-        stage = fid.variables['stage']        
-        i = len(time)
-
-        #Store stage
-        time[i] = self.domain.time
-
-        # Get quantity
-        Q = domain.quantities[name]
-        A,V = Q.get_vertex_values(xy=False,
-                                  precision = self.precision)
-        
-        stage[i,:] = A.astype(self.precision)
-
-        #Flush and close
-        fid.sync()
-        fid.close()
-
-
 
 

inundation/ga/storm_surge/pyvolution/shallow_water.py

@@ -52 +52 @@
         self.reduction = min  #Looks better near steep slopes
         
-
+        self.quantities_to_be_stored = ['level']
+
+        
         #Establish shortcuts to relevant quantities (for efficiency)
         #self.w = self.quantities['level']
@@ -172 +174 @@
             #Store model data, e.g. for subsequent visualisation    
             if self.store is True:
-                self.store_timestep('level')
+                #self.store_timestep(['level', 'xmomentum', 'ymomentum'])
+                self.store_timestep(self.quantities_to_be_stored)
+                
                 #FIXME: Could maybe be taken from specified list 
                 #of 'store every step' quantities       
-
 
             #Pass control on to outer loop for more specific actions