source: anuga_work/production/gippsland_2005/run_gippsland2.py @ 5587

Last change on this file since 5587 was 3514, checked in by duncan, 18 years ago

Hi all,
I'm doing a change in the anuga structure, moving the code to

\anuga_core\source\anuga

After you have done an svn update, the PYTHONPATH has to be changed to;
PYTHONPATH = anuga_core/source/

This is part of changes required to make installation of anuga quicker and reducing the size of our sandpits.

If any imports are broken, try fixing them. With adding anuga. to them for example. If this seems to have really broken things, email/phone me.

Cheers
Duncan

File size: 4.8 KB
Line 
1
2#Convention for strings representing files
3# #_file has the extention
4#           #name does not have the extension
5
6import time
7
8from anuga.pyvolution.data_manager import convert_dem_from_ascii2netcdf,\
9     dem2pts, asc_csiro2sww
10from anuga.pyvolution.pmesh2domain import pmesh_to_domain_instance
11from caching import cache
12from create_mesh2 import create_mesh
13from anuga.pyvolution.shallow_water import Domain, Reflective_boundary,\
14     File_boundary, Dirichlet_boundary, Time_boundary, Transmissive_boundary
15from anuga.pyvolution.least_squares import fit_to_mesh_file, DEFAULT_ALPHA
16
17# Import this last!  It stuffs up the loading of c extensions otherwise.
18import project
19
20#Data preparation
21#Convert ASC 2 DEM 2 PTS using source data and store result in source data
22# just interested in the boundary/ mesh intereaction first off.
23#cache(convert_dem_from_ascii2netcdf, demname, {'verbose': True},
24#      dependencies = [demname + '.asc'],
25#      verbose = True)
26#      #evaluate = True)
27
28cache(dem2pts, project.demname[:-4], {'verbose': True},
29      dependencies = [project.demname]     
30      ,verbose = False
31      )
32
33
34# Convert CSIRO boundary
35# NOTE: This function is dependent on a lot of files, and I haven't listed
36# them, since there is so many and I don't think they will be changing.
37cache(asc_csiro2sww,
38      (project.bath_dir,
39       project.elevation_dir,
40       project.ucur_dir,
41       project.vcur_dir,
42       project.boundaryname),
43       {'verbose': True,
44       'minlat': -38.25,
45       'maxlat': -37.7,
46       'minlon': 147.0,
47       'maxlon': 148.25,
48       'mean_stage': project.tide,
49       'zscale': 1,                 #Enhance storm surge
50       'fail_on_NaN': False}
51      #,evaluate = True                 
52      ,verbose = False
53      )
54
55meshname = project.meshname
56pointname = project.pointname
57mesh_elevname = project.mesh_elevname
58outputname = project.outputname
59
60t0 = time.time()
61#Create the mesh without elevation data
62# 100000 very course - 16,827 triangle mesh
63
64meshname, triagle_count = cache(create_mesh,(1000000),
65                                {'factor':1000000,
66                                 'mesh_file':meshname,
67                                 'triangles_in_name':True}
68                                ,dependencies = ['create_mesh2.py']
69                                #,evaluate = True     
70                                )
71
72mesh_elevname = mesh_elevname[:-9] + '_' + str(triagle_count) +  \
73                mesh_elevname[-9:]
74outputname = outputname[:-4] + '_' + str(triagle_count) + outputname[-4:]
75
76
77#Add elevation data to the mesh
78#fit_to_mesh_file(mesh_file, point_file, mesh_output_file, DEFAULT_ALPHA,
79#                 verbose=True, expand_search=True, precrop=True)     
80cache(fit_to_mesh_file,(meshname,
81                 pointname,
82                 mesh_elevname,
83                 DEFAULT_ALPHA),
84      {'verbose': True,
85       'expand_search': True,
86       'precrop': True}
87      ,dependencies = [meshname, pointname]
88      #,evaluate = True     
89      ,verbose = False
90      )
91 
92print 'Initialising the mesh took %.2f seconds' %(time.time()-t0) 
93
94#Setup domain
95domain = cache(pmesh_to_domain_instance, (mesh_elevname, Domain),
96               dependencies = [mesh_elevname]                   
97               ,verbose = False
98               )               
99
100
101domain.set_name(project.basename + '_%d' %triagle_count)
102domain.set_datadir(project.outputdir)
103domain.store = True
104domain.quantities_to_be_stored = ['stage']
105
106print 'Number of triangles = ', len(domain)
107print 'The extent is ', domain.get_extent()
108
109#Setup Initial Conditions
110domain.set_quantity('friction', 0)
111domain.set_quantity('stage', project.tide)
112
113#domain.max_size = 3000
114#Setup Boundary Conditions
115print domain.get_boundary_tags()
116
117print "****** run  ****"
118print "project.boundaryname",project.boundaryname
119print "**********"
120
121Bf = File_boundary(project.boundaryname, domain, verbose = True)
122#domain.starttime = 3000  #Obtained from MOST
123domain.starttime = 0  #Obtained from MOST
124
125Br = Reflective_boundary(domain)
126Bt = Transmissive_boundary(domain)
127Bd = Dirichlet_boundary([project.tide,0,0])
128#Bw = Time_boundary(domain=domain,
129#                   f=lambda t: [(60<t<660)*4, 0, 0])
130
131domain.set_boundary( {'back': Br,'side': Br, 'ocean': Bf} ) #CSIRO storm surge
132
133#domain.set_boundary( {'back': Br,'side': Br, 'ocean': Br} ) #NO storm surge
134
135#Evolve
136t0 = time.time()
137
138for t in domain.evolve(yieldstep = 900, finaltime = 432000):
139#for t in domain.evolve(yieldstep = 10, finaltime = 30):
140#                       skip_initial_step = True):
141    domain.write_time()
142    domain.write_boundary_statistics(tags = 'ocean') #quantities = 'stage')       
143
144print 'That took %.2f seconds' %(time.time()-t0)
Note: See TracBrowser for help on using the repository browser.