Changeset 389 for inundation/ga

Timestamp:

Oct 12, 2004, 2:45:24 PM (20 years ago)

Author:

duncan

Message:

add vert atts to pmesh2domain/pyvolution

Location:

inundation/ga/storm_surge/pyvolution

Files:

• domain.py (modified) (1 diff)
• interpolate_sww.py (modified) (1 diff)
• least_squares.py (modified) (1 diff)
• pmesh2domain.py (modified) (4 diffs)
• quantity.py (modified) (3 diffs)
• run_tsh.py (modified) (4 diffs)

inundation/ga/storm_surge/pyvolution/domain.py

@@ -125 +125 @@
         return q 
 
+
+    def set_quantity_vertices_dict(self, quantity_dict):
+        """Set values for named quantities.
+        The index is the quantity
+        
+        name: Name of quantity
+        X: Compatible list, Numeric array, const or function (see below)
+        
+        The values will be stored in elements following their
+        internal ordering.
+
+        """
+        for key in quantity_dict.keys():
+            self.set_quantity(key,quantity_dict[key], location='vertices')
+        
     def set_quantity(self, name, X, location='vertices'):
         """Set values for named quantity

inundation/ga/storm_surge/pyvolution/interpolate_sww.py

@@ -116 +116 @@
         volumes = fid.variables['volumes'][:] 
         time = fid.variables['time'][:]
-        stage = fid.variables['stage'][:,:]   # 2D  
+        stage = fid.variables['stage'][:,:]   # 2D 
         
         fid.close()

inundation/ga/storm_surge/pyvolution/least_squares.py

@@ -62 +62 @@
     point_dict = load_xya_file(point_file)
     point_coordinates = point_dict['pointlist']
-    point_attributes = point_dict['pointattributelist']
     point_attributes = point_dict['pointattributelist']
     title_string = point_dict['title']

inundation/ga/storm_surge/pyvolution/pmesh2domain.py

@@ -8 +8 @@
 
 
-def pmesh_to_domain(fileName, tags = None, setting_function = None):
+def pmesh_to_domain_instance(fileName, DomainClass,tags = None, setting_function = None):
     """
     """
-    #FIXME:  The current design doesn't seem to accomadate tags and
-    # setting_functions being passed into the domain at this point.
+    #FIXME:  The current design doesn't seem to accomadate tags
+    #  being passed into the domain at this point.
 
-    #FIXME: Plus the names of the functions are no longer appropriate,
-    #since domain objects aren't being returned.
-    
     import sys
-    from config import pmesh_filename
-    sys.path.append(pmesh_filename)
-    from load_mesh.loadASCII import import_trianglulation
+    from domain import Domain
 
-    try:
-        meshdic = import_trianglulation(fileName)
-    except IOError, e:        
-        msg = 'Could not open file %s ' %fileName
-        raise IOError, msg 
+    vertex_coordinates, volumes, marker_dict, vertex_quantity_dict = \
+                        pmesh_to_domain(fileName,
+                                                   setting_function=setting_function)
+        
+    assert issubclass(DomainClass, Domain), "DomainClass is not a subclass of Domain."
+    domain = DomainClass(vertex_coordinates, volumes, marker_dict)
 
-        
-    return pmesh_dictionary_to_domain(meshdic,
-                                      setting_function = setting_function)
+    # set the water level to be the elevation
+    if vertex_quantity_dict.has_key('elevation') and not vertex_quantity_dict.has_key('level'):
+        vertex_quantity_dict['level'] = vertex_quantity_dict['elevation']
+    domain.set_quantity_vertices_dict(vertex_quantity_dict)
+    #print "vertex_quantity_dict",vertex_quantity_dict
+    return domain
 
 
 
-def pmesh_dictionary_to_domain(meshdic, setting_function = None):
+def pmesh_to_domain(fileName, setting_function = None):
     """
     convert a pmesh dictionary to a list of Volumes.
     Also, return a list of triangles which have boundary tags
-    meshdic structure;
+    mesh_dict structure;
     generated point list: [(x1,y1),(x2,y2),...] (Tuples of doubles)  
     generated point attribute list:[(P1att1,P1attt2, ...),(P2att1,P2attt2,...),...]
@@ -47 +46 @@
     triangle attribute list: [T1att, T2att, ...] (list of floats)
     """
-
-    vertex_coordinates = meshdic['generatedpointlist']
-    volumes = meshdic['generatedtrianglelist']
+    
+    from Numeric import transpose
+    from load_mesh.loadASCII import mesh_file_to_mesh_dictionary
+    
+    mesh_dict = mesh_file_to_mesh_dictionary(fileName)
+    #print "mesh_dict",mesh_dict
+    vertex_coordinates = mesh_dict['generatedpointlist']
+    volumes = mesh_dict['generatedtrianglelist']
 
     #if setting_function:
@@ -55 +59 @@
     #        setting_function = [setting_function]
     #    for funct in setting_function:
-    #        meshdic = funct(meshdic, vertices = mesh_vertices,
+    #        mesh_dict = funct(mesh_dict, vertices = mesh_vertices,
     #                        volumes = volumes)
-        
-    marker_dict = pmesh_dict_to_marker_dict(volumes, meshdic,
+
+
+    vertex_quantity_dict = {}
+    point_atts = transpose(mesh_dict['generatedpointattributelist'])
+    point_titles  = mesh_dict['generatedpointattributetitlelist']
+    #print "point_titles",point_titles
+    for quantity, value_vector in map (None, point_titles, point_atts):
+        vertex_quantity_dict[quantity] = value_vector
+    marker_dict = pmesh_dict_to_marker_dict(volumes, mesh_dict,
                                             vertex_coordinates)
-
-    return vertex_coordinates, volumes, marker_dict
+    
+    return vertex_coordinates, volumes, marker_dict, vertex_quantity_dict
 
 
@@ -75 +86 @@
 
    
-def pmesh_dict_to_marker_dict(triangles,meshdic,Vertices):
-    """ Convert the pmesh dictionary (meshdic) description of boundary tags
+def pmesh_dict_to_marker_dict(triangles,mesh_dict,Vertices):
+    """ Convert the pmesh dictionary (mesh_dict) description of boundary tags
     to a dictionary of markers, indexed with volume id and face number.
     """
     sides = calc_sides(triangles)
     marker_dict = {}
-    for seg, marker in map(None,meshdic['generatedsegmentlist'],
-                           meshdic['generatedsegmentmarkerlist']):
+    for seg, marker in map(None,mesh_dict['generatedsegmentlist'],
+                           mesh_dict['generatedsegmentmarkerlist']):
         v1 = seg[0]
         v2 = seg[1]

inundation/ga/storm_surge/pyvolution/quantity.py

@@ -56 +56 @@
     def __len__(self):
         return self.centroid_values.shape[0]
-        
+    
     def interpolate(self):
         """Compute interpolated values at edges and centroid
@@ -201 +201 @@
             if len(values.shape) == 1:
                 #Values are being specified once for each unique vertex
-                msg = 'Number of values must mact number of vertices'
+                msg = 'Number of values must match number of vertices'
                 assert values.shape[0] == self.domain.coordinates.shape[0], msg
 
@@ -218 +218 @@
                 
 
-
+    # FIXME have a get_vertex_values as well, so the 'level' quantity can be
+    # set, based on the elevation    
     def set_vertex_values(self, A):
         """Set vertex values for all triangles based on input array A

inundation/ga/storm_surge/pyvolution/run_tsh.py

@@ -16 +16 @@
 
 from mesh_factory import rectangular
-from pmesh2domain import pmesh_to_domain
+from pmesh2domain import pmesh_to_domain, pmesh_to_domain_instance
 
 from Numeric import array
@@ -30 +30 @@
 
 print 'Creating domain from', filename
-points, vertices, boundary = pmesh_to_domain(filename)
-domain = Domain(points, vertices, boundary)
+domain = pmesh_to_domain_instance(filename, Domain)
 print "Number of triangles = ", len(domain)
 
@@ -48 +47 @@
 print 'Field values'
 
-domain.set_quantity('elevation', W)
+#domain.set_quantity('elevation', W)
 domain.set_quantity('friction', manning)
 
@@ -73 +72 @@
 
 #Set boundary conditions
-domain.set_boundary({'left': Bw, '0': Br})
+domain.set_boundary({'left': Bw, '0': Br, '1':Bw, 'external':Br})
 
+
+print domain.quantities['elevation'].vertex_values
+#print domain.quantities['level'].vertex_values
+         
 ######################
 #Initial condition
 print 'Initial condition'
-domain.set_quantity('level', Constant_height(W, 0.))
+#domain.set_quantity('level', Constant_height(W, 0.))
 domain.check_integrity()