Changeset 8930
- Timestamp:
- Jun 26, 2013, 9:29:55 AM (12 years ago)
- Location:
- trunk/anuga_core/source
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py
r8871 r8930 239 239 Q.set_values(other) 240 240 241 result = Quantity(self.domain)242 241 243 242 # The maximum of vertex_values, edge_values and centroid_values … … 245 244 # set_values (which calls interpolate). Otherwise 246 245 # edge and centroid values wouldn't be max from q1 and q2 247 result.vertex_values= num.maximum(self.vertex_values, Q.vertex_values)248 result.edge_values= num.maximum(self.edge_values, Q.edge_values)249 result.centroid_values= num.maximum(self.centroid_values, Q.centroid_values)250 251 return result246 self.vertex_values[:] = num.maximum(self.vertex_values, Q.vertex_values) 247 self.edge_values[:] = num.maximum(self.edge_values, Q.edge_values) 248 self.centroid_values[:] = num.maximum(self.centroid_values, Q.centroid_values) 249 250 252 251 253 252 … … 266 265 Q.set_values(other) 267 266 268 result = Quantity(self.domain)269 267 270 268 # The minimum of vertex_values, edge_values and centroid_values 271 269 # are calculated and assigned directly without using 272 270 # set_values (which calls interpolate). Otherwise 273 result.vertex_values= num.minimum(self.vertex_values, Q.vertex_values)274 result.edge_values= num.minimum(self.edge_values, Q.edge_values)275 result.centroid_values= num.minimum(self.centroid_values, Q.centroid_values)276 277 return result 271 self.vertex_values[:]= num.minimum(self.vertex_values, Q.vertex_values) 272 self.edge_values[:] = num.minimum(self.edge_values, Q.edge_values) 273 self.centroid_values[:] = num.minimum(self.centroid_values, Q.centroid_values) 274 275 278 276 279 277 -
trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py
r8486 r8930 2438 2438 2439 2439 2440 def test_maximum(self): 2441 quantity = Quantity(self.mesh4) 2442 2443 # get referece to data arrays 2444 centroid_values = quantity.centroid_values 2445 vertex_values = quantity.vertex_values 2446 edge_values = quantity.edge_values 2447 2448 quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]], 2449 location = 'vertices') 2450 assert num.allclose(quantity.vertex_values, 2451 [[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]]) 2452 2453 assert id(vertex_values) == id(quantity.vertex_values) 2454 assert num.allclose(quantity.centroid_values, [2., 5., 3., 0.]) #Centroid 2455 assert num.allclose(quantity.edge_values, [[2.5, 2.0, 1.5], 2456 [5., 5., 5.], 2457 [4.5, 4.5, 0.], 2458 [3.0, -1.5, -1.5]]) 2459 2460 2461 2462 other_quantity = Quantity(self.mesh4) 2463 2464 other_quantity.set_values([[0,0,0], [1,1,6], [10,10,10], [0,0,4]], 2465 location = 'vertices') 2466 2467 2468 #=============================== 2469 quantity.maximum(other_quantity) 2470 #=============================== 2471 2472 exact_vertex_values = num.array([[ 1., 2., 3.], 2473 [ 5., 5., 6.], 2474 [ 10., 10., 10.], 2475 [ 0., 3., 4.]]) 2476 exact_centroid_values = num.array([ 2., 5., 10., 1.33333333]) 2477 exact_edge_values = num.array([[ 2.5, 2., 1.5], 2478 [ 5., 5., 5., ], 2479 [ 10., 10., 10. ], 2480 [ 3., 2., 0. ]]) 2481 2482 2483 2484 assert num.allclose(quantity.vertex_values, 2485 exact_vertex_values) 2486 assert num.allclose(quantity.centroid_values, exact_centroid_values) #Centroid 2487 assert num.allclose(quantity.edge_values, exact_edge_values) 2488 2489 def test_minimum(self): 2490 quantity = Quantity(self.mesh4) 2491 2492 # get referece to data arrays 2493 centroid_values = quantity.centroid_values 2494 vertex_values = quantity.vertex_values 2495 edge_values = quantity.edge_values 2496 2497 quantity.set_values([[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]], 2498 location = 'vertices') 2499 assert num.allclose(quantity.vertex_values, 2500 [[1,2,3], [5,5,5], [0,0,9], [-6, 3, 3]]) 2501 2502 assert id(vertex_values) == id(quantity.vertex_values) 2503 assert num.allclose(quantity.centroid_values, [2., 5., 3., 0.]) #Centroid 2504 assert num.allclose(quantity.edge_values, [[2.5, 2.0, 1.5], 2505 [5., 5., 5.], 2506 [4.5, 4.5, 0.], 2507 [3.0, -1.5, -1.5]]) 2508 2509 2510 2511 other_quantity = Quantity(self.mesh4) 2512 2513 other_quantity.set_values([[0,0,0], [1,1,6], [10,10,10], [0,0,4]], 2514 location = 'vertices') 2515 2516 2517 #=============================== 2518 quantity.minimum(other_quantity) 2519 #=============================== 2520 2521 exact_vertex_values = num.array([[ 0., 0., 0.], 2522 [ 1., 1., 5.], 2523 [ 0., 0., 9.], 2524 [ -6., 0., 3.]]) 2525 exact_centroid_values = num.array([ 0., 2.66666667, 3., 0.]) 2526 exact_edge_values = num.array([[ 0., 0., 0.], 2527 [ 3.5, 3.5, 1., ], 2528 [ 4.5, 4.5, 0. ], 2529 [ 2., -1.5, -1.5 ]]) 2530 2531 2532 assert num.allclose(quantity.vertex_values, 2533 exact_vertex_values) 2534 assert num.allclose(quantity.centroid_values, exact_centroid_values) #Centroid 2535 assert num.allclose(quantity.edge_values, exact_edge_values) 2440 2536 2441 2537 -
trunk/anuga_core/source/anuga/operators/set_stage_operator.py
r8928 r8930 1 1 """ 2 Set value operators2 Set stage operator 3 3 4 4 Constraints: See GPL license in the user guide … … 17 17 from anuga.geometry.polygon import inside_polygon 18 18 19 from anuga.operators. base_operator import Operator20 from anuga.fit_interpolate.interpolate import Modeltime_too_early, \ 21 Modeltime_too_late 19 from anuga.operators.set_quantity_operator import Set_quantity_operator 20 21 22 22 from anuga import indent 23 23 24 24 25 25 26 class Set_stage_operator( Operator):26 class Set_stage_operator(Set_quantity_operator): 27 27 """ 28 28 Set the stage in a region … … 34 34 """ 35 35 36 get_stage = Set_quantity_operator.get_value 37 set_stage = Set_quantity_operator.set_value 38 36 39 def __init__(self, 37 40 domain, 38 41 stage=None, 39 42 indices=None, 43 polygon=None, 44 center=None, 45 radius=None, 40 46 description = None, 41 47 label = None, … … 44 50 45 51 46 Operator.__init__(self, domain, description, label, logging, verbose) 47 48 #------------------------------------------ 49 # Local variables 50 #------------------------------------------ 51 self.stage = stage 52 self.indices = indices 53 54 55 def __call__(self): 56 """ 57 Apply rate to those triangles defined in indices 58 59 indices == [], then don't apply anywhere 60 indices == None, then apply everywhere 61 otherwise apply for the specific indices 62 """ 63 64 65 if self.indices is []: 66 return 67 68 stage = self.get_stage() 69 70 if stage is None: 71 return 72 73 if self.verbose is True: 74 log.critical('Stage of %s at time = %.2f = %f' 75 % (self.quantity_name, domain.get_time(), stage)) 76 77 if self.indices is None: 78 self.stage_c[:] = stage 79 else: 80 self.stage_c[self.indices] = stage 81 82 83 def get_stage(self, t=None): 84 """Get value of stage at time t. 85 If t not specified, return stage at current domain time 86 """ 87 88 if t is None: 89 t = self.domain.get_time() 90 91 if callable(self.stage): 92 try: 93 stage = self.stage(t) 94 except Modeltime_too_early, e: 95 raise Modeltime_too_early(e) 96 except Modeltime_too_late, e: 97 msg = '%s: ANUGA is trying to run longer than specified data.\n' %str(e) 98 msg += 'You can specify keyword argument default_rate in the ' 99 msg += 'stage function to tell it what to do in the absence of time data.' 100 raise Modeltime_too_late(msg) 101 else: 102 stage = self.stage 103 104 105 if stage is None: 106 msg = ('Attribute stage must be specified in '+self.__name__+ 107 ' before attempting to call it') 108 raise Exception(msg) 109 110 return stage 52 Set_quantity_operator.__init__(self, domain, 53 quantity = 'stage', 54 value = stage, 55 indices = indices, 56 polygon = polygon, 57 center = center, 58 radius = radius, 59 description = description, 60 label = label, 61 logging = logging, 62 verbose = verbose) 111 63 112 64 … … 152 104 assert radius is not None 153 105 154 155 # Determine indices in update region156 N = domain.get_number_of_triangles()157 points = domain.get_centroid_coordinates(absolute=True)158 159 160 indices = []161 162 c = center163 r = radius164 165 self.center = center166 self.radius = radius167 168 intersect = False169 for k in range(N):170 x, y = points[k,:] # Centroid171 172 if ((x-c[0])**2+(y-c[1])**2) < r**2:173 intersect = True174 indices.append(k)175 176 177 msg = 'No centroids intersect circle center'+str(center)+' radius '+str(radius)178 assert intersect, msg179 180 181 182 183 # It is possible that circle doesn't intersect with mesh (as can happen184 # for parallel runs185 186 187 106 Set_stage_operator.__init__(self, 188 107 domain, 189 108 stage=stage, 190 indices=indices, 109 center=center, 110 radius=radius, 191 111 verbose=verbose) 192 193 194 195 196 112 197 113 #=============================================================================== … … 213 129 214 130 215 # Determine indices in update region216 N = domain.get_number_of_triangles()217 points = domain.get_centroid_coordinates(absolute=True)218 219 220 indices = inside_polygon(points, polygon)221 self.polygon = polygon222 223 # It is possible that circle doesn't intersect with mesh (as can happen224 # for parallel runs225 226 131 227 132 Set_stage_operator.__init__(self, 228 133 domain, 229 134 stage=stage, 230 indices=indices,135 polygon=polygon, 231 136 verbose=verbose) 232 137 -
trunk/anuga_core/source/anuga/operators/test_set_elevation_operator.py
r8928 r8930 869 869 870 870 if __name__ == "__main__": 871 suite = unittest.makeSuite(Test_set_elevation_operator s, 'test')871 suite = unittest.makeSuite(Test_set_elevation_operator, 'test') 872 872 runner = unittest.TextTestRunner(verbosity=1) 873 873 runner.run(suite) -
trunk/anuga_core/source/anuga/operators/test_set_stage_operator.py
r8928 r8930 13 13 from anuga.config import time_format 14 14 15 from set_stage_operator simport *15 from set_stage_operator import * 16 16 17 17 import numpy as num … … 253 253 254 254 polygon = [(0.5,0.5), (1.5,0.5), (1.5,1.5), (0.5,1.5)] 255 #operator = Polygonal_set_stage_operator(domain, stage=stage, polygon=polygon) 255 256 operator = Polygonal_set_stage_operator(domain, stage=stage, polygon=polygon) 256 257 -
trunk/anuga_core/source/anuga/utilities/function_utils.py
r8920 r8930 20 20 # or csllable 21 21 #------------------------------------------ 22 msg = "Input argument must be a scalar, or a function "22 msg = "Input argument must be a scalar, or a function or None" 23 23 24 if function is None: 25 return None 26 24 27 assert (isinstance(function, (int, float)) or 25 28 callable(function)), msg -
trunk/anuga_core/source/anuga_validation_tests/case_studies/patong/README.txt
r8908 r8930 28 28 -- Debugged it, and got rid of 'setup_model.py' to simplify the code 29 29 30 -- Go dit running in parallel30 -- Got it running in parallel
Note: See TracChangeset
for help on using the changeset viewer.