- Timestamp:
- Dec 12, 2007, 3:53:24 PM (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_work/production/shark_bay_2007/run_shark_bay_embayment_study.py
r4884 r4885 52 52 53 53 54 def run_model(**kwargs):55 54 56 alpha = kwargs['alpha']57 friction = kwargs['friction']58 time_thinning = kwargs['time_thinning']59 scenario_name = kwargs['aa_scenario_name']60 55 61 62 63 64 56 #----------------------------------------------------------------- 57 # Copy scripts to time stamped output directory and capture screen 58 # output to file 59 #----------------------------------------------------------------- 65 60 66 #copy_code_files(project.output_run_time_dir,67 #__file__,68 #dirname(project.__file__)+sep+project.__name__+'.py' )61 copy_code_files(project.output_run_time_dir, 62 __file__, 63 dirname(project.__file__)+sep+project.__name__+'.py' ) 69 64 70 65 71 #start_screen_catcher(project.output_run_time_dir)66 start_screen_catcher(project.output_run_time_dir) 72 67 73 68 74 75 76 77 78 79 69 #-------------------------------------------------------------------------- 70 # Create the triangular mesh based on overall clipping polygon with a 71 # tagged 72 # boundary and interior regions defined in project.py along with 73 # resolutions (maximal area of per triangle) for each polygon 74 #-------------------------------------------------------------------------- 80 75 81 82 83 84 85 86 87 88 89 76 print 'Create mesh from regions' 77 create_mesh_from_regions(project.small_bounding_polygon, 78 boundary_tags=project.boundary_tags_small, 79 maximum_triangle_area=project.res_bounding_polygon, 80 interior_regions=project.interior_regions, 81 filename=project.mesh_name+'.msh', 82 fail_if_polygons_outside=False, 83 use_cache=False, 84 verbose=True) 90 85 91 92 93 94 95 96 97 98 86 #------------------------------------------------------------------------- 87 # Setup computational domain 88 #------------------------------------------------------------------------- 89 print 'Setup computational domain' 90 domain = Domain(project.mesh_name+'.msh', 91 use_cache=False, verbose=True) 92 93 print domain.statistics() 99 94 100 95 101 102 103 96 #------------------------------------------------------------------------- 97 # Setup initial conditions 98 #------------------------------------------------------------------------- 104 99 105 106 107 domain.set_quantity('friction',friction)108 109 110 111 112 100 print 'Set initial conditions' 101 domain.set_quantity('stage', project.tide) 102 domain.set_quantity('friction', project.friction) 103 domain.set_quantity('elevation', 104 filename = project.combined_dir_name + '.pts', 105 use_cache = True, 106 verbose = True, 107 alpha = project.alpha) 113 108 114 109 115 116 117 118 domain.set_name(scenario_name)119 120 121 122 123 110 #------------------------------------------------------ 111 # Set domain parameters 112 #------------------------------------------------------ 113 domain.set_name(project.scenario_name) 114 domain.set_datadir(project.output_run_time_dir) 115 domain.set_default_order(2) # Apply second order scheme 116 domain.set_minimum_storable_height(0.01) # Don't store anything < 1cm 117 domain.set_store_vertices_uniquely(False) 118 domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum']) 124 119 125 126 120 domain.beta_h = 0 121 domain.tight_slope_limiters = 1 127 122 128 123 129 130 131 132 133 124 #------------------------------------------------------------------------- 125 # Setup boundary conditions 126 #------------------------------------------------------------------------- 127 print 'Available boundary tags', domain.get_boundary_tags() 128 print 'domain id', id(domain) 134 129 135 136 137 138 139 140 141 142 143 130 def waveform(t): 131 delay = 900 132 A = project.amplitude 133 T = project.period 134 if t < delay: 135 # Let things settle 136 return project.tide 137 else: 138 return project.tide + A*sin(2*pi*(t-delay)/T) 144 139 145 146 147 140 Bf = Transmissive_Momentum_Set_Stage_boundary(domain, waveform) 141 Br = Reflective_boundary(domain) 142 Bd = Dirichlet_boundary([project.tide,0,0]) 148 143 149 144 150 151 152 145 print 'Set boundary' 146 domain.set_boundary({'tide': Bd, 147 'ocean': Bf}) 153 148 154 149 155 156 157 158 159 160 161 162 domain.write_boundary_statistics(tags ='ocean')150 #------------------------------ 151 # Evolve system through time 152 #------------------------------ 153 t0 = time.time() 154 for t in domain.evolve(yieldstep = 5, 155 finaltime = 10000): 156 domain.write_time() 157 domain.write_boundary_statistics(tags='ocean') 163 158 164 159 165 166 167 168 169 160 #------------------------------ 161 # Post processing 162 #------------------------------ 163 x, y = domain.get_maximum_inundation_location() 164 q = domain.get_maximum_inundation_elevation() 170 165 171 166 print 'Maximum runup observed at (%.2f, %.2f) with elevation %.2f' %(x,y,q) 172 167 173 168 print 'That took %.2f seconds' %(time.time()-t0) 174 169 175 170 176 #-------------------------------------------------------------177 if __name__ == "__main__":178 171 179 run_model(file_name=project.home+'detail.csv',180 aa_scenario_name=project.scenario_name,181 ab_time=project.time,182 res_factor= project.res_factor,183 tide=project.tide,184 user=project.user,185 alpha = project.alpha,186 friction=project.friction,187 time_thinning = project.time_thinning,188 dir_comment=project.dir_comment)189 190
Note: See TracChangeset
for help on using the changeset viewer.