source: anuga_validation/convergence_study/dam_break.py @ 5300

Last change on this file since 5300 was 5300, checked in by steve, 16 years ago

Working version of wave.py and dam_break.py

File size: 4.0 KB
Line 
1"""Simple water flow example using ANUGA
2
3Water driven up a linear slope and time varying boundary,
4similar to a beach environment
5"""
6
7#------------------------------------------------------------------------------
8# Import necessary modules
9#------------------------------------------------------------------------------
10
11import sys
12from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
13from anuga.shallow_water import Domain
14from anuga.shallow_water import Reflective_boundary
15from anuga.shallow_water import Dirichlet_boundary
16from anuga.shallow_water import Time_boundary
17from anuga.shallow_water import Transmissive_boundary
18from anuga.shallow_water import Transmissive_Momentum_Set_Stage_boundary
19from anuga.shallow_water.data_manager import start_screen_catcher, copy_code_files
20
21from math import cos
22from Numeric import zeros, Float
23from time import localtime, strftime, gmtime
24
25
26
27#-------------------------------------------------------------------------------
28# Copy scripts to time stamped output directory and capture screen
29# output to file
30#-------------------------------------------------------------------------------
31time = strftime('%Y%m%d_%H%M%S',localtime())
32
33output_dir = 'dam_break_'+time
34output_file = 'dam_break'
35
36copy_code_files(output_dir,__file__)
37#start_screen_catcher(output_dir+'_')
38
39
40#------------------------------------------------------------------------------
41# Setup domain
42#------------------------------------------------------------------------------
43dx = 200.
44dy = dx
45L = 100000.
46W = 10*dx
47
48# structured mesh
49points, vertices, boundary = rectangular_cross(int(L/dx), int(W/dy), L, W, (0.0, -W/2))
50
51domain = Domain(points, vertices, boundary) 
52
53domain.set_name(output_file)               
54domain.set_datadir(output_dir) 
55
56#------------------------------------------------------------------------------
57# Setup Algorithm
58#------------------------------------------------------------------------------
59domain.set_timestepping_method('rk2')
60domain.set_default_order(2)
61
62print domain.get_timestepping_method()
63
64domain.use_edge_limiter = True
65domain.tight_slope_limiters = True
66domain.use_centroid_velocities = False
67
68domain.CFL = 1.0
69
70domain.beta_w      = 0.6
71domain.beta_uh     = 0.6
72domain.beta_vh     = 0.6
73domain.beta_h    = 0.0
74
75
76#------------------------------------------------------------------------------
77# Setup initial conditions
78#------------------------------------------------------------------------------
79domain.set_quantity('elevation',0.0)
80domain.set_quantity('friction', 0.0)
81
82h0 = 10.0
83h1 = 0.0
84
85def height(x,y):
86    z = zeros(len(x),Float)
87    for i in range(len(x)):
88        if x[i]<=50000.0:
89            z[i] = h0
90        else:
91            z[i] = h1
92    return z
93
94
95domain.set_quantity('stage', height)
96#-----------------------------------------------------------------------------
97# Setup boundary conditions
98#------------------------------------------------------------------------------
99from math import sin, pi, exp
100Br = Reflective_boundary(domain)      # Solid reflective wall
101Bt = Transmissive_boundary(domain)    # Continue all values on boundary
102Bd = Dirichlet_boundary([1,0.,0.]) # Constant boundary values
103
104
105# Associate boundary tags with boundary objects
106domain.set_boundary({'left': Bt, 'right': Bt, 'top': Br, 'bottom': Br})
107
108
109#===============================================================================
110from anuga.visualiser import RealtimeVisualiser
111vis = RealtimeVisualiser(domain)
112vis.render_quantity_height("stage", zScale =h0*500, dynamic=True)
113vis.colour_height_quantity('stage', (0.0, 0.5, 1.0))
114vis.start()
115#===============================================================================
116
117
118#------------------------------------------------------------------------------
119# Evolve system through time
120#------------------------------------------------------------------------------
121
122for t in domain.evolve(yieldstep = 100.0, finaltime = 60*60.):
123    print domain.timestepping_statistics(track_speeds=True)
124    vis.update()
125   
126vis.evolveFinished()
127
Note: See TracBrowser for help on using the repository browser.