Changeset 288

Timestamp:

Sep 10, 2004, 5:56:19 PM (20 years ago)

Author:

ole

Message:

Started playing with checkpointing

Location:

inundation/ga/storm_surge/pyvolution

Files:

• data_manager.py (modified) (3 diffs)
• mesh.py (modified) (1 diff)
• netherlands.py (modified) (2 diffs)
• quantity.py (modified) (2 diffs)

inundation/ga/storm_surge/pyvolution/data_manager.py

@@ -284 +284 @@
         self.domain = domain
     
-
-    # def store_connectivity(self):
-#         """Writes x,y,z coordinates of triangles constituting
-#         the bed elevation.
-#         """
-
-#         msg = 'Method should be overridden in concrete class'
-#         raise msg
-
-#     def store_timestep(self, t, V0, V1, V2):
-#         """Store time, water heights (and momentums) to file
-#         """
-
-#         msg = 'Method should be overridden in concrete class'
-#         raise msg        
     
     #FIXME: Should we have a general set_precision function?
     
 
-#Function for storing output to e.g. visualisation
+
+#Class for storing output to e.g. visualisation
 class Data_format_sww(Data_format):
     """Interface to native NetCDF format (.sww)
@@ -325 +311 @@
             #Create new file
             fid.institution = 'Geoscience Australia'
+            fid.description = 'Output from pyvolution suitable for plotting'
+            
             if domain.smooth:
                 fid.smoothing = 'Yes'
@@ -454 +442 @@
 
 
+#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/mesh.py

@@ -541 +541 @@
         return vertex_coordinates
   
-
-
-
+    def get_vertices(unique = False):
+        """Get connectivity
+        If unique is True give them only once as stored internally.
+        If unique is False return
+        """
+
+        if unique is True:
+            V = self.domain.vertices
+        else:
+            
+
+
+#FIXME: May get rid of
 def rectangular_mesh(m, n, len1=1.0, len2=1.0, origin = (0.0, 0.0)):
     from mesh_factory import rectangular

inundation/ga/storm_surge/pyvolution/netherlands.py

@@ -92 +92 @@
 
 N = 600 #Size = 720000
-N = 100
+N = 220
 
 
@@ -171 +171 @@
 t0 = time.time()
 
-for t in domain.evolve(yieldstep = 0.1, finaltime = 1.5):
+for t in domain.evolve(yieldstep = 0.01, finaltime = 1.0):
     domain.write_time()
     

inundation/ga/storm_surge/pyvolution/quantity.py

@@ -234 +234 @@
     
     #Method for outputting model results
+    #FIXME: Split up into geometric and numeric stuff.
+    #FIXME: Geometric (X,Y,V) should live in mesh.py
     def get_vertex_values(self,
                           xy=True, 
@@ -282 +284 @@
             N = len(self.domain.vertexlist)
             A = zeros(N, precision)
-            V = self.domain.vertices
+            V = self.domain.get_vertices(smooth=True)
             
             #Smoothing loop