Changeset 5686
- Timestamp:
- Aug 26, 2008, 3:23:31 PM (15 years ago)
- Location:
- anuga_validation/Hinwood_2008
- Files:
-
- 1 added
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_validation/Hinwood_2008/project.py
r5682 r5686 1 """Common filenames and locations for topographic data, meshes and outputs. 2 Also includes origin for slump scenario. 1 """Common filenames and locations. 3 2 """ 4 3 5 from os import sep, getenv, getcwd, mkdir, access, F_OK 6 import sys 7 from time import localtime, strftime 8 from os.path import join 9 from anuga.shallow_water.data_manager import copy_code_files 10 from anuga.abstract_2d_finite_volumes.util import add_directories 11 12 from anuga.utilities.system_tools import get_user_name 13 14 15 class Project: 16 def __init__(self, 17 trunk=['data','flumes','Hinwood_2008'], 18 outputdir_name = None, 19 home = None): 20 21 self.user = get_user_name() 22 if home is None: 23 try: 24 home = getenv('INUNDATIONHOME') #Sandpit's parent dir 25 except: 26 home = '.' 27 self.home = home 28 # Create the structure of where the output directories will go 29 scenariodir=add_directories(home, trunk) 30 31 #Derive subdirectories and filenames 32 if outputdir_name is None: 33 #gets time for dir 34 outputdir_name = strftime('%Y%m%d_%H%M%S',localtime()) 35 36 general_outputdir = scenariodir+sep+'output'+sep 37 self.outputdir = general_outputdir+outputdir_name+sep 38 self.meshdir = scenariodir+sep+'meshes'+sep 39 self.scenariodir = scenariodir+sep 40 self.boundarydir = self.outputdir 41 self.raw_data_dir = join(scenariodir, 'raw_data') 42 self.plots_dir = join(scenariodir, 'plots') 43 self.rmsd_dir = join(scenariodir, 'rmsd') 44 45 # creates copy of output dir structure, if it doesn't exist 46 if not access(self.meshdir,F_OK): 47 mkdir (self.meshdir) 48 if not access(general_outputdir,F_OK): 49 mkdir (general_outputdir) 50 if not access(self.outputdir,F_OK): 51 mkdir (self.outputdir) 52 if not access(self.plots_dir,F_OK): 53 mkdir (self.plots_dir) 54 if not access(self.rmsd_dir,F_OK): 55 mkdir (self.rmsd_dir) 56 57 self.codedir = getcwd()+sep 58 59 60 #------------------------------------------------------------- 61 if __name__ == "__main__": 62 p = Project(['eagle'], 'lego','.') 63 print p.outputdir 64 print p.meshdir 4 boundary_file_tag = '_boundary.tsm' -
anuga_validation/Hinwood_2008/run_dam.py
r5682 r5686 2 2 3 3 Script for running a breaking wave simulation of Jon Hinwoods wave tank. 4 Note: this is based on the frinction_ua_flume_2006 structure.5 6 4 7 5 Duncan Gray, GA - 2007 8 9 10 11 6 """ 12 7 … … 18 13 # Standard modules 19 14 import time 20 from time import localtime, strftime 21 import sys 22 from shutil import copy 23 from os import path, sep 24 from os.path import dirname, join #, basename 25 from Numeric import zeros, size, Float 15 26 16 27 17 # Related major packages 28 from anuga.shallow_water import Domain, Reflective_boundary, \ 29 Dirichlet_boundary, Time_boundary, \ 30 File_boundary, \ 31 Transmissive_Momentum_Set_Stage_boundary 18 from anuga.shallow_water import Domain, Reflective_boundary, Time_boundary 32 19 from anuga.fit_interpolate.interpolate import interpolate_sww2csv 33 from anuga.abstract_2d_finite_volumes.util import start_screen_catcher, \ 34 file_function 35 from anuga.shallow_water.data_manager import copy_code_files 36 from anuga.abstract_2d_finite_volumes.generic_boundary_conditions\ 37 import File_boundary_time 20 from anuga.abstract_2d_finite_volumes.util import file_function 21 from anuga.utilities.interp import interp 38 22 39 23 # Scenario specific imports 40 import project # Definition of file names and polygons41 24 import create_mesh 42 from prepare_time_boundary import prepare_time_boundary 43 from anuga.utilities.interp import interp 44 45 46 class Elevation_function: 47 def __init__(self, slope): 48 self.xslope_position = [slope['xleft'][0],slope['xtoe'][0], 49 slope['xbeach'][0],slope['xright'][0]] 50 self.yslope_height = [slope['xleft'][1],slope['xtoe'][1], 51 slope['xbeach'][1],slope['xright'][1]] 52 53 def __call__(self, x,y): 54 55 z = interp(self.yslope_height, self.xslope_position, x) 56 return z 57 58 def main(boundary_file, 59 metadata_dic, 25 from prepare_time_boundary import csv2tms 26 27 28 def main(metadata_dic, 60 29 maximum_triangle_area, 61 30 yieldstep, 62 31 boundary_path=None, 63 friction=0. 012, # planed wood.http://www.lmnoeng.com/manningn.htm32 friction=0., # planed wood 0.012 http://www.lmnoeng.com/manningn.htm 64 33 outputdir_name=None, 65 run_type=1,34 isTest=False, 66 35 width=1.0, 67 36 use_limits=True, 68 end_tag = '_limiterD'): 69 70 71 basename = 'zz_' + metadata_dic['scenario_id'] 72 73 if run_type == 0: 37 end_tag = ''): 38 39 id = metadata_dic['scenario_id'] 40 41 if isTest is True: 74 42 yieldstep = 1.0 75 43 finaltime = 15. 76 44 maximum_triangle_area=0.1 77 outputdir_name += '_test' 78 79 elif run_type == 1: 80 finaltime = None 81 82 elif run_type == 2: 83 yieldstep = 0.5 45 outputdir_name += '_test' 46 else: 84 47 finaltime = None 85 maximum_triangle_area=0.01 86 outputdir_name += '_test_long_time' 87 88 elif run_type == 3: 89 yieldstep = 0.1 90 finaltime = None 91 maximum_triangle_area=0.01 92 #outputdir_name += '_yieldstep_0.1' 93 94 elif run_type == 4: 95 # this is not a test 96 # Output will go to a file 97 # The sww file will be interpolated 98 yieldstep = 0.01 99 finaltime = None 100 maximum_triangle_area=0.01 101 102 elif run_type == 5: 103 # this is not a test 104 # Output will go to a file 105 # The sww file will be interpolated 106 yieldstep = 0.01 107 finaltime = None 108 maximum_triangle_area=0.001 109 #outputdir_name += '_good' 110 111 elif run_type == 6: 112 # this is not a test 113 # Output will go to a file 114 # The sww file will be interpolated 115 yieldstep = 0.01 116 finaltime = None 117 maximum_triangle_area=0.0001 118 #outputdir_name += '_good' 119 120 elif run_type == 7: 121 # this is not a test 122 # Output will go to a file 123 # The sww file will be interpolated 124 yieldstep = 0.01 125 finaltime = None 126 maximum_triangle_area=0.00001 127 #outputdir_name += '_good' 128 129 elif run_type == 8: 130 # this is not a test 131 # Output will go to a file 132 # The sww file will be interpolated 133 yieldstep = 0.05 134 finaltime = None 135 maximum_triangle_area=0.00001 136 #outputdir_name += '_good' 137 #outputdir_name += '_no_velocity' 138 48 139 49 outputdir_name = paras2outputdir_tag(mta=maximum_triangle_area, 140 50 yieldstep=yieldstep, … … 143 53 friction=friction, 144 54 end_tag=end_tag, 145 outputdir_name=outputdir_name 146 ) 147 55 outputdir_name=outputdir_name) 56 basename = outputdir_name 148 57 metadata_dic = set_z_origin_to_water_depth(metadata_dic) 149 150 pro_instance = project.Project(['data','flumes','Hinwood_2008'], 151 outputdir_name=outputdir_name) 152 print "The output dir is", pro_instance.outputdir 153 copy_code_files(pro_instance.outputdir,__file__, 154 dirname(project.__file__) \ 155 + sep + project.__name__+'.py') 156 copy (pro_instance.codedir + 'run_dam.py', 157 pro_instance.outputdir + 'run_dam.py') 158 copy (pro_instance.codedir + 'create_mesh.py', 159 pro_instance.outputdir + 'create_mesh.py') 160 161 boundary_final_time = prepare_time_boundary(metadata_dic, 162 pro_instance.raw_data_dir, 163 pro_instance.boundarydir) 164 #return pro_instance 58 165 59 if finaltime is None: 166 finaltime = boundary_final_time - 0.1 # Edge boundary problems 167 # Boundary file manipulation 168 if boundary_path is None: 169 boundary_path = pro_instance.boundarydir 170 boundary_file_path = join(boundary_path, boundary_file) 171 # # Convert the boundary file, .csv to .tsm 172 # try: 173 # temp = open(boundary_file_path) 174 # temp.close() 175 # except IOError: 176 # prepare_time_boundary(boundary_file_path) 177 178 mesh_filename = pro_instance.meshdir + basename + '.msh' 179 180 #-------------------------------------------------------------------------- 181 # Copy scripts to output directory and capture screen 182 # output to file 183 #-------------------------------------------------------------------------- 184 185 # creates copy of code in output dir 186 if run_type >= 1: 187 #start_screen_catcher(pro_instance.outputdir, rank, pypar.size()) 188 start_screen_catcher(pro_instance.outputdir) 189 190 print 'USER: ', pro_instance.user 60 finaltime = metadata_dic['wave_times'][1] + 0.1 61 62 boundary_file_path = id + '_boundary.tsm' 63 64 # Convert the boundary file, .csv to .tsm 65 csv2tms(boundary_file_path, metadata_dic['xleft'][1]) 66 67 mesh_filename = basename + '.msh' 68 191 69 #------------------------------------------------------------------------- 192 70 # Create the triangular mesh 193 71 #------------------------------------------------------------------------- 194 72 195 # this creates the mesh196 #gate_position = 12.0197 73 create_mesh.generate(mesh_filename, metadata_dic, width=width, 198 74 maximum_triangle_area=maximum_triangle_area) 199 200 head,tail = path.split(mesh_filename) 201 copy (mesh_filename, 202 pro_instance.outputdir + tail ) 75 203 76 #------------------------------------------------------------------------- 204 77 # Setup computational domain 205 78 #------------------------------------------------------------------------- 206 domain = Domain(mesh_filename, use_cache = False, verbose = True) 207 79 domain = Domain(mesh_filename, use_cache = False, verbose = True) 208 80 209 81 print 'Number of triangles = ', len(domain) 210 82 print 'The extent is ', domain.get_extent() 211 83 print domain.statistics() 212 213 84 214 85 domain.set_name(basename) 215 domain.set_datadir( pro_instance.outputdir)86 domain.set_datadir('.') 216 87 domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum']) 217 88 domain.set_minimum_storable_height(0.0001) … … 237 108 domain.set_quantity('elevation', elevation_function) 238 109 239 240 print 'Available boundary tags', domain.get_boundary_tags() 241 242 # Create boundary function from timeseries provided in file 243 #function = file_function(project.boundary_file, domain, verbose=True) 244 #Bts = Transmissive_Momentum_Set_Stage_boundary(domain, function) 245 try: 246 function = file_function(boundary_file_path, domain, 247 verbose=True) 248 print "boundary_file_path", boundary_file_path 249 except IOError: 250 msg = 'Run prepare_time_boundary.py. File "%s" could not be opened.'\ 251 %(pro_instance.boundary_file) 252 raise msg 110 function = file_function(boundary_file_path, domain, verbose=True) 253 111 254 112 Br = Reflective_boundary(domain) 255 Bd = Dirichlet_boundary([0.3,0,0])256 113 Bts = Time_boundary(domain, function) 257 #Bts = Transmissive_Momentum_Set_Stage_boundary(domain, function)258 114 domain.set_boundary( {'wall': Br, 'wave': Bts} ) 259 #domain.set_boundary( {'wall': Br, 'wave': Bd} )260 115 261 116 #------------------------------------------------------------------------- … … 264 119 t0 = time.time() 265 120 266 # It seems that ANUGA can't handle a starttime that is >0.267 #domain.starttime = 1.0 #!!! what was this doing?268 121 for t in domain.evolve(yieldstep, finaltime): 269 122 domain.write_time() … … 271 124 print 'finished' 272 125 273 flume_y_middle = 0.0274 points = []275 for gauge_x in metadata_dic['gauge_x']:276 points.append([gauge_x, flume_y_middle])277 278 126 279 127 #------------------------------------------------------------------------- … … 281 129 #------------------------------------------------------------------------- 282 130 283 if run_type >= 1: 131 if isTest is not True: 132 points = [] 133 for gauge_x in metadata_dic['gauge_x']: 134 points.append([gauge_x, 0.0]) 135 284 136 id = metadata_dic['scenario_id'] + ".csv" 285 interpolate_sww2csv( pro_instance.outputdir +basename +".sww",137 interpolate_sww2csv(basename +".sww", 286 138 points, 287 pro_instance.outputdir + "depth_" + id, 288 pro_instance.outputdir + "velocity_x_" + id, 289 pro_instance.outputdir + "velocity_y_" + id, 290 pro_instance.outputdir + "stage_" + id, 291 pro_instance.outputdir + "froude_" + id) 139 basename + "_depth_" + id, 140 basename + "_velocity_x_" + id, 141 basename + "_velocity_y_" + id, 142 basename + "_stage_" + id) 292 143 293 return pro_instance294 295 296 144 297 145 def paras2outputdir_tag(mta, … … 301 149 friction, 302 150 end_tag, 303 outputdir_name=None 304 ): 305 151 outputdir_name=None): 152 """ 153 154 """ 306 155 if outputdir_name is None: 307 156 outputdir_name = '' … … 319 168 return outputdir_name 320 169 321 322 def set_z_origin_to_water_depth(seabed_coords): 323 offset = seabed_coords['offshore_water_depth'] 170 171 class Elevation_function: 172 """ 173 """ 174 def __init__(self, slope): 175 self.xslope_position = [slope['xleft'][0],slope['xtoe'][0], 176 slope['xbeach'][0],slope['xright'][0]] 177 self.yslope_height = [slope['xleft'][1],slope['xtoe'][1], 178 slope['xbeach'][1],slope['xright'][1]] 179 180 def __call__(self, x,y): 181 z = interp(self.yslope_height, self.xslope_position, x) 182 return z 183 184 def set_z_origin_to_water_depth(metadata_dic): 185 """ 186 """ 187 offset = metadata_dic['offshore_water_depth'] 324 188 keys = ['xleft', 'xtoe', 'xbeach', 'xright'] 325 189 for x in keys: 326 seabed_coords[x][1] -= offset 327 return seabed_coords 190 metadata_dic[x][1] -= offset 191 return metadata_dic 192 328 193 #------------------------------------------------------------- 329 194 if __name__ == "__main__": 330 195 331 196 from scenarios import scenarios 332 from slope import gauges_for_slope 333 #from plot import plot 334 335 # 1 is fast and dirty 336 # 4 is 0.01 337 # 5 is 0.001 338 # 6 is 0.0001 339 # 7 is 0.00001 340 # 8 is 0.00001 yieldstep = 0.05 341 342 run_type = 1 343 run_type = 0 # for testing 344 345 #for run_data in [scenarios[5]]: 346 #scenarios = scenarios[3:] 197 347 198 scenarios = [scenarios[0]] 348 199 349 width = 1.0350 200 width = 0.1 351 #width = 0.01 352 353 maximum_triangle_area=0.0001 201 maximum_triangle_area=0.01 354 202 yieldstep = 0.01 203 yieldstep = 0.5 355 204 friction=0.0 205 isTest=True 206 #isTest=False 356 207 357 208 for run_data in scenarios: 358 pro_instance = main( run_data['scenario_id'] + '_boundary.tsm' , 359 run_data, 360 maximum_triangle_area=maximum_triangle_area, 361 yieldstep=yieldstep, 362 width=width, 363 run_type=run_type, 364 outputdir_name=run_data['scenario_id'], 365 use_limits=False, 366 friction=friction, 367 end_tag='_I') 368 #gauges_for_slope(pro_instance.outputdir,[run_data]) 209 main(run_data, 210 maximum_triangle_area=maximum_triangle_area, 211 yieldstep=yieldstep, 212 width=width, 213 isTest=isTest, 214 outputdir_name=run_data['scenario_id'], 215 use_limits=False, 216 friction=friction, 217 end_tag='') -
anuga_validation/Hinwood_2008/scenarios.py
r5683 r5686 19 19 'start_slope_x':4, 20 20 'finish_slope_x':6, 21 'break_times':[34.15, 38.35, 42.95, 47.85, 53.15, 58.15],22 21 'break_xs':[5.45, 5.45, 5.28, 5.28, 5.29, 5.29], 23 22 'break_type':['none', 'plunge', 'plunge', 'plunge', 'plunge', \ … … 45 44 'start_slope_x':4, 46 45 'finish_slope_x':6, 47 'break_times':[28.68, 33.18, 37.68, 42.18, 47.38, 52.18],48 46 'break_xs':[5.45, 5.58, 5.29, 5.245, 5.21, 5.21], 49 47 'break_type':['none', 'plunge', 'plunge', 'plunge', 'plunge', \ … … 72 70 'start_slope_x':5, 73 71 'finish_slope_x':6.5, 74 'break_times':[33.95, 41.25, 49.55, 56.75], # last 2 removed75 72 'break_xs':[5.93,5.93,5.907,5.78], # since depth finishes at 65 sec 76 73 'break_type':['none', 'front steepened', 'front steepened', 'spill'], … … 96 93 'start_slope_x':5, 97 94 'finish_slope_x':6.5, 98 'break_times':[47.00,55.00,62.00,70.00], # last 2 removed99 95 'break_xs':[6.32,6.32,6.033,5.935], # since depth finishes at 75 sec 100 96 'break_type':['none', 'front steepened', 'front steepened', … … 128 124 'start_slope_x':8., 129 125 'finish_slope_x':9.5, 130 'break_times':[61.21,68.51,76.11,84.11],131 126 'break_xs':[9.17,9.135,9.135,9.105], 132 127 'break_type':['collapse', 'collapse', 'collapse', … … 161 156 'start_slope_x':8., 162 157 'finish_slope_x':9.5, 163 'break_times':[61.61,66.61,74.31,81.71],164 158 'break_xs':[9.063,9.063,9.043,9.043], 165 159 'break_type':['none','collapse', 'collapse', 'collapse'], … … 193 187 'start_slope_x':6, 194 188 'finish_slope_x':9, 195 'break_times':[65.99,69.89,75.49,81.49,86.19],196 189 'break_xs':[8.53,7.46,7.492,7.492,7.444], 197 190 'break_type':['collapse', 'spill', 'spill', 'spill', 'spill'], … … 212 205 'gauge_names':['B','1','2','3','5','7', 213 206 '9','10','11','12'], 214 215 207 'gauge_x':[-0.325, 1.572, 2.571, 3.571, 5.57, 216 208 7.57, 9.569, 10.569, 11.569,
Note: See TracChangeset
for help on using the changeset viewer.