Changeset 6264 for anuga_work/production/patong/run_patong.py

Timestamp:

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

Author:

ole

Message:

Patong work, buildings, transmissive bdry

File:

• anuga_work/production/patong/run_patong.py (modified) (5 diffs)

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')