Changeset 4147 for anuga_work/production/dampier_2006/run_dampier.py

Timestamp:

Jan 8, 2007, 5:59:47 PM (18 years ago)

Author:

sexton

Message:

(1) updates to Dampier script based on Perth script (2) minor updates to Onslow report

File:

• anuga_work/production/dampier_2006/run_dampier.py (modified) (12 diffs)

anuga_work/production/dampier_2006/run_dampier.py

@@ -6 +6 @@
 
 The scenario is defined by a triangular mesh created from project.polygon,
-the elevation data and a simulated submarine landslide.
+the elevation data and a simulated tsunami generated with URS code.
 
 Ole Nielsen and Duncan Gray, GA - 2005 and Jane Sexton, Nick Bartzis, GA - 2006
@@ -23 +23 @@
 import sys
 
-
 # Related major packages
 from anuga.shallow_water import Domain
@@ -32 +31 @@
 
 from anuga.pmesh.mesh_interface import create_mesh_from_regions
-
-from anuga.geospatial_data.geospatial_data import *
 from anuga.abstract_2d_finite_volumes.util import start_screen_catcher, copy_code_files
 from anuga_parallel.parallel_api import distribute, numprocs, myid, barrier
@@ -45 +42 @@
 #------------------------------------------------------------------------------
 
+start_screen_catcher(project.output_run_time_dir, myid, numprocs)
 
 # filenames
-
-#build_time = '20061029_231935_build_tide_24'
-#build_time = '20061030_165746_build_tide_24'
-#build_time = '20061102_215532_build_plus_old_data'
-#build_time = '20061103_055258_build'
 build_time = '20061107_063840_build'
-#build_time = '20061025_153643_build_basic'
 bound_time = '20061102_221245_build'
 
 boundaries_name = project.boundaries_name
-#meshes_time_dir_name = project.meshes_time_dir_name+'.msh'
 meshes_dir_name = project.meshes_dir_name+'.msh'
-#source_dir = project.boundarydir
-#boundaries_time_dir_name = project.boundaries_dir + build_time + sep + boundaries_name
-boundaries_time_dir_name = project.boundaries_time_dir_name
 boundaries_dir_name = project.boundaries_dir_name
+
 tide = project.tide
 
@@ -70 +59 @@
                  dirname(project.__file__)+sep+ project.__name__+'.py' )
 barrier()
-#start_screen_catcher(project.output_run_time_dir, myid, numprocs)
 
 print 'USER: ', project.user
-#sys.exit()
+print 'min triangles', project.trigs_min,
+print 'Note: This is generally about 20% less than the final amount'
+
 #--------------------------------------------------------------------------
 # Create the triangular mesh based on overall clipping polygon with a
@@ -84 +74 @@
     
     print 'start create mesh from regions'
-    interior_regions = [#[project.karratha_polygon, 25000],
-                    [project.poly_coast, 10000],
-                    [project.poly_pipeline, 2000],
-                    [project.poly_facility, 500]]    
-#                    [project.poly_interior, 1000]]       
-    from anuga.utilities.polygon import plot_polygons
-    figname = project.output_run_time_dir + 'poly_pic'
-    plot_polygons([project.poly_coast,project.poly_pipeline,project.poly_facility,
-               project.bounding_polygon],
-               figname,
-               verbose = True)
-#    if access(project.meshes_time_dir,F_OK) == 0:
-#        mkdir(project.meshes_time_dir)
-#    if access(project.meshes_dir,F_OK) == 0:
-#        mkdir(project.meshes_dir)
-    import sys; sys.exit()
-    print 'start create mesh from regions'
-    interior_regions = [#[project.karratha_polygon, 25000],
-#                    [project.cipma_polygon, 1000],
-#                    [project.poly_pipeline, 5000],
-#                    [project.poly_facility, 500]]    
-                    [project.poly_interior, 1000]]    
-#    meshes_dir_name = project.meshes_dir_name + '.msh'
-
-    create_mesh_from_regions(project.bounding_polygon,
-                         boundary_tags={'back': [7, 8], 'side': [0, 6],
-                                        'ocean': [1, 2, 3, 4, 5]},
-                         maximum_triangle_area=100000,
-                         interior_regions=interior_regions,
-#                         filename=meshes_time_dir_name,
-                         filename=meshes_dir_name,
-                         use_cache=True,
-                         verbose=True)
+    create_mesh_from_regions(project.poly_all,
+                             boundary_tags={'back': [7, 8], 'side': [0, 6],
+                                            'ocean': [1, 2, 3, 4, 5]},
+                             maximum_triangle_area=project.res_poly_all,
+                             interior_regions=project.interior_regions,
+                             filename=meshes_dir_name,
+                             use_cache=True,
+                             verbose=True)
 
 # to sync all processors are ready 
@@ -125 +90 @@
 #-------------------------------------------------------------------------
 print 'Setup computational domain'
-#domain = Domain(meshes_time_dir_name, use_cache=True, verbose=True)
 domain = Domain(meshes_dir_name, use_cache=True, verbose=True)
 print domain.statistics()
 
+"""
 print 'starting to create boundary conditions'
 boundaries_in_dir_name = project.boundaries_in_dir_name
 
-from anuga.shallow_water.data_manager import urs2sww, ferret2sww
-
-print 'maxlat=project.south_boundary, minlat=project.north_boundary', project.south_boundary,project.north_boundary
-print 'minlon= project.west_boundary, maxlon=project.east_boundary',project.west_boundary, project.east_boundary
-print ' maxlon=project.east',project.east
-
-print 'origin: domain.geo_reference.get_origin()',domain.geo_reference.get_origin()
-
-#import sys; sys.exit()
-
-#if access(project.boundaries_time_dir,F_OK) == 0:
-#    mkdir (project.boundaries_time_dir)
+from anuga.shallow_water.data_manager import urs2sw
+
 # put above distribute
 print 'boundary file is: ',boundaries_dir_name
 from caching import cache
 if myid == 0:
-    cache(ferret2sww,
+    cache(urs2sww,
           (boundaries_in_dir_name,
     #       boundaries_time_dir_name), 
@@ -169 +124 @@
            )
 barrier()
-
+"""
 
 #-------------------------------------------------------------------------
@@ -179 +134 @@
 
     domain.set_quantity('stage', tide)
-    domain.set_quantity('friction', 0.0) 
+    domain.set_quantity('friction', 0.01) 
     #combined_time_dir_name = project.topographies_dir+build_time+project.combined_name
     print 'Start Set quantity'
 
-
     domain.set_quantity('elevation', 
-                    filename = project.topographies_dir + build_time + sep + project.combined_final_name + '.pts',
+                    filename = project.combined_dir_name + '.txt',
                     use_cache = True,
                     verbose = True,
                     alpha = 0.1)
-    #domain.set_quantity('elevation', -50)
     print 'Finished Set quantity'
 barrier()
@@ -211 +164 @@
 domain.set_store_vertices_uniquely(False)
 domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
-domain.set_maximum_allowed_speed(0.0) # Allow a little runoff (0.1 is OK)
+domain.set_maximum_allowed_speed(0.1) # Allow a little runoff (0.1 is OK)
 
 #-------------------------------------------------------------------------
@@ -217 +170 @@
 #-------------------------------------------------------------------------
 print 'Available boundary tags', domain.get_boundary_tags()
-
+print 'domain id', id(domain)
 print 'Reading Boundary file'
-#boundariesname = project.boundaries_dir + '20061101_003322_build'+sep+boundaries_name
-#print'boundariesname',boundariesname
-#Bf = File_boundary(boundaries_time_dir_name + '.sww',
-
-#Bf = File_boundary(boundariesname + '.sww',
-
-print 'domain id', id(domain)
-Bf = File_boundary(boundaries_dir_name + '.sww',
-                  domain, time_thinning=5, use_cache=True, verbose=True)
+#Bf = File_boundary(boundaries_dir_name + '.sww',
+#                  domain, time_thinning=5, use_cache=True, verbose=True)
 
 print 'finished reading boundary file'
@@ -246 +192 @@
 t0 = time.time()
 
-#for t in domain.evolve(yieldstep = 60, finaltime = 34000):
+for t in domain.evolve(yieldstep = 60, finaltime = 34000):
+    domain.write_time()
+    domain.write_boundary_statistics(tags = 'ocean')
+
+#for t in domain.evolve(yieldstep = 120, finaltime = 9000):
 #    domain.write_time()
 #    domain.write_boundary_statistics(tags = 'ocean')
-
-for t in domain.evolve(yieldstep = 120, finaltime = 9000):
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'ocean')
-#    print 'time: ',time.time()
- 
-    if allclose(t, 6000):
-        domain.set_quantity('xmomentum', 0)
-        domain.set_quantity('ymomentum', 0)
-#import sys; sys.exit()
      
-for t in domain.evolve(yieldstep = 60, finaltime = 28800
-                       ,skip_initial_step = True):
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'ocean')    
-
-for t in domain.evolve(yieldstep = 120, finaltime = 34800
-                       ,skip_initial_step = True): 
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'ocean')    
+#for t in domain.evolve(yieldstep = 60, finaltime = 28800
+#                       ,skip_initial_step = True):
+#    domain.write_time()
+#    domain.write_boundary_statistics(tags = 'ocean')    
+
+#for t in domain.evolve(yieldstep = 120, finaltime = 34800
+#                       ,skip_initial_step = True): 
+#    domain.write_time()
+#    domain.write_boundary_statistics(tags = 'ocean')    
     
 print 'That took %.2f seconds' %(time.time()-t0)