# Contents of /trunk/speckley/test/python/run_specialOnSpeckley.py

Revision 5707 - (show annotations)
Mon Jun 29 03:59:06 2015 UTC (3 years, 6 months ago) by sshaw
File MIME type: text/x-python
File size: 31235 byte(s)
```adding copyright headers to files without copyright info, moved header to top of file in some cases where it wasn't
```
 1 2 ############################################################################## 3 # 4 # Copyright (c) 2003-2015 by The University of Queensland 5 6 # 7 # Primary Business: Queensland, Australia 8 # Licensed under the Open Software License version 3.0 9 10 # 11 # Development until 2012 by Earth Systems Science Computational Center (ESSCC) 12 # Development 2012-2013 by School of Earth Sciences 13 # Development from 2014 by Centre for Geoscience Computing (GeoComp) 14 # 15 ############################################################################## 16 17 from __future__ import division, print_function 18 19 __copyright__="""Copyright (c) 2003-2015 by The University of Queensland 20 http://www.uq.edu.au 21 Primary Business: Queensland, Australia""" 22 __license__="""Licensed under the Open Software License version 3.0 23 24 __url__= 25 26 import os 27 import sys 28 import esys.escriptcore.utestselect as unittest 29 from esys.escriptcore.testing import * 30 from esys.escript import * 31 from esys.speckley import Rectangle, Brick, speckleycpp 32 33 class Test_Speckley_Assemblers(unittest.TestCase): 34 TOLERANCE = 1e-10 35 36 def test_Brick_XY_single(self): 37 ranks = getMPISizeWorld() 38 for expanded in (True, False): 39 for order in range(2,11): 40 dom = Brick(order, 3, 3*ranks, 3, l0=6, l1=6, l2=6, d1=ranks) 41 Y = Data(3, Function(dom), expanded) 42 X = Data(0, (3,), Function(dom), expanded) 43 X[0] = dom.getX()[0] 44 X[1] = dom.getX()[1] 45 X[2] = dom.getX()[2] 46 f = Data(0, Solution(dom), expanded) 47 48 dom.addToRHS(f, [("Y",Y)], 49 dom.createAssembler("DefaultAssembler", [])) 50 dom.addToRHS(f, [("X", X)], 51 dom.createAssembler("DefaultAssembler", [])) 52 #nuke the boundary back to zero since it's not dealt with here 53 for dim in range(3): 54 f *= whereNegative(dom.getX()[dim]-6) 55 res = Lsup(f) 56 self.assertLess(res, self.TOLERANCE, 57 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 58 res, self.TOLERANCE)).format("" if expanded else "un-")) 59 60 def test_Rectangle_XY_single(self): 61 ranks = getMPISizeWorld() 62 for expanded in (True, False): 63 for order in range(2,11): 64 dom = Rectangle(order, 3, 3*ranks, d1=ranks, l0=6, l1=6) 65 Y = Data(2, Function(dom), expanded) 66 X = Data(0, (2,), Function(dom), expanded) 67 X[0] = dom.getX()[0] 68 X[1] = dom.getX()[1] 69 f = Data(0, Solution(dom), expanded) 70 71 dom.addToRHS(f, [("Y",Y)], 72 dom.createAssembler("DefaultAssembler", [])) 73 dom.addToRHS(f, [("X", X)], 74 dom.createAssembler("DefaultAssembler", [])) 75 #nuke the boundary back to zero since it's not dealt with here 76 for dim in range(2): 77 f *= whereNegative(dom.getX()[dim]-6) 78 res = Lsup(f) 79 self.assertLess(res, self.TOLERANCE, 80 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 81 res, self.TOLERANCE)).format("" if expanded else "un-")) 82 83 def test_Brick_XY_system(self): 84 ranks = getMPISizeWorld() 85 for expanded in (True, False): 86 for order in range(2,11): 87 dom = Brick(order, 3, 3*ranks, 3, l0=6, l1=6, l2=6, d1=ranks) 88 Y = Data(1, (3,), Function(dom), expanded) 89 X = Data(0, (3,3), Function(dom), expanded) 90 X[0,0] = dom.getX()[0] 91 X[1,1] = dom.getX()[1] 92 X[2,2] = dom.getX()[2] 93 94 f = Data(0, (3,), Solution(dom), expanded) 95 96 dom.addToRHS(f, [("Y",Y)], 97 dom.createAssembler("DefaultAssembler", [])) 98 dom.addToRHS(f, [("X", X)], 99 dom.createAssembler("DefaultAssembler", [])) 100 #nuke the boundary back to zero since it's not dealt with here 101 for dim in range(3): 102 f *= whereNegative(dom.getX()[dim]-6) 103 res = Lsup(f) 104 self.assertLess(res, self.TOLERANCE, 105 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 106 res, self.TOLERANCE)).format("" if expanded else "un-")) 107 108 def test_Rectangle_XY_system(self): 109 ranks = getMPISizeWorld() 110 for expanded in (True, False): 111 for order in range(2,11): 112 dom = Rectangle(order, 3, 3*ranks, l0=6, l1=6, d1=ranks) 113 Y = Data(1, (2,), Function(dom), expanded) 114 X = Data(0, (2,2), Function(dom), expanded) 115 X[0,0] = dom.getX()[0] 116 X[1,1] = dom.getX()[1] 117 118 f = Data(0, (2,), Solution(dom), expanded) 119 120 dom.addToRHS(f, [("Y",Y)], 121 dom.createAssembler("DefaultAssembler", [])) 122 dom.addToRHS(f, [("X", X)], 123 dom.createAssembler("DefaultAssembler", [])) 124 #nuke the boundary back to zero since it's not dealt with here 125 f *= whereNegative(dom.getX()[0]-6)*whereNegative(dom.getX()[1]-6) 126 res = Lsup(f) 127 self.assertLess(res, self.TOLERANCE, 128 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 129 res, self.TOLERANCE)).format("" if expanded else "un-")) 130 131 def test_Brick_Du_Y_single(self): 132 #solves for u in Du = Y, where D = 1, Y = 2 133 ranks = getMPISizeWorld() 134 for expanded in (True, False): 135 for order in range(2,11): 136 dom = Brick(order, 3, 3*ranks, 3, d1=ranks) 137 D = Data(1, ContinuousFunction(dom), expanded) 138 Y = Data(2, ContinuousFunction(dom), expanded) 139 140 Y = interpolate(Y, Function(dom)) 141 D = interpolate(D, Function(dom)) 142 143 rhs = Data(0, Solution(dom), expanded) 144 lhs = Data(0, Solution(dom), expanded) 145 146 dom.addToRHS(rhs, [("Y",Y)], 147 dom.createAssembler("DefaultAssembler", [])) 148 dom.addToRHS(lhs, [("D",D)], 149 dom.createAssembler("DefaultAssembler", [])) 150 151 res = Lsup((rhs/lhs)-2) 152 self.assertLess(res, self.TOLERANCE, 153 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 154 res, self.TOLERANCE)).format("" if expanded else "un-")) 155 156 def test_Rectangle_Du_Y_single(self): 157 #solves for u in Du = Y, where D = 1, Y = 2 158 ranks = getMPISizeWorld() 159 for expanded in (True, False): 160 for order in range(2,11): 161 dom = Rectangle(order, 3, 3*ranks, d1=ranks) 162 D = Data(1, ContinuousFunction(dom), expanded) 163 Y = Data(2, ContinuousFunction(dom), expanded) 164 165 Y = interpolate(Y, Function(dom)) 166 D = interpolate(D, Function(dom)) 167 168 rhs = Data(0, Solution(dom), expanded) 169 lhs = Data(0, Solution(dom), expanded) 170 171 dom.addToRHS(rhs, [("Y",Y)], 172 dom.createAssembler("DefaultAssembler", [])) 173 dom.addToRHS(lhs, [("D",D)], 174 dom.createAssembler("DefaultAssembler", [])) 175 176 res = Lsup((rhs/lhs)-2) 177 self.assertLess(res, self.TOLERANCE, 178 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 179 res, self.TOLERANCE)).format("" if expanded else "un-")) 180 181 def test_Brick_Du_Y_system(self): 182 #solves for u in Du = Y, where D = [1,2], Y = [2,4] 183 ranks = getMPISizeWorld() 184 for expanded in (True, False): 185 for order in range(2,11): 186 dom = Brick(order, 3, 3*ranks, 3, d1=ranks) 187 D = Data(1, (2,), ContinuousFunction(dom), expanded) 188 D[1] = 2 189 Y = Data(2, (2,), ContinuousFunction(dom), expanded) 190 Y[1] = 4 191 192 Y = interpolate(Y, Function(dom)) 193 D = interpolate(D, Function(dom)) 194 195 rhs = Data(0, (2,), Solution(dom), expanded) 196 lhs = Data(0, (2,), Solution(dom), expanded) 197 198 dom.addToRHS(rhs, [("Y",Y)], 199 dom.createAssembler("DefaultAssembler", [])) 200 dom.addToRHS(lhs, [("D",D)], 201 dom.createAssembler("DefaultAssembler", [])) 202 203 res = Lsup((rhs/lhs)-2) 204 self.assertLess(res, self.TOLERANCE, 205 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 206 res, self.TOLERANCE)).format("" if expanded else "un-")) 207 208 def test_Rectangle_Du_Y_system(self): 209 #solves for u in Du = Y, where D = [1,2], Y = [2,4] 210 ranks = getMPISizeWorld() 211 for expanded in (True, False): 212 for order in range(2,11): 213 dom = Rectangle(order, 3, 3*ranks, d1=ranks) 214 D = Data(1, (2,), ContinuousFunction(dom), expanded) 215 D[1] = 2 216 Y = Data(2, (2,), ContinuousFunction(dom), expanded) 217 Y[1] = 4 218 219 Y = interpolate(Y, Function(dom)) 220 D = interpolate(D, Function(dom)) 221 222 rhs = Data(0, (2,), Solution(dom), expanded) 223 lhs = Data(0, (2,), Solution(dom), expanded) 224 225 dom.addToRHS(rhs, [("Y",Y)], 226 dom.createAssembler("DefaultAssembler", [])) 227 dom.addToRHS(lhs, [("D",D)], 228 dom.createAssembler("DefaultAssembler", [])) 229 230 res = Lsup((rhs/lhs)-2) 231 self.assertLess(res, self.TOLERANCE, 232 ("assembly for {0}expanded order %d failed with %g >= %g"%(order, 233 res, self.TOLERANCE)).format("" if expanded else "un-")) 234 235 class Test_Speckley(unittest.TestCase): 236 TOLERANCE = 1e-10 237 def test_Rectangle_ReducedFunction(self): 238 ranks = getMPISizeWorld() 239 for order in range(2, 11): 240 dom = Rectangle(order, 3, 3*ranks, l0=3, l1=3*ranks, d1=ranks) 241 X = dom.getX() 242 redData = interpolate(X, ReducedFunction(dom)) 243 data = [(interpolate(redData, ReducedFunction(dom)), "ReducedFunction"), 244 (interpolate(redData, Function(dom)), "Function"), 245 (interpolate(redData, ContinuousFunction(dom)), "ContinuousFunction")] 246 for d, fs in data: 247 self.assertLess(inf(d-[0.5]*2), self.TOLERANCE, 248 "reduced->%s failure with order %d: %g != 0"%(fs, order, inf(d-[0.5]*2))) 249 self.assertLess(sup(d[0]+0.5) - 3, self.TOLERANCE, 250 "reduced->%s failure with order %d: %g != 3"%(fs, order, sup(d[0]+0.5))) 251 self.assertLess(sup(d[1]+0.5) - 3*ranks, self.TOLERANCE, 252 "reduced->%s failure with order %d: %g != %g"%(fs, order, sup(d[1]+0.5), 3*ranks)) 253 254 def test_Brick_ReducedFunction(self): 255 ranks = getMPISizeWorld() 256 for order in range(2, 11): 257 dom = Brick(order, 3, 3*ranks, 3, l0=3, l1=3*ranks, l2=3, d1=ranks) 258 X = dom.getX() 259 redData = interpolate(X, ReducedFunction(dom)) 260 data = [(interpolate(redData, ReducedFunction(dom)), "ReducedFunction"), 261 (interpolate(redData, Function(dom)), "Function"), 262 (interpolate(redData, ContinuousFunction(dom)), "ContinuousFunction")] 263 for d, fs in data: 264 self.assertLess(inf(d-[0.5]*3), self.TOLERANCE, 265 "reduced->%s failure with order %d: %g != 0"%(fs, order, inf(d-[0.5]*3))) 266 self.assertLess(sup(d[0]+0.5) - 3, self.TOLERANCE, 267 "reduced->%s failure with order %d: %g != 3"%(fs, order, sup(d[0]+0.5))) 268 self.assertLess(sup(d[1]+0.5) - 3*ranks, self.TOLERANCE, 269 "reduced->%s failure with order %d: %g != %g"%(fs, order, sup(d[1]+0.5), 3*ranks)) 270 self.assertLess(sup(d[2]+0.5) - 3, self.TOLERANCE, 271 "reduced->%s failure with order %d: %g != 3"%(fs, order, sup(d[2]+0.5))) 272 273 def test_Rectangle_Function_gradient(self): #expanded and non-expanded 274 ranks = getMPISizeWorld() 275 for expanded in [True, False]: 276 for order in range(2,11): 277 dom = Rectangle(order, 3, 3*ranks, d1=ranks) 278 x = Data(5, Function(dom), True) 279 self.assertLess(Lsup(grad(x)), 1e-10, 280 "single component failure, order %d%s, %g >= 1e-10"%(order, 281 (" expanded" if expanded else ""), Lsup(grad(x)))) 282 for data in [[5,1], [-5,-1], [5,1,1e-5]]: 283 x = Data(data, Function(dom), True) 284 g = grad(x) 285 for n,d in enumerate(data): 286 self.assertLess(Lsup(g[n]), 1e-10, 287 "%d-component failure, order %d %sexpanded, %g >= 1e-10"%(len(data), 288 order, ("" if expanded else "un-"), Lsup(g[n]))) 289 290 def test_Brick_Function_gradient(self): 291 ranks = getMPISizeWorld() 292 for expanded in [True, False]: 293 for order in range(2,11): 294 dom = Brick(order, 3, 3*ranks, 3, d1=ranks) 295 x = Data(5, Function(dom), True) 296 self.assertLess(Lsup(grad(x)), 1e-10, 297 "single component failure, order %d%s, %g >= 1e-10"%(order, 298 (" expanded" if expanded else ""), Lsup(grad(x)))) 299 for data in [[5,1], [-5,-1], [5,1,1e-5]]: 300 x = Data(data, Function(dom), True) 301 g = grad(x) 302 for n,d in enumerate(data): 303 self.assertLess(Lsup(g[n]), 1e-10, 304 "%d-component failure, order %d %sexpanded, %g >= 1e-10"%(len(data), 305 order, ("" if expanded else "un-"), Lsup(g[n]))) 306 307 308 def test_Rectangle_ContinuousFunction_gradient(self): 309 ranks = getMPISizeWorld() 310 for order in range(2,11): 311 dom = Rectangle(order, 3, 3*ranks, d1=ranks, l0=100, l1=100) 312 X = dom.getX() 313 u = X[0] + X[1] + 1 314 v = Lsup(grad(u) - 1) 315 self.assertLess(v, 1e-10, "order %d, %g >= 1e-10, %s"%(order, v, str(grad(u)-1))) 316 for power in range(1, order+1): 317 for power2 in range(1, order+1): 318 a = X[0]**power * X[1]**power2 319 da = grad(a) 320 first = Lsup(da[0] - power*X[0]**(power-1) * X[1]**power2) \ 321 /Lsup(power*X[0]**(power-1) * X[1]**power2) 322 second = Lsup(da[1] - power2*X[1]**(power2-1) * X[0]**power) \ 323 /Lsup(power2*X[1]**(power2-1) * X[0]**power) 324 self.assertLess(first, 1e-10, 325 "order %d and degree %d,%d, %g >= 1e-9"%(order, 326 power, power2, first)) 327 self.assertLess(second, 1e-10, 328 "order %d and degree %d,%d, %g >= 1e-9"%(order, 329 power, power2, second)) 330 331 def test_Brick_ContinuousFunction_gradient(self): 332 ranks = getMPISizeWorld() 333 for order in range(2,11): 334 dom = Brick(order, 3, 3*ranks, 3, l0=100, l1=100, l2=100, d1=ranks) 335 X = dom.getX() 336 u = X[0] + X[1] + X[2] + 1 337 v = Lsup(grad(u) - 1) 338 self.assertLess(v, 1e-10, "order %d, %g >= 1e-10, %s"%(order, v, 339 str(grad(u)-1))) 340 for power1 in range(1, order+1, order//2): 341 for power2 in range(1, order+1, order//2): 342 for power3 in range(1, order+1, order//2): 343 a = X[0]**power1 * X[1]**power2 * X[2]**power3 344 da = grad(a) 345 346 temp = power1*X[0]**(power1-1) * X[1]**power2*X[2]**power3 347 first = Lsup(da[0] - temp) / Lsup(temp) 348 temp = power2*X[1]**(power2-1) * X[0]**power1*X[2]**power3 349 second = Lsup(da[1] - temp) / Lsup(temp) 350 temp = power3*X[2]**(power3-1) * X[0]**power1*X[1]**power2 351 third = Lsup(da[2] - temp) / Lsup(temp) 352 353 self.assertLess(first, 1e-10, 354 "order %d and degree %d,%d,%d, %g >= 1e-9"%(order, 355 power1, power2, power3, first)) 356 self.assertLess(second, 1e-10, 357 "order %d and degree %d,%d,%d, %g >= 1e-9"%(order, 358 power1, power2, power3, second)) 359 self.assertLess(third, 1e-10, 360 "order %d and degree %d,%d,%d, %g >= 1e-9"%(order, 361 power1, power2, power3, third)) 362 363 def test_Rectangle_interpolation_continuous_noncontinuous_and_back(self): 364 ranks = getMPISizeWorld() 365 for order in range(2,11): 366 dom = Rectangle(order, 3, 3*ranks, d1=ranks, l0=6, l1=6) 367 original = Data(5, Function(dom), True) 368 cont = interpolate(original, ContinuousFunction(dom)) 369 func = interpolate(cont, Function(dom)) 370 self.assertEqual(Lsup(original-func), 0, 371 "interpolation of constant, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 372 x = dom.getX() 373 original = x[0] + x[1] + 1 374 cont = interpolate(original, ContinuousFunction(dom)) 375 func = interpolate(cont, Function(dom)) 376 self.assertEqual(Lsup(original-func), 0, 377 "interpolation of expanded, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 378 original = whereZero(x[0]-2) + whereZero(x[1]-2) 379 cont = interpolate(original, ContinuousFunction(dom)) 380 func = interpolate(cont, Function(dom)) 381 self.assertEqual(Lsup(original-func), 0, 382 "interpolation of point, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 383 384 def test_Brick_interpolation_continuous_noncontinuous_and_back(self): 385 ranks = getMPISizeWorld() 386 for order in range(2,11): 387 dom = Brick(order, 3, 3*ranks, 3, l0=6, l1=ranks, l2=6, d1=ranks) 388 original = Data(5, Function(dom), True) 389 cont = interpolate(original, ContinuousFunction(dom)) 390 func = interpolate(cont, Function(dom)) 391 self.assertEqual(Lsup(original-func), 0, 392 "interpolation of constant, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 393 x = dom.getX() 394 original = x[0] + x[1] + x[2] + 1 395 cont = interpolate(original, ContinuousFunction(dom)) 396 func = interpolate(cont, Function(dom)) 397 self.assertEqual(Lsup(original-func), 0, 398 "interpolation of expanded, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 399 original = whereZero(x[0]-2) + whereZero(x[1]-2) + whereZero(x[2] - 2) 400 cont = interpolate(original, ContinuousFunction(dom)) 401 func = interpolate(cont, Function(dom)) 402 self.assertEqual(Lsup(original-func), 0, 403 "interpolation of point, order %d: original and final not equal, %g != 0"%(order, Lsup(original-func))) 404 405 def test_Rectangle_integration(self): 406 ranks = getMPISizeWorld() 407 for order in range(2,11): 408 size = 6 409 dom = Rectangle(order, 3, 3*ranks, d1=ranks, l0=size, l1=size) 410 X = dom.getX() 411 for k in range(1, order * 2): 412 for l in range(1, order * 2): 413 powered = X[0]**k * X[1]**l 414 integral = integrate(powered) 415 actual = size**(k+1) * size**(l+1) / ((k+1.)*(l+1.)) 416 self.assertLess(abs(integral - actual)/actual, 1e-11, 417 "too much variance in integral result (order %d, degrees %d %d)"%(order, k, l)) 418 419 def test_Brick_integration(self): 420 ranks = getMPISizeWorld() 421 for order in range(2,11): 422 size = 6 423 dom = Brick(order, 3, 3*ranks, 3, l0=6, l1=6, l2=6, d1=ranks) 424 X = dom.getX() 425 for k in [1, order, order*2 - 1]: 426 for l in [1, order, order*2 - 1]: 427 for m in [1, order, order*2 - 1]: 428 powered = X[0]**k * X[1]**l * X[2]**m 429 integral = integrate(powered) 430 actual = size**(k+1) * size**(l+1) * size**(m+1) \ 431 /((k+1.)*(l+1.)*(m+1.)) 432 res = abs(integral - actual)/actual 433 self.assertLess(res, 1e-11, 434 "too much variance in integral result (order %d, degrees %d %d, %g >= 1e-11)"%(order, 435 k, l, res)) 436 437 @unittest.skipIf(getMPISizeWorld() == 1, 438 "only works with more than one rank") 439 def test_Rectangle_MPI_construction_single_dimensional(self): 440 ranks = getMPISizeWorld() 441 for order in range(2,11): 442 dom = Rectangle(order, 2, 2*ranks, l1=ranks, d1=ranks) 443 self.assertEqual(Lsup(dom.getX()[1]+dom.getX()[0]), ranks+1, 444 "invalid getX() for y-splits order %d"%order) 445 446 filt = whereZero(dom.getX()[1] - 1) 447 for i in range(2,ranks): 448 filt += whereZero(dom.getX()[1] - i) 449 if getMPIRankWorld() % 2: 450 filt *= -1 451 d = Vector(0, Function(dom)) 452 d[0] = 1 * filt 453 d[1] = 10 * filt 454 X = interpolate(d, Function(dom)) 455 res = interpolate(X, ContinuousFunction(dom)) 456 val = Lsup(res) 457 self.assertEqual(val, 0, 458 "summation stage failure for y-splits in order %d"%order) 459 460 dom = Rectangle(2, 2*ranks, 2, l0=ranks, d0=ranks) 461 self.assertEqual(Lsup(dom.getX()[1]+dom.getX()[0]), ranks+1, 462 "invalid getX() for x-splits order %d"%order) 463 filt = whereZero(dom.getX()[0] - 1) 464 for i in range(2,getMPISizeWorld()): 465 filt += whereZero(dom.getX()[0] - i) 466 if getMPIRankWorld() % 2: 467 filt *= -1 468 d = Vector(0, Function(dom)) 469 d[0] = 1 * filt 470 d[1] = 10 * filt 471 X = interpolate(d, Function(dom)) 472 res = interpolate(X, ContinuousFunction(dom)) 473 val = Lsup(res) 474 self.assertEqual(val, 0, 475 "summation stage failure for x-splits in order %d"%order) 476 477 X = interpolate(d+2, Function(dom)) 478 res = interpolate(X, ContinuousFunction(dom)) 479 val = Lsup(res-2) 480 self.assertEqual(val, 0, 481 "averaging stage failure for x-splits in order %d"%order) 482 483 @unittest.skipIf(getMPISizeWorld() != 4, "requires 4 ranks") 484 def test_Rectangle_MPI_construction_multi_dimensional(self): 485 ranks = getMPISizeWorld() 486 half = 2 #change if ranks != 4 (sqrt(ranks)) 487 for order in range(2, 11): 488 dom = Rectangle(order, ranks, ranks, l0=half, l1=half, 489 d0=half, d1=half) 490 self.assertEqual(Lsup(dom.getX()[0]), half, 491 "invalid getX() for multidimensional splits in order %d"%order) 492 self.assertEqual(Lsup(dom.getX()[1]), half, 493 "invalid getX() for multidimensional splits in order %d"%order) 494 xfilt = whereZero(dom.getX()[0] - 1) + whereZero(dom.getX()[1] - 1) 495 for i in range(2,half): 496 xfilt += whereZero(dom.getX()[0] - i) 497 xfilt += whereZero(dom.getX()[1] - i) 498 xfilt = whereNonZero(xfilt) 499 if getMPIRankWorld() in [1,2]: #change if ranks != 4 500 xfilt *= -1 501 X = interpolate(xfilt, Function(dom)) 502 res = interpolate(X, ContinuousFunction(dom)) 503 val = Lsup(res) 504 self.assertEqual(val, 0, 505 "summation failure for mixed-splits in order %d"%order) 506 X = interpolate(xfilt+2, Function(dom)) 507 res = interpolate(X, ContinuousFunction(dom)) 508 val = Lsup(res-2) 509 self.assertEqual(val, 0, 510 "averaging failure for mixed-splits in order %d"%order) 511 512 @unittest.skipIf(getMPISizeWorld() == 1, 513 "only works with more than one rank") 514 def test_Brick_MPI_construction(self): 515 for order in range(2,11): 516 dom = Brick(order, 2*getMPISizeWorld(), 2, 2, d0=getMPISizeWorld()) 517 self.assertEqual(Lsup(dom.getX()[0] + dom.getX()[1] + dom.getX()[2]), 518 3, "invalid Lsup(getX()) for x-split order %d"%order) 519 X = interpolate(dom.getX()[0], Function(dom)) \ 520 + interpolate(dom.getX()[1], Function(dom)) 521 val = Lsup(interpolate(X, ContinuousFunction(dom)) \ 522 - (dom.getX()[0] + dom.getX()[1])) 523 self.assertLess(val, 1e-10, 524 "interpolation failure for x-split order %d"%order) 525 526 dom = Brick(order, 2, 2*getMPISizeWorld(), 2, d1=getMPISizeWorld()) 527 self.assertEqual(Lsup(dom.getX()[0] + dom.getX()[1] + dom.getX()[2]), 528 3, "invalid Lsup(getX()) for y-split order %d"%order) 529 X = interpolate(dom.getX()[0], Function(dom)) \ 530 + interpolate(dom.getX()[1], Function(dom)) 531 val = Lsup(interpolate(X, ContinuousFunction(dom)) \ 532 - (dom.getX()[0] + dom.getX()[1])) 533 self.assertLess(val, 1e-10, 534 "interpolation failure for y-split order %d"%order) 535 536 dom = Brick(order, 2, 2, 2*getMPISizeWorld(), d2=getMPISizeWorld()) 537 self.assertEqual(Lsup(dom.getX()[0] + dom.getX()[1] + dom.getX()[2]), 538 3, "invalid Lsup(getX()) for z-split order %d"%order) 539 X = interpolate(dom.getX()[0], Function(dom)) \ 540 + interpolate(dom.getX()[1], Function(dom)) 541 val = Lsup(interpolate(X, ContinuousFunction(dom)) \ 542 - (dom.getX()[0] + dom.getX()[1])) 543 self.assertLess(val, 1e-10, 544 "interpolation failure for z-split order %d"%order) 545 546 @unittest.skipIf(getMPISizeWorld() == 1, "requires multiple MPI processes") 547 def test_Brick_singledim_subdivision(self): 548 ranks = getMPISizeWorld() 549 for dim in range(0,3): 550 label = ["x","y","z"][dim] 551 size = [2,2,2] 552 size[dim] *= ranks 553 lengths = [1,1,1] 554 lengths[dim] *= ranks 555 splits = [1,1,1] 556 splits[dim] *= ranks 557 558 for order in range(2, 11): 559 dom = Brick(order, size[0], size[1], size[2], 560 l0=lengths[0], l1=lengths[1], l2=lengths[2], 561 d0=splits[0], d1=splits[1], d2=splits[2]) 562 self.assertEqual(Lsup(dom.getX()[1]+dom.getX()[0]+dom.getX()[2]), 563 ranks+2, "invalid getX() for %s-splits order %d"%(\ 564 label, order)) 565 566 filt = whereZero(dom.getX()[dim] - 1) 567 for i in range(2,ranks): 568 filt += whereZero(dom.getX()[0] - i) 569 if getMPIRankWorld() % 2: 570 filt *= -1 571 d = Vector(0, Function(dom)) 572 d[0] = 1 * filt 573 d[1] = 10 * filt 574 d[2] = 100 * filt 575 X = interpolate(d, Function(dom)) 576 res = interpolate(X, ContinuousFunction(dom)) 577 val = Lsup(res) 578 self.assertEqual(val, 0, 579 "summation stage failure for %s-splits in order %d,"%(\ 580 label, order)) 581 X = interpolate(d+2, Function(dom)) 582 res = interpolate(X, ContinuousFunction(dom)) 583 val = Lsup(res-2) 584 self.assertEqual(val, 0, 585 "averaging stage failure for %s-splits in order %d,"%(\ 586 label, order)) 587 588 @unittest.skipIf(getMPISizeWorld() != 4, "requires 4 ranks exactly") 589 def test_Brick_multidim_subdivision(self): 590 ranks = getMPISizeWorld() 591 half = 2 #change if ranks != 4 (sqrt(ranks)) 592 for order in range(2, 11): 593 for dom,dim1,dim2 in [ 594 (Brick(order, 2*half, 2*half, 2, l0=half, l1=half, 595 d0=half, d1=half), 0,1), 596 (Brick(order, 2*half, 2, 2*half, l0=half, l2=half, 597 d0=half, d2=half), 0,2), 598 (Brick(order, 2, 2*half, 2*half, l1=half, l2=half, 599 d1=half, d2=half), 1,2)]: 600 self.assertEqual(ranks + 1, 601 Lsup(dom.getX()[0] + dom.getX()[1] + dom.getX()[2]), 602 "invalid getX() for multidimensional split " + \ 603 "(dims %d,%d) in order %d"%(dim1,dim2,order)) 604 xfilt = whereZero(dom.getX()[dim1] - 1) \ 605 + whereZero(dom.getX()[dim2] - 1) 606 for i in range(2,half): 607 xfilt += whereZero(dom.getX()[dim1] - i) 608 xfilt += whereZero(dom.getX()[dim2] - i) 609 xfilt = whereNonZero(xfilt) 610 if getMPIRankWorld() in [1,2]: #change if ranks != 4 611 xfilt *= -1 612 X = interpolate(xfilt, Function(dom)) 613 res = interpolate(X, ContinuousFunction(dom)) 614 val = Lsup(res) 615 self.assertEqual(val, 0, 616 "summation failure for mixed-splits " \ 617 + "(dims %d,%d) in order %d"%(dim1,dim2,order)) 618 X = interpolate(xfilt+2, Function(dom)) 619 res = interpolate(X, ContinuousFunction(dom)) 620 val = Lsup(res-2) 621 self.assertEqual(val, 0, 622 "averaging failure for mixed-splits "\ 623 + "(dims %d,%d) in order %d"%(dim1,dim2,order)) 624 625 if __name__ == '__main__': 626 run_tests(__name__, exit_on_failure=True) 627

 ViewVC Help Powered by ViewVC 1.1.26