source: development/dam_2006/run_dam.py @ 3226

Last change on this file since 3226 was 3226, checked in by duncan, 19 years ago

tweaking dam simulation

File size: 3.6 KB
RevLine 
[3125]1"""Script for running a tsunami inundation scenario for Onslow, WA, Australia.
2
3Source data such as elevation and boundary data is assumed to be available in
4directories specified by project.py
5The output sww file is stored in project.outputtimedir
6
7The scenario is defined by a triangular mesh created from project.polygon,
8the elevation data and a simulated submarine landslide.
9
10Ole Nielsen and Duncan Gray, GA - 2005 and Nick Bartzis, GA - 2006
11"""
12
13
[3178]14#----------------------------------------------------------------------------
15# Import necessary modules
16#----------------------------------------------------------------------------
[3125]17
18# Standard modules
19import time
20import sys
21from shutil import copy
[3226]22from os import mkdir, access, F_OK, path
[3125]23
24# Related major packages
25from pyvolution.shallow_water import Domain, Reflective_boundary, \
26                            Dirichlet_boundary, Time_boundary, File_boundary
27from pyvolution.util import Screen_Catcher
28from pyvolution.region import Set_region
29
30# Application specific imports
31import project                 # Definition of file names and polygons
32import create_mesh
33
[3145]34
[3178]35#-----------------------------------------------------------------------------
[3145]36# Setup archiving of simulations
[3178]37#-----------------------------------------------------------------------------
[3145]38
[3125]39copy (project.codedirname, project.outputtimedir + 'project.py')
40copy (project.codedir + 'run_dam.py', project.outputtimedir + 'run_dam.py')
[3178]41copy (project.codedir + 'create_mesh.py',
42      project.outputtimedir + 'create_mesh.py')
[3125]43print'output dir', project.outputtimedir
44
45#normal screen output is stored in
46screen_output_name = project.outputtimedir + "screen_output.txt"
47screen_error_name = project.outputtimedir + "screen_error.txt"
48
[3178]49#-----------------------------------------------------------------------------
[3125]50# Create the triangular mesh
[3178]51#-----------------------------------------------------------------------------
[3125]52
53create_mesh.generate() # this creates the mesh
54
[3226]55head,tail = path.split(project.mesh_filename)
56copy (project.mesh_filename,
57      project.outputtimedir + tail )
[3178]58#-----------------------------------------------------------------------------
[3125]59# Setup computational domain
[3178]60#-----------------------------------------------------------------------------
[3125]61domain = Domain(project.mesh_filename, use_cache = False, verbose = True)
62   
63
64print 'Number of triangles = ', len(domain)
65print 'The extent is ', domain.get_extent()
66print domain.statistics()
67
68domain.set_name(project.basename)
69domain.set_datadir(project.outputtimedir)
70domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
71
[3178]72#-----------------------------------------------------------------------------
[3125]73# Setup initial conditions
[3178]74#-----------------------------------------------------------------------------
[3125]75
76tide = 0.0
77
78def elevation_tilt(x, y):
79    return -x*0.15
80       
81domain.set_quantity('stage', elevation_tilt)
82domain.set_quantity('friction', 0.03) 
83domain.set_quantity('elevation',elevation_tilt)
84
85print 'Available boundary tags', domain.get_boundary_tags()
86
87domain.set_region(Set_region('dam','stage',0.4,location = 'unique vertices'))
88
89Br = Reflective_boundary(domain)
90Bd = Dirichlet_boundary([tide,0,0])
91
92#domain.set_boundary( {'wall': Br, 'wave': Bw} )
93domain.set_boundary( {'wall': Br, 'wave': Br} )
94
[3178]95#-----------------------------------------------------------------------------
[3125]96# Evolve system through time
[3178]97#-----------------------------------------------------------------------------
[3125]98import time
99t0 = time.time()
100
[3145]101for t in domain.evolve(yieldstep = 0.1, finaltime = 5):
[3125]102    domain.write_time()
103   
104print 'That took %.2f seconds' %(time.time()-t0)
105
106print 'finished'
Note: See TracBrowser for help on using the repository browser.