Changeset 2056

Timestamp:

Nov 24, 2005, 10:27:59 AM (19 years ago)

Author:

duncan

Message:

create a domain instance from a pmesh.mesh instance

Location:

inundation

Files:

• examples/run_gippsland.py (modified) (1 diff)
• pyvolution/pmesh2domain.py (modified) (5 diffs)
• pyvolution/test_pmesh2domain.py (modified) (3 diffs)

inundation/examples/run_gippsland.py

@@ -1 +1 @@
 from pyvolution.data_manager import asc_csiro2sww
+import time
+
+t0 = time.time()
+
 
 asc_csiro2sww('bathymetry_expanded','elev_expanded','ucur_expanded',
-              'vcur_expanded', 'test.sww',zscale=1000, mean_stage = 100,
+              'vcur_expanded', 'test_extent_true_38.25__37.7__147__148.25.sww',zscale=1,
+              mean_stage = 0,
               fail_on_NaN = False,
-              elevation_NaN_filler = 0)
+              elevation_NaN_filler = 0
+              ,
+              minlat = -38.25, maxlat = -37.7, ###-37.75,
+              minlon = 147.0, maxlon = 148.25
+              , verbose = True
+              )
+
+print 'That took %.2f seconds' %(time.time()-t0)

inundation/pyvolution/pmesh2domain.py

@@ -7 +7 @@
 """
 
+# FXME (DSG-DSG): Combine with function pmesh_to_domain_instance
+def pmesh_instance_to_domain_instance(mesh,
+                                      DomainClass,
+                                      setting_function=None):
+    """
+    """
+    import sys
+    from domain import Domain
+
+    vertex_coordinates, vertices, tag_dict, vertex_quantity_dict \
+                        ,tagged_elements_dict, geo_reference = \
+                pmesh_to_domain(mesh_instance=mesh,
+                                setting_function=setting_function)
+
+
+    assert issubclass(DomainClass, Domain),"DomainClass is not a subclass of Domain."
+
+
+    domain = DomainClass(coordinates = vertex_coordinates,
+                         vertices = vertices,
+                         boundary = tag_dict,
+                         tagged_elements = tagged_elements_dict,
+                         geo_reference = geo_reference )
+
+    # set the water stage to be the elevation
+    if vertex_quantity_dict.has_key('elevation') and not vertex_quantity_dict.has_key('stage'):
+        vertex_quantity_dict['stage'] = vertex_quantity_dict['elevation']
+
+    domain.set_quantity_vertices_dict(vertex_quantity_dict)
+    #print "vertex_quantity_dict",vertex_quantity_dict
+    return domain
 
 def pmesh_to_domain_instance(fileName, DomainClass, setting_function = None):
@@ -16 +47 @@
     vertex_coordinates, vertices, tag_dict, vertex_quantity_dict \
                         ,tagged_elements_dict, geo_reference = \
-                pmesh_to_domain(fileName, setting_function=setting_function)
+                pmesh_to_domain(fileName=fileName,
+                                setting_function=setting_function)
 
 
@@ -52 +84 @@
 
 
-def pmesh_to_domain(fileName, setting_function = None):
+def pmesh_to_domain(fileName=None,
+                    mesh_instance=None,
+                    setting_function=None):
     """
     convert a pmesh dictionary to a list of Volumes.
@@ -69 +103 @@
     from load_mesh.loadASCII import import_mesh_file
 
-    mesh_dict = import_mesh_file(fileName)
+    if fileName is None:
+        mesh_dict = mesh_instance.Mesh2IODict()
+    else:
+        mesh_dict = import_mesh_file(fileName)
     #print "mesh_dict",mesh_dict
     vertex_coordinates = mesh_dict['vertices']
@@ -98 +135 @@
     #print "DSG pm2do tagged_elements", tagged_elements
     return tagged_elements
-
-#FIXME: The issue is whether this format should be stored in the tsh file
-#instead of having to be created here?
-
-#This information is pyvolution focused, not mesh generation focused.
-#This is an appropriate place for the info to be created. -DSG
-
-#FIXME: Another issue is that the tsh file stores consecutive
-#indices explicitly. This is really redundant.
-#Suggest looking at obj and our own sww format and also consider
-#using netCDF.
-
-# It is redundant.  It was the format originally decided on.
-# I'm happy for it to change. -DSG
 
 def pmesh_dict_to_tag_dict(mesh_dict):

inundation/pyvolution/test_pmesh2domain.py

@@ -3 +3 @@
 
 import unittest
-#from math import sqrt
-#from Numeric import array, allclose, minimum, maximum, Float
-
-#from pytools.stats import mean
-#from domain import *
-#from boundary import *
-#from python_versions import compute_gradient
-#from config import epsilon
-#from shallow_water import attributes_2_field_values
-
 
 from Numeric import allclose, array
@@ -25 +15 @@
 from coordinate_transforms.geo_reference import Geo_reference
 
+#This is making pyvolution dependent on pmesh.
+# not good.  this should be in a seperate package.(directory)
+from pmesh.mesh import importMeshFromFile
 
 class Test_pmesh2domain(unittest.TestCase):
@@ -123 +116 @@
          self.failUnless(domain.geo_reference.xllcorner  == 140.0,
                           "bad geo_referece")
+    #************
+    
+    def test_pmesh2Domain_instance(self):
+         import os
+         import tempfile
+
+         fileName = tempfile.mktemp(".tsh")
+         file = open(fileName,"w")
+         file.write("4 3 # <vertex #> <x> <y> [attributes]\n \
+0 0.0 0.0 0.0 0.0 0.01 \n \
+1 1.0 0.0 10.0 10.0 0.02  \n \
+2 0.0 1.0 0.0 10.0 0.03  \n \
+3 0.5 0.25 8.0 12.0 0.04  \n \
+# Vert att title  \n \
+elevation  \n \
+stage  \n \
+friction  \n \
+2 # <triangle #> [<vertex #>] [<neigbouring triangle #>]  \n\
+0 0 3 2 -1  -1  1 dsg\n\
+1 0 1 3 -1  0 -1   ole nielsen\n\
+4 # <segment #> <vertex #>  <vertex #> [boundary tag] \n\
+0 1 0 2 \n\
+1 0 2 3 \n\
+2 2 3 \n\
+3 3 1 1 \n\
+3 0 # <x> <y> [attributes] ...Mesh Vertices... \n \
+0 216.0 -86.0 \n \
+1 160.0 -167.0 \n \
+2 114.0 -91.0 \n \
+3 # <vertex #>  <vertex #> [boundary tag] ...Mesh Segments... \n \
+0 0 1 0 \n \
+1 1 2 0 \n \
+2 2 0 0 \n \
+0 # <x> <y> ...Mesh Holes... \n \
+0 # <x> <y> <attribute>...Mesh Regions... \n \
+0 # <x> <y> <attribute>...Mesh Regions, area... \n\
+#Geo reference \n \
+56 \n \
+140 \n \
+120 \n")
+         file.close()
+
+         mesh_instance = importMeshFromFile(fileName)
+        
+         tags = {}
+         b1 =  Dirichlet_boundary(conserved_quantities = array([0.0]))
+         b2 =  Dirichlet_boundary(conserved_quantities = array([1.0]))
+         b3 =  Dirichlet_boundary(conserved_quantities = array([2.0]))
+         tags["1"] = b1
+         tags["2"] = b2
+         tags["3"] = b3
+
+         domain = pmesh_instance_to_domain_instance(mesh_instance, Domain)
+         os.remove(fileName)
+         #print "domain.tagged_elements", domain.tagged_elements
+         ## check the quantities
+         #print domain.quantities['elevation'].vertex_values
+         answer = [[0., 8., 0.],
+                   [0., 10., 8.]]
+         assert allclose(domain.quantities['elevation'].vertex_values,
+                        answer)
+
+         #print domain.quantities['stage'].vertex_values
+         answer = [[0., 12., 10.],
+                   [0., 10., 12.]]
+         assert allclose(domain.quantities['stage'].vertex_values,
+                        answer)
+
+         #print domain.quantities['friction'].vertex_values
+         answer = [[0.01, 0.04, 0.03],
+                   [0.01, 0.02, 0.04]]
+         assert allclose(domain.quantities['friction'].vertex_values,
+                        answer)
+
+         #print domain.quantities['friction'].vertex_values
+         assert allclose(domain.tagged_elements['dsg'][0],0)
+         assert allclose(domain.tagged_elements['ole nielsen'][0],1)
+
+         self.failUnless( domain.boundary[(1, 0)]  == '1',
+                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
+         self.failUnless( domain.boundary[(1, 2)]  == '2',
+                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
+         self.failUnless( domain.boundary[(0, 1)]  == '3',
+                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
+         self.failUnless( domain.boundary[(0, 0)]  == 'exterior',
+                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
+         #print "domain.boundary",domain.boundary
+         self.failUnless( len(domain.boundary)  == 4,
+                          "test_pmesh2Domain Too many boundaries")
+         #FIXME change to use get_xllcorner
+         #print "d.geo_reference.xllcorner",domain.geo_reference.xllcorner 
+         self.failUnless(domain.geo_reference.xllcorner  == 140.0,
+                          "bad geo_referece")
          
+    #***********
     def old_test_tags_to_boundaries (self):
          meshDict = {}