Changeset 5575 for anuga_work/production/busselton/run_busselton.py

Timestamp:

Jul 25, 2008, 4:17:26 PM (15 years ago)

Author:

kristy

Message:

Addition of Polyline information with sts files

File:

• anuga_work/production/busselton/run_busselton.py (modified) (5 diffs)

anuga_work/production/busselton/run_busselton.py

@@ -39 +39 @@
 from anuga.damage_modelling.inundation_damage import add_depth_and_momentum2csv, inundation_damage
 from anuga.fit_interpolate.benchmark_least_squares import mem_usage 
+from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
 
 # Application specific imports
@@ -71 +72 @@
     print "Processor Name:",get_processor_name()
 
-    # filenames
-#    meshes_dir_name = project.meshes_dir_name+'.msh'
-
-    # creates copy of code in output dir
-    print 'min triangles', project.trigs_min,
-    print 'Note: This is generally about 20% less than the final amount'
-
+    #-----------------------------------------------------------------------
+    # Domain definitions
+    #-----------------------------------------------------------------------
+
+    # Read in boundary from ordered sts file
+    urs_bounding_polygon=create_sts_boundary(project.scenario_name)
+
+    # Reading the landward defined points, this incorporates the original clipping
+    # polygon minus the 100m contour
+    landward_bounding_polygon = read_polygon(project.polygons_dir+'landward_bounding_polygon.txt')
+
+    # Combine sts polyline with landward points
+    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
+    
+    # counting segments
+    N = len(urs_bounding_polygon)-1
+    boundary_tags={'back': [N+1,N+2,N+3,N+4, N+5], 'side': [N,N+6],'ocean': range(N)}
+
+    
     #--------------------------------------------------------------------------
     # Create the triangular mesh based on overall clipping polygon with a
@@ -106 +119 @@
     print 'Setup computational domain'
 
-    #domain = cache(Domain, (meshes_dir_name), {'use_cache':True, 'verbose':True}, verbose=True)
-    #above don't work
     domain = Domain(project.meshes_dir_name+'.msh', use_cache=False, verbose=True)
     print 'memory usage before del domain',mem_usage()
@@ -168 +179 @@
     print 'Available boundary tags', domain.get_boundary_tags()
     print 'domain id', id(domain)
-    #print 'Reading Boundary file',project.boundaries_dir_namea + '.sww'
-
+   
     print'set_boundary'
 
+    boundary_urs_out=project.boundaries_dir_name
+    
+    print 'Available boundary tags', domain.get_boundary_tags()
+    Bf = File_boundary(boundary_urs_out+'.sts', 
+                   domain, time_thinning=1,
+                   use_cache=True,
+                   verbose = True,
+                   boundary_polygon=bounding_polygon)
+    
     Br = Reflective_boundary(domain)
     Bd = Dirichlet_boundary([kwargs['tide'],0,0])
-    Bw = Dirichlet_boundary([kwargs['tide']+10.0,0,0]) # To be deleted for FESA runs
-
-    Bf = Field_boundary(kwargs['boundary_file'],
-                    domain, time_thinning=kwargs['time_thinning'], mean_stage=kwargs['tide'], 
-                    use_cache=True, verbose=True)
+    #Bw = Dirichlet_boundary([kwargs['tide']+10.0,0,0]) # To be deleted for FESA runs
+
     print 'finished reading boundary file'
+
     domain.set_boundary({'back': Bd,
-                             'side': Bd,
-                             'ocean': Bd}) #changed from Bf to Bd for large wave
+                         'side': Bd,
+                         'ocean': Bf}) #changed from Bf to Bd for large wave
 
     kwargs['input_start_time']=domain.starttime
@@ -198 +215 @@
         domain.write_boundary_statistics(tags = 'ocean')
 
-        # To be deleted for FESA runs
-        if allclose(t, 240):
-            domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bw})
-
-        if allclose(t, 1440):
-            domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bd})
+##        # To be deleted for FESA runs
+##        if allclose(t, 240):
+##            domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bw})
+##
+##        if allclose(t, 1440):
+##            domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bd})