Changeset 7064 for anuga_core/documentation/user_manual/demos
- Timestamp:
- May 22, 2009, 4:40:11 PM (15 years ago)
- Location:
- anuga_core/documentation/user_manual/demos
- Files:
-
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/documentation/user_manual/demos/cairns/ExportResults.py
r5408 r7064 1 import project,os1 import os 2 2 import sys 3 import project 3 4 4 5 from anuga.shallow_water.data_manager import sww2dem … … 7 8 8 9 scenario = 'slide' 9 #scenario = 'fixed_wave'10 10 11 11 name = scenario + sep + scenario + 'source' … … 13 13 print 'output dir:', name 14 14 which_var = 4 15 if which_var == 0: # Stage 15 16 if which_var == 0: # Stage 16 17 outname = name + '_stage' 17 18 quantityname = 'stage' 18 19 19 if which_var == 1: # Absolute Momentum20 if which_var == 1: # Absolute Momentum 20 21 outname = name + '_momentum' 21 quantityname = '(xmomentum**2 + ymomentum**2)**0.5' #Absolute momentum22 quantityname = '(xmomentum**2 + ymomentum**2)**0.5' #Absolute momentum 22 23 23 if which_var == 2: # Depth24 if which_var == 2: # Depth 24 25 outname = name + '_depth' 25 26 quantityname = 'stage-elevation' #Depth 26 27 27 if which_var == 3: # Speed28 if which_var == 3: # Speed 28 29 outname = name + '_speed' 29 30 quantityname = '(xmomentum**2 + ymomentum**2)**0.5/(stage-elevation+1.e-30)' #Speed 30 31 31 if which_var == 4: # Elevation32 if which_var == 4: # Elevation 32 33 outname = name + '_elevation' 33 34 quantityname = 'elevation' #Elevation … … 35 36 print 'start sww2dem' 36 37 37 sww2dem(name, basename_out = outname,38 quantity = quantityname,39 cellsize = 100,40 easting_min = project.eastingmin,41 easting_max = project.eastingmax,42 northing_min = project.northingmin,43 northing_max = project.northingmax,44 reduction = max,45 verbose = True,46 format = 'ers')47 38 sww2dem(name, 39 basename_out=outname, 40 quantity=quantityname, 41 cellsize=100, 42 easting_min=project.eastingmin, 43 easting_max=project.eastingmax, 44 northing_min=project.northingmin, 45 northing_max=project.northingmax, 46 reduction=max, 47 verbose=True, 48 format='ers') -
anuga_core/documentation/user_manual/demos/cairns/GetTimeseries.py
r4966 r7064 7 7 8 8 Note, this script will only work if pylab is installed on the platform 9 10 9 """ 11 10 … … 16 15 sww2csv_gauges('slide'+sep+'slidesource.sww', 17 16 project.gauge_filename, 18 quantities =['stage','speed','depth','elevation'],17 quantities=['stage','speed','depth','elevation'], 19 18 verbose=True) 20 19 21 20 sww2csv_gauges('fixed_wave'+sep+'fixed_wavesource.sww', 22 21 project.gauge_filename, 23 quantities =['stage', 'speed','depth','elevation'],22 quantities=['stage', 'speed','depth','elevation'], 24 23 verbose=True) 25 24 26 25 try: 27 26 import pylab 28 csv2timeseries_graphs(directories_dic={'slide'+sep: ['Slide',0, 0],29 'fixed_wave'+sep: ['Fixed Wave',0,0]},30 31 32 33 34 35 36 37 27 csv2timeseries_graphs(directories_dic={'slide'+sep: ['Slide',0, 0], 28 'fixed_wave'+sep: ['Fixed Wave',0,0]}, 29 output_dir='fixed_wave'+sep, 30 base_name='gauge_', 31 plot_numbers='', 32 quantities=['stage','speed','depth'], 33 extra_plot_name='', 34 assess_all_csv_files=True, 35 create_latex=False, 36 verbose=True) 38 37 except ImportError: 39 38 #ANUGA does not rely on pylab to work 40 39 print 'must have pylab installed to generate plots' 41 42 -
anuga_core/documentation/user_manual/demos/cairns/project.py
r6889 r7064 1 """Common filenames and locations for topographic data, meshes and outputs. 2 """ 1 """Common filenames and locations for topographic data, meshes and outputs.""" 3 2 4 3 from anuga.utilities.polygon import read_polygon, plot_polygons, \ 5 4 polygon_area, is_inside_polygon 6 5 7 8 6 #------------------------------------------------------------------------------ 9 7 # Define scenario as either slide or fixed_wave. 10 8 #------------------------------------------------------------------------------ 11 #scenario = 'slide'12 9 scenario = 'fixed_wave' 13 14 10 #scenario = 'slide' 11 15 12 #------------------------------------------------------------------------------ 16 13 # Filenames 17 14 #------------------------------------------------------------------------------ 18 demname = 'cairns' 15 demname = 'cairns' 19 16 meshname = demname + '.msh' 20 17 21 18 # Filename for locations where timeseries are to be produced 22 19 gauge_filename = 'gauges.csv' 23 24 20 25 21 #------------------------------------------------------------------------------ 26 22 # Domain definitions 27 23 #------------------------------------------------------------------------------ 28 29 24 # bounding polygon for study area 30 25 bounding_polygon = read_polygon('extent.csv') 31 26 32 A = polygon_area(bounding_polygon) /1000000.027 A = polygon_area(bounding_polygon) / 1000000.0 33 28 print 'Area of bounding polygon = %.2f km^2' % A 34 35 29 36 30 #------------------------------------------------------------------------------ 37 31 # Interior region definitions 38 32 #------------------------------------------------------------------------------ 39 40 33 # Read interior polygons 41 34 poly_cairns = read_polygon('cairns.csv') … … 47 40 48 41 # Optionally plot points making up these polygons 49 #plot_polygons([bounding_polygon,poly_cairns,poly_island0,poly_island1,\ 50 # poly_island2,poly_island3,poly_shallow],\ 51 # style='boundingpoly',verbose=False) 52 53 42 #plot_polygons([bounding_polygon, poly_cairns, poly_island0, poly_island1, 43 # poly_island2, poly_island3, poly_shallow], 44 # style='boundingpoly', verbose=False) 54 45 55 46 # Define resolutions (max area per triangle) for each polygon 56 default_res = 10000000 # Background resolution47 default_res = 10000000 # Background resolution 57 48 islands_res = 100000 58 49 cairns_res = 100000 … … 60 51 61 52 # Define list of interior regions with associated resolutions 62 interior_regions = [[poly_cairns, cairns_res],53 interior_regions = [[poly_cairns, cairns_res], 63 54 [poly_island0, islands_res], 64 55 [poly_island1, islands_res], … … 66 57 [poly_island3, islands_res], 67 58 [poly_shallow, shallow_res]] 68 69 70 59 71 60 #------------------------------------------------------------------------------ … … 80 69 # Data for landslide 81 70 #------------------------------------------------------------------------------ 82 slide_origin = [451871, 8128376] # Assume to be on continental shelf71 slide_origin = [451871, 8128376] # Assume to be on continental shelf 83 72 slide_depth = 500. 84 85 86 87 88 -
anuga_core/documentation/user_manual/demos/cairns/runcairns.py
r6889 r7064 15 15 # Import necessary modules 16 16 #------------------------------------------------------------------------------ 17 18 17 # Standard modules 19 18 import os … … 37 36 import project # Definition of file names and polygons 38 37 39 40 38 #------------------------------------------------------------------------------ 41 39 # Preparation of topographic data 42 40 # Convert ASC 2 DEM 2 PTS using source data and store result in source data 43 41 #------------------------------------------------------------------------------ 44 45 42 # Create DEM from asc data 46 43 convert_dem_from_ascii2netcdf(project.demname, use_cache=True, verbose=True) … … 49 46 dem2pts(project.demname, use_cache=True, verbose=True) 50 47 51 52 48 #------------------------------------------------------------------------------ 53 49 # Create the triangular mesh and domain based on … … 55 51 # boundary and interior regions as defined in project.py 56 52 #------------------------------------------------------------------------------ 57 58 53 domain = create_domain_from_regions(project.bounding_polygon, 59 54 boundary_tags={'top': [0], … … 67 62 verbose=True) 68 63 69 70 64 # Print some stats about mesh and domain 71 65 print 'Number of triangles = ', len(domain) … … 76 70 # Setup parameters of computational domain 77 71 #------------------------------------------------------------------------------ 78 79 80 72 domain.set_name('cairns_' + project.scenario) # Name of sww file 81 73 domain.set_datadir('.') # Store sww output here 82 74 domain.set_minimum_storable_height(0.01) # Store only depth > 1cm 83 75 84 85 76 #------------------------------------------------------------------------------ 86 77 # Setup initial conditions 87 78 #------------------------------------------------------------------------------ 88 89 79 tide = 0.0 90 80 domain.set_quantity('stage', tide) … … 96 86 alpha=0.1) 97 87 98 99 88 #------------------------------------------------------------------------------ 100 89 # Setup information for slide scenario (to be applied 1 min into simulation 101 90 #------------------------------------------------------------------------------ 102 103 91 if project.scenario == 'slide': 104 92 # Function for submarine slide … … 113 101 verbose=True) 114 102 115 116 103 #------------------------------------------------------------------------------ 117 104 # Setup boundary conditions 118 105 #------------------------------------------------------------------------------ 119 120 106 print 'Available boundary tags', domain.get_boundary_tags() 121 122 107 123 108 Bd = Dirichlet_boundary([tide,0,0]) # Mean water level 124 109 Bs = Transmissive_stage_zero_momentum_boundary(domain) # Neutral boundary 125 126 110 127 111 if project.scenario == 'fixed_wave': … … 129 113 Bw = Time_boundary(domain=domain, 130 114 function=lambda t: [(60<t<3660)*50, 0, 0]) 131 132 115 domain.set_boundary({'ocean_east': Bw, 133 116 'bottom': Bs, … … 141 124 'onshore': Bd, 142 125 'top': Bd}) 143 144 126 145 127 #------------------------------------------------------------------------------ 146 128 # Evolve system through time 147 129 #------------------------------------------------------------------------------ 148 149 130 import time 150 131 t0 = time.time() … … 154 135 155 136 if project.scenario == 'slide': 156 157 137 for t in domain.evolve(yieldstep=10, finaltime=60): 158 138 print domain.timestepping_statistics() … … 171 151 print domain.boundary_statistics(tags='ocean_east') 172 152 173 174 153 if project.scenario == 'fixed_wave': 175 176 154 # Save every two mins leading up to wave approaching land 177 155 for t in domain.evolve(yieldstep=120, finaltime=5000): -
anuga_core/documentation/user_manual/demos/channel1.py
r5173 r7064 12 12 from anuga.shallow_water import Dirichlet_boundary 13 13 14 15 14 #------------------------------------------------------------------------------ 16 15 # Setup computational domain … … 21 20 domain = Domain(points, vertices, boundary) # Create domain 22 21 domain.set_name('channel1') # Output name 23 24 22 25 23 #------------------------------------------------------------------------------ … … 34 32 expression='elevation + 0.0') 35 33 36 37 34 #------------------------------------------------------------------------------ 38 35 # Setup boundary conditions … … 43 40 domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br}) 44 41 45 46 42 #------------------------------------------------------------------------------ 47 43 # Evolve system through time 48 44 #------------------------------------------------------------------------------ 49 for t in domain.evolve(yieldstep = 0.2, finaltime =40.0):45 for t in domain.evolve(yieldstep=0.2, finaltime=40.0): 50 46 print domain.timestepping_statistics() 51 52 -
anuga_core/documentation/user_manual/demos/channel2.py
r6889 r7064 13 13 from anuga.shallow_water import Time_boundary 14 14 15 16 15 #------------------------------------------------------------------------------ 17 16 # Setup computational domain … … 26 25 domain.set_name('channel2') # Output name 27 26 28 29 27 #------------------------------------------------------------------------------ 30 28 # Setup initial conditions … … 38 36 expression='elevation') # Dry initial condition 39 37 40 41 38 #------------------------------------------------------------------------------ 42 39 # Setup boundary conditions … … 48 45 domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br}) 49 46 50 51 47 #------------------------------------------------------------------------------ 52 48 # Evolve system through time 53 49 #------------------------------------------------------------------------------ 54 for t in domain.evolve(yieldstep = 0.2, finaltime =40.0):50 for t in domain.evolve(yieldstep=0.2, finaltime=40.0): 55 51 print domain.timestepping_statistics() 56 52 … … 59 55 print 'Stage > 0: Changing to outflow boundary' 60 56 domain.set_boundary({'right': Bo}) 61 62 63 64 -
anuga_core/documentation/user_manual/demos/channel3.py
r6889 r7064 13 13 from anuga.shallow_water import Time_boundary 14 14 15 16 15 #------------------------------------------------------------------------------ 17 16 # Setup computational domain … … 23 22 points, vertices, boundary = rectangular_cross(int(length/dx), int(width/dy), 24 23 len1=length, len2=width) 25 domain = Domain(points, vertices, boundary) 24 domain = Domain(points, vertices, boundary) 26 25 domain.set_name('channel3') # Output name 27 26 print domain.statistics() 28 29 27 30 28 #------------------------------------------------------------------------------ … … 32 30 #------------------------------------------------------------------------------ 33 31 def topography(x,y): 34 """Complex topography defined by a function of vectors x and y 35 """ 36 37 z = -x/10 32 """Complex topography defined by a function of vectors x and y.""" 33 34 z = -x/10 38 35 39 36 N = len(x) 40 37 for i in range(N): 41 42 38 # Step 43 39 if 10 < x[i] < 12: 44 z[i] += 0.4 - 0.05*y[i] 45 40 z[i] += 0.4 - 0.05*y[i] 41 46 42 # Constriction 47 43 if 27 < x[i] < 29 and y[i] > 3: 48 z[i] += 2 49 44 z[i] += 2 45 50 46 # Pole 51 47 if (x[i] - 34)**2 + (y[i] - 2)**2 < 0.4**2: … … 54 50 return z 55 51 56 57 domain.set_quantity('elevation', topography) # Use function for elevation 58 domain.set_quantity('friction', 0.01) # Constant friction 59 domain.set_quantity('stage', 60 expression='elevation') # Dry initial condition 61 52 domain.set_quantity('elevation', topography) # elevation is a function 53 domain.set_quantity('friction', 0.01) # Constant friction 54 domain.set_quantity('stage', expression='elevation') # Dry initial condition 62 55 63 56 #------------------------------------------------------------------------------ … … 70 63 domain.set_boundary({'left': Bi, 'right': Bo, 'top': Br, 'bottom': Br}) 71 64 72 73 65 #------------------------------------------------------------------------------ 74 66 # Evolve system through time 75 67 #------------------------------------------------------------------------------ 76 for t in domain.evolve(yieldstep = 0.1, finaltime =16.0):68 for t in domain.evolve(yieldstep=0.1, finaltime=16.0): 77 69 print domain.timestepping_statistics() 78 70 79 80 71 if domain.get_quantity('stage').\ 81 get_values(interpolation_points=[[10, 2.5]]) > 0: 72 get_values(interpolation_points=[[10, 2.5]]) > 0: 82 73 print 'Stage > 0: Changing to outflow boundary' 83 74 domain.set_boundary({'right': Bo}) -
anuga_core/documentation/user_manual/demos/runup.py
r6889 r7064 5 5 """ 6 6 7 8 7 #------------------------------------------------------------------------------ 9 8 # Import necessary modules 10 9 #------------------------------------------------------------------------------ 11 12 10 from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross 13 11 from anuga.shallow_water import Domain … … 22 20 # Setup computational domain 23 21 #------------------------------------------------------------------------------ 24 25 22 points, vertices, boundary = rectangular_cross(10, 10) # Basic mesh 26 23 … … 29 26 domain.set_datadir('.') # Use current directory for output 30 27 31 32 28 #------------------------------------------------------------------------------ 33 29 # Setup initial conditions 34 30 #------------------------------------------------------------------------------ 35 36 31 def topography(x,y): 37 return -x/2 # linear bed slope32 return -x/2 # linear bed slope 38 33 #return x*(-(2.0-x)*.5) # curved bed slope 39 34 … … 42 37 domain.set_quantity('stage', -.4) # Constant negative initial stage 43 38 44 45 39 #------------------------------------------------------------------------------ 46 40 # Setup boundary conditions 47 41 #------------------------------------------------------------------------------ 48 49 42 Br = Reflective_boundary(domain) # Solid reflective wall 50 43 Bt = Transmissive_boundary(domain) # Continue all values on boundary … … 56 49 domain.set_boundary({'left': Br, 'right': Bw, 'top': Br, 'bottom': Br}) 57 50 58 59 51 #------------------------------------------------------------------------------ 60 52 # Evolve system through time 61 53 #------------------------------------------------------------------------------ 62 63 54 for t in domain.evolve(yieldstep = 0.1, finaltime = 10.0): 64 55 print domain.timestepping_statistics() 65 66 67
Note: See TracChangeset
for help on using the changeset viewer.