Changeset 4427

Timestamp:

May 15, 2007, 4:12:27 PM (18 years ago)

Author:

nick

Message:

updated alpha stuff

File:

• anuga_work/development/alpha_validation/find_alpha.py (modified) (1 diff)

anuga_work/development/alpha_validation/find_alpha.py

@@ -1 +1 @@
 
-from abstract_2d_finite_volumes import domain
+from anuga.abstract_2d_finite_volumes import domain
+from anuga.shallow_water import Domain
+from anuga.geospatial_data.geospatial_data import Geospatial_data
+from anuga.coordinate_transforms import Geo_reference
+from anuga.pmesh.mesh_interface import create_mesh_from_regions
+from anuga.shallow_water import Reflective_boundary
 
 
 
-create_mesh_from_regions(project.bounding_polygon,
-                         boundary_tags={'back': [7, 8], 'side': [0, 6],
-                                        'ocean': [1, 2, 3, 4, 5]},
+poly_topo = [[664000,7752000],[664000,7754000],
+             [666000,7754000],[666000,7752000]]
+mesh_dir_name = 'temp.msh'
+                  
+create_mesh_from_regions(poly_topo,
+                         boundary_tags={'back': [2], 'side': [1,3],
+                                        'ocean': [0]},
                          maximum_triangle_area=100000,
-                         interior_regions=interior_regions,
-#                         filename=meshes_time_dir_name,
-                         filename=meshes_dir_name,
-                         use_cache=True,
+                         filename=mesh_dir_name,
+                         use_cache=False,
                          verbose=True)
 
 
+topo_dir_name = 'pt_hedland_export.txt'
+sample = 'sample.txt'
+remainder ='remainder.txt'
+#split topo data
+G = Geospatial_data(file_name = topo_dir_name)
+print 'start split'
+G_small, G_other = G.split(0.1,True)
 
-alphas = [0.001, 0.01, 0.1, 1.0, 10.0, 100.0]
+#print 'start export'
+#G_small.export_points_file(sample)
+#G_other.export_points_file(remainder)
+
+
+#alphas = [0.001, 0.01, 0.1, 1.0, 10.0, 100.0]
+alphas = [1]
 domains = {}
+
 print 'Setup computational domain'
 for alpha in alphas:
     
-    domain = Domain(meshes_dir_name, use_cache=True, verbose=True)
+    #think about making domain absolute when it is created
+    geo_ref=Geo_reference(zone = 0, xllcorner = 0.0, yllcorner = 0.0)
+    
+#    domain = Domain(mesh_dir_name,geo_reference=geo_ref, use_cache=False, verbose=True)
+    domain = Domain(mesh_dir_name, use_cache=False, verbose=True)
     print domain.statistics()
     domain.set_quantity('elevation', 
-                    filename = project.topographies_dir + build_time + sep + project.combined_final_name + '.pts',
+                    geospatial_data = G_other,
                     use_cache = True,
                     verbose = True,
                     alpha = alpha)
+    domains[alpha]=domain
+    domain.set_name('temp'+str(alpha))
+    domain.set_datadir('.')
+    print'set_boundary'
+#    points=[[665000,7753000],[665100,7753000]]
+    points=G_small.get_data_points()
+    print'geo',domain.geo_reference,points[0:10]
+    #gets relative point from sample
+    points=domain.geo_reference.change_points_geo_ref(points)
+    
+    point_elevation=G_small.get_attributes(attribute_name='elevation')
+
+    print'poin',points[0:10],point_elevation[0:10]
+    print'quantities',domain.get_quantity_names()
+    hmm=domain.quantities['elevation'].get_values(interpolation_points=points)
+    print'hmm',hmm[0:10],len(points),len(hmm)
+    
+#    for i in si
+    
+#    Br = Reflective_boundar[y(domain)
+    
+#    domain.set_boundary({'back': Br,
+#                     'side': Br,
+#                     'ocean': Br}) 
+
+#    for t in domain.evolve(yieldstep = 1, finaltime = 1):
+#        domain.write_time()
+        
+
+
+
+# from file load x,y,z data
+# use file xy data to find predicted z data from domain
+# determine error and repeat