source: trunk/anuga_validation/okushiri_2005/run_okushiri.py @ 7877

Last change on this file since 7877 was 7492, checked in by steve, 15 years ago

Changed the order of update_timestep and compute_forcing_terms.
Can now update the parameter domain.flux_timestep within
the computation of the forcing terms, so as to control stability of
the evolver.

Hopefully this will help in reducing the oscillations seen within the
culvert computation.

File size: 3.2 KB
Line 
1"""Validation of the AnuGA implementation of the shallow water wave equation.
2
3This script sets up Okushiri Island benchmark as published at the
4
5THE THIRD INTERNATIONAL WORKSHOP ON LONG-WAVE RUNUP MODELS
6June 17-18 2004
7Wrigley Marine Science Center
8Catalina Island, California
9http://www.cee.cornell.edu/longwave/
10
11
12The validation data was downloaded and made available in this directory
13for convenience but the original data is available at
14http://www.cee.cornell.edu/longwave/index.cfm?page=benchmark&problem=2
15where a detailed description of the problem is also available.
16
17
18Run create_okushiri.py to process the boundary condition and build a the
19mesh before running this script.
20
21"""
22
23# Module imports
24from anuga.interface import Domain
25from anuga.interface import Reflective_boundary
26from anuga.interface import Time_boundary
27from anuga.interface import Transmissive_momentum_set_stage_boundary
28from anuga.interface import Transmissive_n_momentum_zero_t_momentum_set_stage_boundary
29from anuga.interface import file_function
30
31import project
32
33
34#-------------------------
35# Create Domain from mesh
36#-------------------------
37domain = Domain(project.mesh_filename, use_cache=True, verbose=True)
38print domain.statistics()
39
40
41#-------------------------
42# Initial Conditions
43#-------------------------
44domain.set_quantity('friction', 0.0)
45domain.set_quantity('stage', 0.0)
46domain.set_quantity('elevation',
47                    filename=project.bathymetry_filename,
48                    alpha=0.02,                   
49                    verbose=True,
50                    use_cache=True)
51
52
53#-------------------------
54# Set simulation parameters
55#-------------------------
56domain.set_name(project.output_filename)  # Name of output sww file
57domain.set_default_order(2)               # Apply second order scheme
58domain.set_minimum_storable_height(0.001) # Don't store h < 0.001m
59
60#Timings on AMD64-242 (beta_h=0)
61#  tight_slope_limiters = 0:
62#    3035s - 3110s
63#  tight_slope_limiters = 1:
64#    3000s - 3008s
65#
66# beta_h>0: In the order of 3200s
67
68#-------------------------
69# Boundary Conditions
70#-------------------------
71
72# Create boundary function from timeseries provided in file
73function = file_function(project.boundary_filename,
74                         domain, verbose=True)
75
76# Create and assign boundary objects
77Bts = Transmissive_momentum_set_stage_boundary(domain, function)
78Br = Reflective_boundary(domain)
79domain.set_boundary({'wave': Bts, 'wall': Br})
80
81
82# Select triangle containing ch5 for diagnostic output
83# around known gauge
84triangle_id = domain.get_triangle_containing_point([4.521, 1.196])
85# This should get triangle id 32833 with centroid (4.5244, 1.1972)
86bdry_id = domain.get_triangle_containing_point([0.000, 1.696])
87print 'bdry_id = ',bdry_id
88
89
90#-------------------------
91# Evolve through time
92#-------------------------
93import time
94t0 = time.time()
95
96for t in domain.evolve(yieldstep = 0.05, finaltime = 22.5):
97    domain.write_time()
98    print  function(domain.get_time())[0],' bdry_d = ',bdry_id,' ',\
99          domain.get_conserved_quantities(bdry_id, edge=0)[0],' ',\
100          domain.get_conserved_quantities(bdry_id, edge=1)[0],' ',\
101          domain.get_conserved_quantities(bdry_id, edge=2)[0]   
102   
103
104print 'That took %.2f seconds' %(time.time()-t0)
Note: See TracBrowser for help on using the repository browser.