1 | |
---|
2 | |
---|
3 | """test_inflow_using_flowline |
---|
4 | Test the ability of a flowline to match inflow above the flowline by |
---|
5 | creating constant inflow onto a circle at the head of a 20m |
---|
6 | wide by 400m long plane dipping at various slopes with a perpendicular flowline and gauge |
---|
7 | downstream of the inflow and a 45 degree flowlines at 200m downstream |
---|
8 | """ |
---|
9 | |
---|
10 | verbose = True |
---|
11 | import Numeric as num |
---|
12 | |
---|
13 | #------------------------------------------------------------------------------ |
---|
14 | # Import necessary modules |
---|
15 | #------------------------------------------------------------------------------ |
---|
16 | from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross |
---|
17 | from anuga.shallow_water import Domain |
---|
18 | from anuga.shallow_water.shallow_water_domain import Reflective_boundary |
---|
19 | from anuga.shallow_water.shallow_water_domain import Dirichlet_boundary |
---|
20 | from anuga.shallow_water.shallow_water_domain import Transmissive_momentum_set_stage_boundary |
---|
21 | from anuga.shallow_water.shallow_water_domain import Inflow |
---|
22 | from anuga.shallow_water.data_manager import get_flow_through_cross_section |
---|
23 | from anuga.abstract_2d_finite_volumes.util import sww2csv_gauges, csv2timeseries_graphs |
---|
24 | |
---|
25 | |
---|
26 | #------------------------------------------------------------------------------ |
---|
27 | # Setup computational domain |
---|
28 | #------------------------------------------------------------------------------ |
---|
29 | |
---|
30 | number_of_inflows = 2 # Number of inflows on top of each other |
---|
31 | finaltime = 1000.0 |
---|
32 | |
---|
33 | length = 400. |
---|
34 | width = 20. |
---|
35 | dx = dy = 2 # Resolution: of grid on both axes |
---|
36 | |
---|
37 | points, vertices, boundary = rectangular_cross(int(length/dx), int(width/dy), |
---|
38 | len1=length, len2=width) |
---|
39 | |
---|
40 | for mannings_n in [0.012, 0.035, 0.070, 0.150]: |
---|
41 | # Loop over a range of roughnesses |
---|
42 | |
---|
43 | for slope in [1.0/300, 1.0/150, 1.0/75]: |
---|
44 | # Loop over a range of bedslopes |
---|
45 | |
---|
46 | |
---|
47 | domain = Domain(points, vertices, boundary) |
---|
48 | domain.set_name('inflow_flowline_test') |
---|
49 | |
---|
50 | |
---|
51 | #---------------------------------------------------------------------- |
---|
52 | # Setup initial conditions |
---|
53 | #---------------------------------------------------------------------- |
---|
54 | |
---|
55 | def topography(x, y): |
---|
56 | z=-x * slope |
---|
57 | return z |
---|
58 | |
---|
59 | domain.set_quantity('elevation', topography) |
---|
60 | domain.set_quantity('friction', mannings_n) |
---|
61 | domain.set_quantity('stage', |
---|
62 | expression='elevation') # Dry initial condition |
---|
63 | |
---|
64 | |
---|
65 | #---------------------------------------------------------------------- |
---|
66 | # Seup Inflow and compute flow characteristics |
---|
67 | #---------------------------------------------------------------------- |
---|
68 | |
---|
69 | # Fixed Flowrate onto Area |
---|
70 | fixed_inflow = Inflow(domain, |
---|
71 | center=(10.0, 10.0), |
---|
72 | radius=5.00, |
---|
73 | rate=10.00) |
---|
74 | |
---|
75 | # Stack this flow |
---|
76 | for i in range(number_of_inflows): |
---|
77 | domain.forcing_terms.append(fixed_inflow) |
---|
78 | |
---|
79 | ref_flow = fixed_inflow.rate*number_of_inflows |
---|
80 | |
---|
81 | # Compute normal depth on plane using Mannings equation |
---|
82 | # v=1/n*(r^2/3)*(s^0.5) or r=(Q*n/(s^0.5*W))^0.6 |
---|
83 | normal_depth=(ref_flow*mannings_n/(slope**0.5*width))**0.6 |
---|
84 | flow_velocity = ref_flow/(width*normal_depth) |
---|
85 | froude_no = flow_velocity/(9.8*normal_depth)**0.5 |
---|
86 | if verbose: |
---|
87 | print |
---|
88 | print 'Slope:', slope, 'Mannings n:', mannings_n, 'Velocity:',flow_velocity,'Froude:',froude_no |
---|
89 | |
---|
90 | |
---|
91 | #---------------------------------------------------------------------- |
---|
92 | # Setup boundary conditions |
---|
93 | #---------------------------------------------------------------------- |
---|
94 | |
---|
95 | Br = Reflective_boundary(domain) # Solid reflective wall |
---|
96 | |
---|
97 | # Define downstream boundary based on predicted depth |
---|
98 | def normal_depth_stage_downstream(t): |
---|
99 | return (-slope*length) + normal_depth |
---|
100 | |
---|
101 | Bt = Transmissive_momentum_set_stage_boundary(domain=domain, |
---|
102 | function=normal_depth_stage_downstream) |
---|
103 | |
---|
104 | domain.set_boundary({'left': Br, 'right': Bt, 'top': Br, 'bottom': Br}) |
---|
105 | |
---|
106 | #---------------------------------------------------------------------- |
---|
107 | # Evolve system through time |
---|
108 | #---------------------------------------------------------------------- |
---|
109 | |
---|
110 | x=200.0 # Reference location for flowlines and gauges |
---|
111 | y=10.00 |
---|
112 | |
---|
113 | for t in domain.evolve(yieldstep=100.0, finaltime=finaltime): |
---|
114 | if verbose : |
---|
115 | print domain.timestepping_statistics() |
---|
116 | print domain.volumetric_balance_statistics() |
---|
117 | |
---|
118 | |
---|
119 | #---------------------------------------------------------------------- |
---|
120 | # Compute flow thru flowlines ds of inflow and check against normal depth |
---|
121 | #---------------------------------------------------------------------- |
---|
122 | |
---|
123 | # Square on flowline at 200m |
---|
124 | q=domain.get_flow_through_cross_section([[200.0,0.0],[200.0,20.0]]) |
---|
125 | msg = 'Predicted square on flow was %f, should have been %f' % (q, ref_flow) |
---|
126 | if verbose: |
---|
127 | print '90 degree flowline: ANUGA = %f, Ref = %f' % (q, ref_flow) |
---|
128 | assert num.allclose(q, ref_flow, rtol=1.0e-2), msg |
---|
129 | |
---|
130 | # 45 degree flowline at 200m |
---|
131 | q=domain.get_flow_through_cross_section([[200.0,0.0],[220.0,20.0]]) |
---|
132 | msg = 'Predicted 45 flow was %f, should have been %f' % (q, ref_flow) |
---|
133 | if verbose: |
---|
134 | print '45 degree flowline: ANUGA = %f, Ref = %f' % (q, ref_flow) |
---|
135 | assert num.allclose(q, ref_flow, rtol=1.0e-2), msg |
---|
136 | |
---|
137 | |
---|
138 | # Stage recorder (gauge) in middle of plane at 200m |
---|
139 | w = domain.get_quantity('stage').get_values(interpolation_points=[[x, y]])[0] |
---|
140 | z = domain.get_quantity('elevation').get_values(interpolation_points=[[x, y]])[0] |
---|
141 | domain_depth = w-z |
---|
142 | |
---|
143 | |
---|
144 | msg = 'Predicted depth of flow at 200m was %f, should have been %f' % (normal_depth, domain_depth) |
---|
145 | msg = 'Mannings n = %s, slope = %s' % (mannings_n, slope) |
---|
146 | if verbose: |
---|
147 | print 'Depth at 200m: ANUGA = %f, Mannings = %f' % (domain_depth, normal_depth) |
---|
148 | assert num.allclose(domain_depth,normal_depth, rtol=1.0e-2), msg |
---|
149 | |
---|
150 | |
---|
151 | |
---|