source: trunk/anuga_core/documentation/user_manual/examples/runsydney.py @ 7808

Last change on this file since 7808 was 7808, checked in by hudson, 12 years ago

Examples work with new API.

File size: 5.4 KB
Line 
1"""Script for running a tsunami inundation scenario for Sydney, NSW, 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.outputdir
6
7The scenario is defined by a triangular mesh created from project.polygon,
8the elevation data and a simulated submarine landslide.
9
10"""
11
12
13#------------------------------------------------------------------------------
14# Import necessary modules
15#------------------------------------------------------------------------------
16
17# Standard modules
18import os
19import time
20
21# Related major packages
22import anuga
23
24# Application specific imports
25import project                           # Define file names and polygons
26
27
28
29#------------------------------------------------------------------------------
30# Prepare bathymetric and topographic data
31#------------------------------------------------------------------------------
32
33# filenames
34#coarsedemname = project.coarsedemname + '.pts'
35#finedemname = project.finedemname + '.pts'
36demname = project.demname
37meshname = project.meshname+'.msh'
38
39anuga.asc2dem(demname, use_cache=True, verbose=True)
40
41# creates pts file
42anuga.dem2pts(demname, use_cache=True, verbose=True)
43
44
45#------------------------------------------------------------------------------
46# Setup computational domain
47#------------------------------------------------------------------------------
48
49# Interior regions and mesh resolutions
50interior_res = 5000
51shallow_res = 15000
52##interior_regions = [[project.northern_polygon, interior_res],
53##                    [project.harbour_polygon, interior_res],
54##                    [project.southern_polygon, interior_res],
55##                    [project.manly_polygon, interior_res],
56##                    [project.top_polygon, interior_res]]
57interior_regions = [[project.coastal_polygon, interior_res],
58                    [project.shallow_polygon, shallow_res]]
59
60anuga.create_mesh_from_regions(project.demopoly,
61                         boundary_tags= {'oceannorth': [0], 
62                                         'onshorenorth1': [1],
63                                         'onshorenorthwest1': [2],
64                                         'onshorenorthwest2': [3],
65                                         'onshorenorth2': [4],
66                                         'onshorewest': [5],
67                                         'onshoresouth': [6],
68                                         'oceansouth': [7],
69                                         'oceaneast': [8]},
70                         maximum_triangle_area=100000,
71                         filename=meshname,           
72                         interior_regions=interior_regions)
73
74
75#Create shallow water domain
76
77domain = anuga.Domain(meshname,
78                use_cache = True,
79                verbose = True)                                 
80
81
82print 'Number of triangles = ', len(domain)
83print 'The extent is ', domain.get_extent()
84
85domain.set_name(project.basename)
86domain.set_datadir(project.outputdir)
87domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
88
89
90#------------------------------------------------------------------------------
91# Set up scenario (tsunami_source is a callable object used with set_quantity)
92#------------------------------------------------------------------------------
93
94tsunami_source = anuga.slump_tsunami(length=30000.0,
95                               depth=400.0,
96                               slope=6.0,
97                               thickness=176.0, 
98                               radius=3330,
99                               dphi=0.23,
100                               x0=project.slump_origin[0], 
101                               y0=project.slump_origin[1], 
102                               alpha=0.0, 
103                               domain=domain)
104
105
106#------------------------------------------------------------------------------
107# Setup initial conditions
108#------------------------------------------------------------------------------
109
110domain.set_quantity('stage', tsunami_source)
111domain.set_quantity('friction', 0.0) 
112domain.set_quantity('elevation',
113                    filename = demname + '.pts',
114                    use_cache = True,
115                    verbose = True)
116
117#------------------------------------------------------------------------------
118# Setup boundary conditions (all Dirichlet)
119#------------------------------------------------------------------------------
120
121print 'Available boundary tags', domain.get_boundary_tags()
122
123Bd = anuga.Dirichlet_boundary([0.0,0.0,0.0])
124domain.set_boundary( {'oceannorth': Bd,
125                      'onshorenorth1': Bd,
126                      'onshorenorthwest1': Bd,
127                      'onshorenorthwest2': Bd,
128                      'onshorenorth2': Bd,
129                      'onshorewest': Bd,
130                      'onshoresouth': Bd,
131                      'oceansouth': Bd,
132                      'oceaneast': Bd} )
133
134#------------------------------------------------------------------------------
135# Evolve system through time
136#------------------------------------------------------------------------------
137
138import time
139t0 = time.time()
140
141for t in domain.evolve(yieldstep = 60, finaltime = 7200):
142    print domain.timestepping_statistics()
143    print domain.boundary_statistics(tags = 'oceaneast')   
144   
145print 'That took %.2f seconds' %(time.time()-t0)
Note: See TracBrowser for help on using the repository browser.