Changeset 7317 for anuga_core/source/anuga/culvert_flows
- Timestamp:
- Jul 22, 2009, 9:22:11 AM (16 years ago)
- Location:
- anuga_core/source/anuga/culvert_flows
- Files:
-
- 4 edited
-
Test_Culvert_Flat_Water_Lev.py (modified) (1 diff)
-
culvert_class.py (modified) (12 diffs)
-
culvert_routines.py (modified) (10 diffs)
-
test_culvert_class.py (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/culvert_flows/Test_Culvert_Flat_Water_Lev.py
r7276 r7317 11 11 12 12 """ 13 13 14 print 'Starting.... Importing Modules...' 15 14 16 #------------------------------------------------------------------------------ 15 17 # Import necessary modules -
anuga_core/source/anuga/culvert_flows/culvert_class.py
r7276 r7317 13 13 from anuga.config import g, epsilon 14 14 from anuga.config import minimum_allowed_height, velocity_protection 15 import anuga.utilities.log as log 15 16 16 17 import numpy as num … … 287 288 msg += ' does not match distance between specified' 288 289 msg += ' end points (%.2f m)' %self.length 289 print msg290 log.critical(msg) 290 291 291 292 self.verbose = verbose … … 441 442 msg += 'Q will be reduced from %.2f m^3/s to %.2f m^3/s.' % (Q, Q_reduced) 442 443 if self.verbose is True: 443 print msg444 log.critical(msg) 444 445 445 446 if self.log_filename is not None: … … 506 507 delta_total_energy = openings[0].total_energy - openings[1].total_energy 507 508 if delta_total_energy > 0: 508 #print 'Flow U/S ---> D/S'509 509 inlet = openings[0] 510 510 outlet = openings[1] 511 511 else: 512 #print 'Flow D/S ---> U/S'513 512 inlet = openings[1] 514 513 outlet = openings[0] … … 530 529 # Flow will be purely controlled by uphill outlet face 531 530 if self.verbose is True: 532 print '%.2fs - WARNING: Flow is running uphill.' % time531 log.critical('%.2fs - WARNING: Flow is running uphill.' % time) 533 532 534 533 if self.log_filename is not None: … … 789 788 msg += ' does not match distance between specified' 790 789 msg += ' end points (%.2f m)' %self.length 791 print msg790 log.critical(msg) 792 791 793 792 self.verbose = verbose … … 906 905 V += d * domain.areas[i] 907 906 908 #Vsimple = mean(depth)*self.inlet.exchange_area # Current volume in exchange area909 #print 'Q', Q, 'dt', delta_t, 'Q*dt', Q*delta_t, 'V', V, 'Vsimple', Vsimple910 911 907 dt = delta_t 912 908 if Q*dt > V: … … 918 914 msg += ' Q will be reduced from %.2f m^3/s to %.2f m^3/s.' % (Q, Q_reduced) 919 915 920 #print msg921 922 916 if self.verbose is True: 923 print msg917 log.critical(msg) 924 918 if hasattr(self, 'log_filename'): 925 919 log_to_file(self.log_filename, msg) … … 1313 1307 1314 1308 opening.total_energy = 0.5*(u*u + v*v)/g + stage 1315 #print 'Et = %.3f m' %opening.total_energy1316 1309 1317 1310 # Store current average stage and depth with each opening object … … 1330 1323 1331 1324 if delta_Et > 0: 1332 #print 'Flow U/S ---> D/S'1333 1325 inlet = openings[0] 1334 1326 outlet = openings[1] … … 1338 1330 1339 1331 else: 1340 #print 'Flow D/S ---> U/S'1341 1332 inlet = openings[1] 1342 1333 outlet = openings[0] … … 1361 1352 # Flow will be purely controlled by uphill outlet face 1362 1353 if self.verbose is True: 1363 print 'WARNING: Flow is running uphill. Watch Out!', inlet.elevation, outlet.elevation 1354 log.critical('WARNING: Flow is running uphill. Watch Out! ' 1355 'inlet.elevation=%s, outlet.elevation%s' 1356 % (str(inlet.elevation), str(outlet.elevation))) 1364 1357 1365 1358 -
anuga_core/source/anuga/culvert_flows/culvert_routines.py
r7276 r7317 86 86 if inlet_depth > 0.1: #this value was 0.01: 87 87 if local_debug =='true': 88 print 'Specific E & Deltat Tot E = ',inlet_specific_energy,delta_total_energy 89 print 'culvert typ = ',culvert_type 88 log.critical('Specific E & Deltat Tot E = %s, %s' 89 % (str(inlet_specific_energy), 90 str(delta_total_energy))) 91 log.critical('culvert type = %s' % str(culvert_type)) 90 92 # Water has risen above inlet 91 93 … … 138 140 flow_width= diameter 139 141 case = 'Inlet CTRL Outlet submerged Circular PIPE FULL' 140 if local_debug =='true': 141 print 'Inlet CTRL Outlet submerged Circular PIPE FULL' 142 if local_debug == 'true': 143 log.critical('Inlet CTRL Outlet submerged Circular ' 144 'PIPE FULL') 142 145 else: 143 146 #alpha = acos(1 - outlet_culvert_depth/diameter) # Where did this Come From ????/ … … 149 152 case = 'INLET CTRL Culvert is open channel flow we will for now assume critical depth' 150 153 if local_debug =='true': 151 print 'INLET CTRL Culvert is open channel flow we will for now assume critical depth' 152 print 'Q Outlet Depth and ALPHA = ',Q,' ',outlet_culvert_depth,' ',alpha 154 log.critical('INLET CTRL Culvert is open channel flow ' 155 'we will for now assume critical depth') 156 log.critical('Q Outlet Depth and ALPHA = %s, %s, %s' 157 % (str(Q), str(outlet_culvert_depth), 158 str(alpha))) 153 159 if delta_total_energy < inlet_specific_energy: # OUTLET CONTROL !!!! 154 160 # Calculate flows for outlet control … … 162 168 case = 'Outlet submerged' 163 169 if local_debug =='true': 164 print 'Outlet submerged'170 log.critical('Outlet submerged') 165 171 else: # Culvert running PART FULL for PART OF ITS LENGTH Here really should use the Culvert Slope to calculate Actual Culvert Depth & Velocity 166 172 # IF outlet_depth < diameter … … 178 184 case = 'Outlet unsubmerged PIPE FULL' 179 185 if local_debug =='true': 180 print 'Outlet unsubmerged PIPE FULL'186 log.critical('Outlet unsubmerged PIPE FULL') 181 187 else: 182 188 alpha = acos(1-2*outlet_culvert_depth/diameter)*2 … … 185 191 perimeter = alpha*diameter/2.0 186 192 case = 'Outlet is open channel flow we will for now assume critical depth' 187 if local_debug =='true': 188 print 'Q Outlet Depth and ALPHA = ',Q,' ',outlet_culvert_depth,' ',alpha 189 print 'Outlet is open channel flow we will for now assume critical depth' 190 if local_debug =='true': 191 print 'FLOW AREA = ',flow_area 192 print 'PERIMETER = ',perimeter 193 print 'Q Interim = ',Q 193 if local_debug == 'true': 194 log.critical('Q Outlet Depth and ALPHA = %s, %s, %s' 195 % (str(Q), str(outlet_culvert_depth), 196 str(alpha))) 197 log.critical('Outlet is open channel flow we ' 198 'will for now assume critical depth') 199 if local_debug == 'true': 200 log.critical('FLOW AREA = %s' % str(flow_area)) 201 log.critical('PERIMETER = %s' % str(perimeter)) 202 log.critical('Q Interim = %s' % str(Q)) 194 203 hyd_rad = flow_area/perimeter 195 204 … … 200 209 # Outlet control velocity using tail water 201 210 if local_debug =='true': 202 print 'GOT IT ALL CALCULATING Velocity'203 print 'HydRad = ',hyd_rad211 log.critical('GOT IT ALL CALCULATING Velocity') 212 log.critical('HydRad = %s' % str(hyd_rad)) 204 213 culvert_velocity = sqrt(delta_total_energy/((sum_loss/2/g)+(manning**2*culvert_length)/hyd_rad**1.33333)) 205 214 Q_outlet_tailwater = flow_area * culvert_velocity 206 215 if local_debug =='true': 207 print 'VELOCITY = ',culvert_velocity208 print 'Outlet Ctrl Q = ',Q_outlet_tailwater216 log.critical('VELOCITY = %s' % str(culvert_velocity)) 217 log.critical('Outlet Ctrl Q = %s' % str(Q_outlet_tailwater)) 209 218 if log_filename is not None: 210 219 s = 'Q_outlet_tailwater = %.6f' %Q_outlet_tailwater … … 212 221 Q = min(Q, Q_outlet_tailwater) 213 222 if local_debug =='true': 214 print ('%s,%.3f,%.3f' %('dcrit 1 , dcit2 =',dcrit1,dcrit2)) 215 print ('%s,%.3f,%.3f,%.3f' %('Q and Velocity and Depth=',Q,culvert_velocity,outlet_culvert_depth)) 223 log.critical('%s,%.3f,%.3f' 224 % ('dcrit 1 , dcit2 =',dcrit1,dcrit2)) 225 log.critical('%s,%.3f,%.3f,%.3f' 226 % ('Q and Velocity and Depth=', Q, 227 culvert_velocity, outlet_culvert_depth)) 216 228 217 229 else: … … 221 233 222 234 #else.... 223 if culvert_type == 'box':224 if local_debug == 'true':225 print 'BOX CULVERT'235 if culvert_type == 'box': 236 if local_debug == 'true': 237 log.critical('BOX CULVERT') 226 238 # Box culvert (rectangle or square) ======================================================================================================================== 227 239 … … 318 330 culv_froude=sqrt(Q**2*flow_width/(g*flow_area**3)) 319 331 if local_debug =='true': 320 print 'FLOW AREA = ',flow_area321 print 'PERIMETER = ',perimeter322 print 'Q final = ',Q323 print 'FROUDE = ',culv_froude332 log.critical('FLOW AREA = %s' % str(flow_area)) 333 log.critical('PERIMETER = %s' % str(perimeter)) 334 log.critical('Q final = %s' % str(Q)) 335 log.critical('FROUDE = %s' % str(culv_froude)) 324 336 if log_filename is not None: 325 337 s = 'Froude in Culvert = %f' % culv_froude -
anuga_core/source/anuga/culvert_flows/test_culvert_class.py
r7276 r7317 135 135 if delta_w < min_delta_w: min_delta_w = delta_w 136 136 137 #print domain.timestepping_statistics()138 137 pass 139 138 … … 236 235 ref_volume = domain.get_quantity('stage').get_integral() 237 236 for t in domain.evolve(yieldstep = 1, finaltime = 25): 238 #print domain.timestepping_statistics()239 237 new_volume = domain.get_quantity('stage').get_integral() 240 238 … … 339 337 ref_volume = domain.get_quantity('stage').get_integral() 340 338 for t in domain.evolve(yieldstep = 0.1, finaltime = 25): 341 #print domain.timestepping_statistics()342 339 new_volume = domain.get_quantity('stage').get_integral() 343 340 344 msg = 'Total volume has changed: Is %.2f m^3 should have been %.2f m^3'\ 345 % (new_volume, ref_volume) 346 if not num.allclose(new_volume, ref_volume): 347 print msg 341 msg = ('Total volume has changed: Is %.2f m^3 should have been %.2f m^3' 342 % (new_volume, ref_volume)) 348 343 assert num.allclose(new_volume, ref_volume), msg 349 344 … … 360 355 ref_volume = domain.get_quantity('stage').get_integral() 361 356 for t in domain.evolve(yieldstep = 0.1, finaltime = 25): 362 #print domain.timestepping_statistics()363 357 new_volume = domain.get_quantity('stage').get_integral() 364 358 365 #print new_volume, ref_volume, new_volume-ref_volume366 359 msg = 'Total volume has changed: Is %.2f m^3 should have been %.2f m^3'\ 367 360 % (new_volume, ref_volume) … … 471 464 for t in domain.evolve(yieldstep = 1, finaltime = 25): 472 465 473 #print domain.timestepping_statistics()474 466 new_volume = domain.get_quantity('stage').get_integral() 475 467 … … 558 550 culvert(domain) 559 551 560 #print 'Inlet flow', culvert.inlet.rate561 #print 'Outlet flow', culvert.outlet.rate562 563 564 552 565 553 #------------------------------------------------------------- 554 566 555 if __name__ == "__main__": 567 556 suite = unittest.makeSuite(Test_Culvert, 'test')
Note: See TracChangeset
for help on using the changeset viewer.
