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