[1845] | 1 | """Script for running a tsunami inundation scenario for Karratha, WA, Australia. |
---|
| 2 | |
---|
| 3 | Source data such as elevation and boundary data is assumed to be available in |
---|
| 4 | directories specified by project.py |
---|
| 5 | The output sww file is stored in project.outputdir |
---|
| 6 | |
---|
| 7 | The scenario is defined by a triangular mesh created from project.polygon, |
---|
| 8 | the elevation data and boundary data obtained from a tsunami simulation done with MOST. |
---|
| 9 | |
---|
| 10 | Ole Nielsen, GA - 2005 |
---|
[1781] | 11 | """ |
---|
| 12 | |
---|
[1855] | 13 | tide = 0.75 #HMWS estimate by Colin French, GA |
---|
| 14 | |
---|
| 15 | |
---|
[1781] | 16 | import os |
---|
| 17 | import time |
---|
| 18 | |
---|
[1784] | 19 | |
---|
[1845] | 20 | from pyvolution.shallow_water import Domain, Reflective_boundary, File_boundary,\ |
---|
| 21 | Dirichlet_boundary, Time_boundary, Transmissive_boundary |
---|
[1796] | 22 | from pyvolution.data_manager import convert_dem_from_ascii2netcdf,\ |
---|
| 23 | dem2pts, ferret2sww |
---|
[1784] | 24 | from pyvolution.pmesh2domain import pmesh_to_domain_instance |
---|
| 25 | from caching import cache |
---|
[1781] | 26 | import project |
---|
| 27 | |
---|
[1796] | 28 | #Data preparation |
---|
[1781] | 29 | #Convert ASC 2 DEM 2 PTS using source data and store result in source data |
---|
[1784] | 30 | demname = project.demname |
---|
| 31 | |
---|
| 32 | cache(convert_dem_from_ascii2netcdf, demname, {'verbose': True}, |
---|
| 33 | dependencies = [demname + '.asc'], |
---|
| 34 | verbose = True) |
---|
[1796] | 35 | #evaluate = True) |
---|
[1784] | 36 | |
---|
| 37 | cache(dem2pts, demname, {'verbose': True}, |
---|
| 38 | dependencies = [demname + '.dem'], |
---|
| 39 | verbose = True) |
---|
| 40 | |
---|
[1801] | 41 | |
---|
[1796] | 42 | #Convert MOST boundary |
---|
[1801] | 43 | source_dir = project.boundarydir |
---|
[1784] | 44 | |
---|
[1801] | 45 | from pyvolution.data_manager import ferret2sww |
---|
[1796] | 46 | |
---|
[1801] | 47 | south = project.south |
---|
| 48 | north = project.north |
---|
| 49 | west = project.west |
---|
| 50 | east = project.east |
---|
| 51 | |
---|
| 52 | |
---|
[1855] | 53 | cache(ferret2sww, |
---|
| 54 | (source_dir+project.boundary_basename, |
---|
| 55 | project.boundary_basename), |
---|
| 56 | {'verbose': True, |
---|
| 57 | 'minlat': south-1, |
---|
| 58 | 'maxlat': north+1, |
---|
| 59 | 'minlon': west-1, |
---|
| 60 | 'maxlon': east+1, |
---|
| 61 | 'origin': project.mesh_origin, |
---|
| 62 | 'mean_stage': tide, |
---|
| 63 | 'zscale': 1, |
---|
| 64 | 'fail_on_NaN': False, |
---|
| 65 | 'inverted_bathymetry': True}, |
---|
| 66 | verbose = True) |
---|
| 67 | #FIXME: Dependencies |
---|
[1845] | 68 | |
---|
[1855] | 69 | |
---|
[1843] | 70 | #ferret2sww(source_dir+project.boundary_basename, |
---|
| 71 | # project.boundary_basename, |
---|
| 72 | # verbose=True, |
---|
| 73 | # minlat=south-1, maxlat=north+1, |
---|
| 74 | # minlon=west-1, maxlon=east+1, |
---|
| 75 | # origin = project.mesh_origin, |
---|
[1855] | 76 | # mean_stage = tide, |
---|
[1843] | 77 | # zscale = 1, |
---|
| 78 | # fail_on_NaN = False, |
---|
[1855] | 79 | # inverted_bathymetry = True) |
---|
[1832] | 80 | |
---|
[1801] | 81 | |
---|
[1855] | 82 | #Create Triangular Mesh |
---|
| 83 | from create_mesh import create_mesh |
---|
[1801] | 84 | |
---|
[1855] | 85 | interior_regions = [[project.karratha_polygon, 10000]] |
---|
| 86 | m = cache(create_mesh, |
---|
| 87 | project.polygon, |
---|
| 88 | {'boundary_tags': {'back': [7, 8], 'side': [0, 6], 'ocean': [1, 2, 3, 4, 5]}, |
---|
| 89 | 'resolution': 80000, |
---|
| 90 | 'filename': project.meshname + '.msh', |
---|
| 91 | 'interior_regions': interior_regions}, |
---|
| 92 | verbose = True) |
---|
| 93 | #verbose = True, |
---|
| 94 | #evaluate = True ) |
---|
| 95 | |
---|
| 96 | |
---|
[1801] | 97 | #Setup domain |
---|
[1784] | 98 | mesh = project.meshname + '.msh' |
---|
| 99 | domain = cache(pmesh_to_domain_instance, (mesh, Domain), |
---|
| 100 | dependencies = [mesh], |
---|
| 101 | verbose = True) |
---|
| 102 | |
---|
| 103 | |
---|
| 104 | domain.set_name(project.basename) |
---|
[1786] | 105 | domain.set_datadir(project.outputdir) |
---|
| 106 | domain.store = True |
---|
| 107 | |
---|
[1855] | 108 | #domain.quantities_to_be_stored = ['stage', 'xmomentum', 'ymomentum'] |
---|
| 109 | domain.quantities_to_be_stored = ['stage'] |
---|
| 110 | |
---|
[1832] | 111 | print 'Number of triangles = ', len(domain) |
---|
[1785] | 112 | print 'The extent is ', domain.get_extent() |
---|
[1784] | 113 | |
---|
[1796] | 114 | |
---|
[1845] | 115 | |
---|
| 116 | #Setup Initial Conditions |
---|
[1843] | 117 | domain.set_quantity('friction', 0) |
---|
[1832] | 118 | domain.set_quantity('stage', tide) |
---|
[1786] | 119 | domain.set_quantity('elevation', |
---|
| 120 | filename = demname + '.pts', |
---|
| 121 | use_cache = True, |
---|
| 122 | verbose = True) |
---|
[1784] | 123 | |
---|
[1801] | 124 | |
---|
[1845] | 125 | |
---|
| 126 | #Setup Boundary Conditions |
---|
[1830] | 127 | print domain.get_boundary_tags() |
---|
| 128 | |
---|
[1832] | 129 | Bf = File_boundary(project.boundary_basename + '.sww', domain, verbose = True) |
---|
[1855] | 130 | #domain.starttime = 3000 #Obtained from MOST |
---|
| 131 | domain.starttime = 0 #Obtained from MOST |
---|
[1830] | 132 | |
---|
[1784] | 133 | Br = Reflective_boundary(domain) |
---|
[1832] | 134 | Bt = Transmissive_boundary(domain) |
---|
| 135 | Bd = Dirichlet_boundary([tide,0,0]) |
---|
[1830] | 136 | Bw = Time_boundary(domain=domain, |
---|
[1832] | 137 | f=lambda t: [(60<t<660)*4, 0, 0]) |
---|
[1784] | 138 | |
---|
[1832] | 139 | domain.set_boundary( {'back': Br,'side': Bd, 'ocean': Bf} ) |
---|
[1784] | 140 | |
---|
[1843] | 141 | #domain.set_boundary( {'back': Br,'side': Bd, 'ocean': Bd} ) #Nothing |
---|
[1832] | 142 | |
---|
[1843] | 143 | |
---|
[1784] | 144 | #Run |
---|
[1855] | 145 | #for t in domain.evolve(yieldstep = 600, finaltime = 15000): |
---|
| 146 | # domain.write_time() |
---|
| 147 | # domain.write_boundary_statistics(tags = 'ocean') #quantities = 'stage') |
---|
| 148 | # |
---|
| 149 | #for t in domain.evolve(yieldstep = 10, finaltime = 35000): |
---|
| 150 | # domain.write_time() |
---|
| 151 | # domain.write_boundary_statistics(tags = 'ocean') #quantities = 'stage') |
---|
| 152 | |
---|
[1843] | 153 | for t in domain.evolve(yieldstep = 60, finaltime = 40000): |
---|
[1784] | 154 | domain.write_time() |
---|
[1855] | 155 | domain.write_boundary_statistics(tags = 'ocean') #quantities = 'stage') |
---|