Changeset 8994
- Timestamp:
- Sep 27, 2013, 10:36:23 AM (11 years ago)
- Location:
- trunk/anuga_core/source
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/anuga_core/source/anuga/shallow_water/shallow_water_ext.c
r8813 r8994 2757 2757 double dk, dv0, dv1, dv2; 2758 2758 2759 double *xmom_centroid_store, *ymom_centroid_store, *stage_centroid_store; 2760 2761 2762 // Associate memory location of Domain varialbe with local aliases 2759 double *xmom_centroid_store; 2760 double *ymom_centroid_store; 2761 //double *stage_centroid_store; 2762 2763 2764 // Associate memory location of Domain varibles with local aliases 2763 2765 number_of_elements = D->number_of_elements; 2764 2766 epsilon = D->epsilon; … … 2836 2838 xmom_centroid_store = malloc(number_of_elements*sizeof(double)); 2837 2839 ymom_centroid_store = malloc(number_of_elements*sizeof(double)); 2838 stage_centroid_store = malloc(number_of_elements*sizeof(double));2840 // stage_centroid_store = malloc(number_of_elements*sizeof(double)); 2839 2841 2840 2842 if (extrapolate_velocity_second_order == 1) { … … 3318 3320 free(xmom_centroid_store); 3319 3321 free(ymom_centroid_store); 3320 free(stage_centroid_store);3322 //free(stage_centroid_store); 3321 3323 3322 3324 -
trunk/anuga_core/source/anuga/structures/boyd_box_operator.py
r8985 r8994 155 155 print 'delta_total_energy ',self.delta_total_energy 156 156 print 'outlet_enquiry_depth ',self.outflow.get_enquiry_depth() 157 print 'inflow_enquiry_depth ',self.inflow.get_enquiry_depth() 158 print 'outlet_enquiry_speed ',self.outflow.get_enquiry_speed() 159 print 'inflow_enquiry_speed ',self.inflow.get_enquiry_speed() 157 160 print 'sum_loss ',self.sum_loss 158 161 print 'manning ',self.manning -
trunk/anuga_core/source/anuga/structures/inlet.py
r8873 r8994 45 45 else: # poly is a polygon 46 46 47 self.triangle_indices = inside_polygon(domain_centroids, self.poly) 47 tris_0 = line_intersect(vertex_coordinates, [self.poly[0],self.poly[1]]) 48 tris_1 = inside_polygon(domain_centroids, self.poly) 49 #print 40*"=" 50 #print tris_0 51 #print tris_1 52 self.triangle_indices = num.union1d(tris_0, tris_1) 53 #print self.triangle_indices 48 54 49 55 if len(self.triangle_indices) == 0: 50 56 msg = 'Inlet poly=%s ' % (self.poly) 51 msg += 'No triangle s intersecting line'57 msg += 'No triangle centroids intersecting poly ' 52 58 raise Exception, msg 53 59 -
trunk/anuga_core/source/anuga/structures/inlet_enquiry.py
r8876 r8994 54 54 if self.enquiry_index in self.triangle_indices: 55 55 msg = 'Enquiry point %s' % (self.enquiry_pt) 56 msg += ' is in an inlet triangle'56 msg += ' is in an inlet triangle' 57 57 import warnings 58 58 warnings.warn(msg) -
trunk/anuga_core/source/anuga/structures/run_culvert.py
r8361 r8994 42 42 domain = anuga.Domain(points, vertices, boundary) 43 43 domain.set_starttime(10) 44 domain.set_name( 'Test_culvert') # Output name44 domain.set_name() # Output name 45 45 domain.set_default_order(2) 46 46 #domain.set_beta(1.5) -
trunk/anuga_core/source/anuga/structures/structure_operator.py
r8950 r8994 101 101 102 102 # Slots for recording current statistics 103 self.accumulated_flow = 0.0 103 104 self.discharge = 0.0 104 105 self.velocity = 0.0 106 self.outlet_depth = 0.0 105 107 self.delta_total_energy = 0.0 106 108 self.driving_energy = 0.0 … … 116 118 self.inlets = [] 117 119 line0 = self.exchange_lines[0] #self.inlet_lines[0] 120 if self.apron is None: 121 poly0 = line0 122 else: 123 offset = -self.apron*self.outward_vector_0 124 #print line0 125 #print offset 126 poly0 = num.array([ line0[0], line0[1], line0[1]+offset, line0[0]+offset]) 127 #print poly0 118 128 if self.invert_elevations is None: 119 129 invert_elevation0 = None … … 126 136 self.inlets.append(inlet_enquiry.Inlet_enquiry( 127 137 self.domain, 128 line0,138 poly0, 129 139 enquiry_point0, 130 140 invert_elevation = invert_elevation0, … … 132 142 verbose = self.verbose)) 133 143 134 144 tris_0 = self.inlets[0].triangle_indices 145 #print tris_0 146 #print self.domain.centroid_coordinates[tris_0] 135 147 136 148 line1 = self.exchange_lines[1] 149 if self.apron is None: 150 poly1 = line1 151 else: 152 offset = -self.apron*self.outward_vector_1 153 #print line1 154 #print offset 155 poly1 = num.array([ line1[0], line1[1], line1[1]+offset, line1[0]+offset]) 156 #print poly1 157 137 158 if self.invert_elevations is None: 138 159 invert_elevation1 = None … … 144 165 self.inlets.append(inlet_enquiry.Inlet_enquiry( 145 166 self.domain, 146 line1,167 poly1, 147 168 enquiry_point1, 148 169 invert_elevation = invert_elevation1, … … 150 171 verbose = self.verbose)) 151 172 173 174 tris_1 = self.inlets[1].triangle_indices 175 #print tris_1 176 #print self.domain.centroid_coordinates[tris_1] 177 152 178 self.set_logging(logging) 153 179 … … 166 192 old_inflow_ymom = self.inflow.get_average_ymom() 167 193 168 169 # Implement the update of flow over a timestep by 170 # using a semi-implict update. This ensures that 171 # the update does not create a negative depth 172 if old_inflow_depth > 0.0 : 173 Q_star = Q/old_inflow_depth 174 else: 175 Q_star = 0.0 176 177 factor = 1.0/(1.0 + Q_star*timestep/self.inflow.get_area()) 178 179 new_inflow_depth = old_inflow_depth*factor 180 new_inflow_xmom = old_inflow_xmom*factor 181 new_inflow_ymom = old_inflow_ymom*factor 182 183 self.inflow.set_depths(new_inflow_depth) 184 185 #inflow.set_xmoms(Q/inflow.get_area()) 186 #inflow.set_ymoms(0.0) 187 188 self.inflow.set_xmoms(new_inflow_xmom) 189 self.inflow.set_ymoms(new_inflow_ymom) 190 191 loss = (old_inflow_depth - new_inflow_depth)*self.inflow.get_area() 192 193 # set outflow 194 if old_inflow_depth > 0.0 : 195 timestep_star = timestep*new_inflow_depth/old_inflow_depth 196 else: 197 timestep_star = 0.0 198 199 200 201 202 outflow_extra_depth = Q*timestep_star/self.outflow.get_area() 203 outflow_direction = - self.outflow.outward_culvert_vector 204 outflow_extra_momentum = outflow_extra_depth*barrel_speed*outflow_direction 205 206 gain = outflow_extra_depth*self.outflow.get_area() 207 208 #print Q, Q*timestep, barrel_speed, outlet_depth, Qstar, factor, timestep_star 209 #print ' ', loss, gain 210 194 semi_implicit = True 195 if semi_implicit: 196 # Implement the update of flow over a timestep by 197 # using a semi-implict update. This ensures that 198 # the update does not create a negative depth 199 if old_inflow_depth > 0.0 : 200 Q_star = Q/old_inflow_depth 201 else: 202 Q_star = 0.0 203 204 factor = 1.0/(1.0 + Q_star*timestep/self.inflow.get_area()) 205 206 new_inflow_depth = old_inflow_depth*factor 207 new_inflow_xmom = old_inflow_xmom*factor 208 new_inflow_ymom = old_inflow_ymom*factor 209 210 self.inflow.set_depths(new_inflow_depth) 211 212 #inflow.set_xmoms(Q/inflow.get_area()) 213 #inflow.set_ymoms(0.0) 214 215 self.inflow.set_xmoms(new_inflow_xmom) 216 self.inflow.set_ymoms(new_inflow_ymom) 217 218 loss = (old_inflow_depth - new_inflow_depth)*self.inflow.get_area() 219 220 # set outflow 221 if old_inflow_depth > 0.0 : 222 timestep_star = timestep*new_inflow_depth/old_inflow_depth 223 else: 224 timestep_star = 0.0 225 226 227 228 229 outflow_extra_depth = Q*timestep_star/self.outflow.get_area() 230 outflow_direction = - self.outflow.outward_culvert_vector 231 outflow_extra_momentum = outflow_extra_depth*barrel_speed*outflow_direction 232 233 gain = outflow_extra_depth*self.outflow.get_area() 234 235 #print gain, loss 236 assert num.allclose(gain-loss, 0.0) 237 238 #print Q, Q*timestep, barrel_speed, outlet_depth, Qstar, factor, timestep_star 239 #print ' ', loss, gain 240 else: 241 242 243 factor = Q*timestep/self.inflow.get_area() 244 245 new_inflow_depth = old_inflow_depth - factor 246 new_inflow_xmom = old_inflow_xmom/old_inflow_depth*new_inflow_depth 247 new_inflow_ymom = old_inflow_ymom/old_inflow_depth*new_inflow_depth 248 249 self.inflow.set_depths(new_inflow_depth) 250 self.inflow.set_xmoms(new_inflow_xmom) 251 self.inflow.set_ymoms(new_inflow_ymom) 252 253 254 loss = (old_inflow_depth - new_inflow_depth)*self.inflow.get_area() 255 256 257 258 outflow_extra_depth = Q*timestep/self.outflow.get_area() 259 outflow_direction = - self.outflow.outward_culvert_vector 260 outflow_extra_momentum = outflow_extra_depth*barrel_speed*outflow_direction 261 262 gain = outflow_extra_depth*self.outflow.get_area() 263 264 assert num.allclose(gain-loss, 0.0) 265 266 #print Q, Q*timestep, barrel_speed, outlet_depth, Qstar, factor, timestep_star 267 #print ' ', loss, gain 268 269 211 270 # Stats 212 self.discharge = Q#outflow_extra_depth*self.outflow.get_area()/timestep 213 self.velocity = barrel_speed#self.discharge/outlet_depth/self.width 271 272 self.accumulated_flow += gain 273 self.discharge = outflow_extra_depth*self.outflow.get_area()/timestep #Q 274 self.velocity = barrel_speed # self.discharge/outlet_depth/self.width 275 self.outlet_depth = outlet_depth 214 276 215 277 new_outflow_depth = self.outflow.get_average_depth() + outflow_extra_depth … … 217 279 if self.use_momentum_jet : 218 280 # FIXME (SR) Review momentum to account for possible hydraulic jumps at outlet 219 #new_outflow_xmom = outflow.get_average_xmom() + outflow_extra_momentum[0]220 #new_outflow_ymom = outflow.get_average_ymom() + outflow_extra_momentum[1]281 #new_outflow_xmom = self.outflow.get_average_xmom() + outflow_extra_momentum[0] 282 #new_outflow_ymom = self.outflow.get_average_ymom() + outflow_extra_momentum[1] 221 283 new_outflow_xmom = barrel_speed*new_outflow_depth*outflow_direction[0] 222 284 new_outflow_ymom = barrel_speed*new_outflow_depth*outflow_direction[1] … … 396 458 message += '--------------------------\n' 397 459 message += 'Type: %s\n' % self.structure_type 460 461 message += 'inlets[0]_enquiry_depth [m]: %.2f\n' %self.inlets[0].get_enquiry_depth() 462 message += 'inlets[0]_enquiry_speed [m/s]: %.2f\n' %self.inlets[0].get_enquiry_speed() 463 message += 'inlets[0]_enquiry_stage [m]: %.2f\n' %self.inlets[0].get_enquiry_stage() 464 message += 'inlets[0]_enquiry_elevation [m]: %.2f\n' %self.inlets[0].get_enquiry_elevation() 465 message += 'inlets[0]_average_depth [m]: %.2f\n' %self.inlets[0].get_average_depth() 466 message += 'inlets[0]_average_speed [m/s]: %.2f\n' %self.inlets[0].get_average_speed() 467 message += 'inlets[0]_average_stage [m]: %.2f\n' %self.inlets[0].get_average_stage() 468 message += 'inlets[0]_average_elevation [m]: %.2f\n' %self.inlets[0].get_average_elevation() 469 470 message += '\n' 471 472 message += 'inlets[1]_enquiry_depth [m]: %.2f\n' %self.inlets[1].get_enquiry_depth() 473 message += 'inlets[1]_enquiry_speed [m/s]: %.2f\n' %self.inlets[1].get_enquiry_speed() 474 message += 'inlets[1]_enquiry_stage [m]: %.2f\n' %self.inlets[1].get_enquiry_stage() 475 message += 'inlets[1]_enquiry_elevation [m]: %.2f\n' %self.inlets[1].get_enquiry_elevation() 476 477 message += 'inlets[1]_average_depth [m]: %.2f\n' %self.inlets[1].get_average_depth() 478 message += 'inlets[1]_average_speed [m/s]: %.2f\n' %self.inlets[1].get_average_speed() 479 message += 'inlets[1]_average_stage [m]: %.2f\n' %self.inlets[1].get_average_stage() 480 message += 'inlets[1]_average_elevation [m]: %.2f\n' %self.inlets[1].get_average_elevation() 481 482 398 483 message += 'Discharge [m^3/s]: %.2f\n' % self.discharge 399 484 message += 'Velocity [m/s]: %.2f\n' % self.velocity 485 message += 'Outlet Depth [m]: %.2f\n' % self.outlet_depth 486 message += 'Accumulated Flow [m^3]: %.2f\n' % self.accumulated_flow 400 487 message += 'Inlet Driving Energy %.2f\n' % self.driving_energy 401 488 message += 'Delta Total Energy %.2f\n' % self.delta_total_energy 402 489 message += 'Control at this instant: %s\n' % self.case 490 491 492 403 493 404 494 print message … … 423 513 message += '%.5f, ' % self.discharge 424 514 message += '%.5f, ' % self.velocity 515 message += '%.5f, ' % self.accumulated_flow 425 516 message += '%.5f, ' % self.driving_energy 426 517 message += '%.5f' % self.delta_total_energy 518 519 427 520 428 521 return message -
trunk/anuga_core/source/anuga_parallel/run_parallel_gate_operator.py
r8877 r8994 2 2 import sys 3 3 4 from anuga.utilities.system_tools import get_pathname_from_package5 from anuga.geometry.polygon_function import Polygon_function4 #from anuga.utilities.system_tools import get_pathname_from_package 5 #from anuga.geometry.polygon_function import Polygon_function 6 6 7 from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross8 from anuga.abstract_2d_finite_volumes.quantity import Quantity7 #from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross 8 #from anuga.abstract_2d_finite_volumes.quantity import Quantity 9 9 10 10 import anuga … … 16 16 import numpy as num 17 17 18 from anuga _parallelimport distribute, myid, numprocs, finalize, barrier18 from anuga import distribute, myid, numprocs, finalize, barrier 19 19 20 from parallel_operator_factoryimport Inlet_operator, Boyd_box_operator20 from anuga import Inlet_operator, Boyd_box_operator 21 21 22 22 … … 56 56 return z 57 57 58 59 """test_that_culvert_runs_rating 60 61 This test exercises the culvert and checks values outside rating curve 62 are dealt with 63 """ 64 65 path = get_pathname_from_package('anuga.culvert_flows') 58 66 59 67 60 length = 40. … … 73 66 if myid == 0: 74 67 75 points, vertices, boundary = rectangular_cross(int(length/dx),76 int(width/dy),77 len1=length,78 len2=width)68 points, vertices, boundary = anuga.rectangular_cross(int(length/dx), 69 int(width/dy), 70 len1=length, 71 len2=width) 79 72 domain = anuga.Domain(points, vertices, boundary) 80 domain.set_name( 'run_gate_operator') # Output name73 domain.set_name() # Output name 81 74 domain.set_flow_algorithm('2_0') 82 75
Note: See TracChangeset
for help on using the changeset viewer.