# Contents of /trunk/finley/test/python/run_models.py

Revision 2264 - (show annotations)
Wed Feb 11 06:48:28 2009 UTC (13 years, 4 months ago) by gross
File MIME type: text/x-python
File size: 30608 byte(s)
```a new darcy flux solver.
```
 1 2 ######################################################## 3 # 4 # Copyright (c) 2003-2008 by University of Queensland 5 # Earth Systems Science Computational Center (ESSCC) 6 7 # 8 # Primary Business: Queensland, Australia 9 # Licensed under the Open Software License version 3.0 10 11 # 12 ######################################################## 13 14 __copyright__="""Copyright (c) 2003-2008 by University of Queensland 15 Earth Systems Science Computational Center (ESSCC) 16 http://www.uq.edu.au/esscc 17 Primary Business: Queensland, Australia""" 18 __license__="""Licensed under the Open Software License version 3.0 19 20 __url__= 21 22 import unittest 23 import tempfile 24 25 from esys.escript import * 26 from esys.finley import Rectangle 27 from esys.escript.models import DarcyFlow 28 import sys 29 import os 30 try: 31 FINLEY_WORKDIR=os.environ['FINLEY_WORKDIR'] 32 except KeyError: 33 FINLEY_WORKDIR='.' 34 35 36 VERBOSE=False # or True 37 DETAIL_VERBOSE=False 38 39 from esys.escript import * 40 from esys.escript.models import StokesProblemCartesian 41 from esys.finley import Rectangle, Brick 42 43 from esys.escript.models import Mountains 44 from math import pi 45 46 #==================================================================================================================== 47 class Test_StokesProblemCartesian2D(unittest.TestCase): 48 def setUp(self): 49 NE=6 50 self.TOL=1e-3 51 self.domain=Rectangle(NE,NE,order=2,useFullElementOrder=True) 52 def tearDown(self): 53 del self.domain 54 def test_PCG_P_0(self): 55 ETA=1. 56 P1=0. 57 58 x=self.domain.getX() 59 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 60 mask=whereZero(x[0]) * [1.,1.] \ 61 +whereZero(x[0]-1) * [1.,1.] \ 62 +whereZero(x[1]) * [1.,0.] \ 63 +whereZero(x[1]-1) * [1.,1.] 64 65 sp=StokesProblemCartesian(self.domain) 66 67 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 68 u0=(1-x[0])*x[0]*[0.,1.] 69 p0=Scalar(P1,ReducedSolution(self.domain)) 70 sp.setTolerance(self.TOL) 71 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=100,usePCG=True) 72 73 error_v0=Lsup(u[0]-u0[0]) 74 error_v1=Lsup(u[1]-u0[1])/0.25 75 error_p=Lsup(p+P1*x[0]*x[1]) 76 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 77 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 78 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 79 80 def test_PCG_P_small(self): 81 ETA=1. 82 P1=1. 83 84 x=self.domain.getX() 85 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 86 mask=whereZero(x[0]) * [1.,1.] \ 87 +whereZero(x[0]-1) * [1.,1.] \ 88 +whereZero(x[1]) * [1.,0.] \ 89 +whereZero(x[1]-1) * [1.,1.] 90 91 sp=StokesProblemCartesian(self.domain) 92 93 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 94 u0=(1-x[0])*x[0]*[0.,1.] 95 p0=Scalar(P1,ReducedSolution(self.domain)) 96 sp.setTolerance(self.TOL*0.2) 97 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=100,usePCG=True) 98 error_v0=Lsup(u[0]-u0[0]) 99 error_v1=Lsup(u[1]-u0[1])/0.25 100 error_p=Lsup(P1*x[0]*x[1]+p) 101 saveVTK("d.vtu",p=p, e=P1*x[0]*x[1]+p, p_ref=P1*x[0]*x[1]) 102 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 103 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 104 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 105 106 def test_PCG_P_large(self): 107 ETA=1. 108 P1=1000. 109 110 x=self.domain.getX() 111 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 112 mask=whereZero(x[0]) * [1.,1.] \ 113 +whereZero(x[0]-1) * [1.,1.] \ 114 +whereZero(x[1]) * [1.,0.] \ 115 +whereZero(x[1]-1) * [1.,1.] 116 117 sp=StokesProblemCartesian(self.domain) 118 119 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 120 u0=(1-x[0])*x[0]*[0.,1.] 121 p0=Scalar(P1,ReducedSolution(self.domain)) 122 sp.setTolerance(self.TOL) 123 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=100,usePCG=True) 124 125 error_v0=Lsup(u[0]-u0[0]) 126 error_v1=Lsup(u[1]-u0[1])/0.25 127 error_p=Lsup(P1*x[0]*x[1]+p)/P1 128 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 129 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 130 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 131 132 def test_GMRES_P_0(self): 133 ETA=1. 134 P1=0. 135 136 x=self.domain.getX() 137 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 138 mask=whereZero(x[0]) * [1.,1.] \ 139 +whereZero(x[0]-1) * [1.,1.] \ 140 +whereZero(x[1]) * [1.,0.] \ 141 +whereZero(x[1]-1) * [1.,1.] 142 143 sp=StokesProblemCartesian(self.domain) 144 145 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 146 u0=(1-x[0])*x[0]*[0.,1.] 147 p0=Scalar(P1,ReducedSolution(self.domain)) 148 sp.setTolerance(self.TOL) 149 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=50,usePCG=False,iter_restart=18) 150 151 error_v0=Lsup(u[0]-u0[0]) 152 error_v1=Lsup(u[1]-u0[1])/0.25 153 error_p=Lsup(P1*x[0]*x[1]+p) 154 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 155 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 156 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 157 158 def test_GMRES_P_small(self): 159 ETA=1. 160 P1=1. 161 162 x=self.domain.getX() 163 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 164 mask=whereZero(x[0]) * [1.,1.] \ 165 +whereZero(x[0]-1) * [1.,1.] \ 166 +whereZero(x[1]) * [1.,0.] \ 167 +whereZero(x[1]-1) * [1.,1.] 168 169 sp=StokesProblemCartesian(self.domain) 170 171 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 172 u0=(1-x[0])*x[0]*[0.,1.] 173 p0=Scalar(P1,ReducedSolution(self.domain)) 174 sp.setTolerance(self.TOL*0.1) 175 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=20,usePCG=False) 176 177 error_v0=Lsup(u[0]-u0[0]) 178 error_v1=Lsup(u[1]-u0[1])/0.25 179 error_p=Lsup(P1*x[0]*x[1]+p) 180 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 181 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 182 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 183 184 def test_GMRES_P_large(self): 185 ETA=1. 186 P1=1000. 187 188 x=self.domain.getX() 189 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.] 190 mask=whereZero(x[0]) * [1.,1.] \ 191 +whereZero(x[0]-1) * [1.,1.] \ 192 +whereZero(x[1]) * [1.,0.] \ 193 +whereZero(x[1]-1) * [1.,1.] 194 195 sp=StokesProblemCartesian(self.domain) 196 197 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 198 u0=(1-x[0])*x[0]*[0.,1.] 199 p0=Scalar(P1,ReducedSolution(self.domain)) 200 sp.setTolerance(self.TOL) 201 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=100,usePCG=False) 202 203 error_v0=Lsup(u[0]-u0[0]) 204 error_v1=Lsup(u[1]-u0[1])/0.25 205 error_p=Lsup(P1*x[0]*x[1]+p)/P1 206 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 207 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 208 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 209 #==================================================================================================================== 210 class Test_StokesProblemCartesian3D(unittest.TestCase): 211 def setUp(self): 212 NE=6 213 self.TOL=1e-4 214 self.domain=Brick(NE,NE,NE,order=2,useFullElementOrder=True) 215 def tearDown(self): 216 del self.domain 217 def test_PCG_P_0(self): 218 ETA=1. 219 P1=0. 220 221 x=self.domain.getX() 222 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 223 x=self.domain.getX() 224 mask=whereZero(x[0]) * [1.,1.,1.] \ 225 +whereZero(x[0]-1) * [1.,1.,1.] \ 226 +whereZero(x[1]) * [1.,0.,1.] \ 227 +whereZero(x[1]-1) * [1.,1.,1.] \ 228 +whereZero(x[2]) * [1.,1.,0.] \ 229 +whereZero(x[2]-1) * [1.,1.,1.] 230 231 232 sp=StokesProblemCartesian(self.domain) 233 234 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 235 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 236 p0=Scalar(P1,ReducedSolution(self.domain)) 237 sp.setTolerance(self.TOL) 238 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE ,max_iter=100,usePCG=True) 239 240 error_v0=Lsup(u[0]-u0[0]) 241 error_v1=Lsup(u[1]-u0[1]) 242 error_v2=Lsup(u[2]-u0[2])/0.25**2 243 error_p=Lsup(P1*x[0]*x[1]*x[2]+p) 244 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 245 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 246 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 247 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 248 249 def test_PCG_P_small(self): 250 ETA=1. 251 P1=1. 252 253 x=self.domain.getX() 254 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 255 mask=whereZero(x[0]) * [1.,1.,1.] \ 256 +whereZero(x[0]-1) * [1.,1.,1.] \ 257 +whereZero(x[1]) * [1.,0.,1.] \ 258 +whereZero(x[1]-1) * [1.,1.,1.] \ 259 +whereZero(x[2]) * [1.,1.,0.] \ 260 +whereZero(x[2]-1) * [1.,1.,1.] 261 262 263 sp=StokesProblemCartesian(self.domain) 264 265 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 266 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 267 p0=Scalar(P1,ReducedSolution(self.domain)) 268 sp.setTolerance(self.TOL*0.1) 269 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE ,max_iter=100,usePCG=True) 270 error_v0=Lsup(u[0]-u0[0]) 271 error_v1=Lsup(u[1]-u0[1]) 272 error_v2=Lsup(u[2]-u0[2])/0.25**2 273 error_p=Lsup(P1*x[0]*x[1]*x[2]+p) 274 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 275 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 276 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 277 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 278 def test_PCG_P_large(self): 279 ETA=1. 280 P1=1000. 281 282 x=self.domain.getX() 283 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 284 mask=whereZero(x[0]) * [1.,1.,1.] \ 285 +whereZero(x[0]-1) * [1.,1.,1.] \ 286 +whereZero(x[1]) * [1.,0.,1.] \ 287 +whereZero(x[1]-1) * [1.,1.,1.] \ 288 +whereZero(x[2]) * [1.,1.,0.] \ 289 +whereZero(x[2]-1) * [1.,1.,1.] 290 291 292 sp=StokesProblemCartesian(self.domain) 293 294 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 295 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 296 p0=Scalar(P1,ReducedSolution(self.domain)) 297 sp.setTolerance(self.TOL) 298 u,p=sp.solve(u0,-p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE ,max_iter=100,usePCG=True) 299 300 error_v0=Lsup(u[0]-u0[0]) 301 error_v1=Lsup(u[1]-u0[1]) 302 error_v2=Lsup(u[2]-u0[2])/0.25**2 303 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1 304 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 305 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 306 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 307 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 308 309 def test_GMRES_P_0(self): 310 ETA=1. 311 P1=0. 312 313 x=self.domain.getX() 314 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 315 x=self.domain.getX() 316 mask=whereZero(x[0]) * [1.,1.,1.] \ 317 +whereZero(x[0]-1) * [1.,1.,1.] \ 318 +whereZero(x[1]) * [1.,1.,1.] \ 319 +whereZero(x[1]-1) * [1.,1.,1.] \ 320 +whereZero(x[2]) * [1.,1.,0.] \ 321 +whereZero(x[2]-1) * [1.,1.,1.] 322 323 324 sp=StokesProblemCartesian(self.domain) 325 326 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 327 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 328 p0=Scalar(P1,ReducedSolution(self.domain)) 329 sp.setTolerance(self.TOL) 330 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE,max_iter=100,usePCG=False,iter_restart=20) 331 332 error_v0=Lsup(u[0]-u0[0]) 333 error_v1=Lsup(u[1]-u0[1]) 334 error_v2=Lsup(u[2]-u0[2])/0.25**2 335 error_p=Lsup(P1*x[0]*x[1]*x[2]+p) 336 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 337 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 338 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 339 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 340 def test_GMRES_P_small(self): 341 ETA=1. 342 P1=1. 343 344 x=self.domain.getX() 345 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 346 mask=whereZero(x[0]) * [1.,1.,1.] \ 347 +whereZero(x[0]-1) * [1.,1.,1.] \ 348 +whereZero(x[1]) * [1.,1.,1.] \ 349 +whereZero(x[1]-1) * [1.,1.,1.] \ 350 +whereZero(x[2]) * [1.,1.,0.] \ 351 +whereZero(x[2]-1) * [1.,1.,1.] 352 353 354 sp=StokesProblemCartesian(self.domain) 355 356 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 357 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 358 p0=Scalar(P1,ReducedSolution(self.domain)) 359 sp.setTolerance(self.TOL*0.1) 360 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE and False, verbose=VERBOSE,max_iter=100,usePCG=False) 361 362 error_v0=Lsup(u[0]-u0[0]) 363 error_v1=Lsup(u[1]-u0[1]) 364 error_v2=Lsup(u[2]-u0[2])/0.25**2 365 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1 366 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 367 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 368 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 369 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 370 def test_GMRES_P_large(self): 371 ETA=1. 372 P1=1000. 373 374 x=self.domain.getX() 375 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.] 376 mask=whereZero(x[0]) * [1.,1.,1.] \ 377 +whereZero(x[0]-1) * [1.,1.,1.] \ 378 +whereZero(x[1]) * [1.,0.,1.] \ 379 +whereZero(x[1]-1) * [1.,1.,1.] \ 380 +whereZero(x[2]) * [1.,1.,0.] \ 381 +whereZero(x[2]-1) * [1.,1.,1.] 382 383 384 sp=StokesProblemCartesian(self.domain) 385 386 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA) 387 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.] 388 p0=Scalar(P1,ReducedSolution(self.domain)) 389 sp.setTolerance(self.TOL) 390 u,p=sp.solve(u0,p0,show_details=DETAIL_VERBOSE, verbose=VERBOSE ,max_iter=100,usePCG=False) 391 392 error_v0=Lsup(u[0]-u0[0]) 393 error_v1=Lsup(u[1]-u0[1]) 394 error_v2=Lsup(u[2]-u0[2])/0.25**2 395 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1 396 self.failUnless(error_p<10*self.TOL, "pressure error too large.") 397 self.failUnless(error_v0<10*self.TOL, "0-velocity error too large.") 398 self.failUnless(error_v1<10*self.TOL, "1-velocity error too large.") 399 self.failUnless(error_v2<10*self.TOL, "2-velocity error too large.") 400 #==================================================================================================================== 401 class Test_Darcy(unittest.TestCase): 402 # this is a simple test for the darcy flux problem 403 # 404 # 405 # p = 1/k * ( 1/2* (fb-f0)/xb* x **2 + f0 * x - ub*x ) + p0 406 # 407 # with f = (fb-f0)/xb* x + f0 408 # 409 # u = f - k * p,x = ub 410 # 411 # we prescribe pressure at x=x0=0 to p0 412 # 413 # if we prescribe the pressure on the bottom x=xb we set 414 # 415 # pb= 1/k * ( 1/2* (fb-f0)* xb + f0 * xb - ub*xb ) + p0 = 1/k * ((fb+f0)/2 - ub ) * xb + p0 416 # 417 # which leads to ub = (fb+f0)/2-k*(pb-p0)/xb 418 # 419 def rescaleDomain(self): 420 x=self.dom.getX().copy() 421 for i in xrange(self.dom.getDim()): 422 x_inf=inf(x[i]) 423 x_sup=sup(x[i]) 424 if i == self.dom.getDim()-1: 425 x[i]=-self.WIDTH*(x[i]-x_sup)/(x_inf-x_sup) 426 else: 427 x[i]=self.WIDTH*(x[i]-x_inf)/(x_sup-x_inf) 428 self.dom.setX(x) 429 def getScalarMask(self,include_bottom=True): 430 x=self.dom.getX().copy() 431 x_inf=inf(x[self.dom.getDim()-1]) 432 x_sup=sup(x[self.dom.getDim()-1]) 433 out=whereZero(x[self.dom.getDim()-1]-x_sup) 434 if include_bottom: out+=whereZero(x[self.dom.getDim()-1]-x_inf) 435 return wherePositive(out) 436 def getVectorMask(self,include_bottom=True): 437 x=self.dom.getX().copy() 438 out=Vector(0.,Solution(self.dom)) 439 for i in xrange(self.dom.getDim()): 440 x_inf=inf(x[i]) 441 x_sup=sup(x[i]) 442 if i != self.dom.getDim()-1: out[i]+=whereZero(x[i]-x_sup) 443 if i != self.dom.getDim()-1 or include_bottom: out[i]+=whereZero(x[i]-x_inf) 444 return wherePositive(out) 445 446 def setSolutionFixedBottom(self, p0, pb, f0, fb, k): 447 d=self.dom.getDim() 448 x=self.dom.getX()[d-1] 449 xb=inf(x) 450 u=Vector(0.,Solution(self.dom))+kronecker(d)[d-1]*((f0+fb)/2.-k*(pb-p0)/xb) 451 p=1./k*((fb-f0)/(xb*2.)* x**2 - (fb-f0)/2.*x)+(pb-p0)/xb*x + p0 452 f= ((fb-f0)/xb* x + f0)*kronecker(Function(self.dom))[d-1] 453 return u,p,f 454 455 def testConstF_FixedBottom_smallK(self): 456 k=1.e-10 457 mp=self.getScalarMask(include_bottom=True) 458 mv=self.getVectorMask(include_bottom=False) 459 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k) 460 p=p_ref*mp 461 u=u_ref*mv 462 df=DarcyFlow(self.dom) 463 df.setValue(g=f, 464 location_of_fixed_pressure=mp, 465 location_of_fixed_flux=mv, 466 permeability=Scalar(k,Function(self.dom))) 467 df.setTolerance(rtol=self.TOL) 468 df.setSubProblemTolerance() 469 v,p=df.solve(u_ref,p, max_iter=100, verbose=VERBOSE) 470 self.failUnless(Lsup(v-u_ref)

 ViewVC Help Powered by ViewVC 1.1.26