"""Script for running a simulation of UQ's wave flume. Ole Nielsen and Duncan Gray, GA - 2006 Modified by Matt Barnes 2012 """ #---------------------------------------------------------------------------- # Import necessary modules #---------------------------------------------------------------------------- # Standard modules import time import sys from shutil import copy from os import path import anuga from anuga.abstract_2d_finite_volumes.region import Set_region from anuga.fit_interpolate.interpolate import interpolate_sww2csv, \ file_function from anuga.utilities.file_utils import copy_code_files from anuga.file.csv_file import load_csv_as_dict from numerical_tools import err from anuga.file_conversion.file_conversion import timefile2netcdf # Application specific imports import create_mesh import project # Definition of file names and polygons def main(): D = 0.2 #initial depth frictions = [0.009] for friction in frictions: scenario(D, friction) def scenario(D, friction): import project #------------------------------------------------------------------------- # Setup archiving of simulations #------------------------------------------------------------------------- id = 'D_'+str(D)+'f_'+str(friction) copy (project.codedirname, project.outputtimedir + 'project.py') run_name = 'run_dam.py' run_name_out = 'run_dam'+id+'.py' copy (project.codedirname, project.outputtimedir + 'project.py') copy (project.codedir + 'run_UQwave.py', project.outputtimedir + 'run_UQwave.py') copy (project.codedir + 'create_mesh.py', project.outputtimedir + 'create_mesh.py') print'output dir', project.outputtimedir #FIXME this isn't working #normal screen output is stored in screen_output_name = project.outputtimedir + "screen_output.txt" screen_error_name = project.outputtimedir + "screen_error.txt" #------------------------------------------------------------------------- # Create the triangular mesh #------------------------------------------------------------------------- flume_length = 7.24 #flume_length = 3.24 create_mesh.generate(project.mesh_filename, flume_length, is_course=True) # this creates the mesh #is_course=False) # this creates the mesh head,tail = path.split(project.mesh_filename) copy (project.mesh_filename, project.outputtimedir + tail ) #------------------------------------------------------------------------- # Setup computational domain #------------------------------------------------------------------------- domain = anuga.Domain(project.mesh_filename, use_cache = False, verbose = True) print 'Number of triangles = ', len(domain) print 'The extent is ', domain.get_extent() print domain.statistics() domain.set_name(project.basename + id) domain.set_datadir(project.outputtimedir) domain.set_quantities_to_be_stored({'stage':2, 'xmomentum': 2, 'ymomentum': 2, 'elevation': 1}) domain.set_minimum_storable_height(0.003) domain.set_store_vertices_uniquely(True) # for writting to sww #------------------------------------------------------------------------- # Setup initial conditions #------------------------------------------------------------------------- slope = 0.084 #beach slope def elevation_tilt(x, y): return x*slope domain.set_quantity('stage', D) #water level domain.set_quantity('friction', friction) #bed friction domain.set_quantity('elevation', 0) domain.set_region('beach', 'elevation', elevation_tilt) print 'Available boundary tags', domain.get_boundary_tags() Br = anuga.Reflective_boundary(domain) Bd = anuga.Dirichlet_boundary([0,0,0]) # to drain the water out. #Bore file input bore_file = 'bore' timefile2netcdf(bore_file + '.txt', quantity_names = ['stage','xmomentum', 'ymomentum'], time_as_seconds=True) F = anuga.File_boundary(bore_file + '.tms', domain) domain.set_boundary( {'wall': Br, 'back':F} ) #------------------------------------------------------------------------- # Evolve system through time #------------------------------------------------------------------------- t0 = time.time() for t in domain.evolve(yieldstep = 0.02, finaltime = 10): #enter timestep and final time domain.write_time() print 'That took %.2f seconds' %(time.time()-t0) print 'finished' points = [] for i in range(41): points.append([float(i)/20 + -2.05, 0.425]) for j in range(120): points.append([float(j)/20, 0.425]) #------------------------------------------------------------------------- # Calculate gauge info #------------------------------------------------------------------------- if True: interpolate_sww2csv(project.outputtimedir + project.basename \ + id+".sww", points, project.depth_filename + id + '.csv', project.velocity_x_filename + id + '.csv', project.velocity_y_filename + id + '.csv') #------------------------------------------------------------- if __name__ == "__main__": main()