# Changeset 5370

Ignore:
Timestamp:
May 27, 2008, 5:39:41 PM (15 years ago)
Message:

Hinwood scenario

Location:
anuga_work/development
Files:
5 edited

Unmodified
Removed
• ## anuga_work/development/Hinwood_2008/create_mesh.py

 r5092 """Create mesh for University of Queensland dam break flume """Create mesh for Jon Hinwoods wave flume """ from anuga.coordinate_transforms.geo_reference import Geo_reference xslope = 5.0 # Distance between the boundary ande the start of the slope xslope = 5.0 # Distance between the boundary and the start of the slope def generate(mesh_filename, maximum_triangle_area=0.01): def generate(mesh_filename, slope, maximum_triangle_area=0.01): """ Generate mesh for University of Aberbeen dam break flume. The gate position is the distance from the right wall of the wave tank to the gate. Generate mesh for Monash University wave generation flume. Origin of x coordinate is the toe of the beach, x measured positive shorewards. """ global xslope xright  = 19.0 xleft = slope['xleft'] xtoe = slope['xtoe'] xbeach = slope['xbeach'] xright = slope['xright'] ybottom = 0 ytop    = 1.00 xleft = 0.0 ###xslope = slope # slope seperation (middle) point_slope_top = [xslope, ytop] point_slope_bottom = [xslope, ybottom] # 1st slope seperation (middle) point_toe_top = [xtoe, ytop] point_toe_bottom = [xtoe, ybottom] # 2nd slope seperation (middle) point_beach_top = [xbeach, ytop] point_beach_bottom = [xbeach, ybottom] m = Mesh() points = [point_sw,   #se point_nw, point_slope_top, point_toe_top, # 2 point_beach_top, # 3 point_ne, point_se, point_slope_bottom point_beach_bottom, # 6 point_toe_bottom # 7 ] [3,4], [4,5], [5,0],  #The outer border [2,5]       # slope separator [5,6], [6,7], [7,0],  #The outer border [2,7],       # toe separator [3,6]       # _beach separator ] segment_tags = {'wall':[1,2,3,4,5], segment_tags = {'wall':[1,2,3,4,5,6,7], 'wave':} # '': m.add_points_and_segments(points, segments, segment_tags) dam = m.add_region(xslope - 0.0000001,(ytop - ybottom)/2) # this is the location of the reservoir region. dam.setTag("flat") #     dam = m.add_region(xslope - 0.0000001,(ytop - ybottom)/2) #     # this is the location of the reservoir region. #     dam.setTag("flat") slope = m.add_region(xslope + 0.0000001,(ytop - ybottom)/2) # this is the location of the slope region. slope.setTag("slope") #     slope = m.add_region(xslope + 0.0000001,(ytop - ybottom)/2) #     # this is the location of the slope region. #     slope.setTag("slope") m.generate_mesh(maximum_triangle_area=maximum_triangle_area)
• ## anuga_work/development/Hinwood_2008/prepare_time_boundary.py

 r5350 def prepare_time_boundary(filename): """Convert benchmark 2 time series to NetCDF tms file. This is a 'throw-away' code taylor made for files like 'Benchmark_2_input.txt' from the LWRU2 benchmark """ # Read the ascii (.csv) version of this file fid = open(filename[:-4] + '.csv') # Skip first line line = fid.readline() # Read remaining lines # Create tms NetCDF file fid = NetCDFFile(filename, 'w') fid.institution = 'Geoscience Australia' """ Convert the rawish velocity and depth values, which have values at Convert the raw velocity and depth values, which have values at different times to a csv file, with values at the same time, with SI units. Set the depth values to be at the same times as the velocity values. The format for the velocity file is; [time, sec], [mm above SWL for sensor A], many other sensors... """ missing = 1e+20 # Read velocity file vfid = open(velocity_file) n_velocity = len(lines) vtimes = zeros(n_velocity, Float)  #Time velocities = zeros(n_velocity, Float)  # x_velocities = zeros(n_velocity, Float)  # y_velocities = zeros(n_velocity, Float)  # for i, line in enumerate(lines): fields = line.split(',') fields = line.split() #(',') vtimes[i] = float(fields) velocities[i] = float(fields) x_velocities[i] = float(fields) y_velocities[i] = float(fields) # Read the depth file dfid.close() n_depth = len(lines) n_sensors = len(lines.split()) dtimes = zeros(n_depth, Float)  #Time depths = zeros(n_depth, Float)  # sensors = zeros((n_depth,n_sensors), Float) for i, line in enumerate(lines): fields = line.split(',') dtimes[i] = float(fields) depths[i] = float(fields) depths_at_vtimes = interp(dtimes, depths, vtimes, missing=1e+20) fields = line.split() #(',') fields = [float(j) for j in fields] dtimes[i] = fields depths[i] = fields sensors[i] = fields #print "dtimes", dtimes #print "depths", depths #print "vtimes", vtimes depths_at_vtimes = interp( depths, dtimes, vtimes, missing=missing) depths_at_vtimes = ensure_numeric(depths_at_vtimes) #print "len(dtimes)", len(vtimes) #print "len(depths_at_vtimes)", len(depths_at_vtimes) #     for i in range(len(depths_at_vtimes)): #         print "i", i #         print "vtimes[i]", vtimes[i] #         print "depths_at_vtimes[i]", depths_at_vtimes[i] #     print "depths_at_vtimes",  depths_at_vtimes depths_at_vtimes = depths_at_vtimes/1000.00 # convert from mm to m velocities = ensure_numeric(velocities) velocities = velocities * -1.0 # Swap axis around missing=missing/1000.00 # Do to missing what is done to depths_at_vtimes x_velocities = ensure_numeric(x_velocities) x_velocities = x_velocities * -1.0 # Swap axis around y_velocities = ensure_numeric(y_velocities) y_velocities = y_velocities * -1.0 # Swap axis around fid = open(out_file,'w') assert len(depths_at_vtimes) == len(vtimes) start_time = None #start_time = 0.0 #for vtime, depth_at_vtime, velocity in map(vtimes, depths_at_vtimes, #                                           velocities): for i in xrange(len(vtimes)): fid.write(str(vtimes[i]) + ',' + str(depths_at_vtimes[i]) \ + ',' + str(velocities[i])+'\n') fid.close() if not depths_at_vtimes[i] == missing: # Make the times start at zero. if start_time is None: start_time = vtimes[i] fid.write(str(vtimes[i]-start_time) \ + ',' + str(depths_at_vtimes[i]) \ + ',' + str(x_velocities[i]) \ + ',' + str(y_velocities[i])+'\n') fid.close() # Since there is a new time reference save the depth info using this # new reference. fid = open(depth_file[:-4] + '_new_time'+depth_file[-4:],'w') sensors[:,0] -= start_time #print "depth_file", depth_file #print "start_time", start_time for i in xrange(len(dtimes)): fid.write(str(sensors[i,0])) for j in xrange(1,len(sensors)): fid.write(' ' + str(sensors[i,j])) fid.write('\n') fid.close() #------------------------------------------------------------------- if __name__ == "__main__": combine_velocity_depth('T2R7velfilt.csv','T2R7pressfilt.csv', 'cyeah.csv') from os import getenv from os.path import join home = getenv('INUNDATIONHOME') Hinwood_dir = join(home,'data','flumes','Hinwood_2008') raw_data_dir = join(Hinwood_dir, 'raw_data') # Test 1 Run 5 combine_velocity_depth(join(raw_data_dir,'T1R5velfilt.txt'), join(raw_data_dir,'T1R5pressfilt.txt'), join(Hinwood_dir, 'T1R5_boundary.csv')) # Create the tsm file prepare_time_boundary(join(Hinwood_dir, 'T1R5_boundary.tsm')) # Test 2 Run 7 combine_velocity_depth(join(raw_data_dir,'T2R7velfilt.txt'), join(raw_data_dir,'T2R7pressfilt.txt'), join(Hinwood_dir, 'T2R7_boundary.csv')) # Create the tsm file prepare_time_boundary(join(Hinwood_dir, 'T2R7_boundary.tsm'))
• ## anuga_work/development/Hinwood_2008/project.py

 r5350 self.meshdir = scenariodir+sep+'meshes'+sep self.scenariodir = scenariodir+sep self.boundary_file = self.scenariodir + 'boundary.tsm' self.boundarydir = self.scenariodir # creates copy of output dir structure, if it doesn't exist
• ## anuga_work/development/Hinwood_2008/run_dam.py

 r5350 from shutil import copy from os import path, sep from os.path import dirname  #, basename from os.path import dirname, join  #, basename from Numeric import zeros, size, Float # Related major packages import create_mesh from prepare_time_boundary import prepare_time_boundary def elevation_function(x,y): from Numeric import zeros, size, Float xslope = create_mesh.xslope  #4 ## Bit of a magic Number print "xslope",xslope z = zeros(size(x), Float) for i in range(len(x)): if x[i] < xslope: z[i] = 0.0  #WARNING: the code in prepare_time_boundary # that calc's momentum assumes this is 0.0 else: z[i] = (x[i]-xslope)*(1./16.) return z def main(friction=0.01, outputdir_name=None, is_trial_run=False): from interp import interp class Elevation_function: def __init__(self, slope): self.xslope_position = [slope['xleft'],slope['xtoe'], slope['xbeach'],slope['xright']] self.yslope_height = [slope['xleft'],slope['xtoe'], slope['xbeach'],slope['xright']] def __call__(self, x,y): z = interp(self.yslope_height, self.xslope_position, x) return z def main(boundary_file, boundary_location, slope, boundary_path=None, friction=0.01, outputdir_name=None, is_trial_run=False): pro_instance.outputdir + 'create_mesh.py') # Boundary file manipulation if boundary_path is None: boundary_path = pro_instance.boundarydir boundary_file_path = join(boundary_path, boundary_file) # Convert the boundary file, .csv to .tsm try: temp = open(pro_instance.boundary_file) temp = open(boundary_file_path) temp.close() except IOError: prepare_time_boundary(pro_instance.boundary_file) prepare_time_boundary(boundary_file_path) mesh_filename = pro_instance.meshdir + basename + '.msh' # this creates the mesh #gate_position = 12.0 create_mesh.generate(mesh_filename, create_mesh.generate(mesh_filename, slope, maximum_triangle_area=maximum_triangle_area) domain.set_quantity('stage', 0.4) domain.set_quantity('friction', friction) domain.set_quantity('friction', friction) elevation_function = Elevation_function(slope) domain.set_quantity('elevation', elevation_function) #Bts = Transmissive_Momentum_Set_Stage_boundary(domain, function) try: function = file_function(pro_instance.boundary_file, domain, function = file_function(boundary_file_path, domain, verbose=True) except IOError: Br = Reflective_boundary(domain) #Bd = Dirichlet_boundary([0.3,0,0]) Bd = Dirichlet_boundary([0.3,0,0]) Bts = Time_boundary(domain, function) domain.set_boundary( {'wall': Br, 'wave': Bts} ) domain.set_boundary( {'wall': Br, 'wave': Bd} ) #------------------------------------------------------------------------- t0 = time.time() # It seems that ANUGA can't handle a starttime that is >0. domain.starttime = 1.0 for t in domain.evolve(yieldstep, finaltime): domain.write_time() #------------------------------------------------------------- if __name__ == "__main__": main( is_trial_run = True, #slopes = [[-4.5,0.0],[0.0,0.0],[1.285,0.090],[16.1,.960]] slope = {'xleft':[-4.5,0.0], 'xtoe':[0.0,0.0], 'xbeach':[1.285,0.090], 'xright':[16.1,.960]} main( 'T1R5_boundary.tsm', 8, slope, is_trial_run = True, outputdir_name='Hinwood_low_stage_low_velocity_draft')
• ## anuga_work/development/friction_UA_flume_2006/parallel_run_dam.py

 r4443 e.g. miprun -c 4 python parallel_run_dam.py or e.g. miprun c4-8 python parallel_run_dam.py mpirun -c 4 python parallel_run_dam.py
Note: See TracChangeset for help on using the changeset viewer.