Changeset 6264

Timestamp:

Feb 3, 2009, 6:22:23 PM (16 years ago)

Author:

ole

Message:

Patong work, buildings, transmissive bdry

Location:

anuga_work/production/patong

Files:

• project.py (modified) (6 diffs)
• run_patong.py (modified) (5 diffs)

anuga_work/production/patong/project.py

@@ -48 +48 @@
 finaltime=15000         # final time for simulation 15000
 
-setup='basic'  # Final can be replaced with trial or basic.
+setup='final'  # Final can be replaced with trial or basic.
                # Either will result in a coarser mesh that will allow a
                # faster, but less accurate, simulation.
@@ -68 +68 @@
     yieldstep=10
 
+use_buildings = True
+
 #------------------------------------------------------------------------------
 # Output Filename
@@ -73 +75 @@
 # Important to distinguish each run - ensure str(user) is included!
 # Note, the user is free to include as many parameters as desired
-dir_comment='_'+setup+'_'+str(tide)+'_buildings_'+str(user)
+if use_buildings:
+    dir_comment='_'+setup+'_'+str(tide)+'_buildings_'+str(user)
+else:
+    dir_comment='_'+setup+'_'+str(tide)+'_nobuildings_'+str(user)    
 
 #------------------------------------------------------------------------------
@@ -170 +175 @@
 # Initial bounding polygon for data clipping 
 poly_all = read_polygon(polygons_dir+'poly_all.csv')
-res_poly_all = 100000*res_factor
+res_poly_all = 150000*res_factor
+
+
+# Inundation area
+poly_ia = read_polygon(polygons_dir+'inundation_area.csv')
+res_ia = 75*res_factor
+poly_saddle = read_polygon(polygons_dir+'saddle_10m.txt') 
 
 # Area of Interest 1 elevation from -10m to 20m
@@ -178 +189 @@
 # Area of Significance 1 elevation from -20m to a 200m buffer of the 20m contour
 poly_aos1 = read_polygon(polygons_dir+'aos.csv')
-res_aos1 = 700*res_factor
+res_aos1 = 900*res_factor
+
 
 # Area of Shallow water and coastal land that needs a finer res than 1000000
-# This is a rectangle.
 poly_sw = read_polygon(polygons_dir+'sw.csv')
 res_sw = 6000*res_factor
+
 
 
@@ -198 +210 @@
 
 # Combined all regions, must check that all are included!
-interior_regions = [[poly_aoi1,res_aoi1],
+interior_regions = [[poly_ia, res_ia],
+                    #[poly_saddle, res_ia],
+                    [poly_aoi1,res_aoi1],
                     [poly_aos1,res_aos1],
                     [poly_sw,res_sw],

anuga_work/production/patong/run_patong.py

@@ -28 +28 @@
 from anuga.interface import create_domain_from_regions
 from anuga.interface import Dirichlet_boundary
+from anuga.interface import Transmissive_stage_zero_momentum_boundary
 from anuga.interface import Reflective_boundary
 from anuga.interface import Field_boundary
@@ -122 +123 @@
                     verbose=True)
 
-
-# Add buildings from file
-print 'Reading building polygons'    
-building_polygons, building_heights = csv2building_polygons(project.building_polygon_file)
-                                                            #clipping_polygons=project.building_area_polygons)
-
-print 'Creating %d building polygons' % len(building_polygons)
-def create_polygon_function(building_polygons, geo_reference=None):
-    L = []
-    for i, key in enumerate(building_polygons):
-        if i%100==0: print i
-        poly = building_polygons[key]
-        elev = building_heights[key]
-        L.append((poly, elev))
-
-        buildings = Polygon_function(L, default=0.0,
-                                     geo_reference=geo_reference)
-    return buildings
-
-buildings = cache(create_polygon_function,
-                  building_polygons,
-                  {'geo_reference': domain.geo_reference},
-                  verbose=True)
-
-print 'Adding buildings'
-domain.add_quantity('elevation',
-                    buildings,
-                    use_cache=True,
-                    verbose=True)
+if project.use_buildings:
+    # Add buildings from file
+    print 'Reading building polygons'    
+    building_polygons, building_heights = csv2building_polygons(project.building_polygon_file)
+    #clipping_polygons=project.building_area_polygons)
+
+    print 'Creating %d building polygons' % len(building_polygons)
+    def create_polygon_function(building_polygons, geo_reference=None):
+        L = []
+        for i, key in enumerate(building_polygons):
+            if i%100==0: print i
+            poly = building_polygons[key]
+            elev = building_heights[key]
+            L.append((poly, elev))
+            
+            buildings = Polygon_function(L, default=0.0,
+                                         geo_reference=geo_reference)
+        return buildings
+
+    buildings = cache(create_polygon_function,
+                      building_polygons,
+                      {'geo_reference': domain.geo_reference},
+                      verbose=True)
+
+    print 'Adding buildings'
+    domain.add_quantity('elevation',
+                        buildings,
+                        use_cache=True,
+                        verbose=True)
 
 
@@ -175 +176 @@
 
 Br = Reflective_boundary(domain)
-Bd = Dirichlet_boundary([project.tide,0,0])
+#Bd = Dirichlet_boundary([project.tide,0,0])
+Bs = Transmissive_stage_zero_momentum_boundary(domain)
 Bf = Field_boundary(urs_boundary_name+'.sts', 
                     domain,
                     mean_stage= project.tide,
                     time_thinning=project.time_thinning,
-                    default_boundary=Bd,
+                    default_boundary=Bs,
                     use_cache=True,
                     verbose=True,
@@ -186 +188 @@
 
 domain.set_boundary({'back': Br,
-                     'side': Bd,
+                     'side': Bs,
                      'ocean': Bf}) 
 
@@ -195 +197 @@
 t0 = time.time()
 
+# Skip over the first 5000 seconds
+for t in domain.evolve(yieldstep=500,
+                       finaltime=5000):
+    print domain.timestepping_statistics()
+    print domain.boundary_statistics(tags='ocean')
+
+# Start detailed model
 for t in domain.evolve(yieldstep=project.yieldstep,
-                       finaltime=project.finaltime):
+                       finaltime=project.finaltime,
+                       skip_initial_step=True):
     print domain.timestepping_statistics()
     print domain.boundary_statistics(tags='ocean')