Changeset 2632
- Timestamp:
- Mar 29, 2006, 4:33:30 PM (19 years ago)
- Location:
- inundation
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
inundation/examples/island.py
r2621 r2632 7 7 8 8 9 #------------------------------------------------------------------------------ -9 #------------------------------------------------------------------------------ 10 10 # Import necessary modules 11 #------------------------------------------------------------------------------ -11 #------------------------------------------------------------------------------ 12 12 13 13 # Standard modules … … 21 21 22 22 23 #------------------------------------------------------------------------------ -23 #------------------------------------------------------------------------------ 24 24 # Setup computational domain 25 #------------------------------------------------------------------------------ -25 #------------------------------------------------------------------------------ 26 26 27 27 #Create basic mesh … … 37 37 maximum_triangle_area = 25, 38 38 filename = 'island.msh', 39 interior_regions=[ ([[5 5,25], [75,25],40 [7 5,75], [55,75]], 3)]39 interior_regions=[ ([[50,25], [70,25], 40 [70,75], [50,75]], 3)] 41 41 ) 42 42 … … 63 63 64 64 65 #------------------------------------------------------------------------------ -65 #------------------------------------------------------------------------------ 66 66 # Setup initial conditions 67 #------------------------------------------------------------------------------ -67 #------------------------------------------------------------------------------ 68 68 69 69 def island(x, y): … … 73 73 return z 74 74 75 domain.set_quantity('friction', 0.0)75 domain.set_quantity('friction', 2.0) 76 76 domain.set_quantity('elevation', island) 77 77 domain.set_quantity('stage', 1) 78 78 79 79 80 #------------------------------------------------------------------------------ -80 #------------------------------------------------------------------------------ 81 81 # Setup boundary conditions (all reflective) 82 #------------------------------------------------------------------------------ -82 #------------------------------------------------------------------------------ 83 83 84 84 Br = Reflective_boundary(domain) … … 88 88 89 89 90 #------------------------------------------------------------------------------ -90 #------------------------------------------------------------------------------ 91 91 # Evolve system through time 92 #------------------------------------------------------------------------------ -92 #------------------------------------------------------------------------------ 93 93 94 94 import time 95 for t in domain.evolve(yieldstep = 0.5, finaltime = 100):95 for t in domain.evolve(yieldstep = 1, finaltime = 100): 96 96 domain.write_time() 97 #for t in domain.evolve(yieldstep = 0.05, finaltime = 5): 98 # domain.write_time() 99 # #print ' Volume: ', domain.get_quantity('stage').get_integral() 97 100 98 101 print 'Done' -
inundation/examples/island_timeseries.py
r2620 r2632 10 10 11 11 swwfile = 'island.sww' 12 #gauges = [[0,0.5], [0.2,0.5], [0.4,0.5], [0.6,0.5], [0.8,0.5], [1,0.5]] 13 #gauges = [[20, 50], [30, 50], [40, 50], [50, 50]] 14 gauges = [[42, 48], [58, 48]] 12 #gauges = [[56, 48], [58, 48], [60, 48], [62, 48], [64, 48]] 13 gauges = [[60, 48], [62, 48], [64, 48]] 15 14 15 16 store = False 16 17 17 18 #Read model output … … 103 104 hold(False) 104 105 105 if elevations[0] < -10:106 if elevations[0] < 0: 106 107 plot(model_time, stages, '-b') 107 108 else: 108 109 plot(model_time, stages, '-b', 109 110 model_time, elevations, '-k') 110 axis([0, 100, z-0.01, z+0.005])111 #axis([0, 100, z-0.01, z+0.005]) 111 112 112 113 #name = 'Gauge_%d: (%.1f, %.1f)' %(k, g[0], g[1]) … … 120 121 shadow=True, 121 122 loc='upper right') 122 #savefig('Gauge_%d_stage' %k) 123 savefig('Gauge_%s_stage' %myloc)123 124 if store is True: savefig('Gauge_%s_stage' %myloc) 124 125 125 126 raw_input('Next') … … 136 137 ylabel('sqrt( uh^2 + vh^2 ) [m^2/s]') 137 138 #savefig('Gauge_%d_momentum' %k) 138 savefig('Gauge_%s_momentum' %myloc)139 if store is True: savefig('Gauge_%s_momentum' %myloc) 139 140 140 141 raw_input('Next') … … 160 161 ylabel(' atan(vh/uh) [degrees from North]') 161 162 #savefig('Gauge_%d_bearing' %k) 162 savefig('Gauge_%s_bearing' %myloc)163 if store is True: savefig('Gauge_%s_bearing' %myloc) 163 164 164 165 raw_input('Next') … … 176 177 ylabel('sqrt( uh^2 + vh^2 ) / depth [m/s]') 177 178 #savefig('Gauge_%d_speed' %k) 178 savefig('Gauge_%s_speed' %myloc)179 if store is True: savefig('Gauge_%s_speed' %myloc) 179 180 180 181 raw_input('Next') -
inundation/pyvolution/shallow_water_ext.c
r2621 r2632 246 246 if (h >= eps) { 247 247 S = -g * eta[k]*eta[k] * sqrt((uh[k]*uh[k] + vh[k]*vh[k])); 248 //S /= pow(h, 7.0/3); //Expensive (on Ole's home computer)249 S /= exp(7.0/3.0*log(h)); //seems to save about 15% over manning_friction248 S /= pow(h, 7.0/3); //Expensive (on Ole's home computer) 249 //S /= exp(7.0/3.0*log(h)); //seems to save about 15% over manning_friction 250 250 //S /= h*h*(1 + h/3.0 - h*h/9.0); //FIXME: Could use a Taylor expansion 251 251 … … 372 372 373 373 if (hc < minimum_allowed_height) { 374 375 //Old code: Set momentum to zero and ensure h is non negative 376 /* 377 xmomc[k] = 0.0; 378 ymomc[k] = 0.0; 379 //wc[k] = zc[k]; 380 if (hc <= 0.0) wc[k] = zc[k]; 381 */ 382 383 //New code: Adjust momentum to guarantee speeds are physical 384 // ensure h is non negative 385 374 386 if (hc <= 0.0) { 375 wc[k] = zc[k]; //Contain 'lost mass' error387 wc[k] = zc[k]; 376 388 xmomc[k] = 0.0; 377 389 ymomc[k] = 0.0; … … 383 395 reduced_speed = maximum_allowed_speed * u/fabs(u); 384 396 //printf("Speed (u) has been reduced from %.3f to %.3f\n", 385 // 397 // u, reduced_speed); 386 398 xmomc[k] = reduced_speed * hc; 387 399 } … … 391 403 reduced_speed = maximum_allowed_speed * v/fabs(v); 392 404 //printf("Speed (v) has been reduced from %.3f to %.3f\n", 393 // 405 // v, reduced_speed); 394 406 ymomc[k] = reduced_speed * hc; 395 407 } 396 408 } 397 398 //Old code399 //wc[k] = zc[k]; //Contain 'lost mass' error400 //printf("Speed has been reduced to %.3f\n", 0.0);401 //xmomc[k] = 0.0;402 //ymomc[k] = 0.0;403 409 } 404 410 }
Note: See TracChangeset
for help on using the changeset viewer.