/[escript]/trunk/finley/test/python/run_amg.py
ViewVC logotype

Contents of /trunk/finley/test/python/run_amg.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3440 - (show annotations)
Fri Jan 14 00:04:53 2011 UTC (8 years, 7 months ago) by gross
File MIME type: text/x-python
File size: 34490 byte(s)
a second AMG interpolation implemented
1 # -*- coding: utf-8 -*-
2
3 ########################################################
4 #
5 # Copyright (c) 2003-2010 by University of Queensland
6 # Earth Systems Science Computational Center (ESSCC)
7 # http://www.uq.edu.au/esscc
8 #
9 # Primary Business: Queensland, Australia
10 # Licensed under the Open Software License version 3.0
11 # http://www.opensource.org/licenses/osl-3.0.php
12 #
13 ########################################################
14
15 __copyright__="""Copyright (c) 2003-2010 by University of Queensland
16 Earth Systems Science Computational Center (ESSCC)
17 http://www.uq.edu.au/esscc
18 Primary Business: Queensland, Australia"""
19 __license__="""Licensed under the Open Software License version 3.0
20 http://www.opensource.org/licenses/osl-3.0.php"""
21 __url__="https://launchpad.net/escript-finley"
22
23 """
24 Test suite for AMG
25
26 :remark:
27
28 :var __author__: name of author
29 :var __licence__: licence agreement
30 :var __url__: url entry point on documentation
31 :var __version__: version
32 :var __date__: date of the version
33 """
34
35 __author__="Lutz Gross, l.gross@uq.edu.au"
36
37 import unittest, sys
38
39 from esys.escript import *
40 from esys.finley import Rectangle,Brick
41 from esys.escript.linearPDEs import LinearPDE, SolverOptions
42 import numpy
43 OPTIMIZE=True
44 SOLVER_VERBOSE=True
45
46 MIN_MATRIX_SIZE=None
47 MIN_SPARSITY=None
48 MAX_LEVEL=None
49 USE_AMG=True or False
50
51 try:
52 FINLEY_TEST_DATA=os.environ['FINLEY_TEST_DATA']
53 except KeyError:
54 FINLEY_TEST_DATA='.'
55 FINLEY_TEST_MESH_PATH=os.path.join(FINLEY_TEST_DATA,"data_meshes")
56
57 # number of elements in the spatial directions
58 NE_TOTAL=4096
59
60 class Test_AMG(unittest.TestCase):
61
62 def test_Poisson(self):
63 x=Solution(self.domain).getX()
64 # --- set exact solution ----
65 u_ex=Scalar(1,Solution(self.domain))
66 g_ex=Vector(0.,Solution(self.domain))
67 for i in xrange(self.domain.getDim()):
68 u_ex+=(i+1)*x[i]
69 g_ex[i]=(i+1)
70
71 # create PDE:
72 pde=LinearPDE(self.domain,numEquations=1)
73 pde.setSymmetryOn()
74 mask=whereZero(x[0])
75 pde.setValue(r=u_ex,q=mask)
76 pde.setValue(A=kronecker(self.domain),y=inner(g_ex,self.domain.getNormal()))
77 # -------- get the solution ---------------------------
78 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
79 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
80 #if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
81 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
82 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
83 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
84 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
85
86 u=pde.getSolution()
87 # -------- test the solution ---------------------------
88 error=Lsup(u-u_ex)/Lsup(u_ex)
89 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
90
91 def test_PoissonWithDirectInterpolation(self):
92 x=Solution(self.domain).getX()
93 # --- set exact solution ----
94 u_ex=Scalar(1,Solution(self.domain))
95 g_ex=Vector(0.,Solution(self.domain))
96 for i in xrange(self.domain.getDim()):
97 u_ex+=(i+1)*x[i]
98 g_ex[i]=(i+1)
99
100 # create PDE:
101 pde=LinearPDE(self.domain,numEquations=1)
102 pde.setSymmetryOn()
103 mask=whereZero(x[0])
104 pde.setValue(r=u_ex,q=mask)
105 pde.setValue(A=kronecker(self.domain),y=inner(g_ex,self.domain.getNormal()))
106 # -------- get the solution ---------------------------
107 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
108 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
109 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
110 pde.getSolverOptions().setNumPreSweeps(3)
111 pde.getSolverOptions().setNumPostSweeps(3)
112 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
113 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
114 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
115 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
116 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().DIRECT_INTERPOLATION)
117
118 u=pde.getSolution()
119 # -------- test the solution ---------------------------
120 error=Lsup(u-u_ex)/Lsup(u_ex)
121 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
122
123 def test_PoissonClassic(self):
124 x=Solution(self.domain).getX()
125 # --- set exact solution ----
126 u_ex=Scalar(1,Solution(self.domain))
127 g_ex=Vector(0.,Solution(self.domain))
128 for i in xrange(self.domain.getDim()):
129 u_ex+=(i+1)*x[i]
130 g_ex[i]=(i+1)
131
132 # create PDE:
133 pde=LinearPDE(self.domain,numEquations=1)
134 pde.setSymmetryOn()
135 mask=whereZero(x[0])
136 pde.setValue(r=u_ex,q=mask)
137 pde.setValue(A=kronecker(self.domain),y=inner(g_ex,self.domain.getNormal()))
138 # -------- get the solution ---------------------------
139 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
140 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
141 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
142 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
143 pde.getSolverOptions().setNumPreSweeps(3)
144 pde.getSolverOptions().setNumPostSweeps(3)
145 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
146 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
147 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
148 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().CLASSIC_INTERPOLATION)
149
150 u=pde.getSolution()
151 # -------- test the solution ---------------------------
152 error=Lsup(u-u_ex)/Lsup(u_ex)
153 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
154
155 def test_PoissonClassicWithFFCoupling(self):
156 x=Solution(self.domain).getX()
157 # --- set exact solution ----
158 u_ex=Scalar(1,Solution(self.domain))
159 g_ex=Vector(0.,Solution(self.domain))
160 for i in xrange(self.domain.getDim()):
161 u_ex+=(i+1)*x[i]
162 g_ex[i]=(i+1)
163
164 # create PDE:
165 pde=LinearPDE(self.domain,numEquations=1)
166 pde.setSymmetryOn()
167 mask=whereZero(x[0])
168 pde.setValue(r=u_ex,q=mask)
169 pde.setValue(A=kronecker(self.domain),y=inner(g_ex,self.domain.getNormal()))
170 # -------- get the solution ---------------------------
171 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
172 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
173 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
174 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
175 pde.getSolverOptions().setNumPreSweeps(3)
176 pde.getSolverOptions().setNumPostSweeps(3)
177 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
178 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
179 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
180 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().CLASSIC_INTERPOLATION_WITH_FF_COUPLING)
181
182 u=pde.getSolution()
183 # -------- test the solution ---------------------------
184 error=Lsup(u-u_ex)/Lsup(u_ex)
185 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
186
187 def test_PoissonSqueezedX(self):
188 x=self.domain.getX().copy()
189 x[0]*=0.5
190 self.domain.setX(x)
191 x=Solution(self.domain).getX()
192 # --- set exact solution ----
193 u_ex=Scalar(1,Solution(self.domain))
194 g_ex=Vector(0.,Solution(self.domain))
195 for i in xrange(self.domain.getDim()):
196 u_ex+=(i+1)*x[i]
197 g_ex[i]=(i+1)
198
199 # create PDE:
200 pde=LinearPDE(self.domain,numEquations=1)
201 pde.setSymmetryOn()
202 mask=whereZero(x[0])
203 pde.setValue(r=u_ex,q=mask)
204 pde.setValue(A=kronecker(self.domain),y=inner(g_ex,self.domain.getNormal()))
205 # -------- get the solution ---------------------------
206 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
207 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
208 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
209 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
210 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
211 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
212 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
213
214 u=pde.getSolution()
215 # -------- test the solution ---------------------------
216 error=Lsup(u-u_ex)/Lsup(u_ex)
217 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
218
219
220 def test_Poisson2(self):
221 x=Solution(self.domain).getX()
222 # --- set exact solution ----
223 u_ex=Data(1.,(2,),Solution(self.domain))
224 g_ex=Data(0.,(2,self.domain.getDim()), Solution(self.domain))
225 A=Data(0.,(2,self.domain.getDim(),2,self.domain.getDim()), Function(self.domain))
226 for i in xrange(self.domain.getDim()):
227 u_ex[0]+= 1*(i+1) *x[i]
228 g_ex[0,i]=1*(i+1)
229 u_ex[1]+= 2*(i+1)*x[i]
230 g_ex[1,i]=2*(i+1)
231 A[0,i,0,i]=1
232 A[1,i,1,i]=1
233
234 # create PDE:
235 pde=LinearPDE(self.domain,numEquations=2)
236 pde.setSymmetryOn()
237 mask=whereZero(x[0])*[1,1]
238 pde.setValue(r=u_ex,q=mask)
239 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
240 # -------- get the solution ---------------------------
241 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
242 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
243 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
244 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
245 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
246 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
247 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
248
249 u=pde.getSolution()
250 # -------- test the solution ---------------------------
251 error=Lsup(u-u_ex)/Lsup(u_ex)
252 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
253
254 def test_Poisson3(self):
255 x=Solution(self.domain).getX()
256 # --- set exact solution ----
257 u_ex=Data(1.,(3,),Solution(self.domain))
258 g_ex=Data(0.,(3,self.domain.getDim()), Solution(self.domain))
259 A=Data(0.,(3,self.domain.getDim(),3,self.domain.getDim()), Function(self.domain))
260 for i in xrange(self.domain.getDim()):
261 u_ex[0]+= 1*(i+1) *x[i]
262 g_ex[0,i]=1*(i+1)
263 u_ex[1]+= 2*(i+1)*x[i]
264 g_ex[1,i]=2*(i+1)
265 u_ex[2]+= 3*(i+1)*x[i]
266 g_ex[2,i]=3*(i+1)
267 A[0,i,0,i]=1
268 A[1,i,1,i]=1
269 A[2,i,2,i]=1
270
271 # create PDE:
272 pde=LinearPDE(self.domain,numEquations=3)
273 pde.setSymmetryOn()
274 mask=whereZero(x[0])*[1,1,1]
275 pde.setValue(r=u_ex,q=mask)
276 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
277 # -------- get the solution ---------------------------
278 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
279 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
280 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
281 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
282 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
283 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
284 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
285
286 u=pde.getSolution()
287 # -------- test the solution ---------------------------
288 error=Lsup(u-u_ex)/Lsup(u_ex)
289 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
290
291 def test_Poisson4(self):
292 x=Solution(self.domain).getX()
293 # --- set exact solution ----
294 u_ex=Data(1.,(4,),Solution(self.domain))
295 g_ex=Data(0.,(4,self.domain.getDim()), Solution(self.domain))
296 A=Data(0.,(4,self.domain.getDim(),4,self.domain.getDim()), Function(self.domain))
297 for i in xrange(self.domain.getDim()):
298 u_ex[0]+= 1*(i+1) *x[i]
299 g_ex[0,i]=1*(i+1)
300 u_ex[1]+= 2*(i+1)*x[i]
301 g_ex[1,i]=2*(i+1)
302 u_ex[2]+= 3*(i+1)*x[i]
303 g_ex[2,i]=3*(i+1)
304 u_ex[3]+= 4*(i+1)*x[i]
305 g_ex[3,i]=4*(i+1)
306 A[0,i,0,i]=1
307 A[1,i,1,i]=1
308 A[2,i,2,i]=1
309 A[3,i,3,i]=1
310
311 # create PDE:
312 pde=LinearPDE(self.domain,numEquations=4)
313 pde.setSymmetryOn()
314 mask=whereZero(x[0])*[1,1,1,1]
315 pde.setValue(r=u_ex,q=mask)
316 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
317 # -------- get the solution ---------------------------
318 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
319 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
320 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
321 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
322 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
323 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
324 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
325
326 u=pde.getSolution()
327 # -------- test the solution ---------------------------
328 error=Lsup(u-u_ex)/Lsup(u_ex)
329 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
330 def test_Poisson2Classic(self):
331 x=Solution(self.domain).getX()
332 # --- set exact solution ----
333 u_ex=Data(1.,(2,),Solution(self.domain))
334 g_ex=Data(0.,(2,self.domain.getDim()), Solution(self.domain))
335 A=Data(0.,(2,self.domain.getDim(),2,self.domain.getDim()), Function(self.domain))
336 for i in xrange(self.domain.getDim()):
337 u_ex[0]+= 1*(i+1) *x[i]
338 g_ex[0,i]=1*(i+1)
339 u_ex[1]+= 2*(i+1)*x[i]
340 g_ex[1,i]=2*(i+1)
341 A[0,i,0,i]=1
342 A[1,i,1,i]=1
343
344 # create PDE:
345 pde=LinearPDE(self.domain,numEquations=2)
346 pde.setSymmetryOn()
347 mask=whereZero(x[0])*[1,1]
348 pde.setValue(r=u_ex,q=mask)
349 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
350 # -------- get the solution ---------------------------
351 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
352 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
353 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
354 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
355 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
356 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
357 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
358 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().CLASSIC_INTERPOLATION)
359
360 u=pde.getSolution()
361 # -------- test the solution ---------------------------
362 error=Lsup(u-u_ex)/Lsup(u_ex)
363 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
364
365 def test_Poisson3Classic(self):
366 x=Solution(self.domain).getX()
367 # --- set exact solution ----
368 u_ex=Data(1.,(3,),Solution(self.domain))
369 g_ex=Data(0.,(3,self.domain.getDim()), Solution(self.domain))
370 A=Data(0.,(3,self.domain.getDim(),3,self.domain.getDim()), Function(self.domain))
371 for i in xrange(self.domain.getDim()):
372 u_ex[0]+= 1*(i+1) *x[i]
373 g_ex[0,i]=1*(i+1)
374 u_ex[1]+= 2*(i+1)*x[i]
375 g_ex[1,i]=2*(i+1)
376 u_ex[2]+= 3*(i+1)*x[i]
377 g_ex[2,i]=3*(i+1)
378 A[0,i,0,i]=1
379 A[1,i,1,i]=1
380 A[2,i,2,i]=1
381
382 # create PDE:
383 pde=LinearPDE(self.domain,numEquations=3)
384 pde.setSymmetryOn()
385 mask=whereZero(x[0])*[1,1,1]
386 pde.setValue(r=u_ex,q=mask)
387 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
388 # -------- get the solution ---------------------------
389 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
390 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
391 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
392 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
393 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
394 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
395 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
396 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().CLASSIC_INTERPOLATION)
397
398 u=pde.getSolution()
399 # -------- test the solution ---------------------------
400 error=Lsup(u-u_ex)/Lsup(u_ex)
401 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
402
403 def test_Poisson4Classic(self):
404 x=Solution(self.domain).getX()
405 # --- set exact solution ----
406 u_ex=Data(1.,(4,),Solution(self.domain))
407 g_ex=Data(0.,(4,self.domain.getDim()), Solution(self.domain))
408 A=Data(0.,(4,self.domain.getDim(),4,self.domain.getDim()), Function(self.domain))
409 for i in xrange(self.domain.getDim()):
410 u_ex[0]+= 1*(i+1) *x[i]
411 g_ex[0,i]=1*(i+1)
412 u_ex[1]+= 2*(i+1)*x[i]
413 g_ex[1,i]=2*(i+1)
414 u_ex[2]+= 3*(i+1)*x[i]
415 g_ex[2,i]=3*(i+1)
416 u_ex[3]+= 4*(i+1)*x[i]
417 g_ex[3,i]=4*(i+1)
418 A[0,i,0,i]=1
419 A[1,i,1,i]=1
420 A[2,i,2,i]=1
421 A[3,i,3,i]=1
422
423 # create PDE:
424 pde=LinearPDE(self.domain,numEquations=4)
425 pde.setSymmetryOn()
426 mask=whereZero(x[0])*[1,1,1,1]
427 pde.setValue(r=u_ex,q=mask)
428 pde.setValue(A=A,y=matrixmult(g_ex,self.domain.getNormal()))
429 # -------- get the solution ---------------------------
430 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
431 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
432 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
433 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
434 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
435 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
436 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
437 pde.getSolverOptions().setAMGInterpolation(pde.getSolverOptions().CLASSIC_INTERPOLATION)
438
439 u=pde.getSolution()
440 # -------- test the solution ---------------------------
441 error=Lsup(u-u_ex)/Lsup(u_ex)
442 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
443
444 def test_WeakCoupled2(self):
445 x=Solution(self.domain).getX()
446 # --- set exact solution ----
447 u_ex=Data(1.,(2,),Solution(self.domain))
448 g_ex=Data(0.,(2,self.domain.getDim()), Solution(self.domain))
449 A=Data(0.,(2,self.domain.getDim(),2,self.domain.getDim()), Function(self.domain))
450 for i in xrange(self.domain.getDim()):
451 u_ex[0]+= 1*(i+1)*x[i]
452 g_ex[0,i]=1*(i+1)
453 u_ex[1]+= 2*(i+1)*x[i]
454 g_ex[1,i]=2*(i+1)
455 A[0,i,0,i]=1
456 A[1,i,1,i]=1
457
458 Y=Data(0.,(2,),Function(self.domain))
459 Y[0]= u_ex[1]
460 Y[1]=u_ex[0]
461 # create PDE:
462 pde=LinearPDE(self.domain,numEquations=2)
463 pde.setSymmetryOn()
464 mask=whereZero(x[0])*[1,1]
465 pde.setValue(r=u_ex,q=mask)
466 pde.setValue(A=A,
467 D=-0.01*(numpy.ones((2,2))-kronecker(2)),
468 Y=-0.01*Y,
469 y=matrixmult(g_ex,self.domain.getNormal()))
470 # -------- get the solution ---------------------------
471 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
472 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
473 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
474 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
475 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
476 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
477 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
478
479 u=pde.getSolution()
480 # -------- test the solution ---------------------------
481 error=Lsup(u-u_ex)/Lsup(u_ex)
482 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
483
484 def test_WeakCoupled3(self):
485 x=Solution(self.domain).getX()
486 # --- set exact solution ----
487 u_ex=Data(1.,(3,),Solution(self.domain))
488 g_ex=Data(0.,(3,self.domain.getDim()), Solution(self.domain))
489 A=Data(0.,(3,self.domain.getDim(),3,self.domain.getDim()), Function(self.domain))
490 for i in xrange(self.domain.getDim()):
491 u_ex[0]+= 1*(i+1)*x[i]
492 g_ex[0,i]=1*(i+1)
493 u_ex[1]+= 2*(i+1)*x[i]
494 g_ex[1,i]=2*(i+1)
495 u_ex[2]+= 3*(i+1)*x[i]
496 g_ex[2,i]=3*(i+1)
497 A[0,i,0,i]=1
498 A[1,i,1,i]=1
499 A[2,i,2,i]=1
500
501 Y=Data(0.,(3,),Function(self.domain))
502 Y[0]= u_ex[1]+u_ex[2]
503 Y[1]=u_ex[0]+ u_ex[2]
504 Y[2]=u_ex[0]+u_ex[1]
505 # create PDE:
506 pde=LinearPDE(self.domain,numEquations=3)
507 pde.setSymmetryOn()
508 mask=whereZero(x[0])*[1,1,1]
509 pde.setValue(r=u_ex,q=mask)
510 pde.setValue(A=A,
511 D=-0.01*(numpy.ones((3,3))-kronecker(3)),
512 Y=-0.01*Y,
513 y=matrixmult(g_ex,self.domain.getNormal()))
514 # -------- get the solution ---------------------------
515 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
516 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
517 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
518 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
519 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
520 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
521 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
522 u=pde.getSolution()
523 # -------- test the solution ---------------------------
524 error=Lsup(u-u_ex)/Lsup(u_ex)
525 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
526
527 def test_WeakCoupled4(self):
528 x=Solution(self.domain).getX()
529 # --- set exact solution ----
530 u_ex=Data(1.,(4,),Solution(self.domain))
531 g_ex=Data(0.,(4,self.domain.getDim()), Solution(self.domain))
532 A=Data(0.,(4,self.domain.getDim(),4,self.domain.getDim()), Function(self.domain))
533 for i in xrange(self.domain.getDim()):
534 u_ex[0]+= 1*(i+1)*x[i]
535 g_ex[0,i]=1*(i+1)
536 u_ex[1]+= 2*(i+1)*x[i]
537 g_ex[1,i]=2*(i+1)
538 u_ex[2]+= 3*(i+1)*x[i]
539 g_ex[2,i]=3*(i+1)
540 u_ex[3]+= 4*(i+1)*x[i]
541 g_ex[3,i]=4*(i+1)
542 A[0,i,0,i]=1
543 A[1,i,1,i]=1
544 A[2,i,2,i]=1
545 A[3,i,3,i]=1
546
547 Y=Data(0.,(4,),Function(self.domain))
548 Y[0]= u_ex[1]+u_ex[2]+u_ex[3]
549 Y[1]=u_ex[0]+ u_ex[2]+u_ex[3]
550 Y[2]=u_ex[0]+u_ex[1]+ u_ex[3]
551 Y[3]=u_ex[0]+u_ex[1]+u_ex[2]
552 # create PDE:
553 pde=LinearPDE(self.domain,numEquations=4)
554 pde.setSymmetryOn()
555 mask=whereZero(x[0])*[1,1,1,1]
556 pde.setValue(r=u_ex,q=mask)
557 pde.setValue(A=A,
558 D=-0.01*(numpy.ones((4,4))-kronecker(4)),
559 Y=-0.01*Y,
560 y=matrixmult(g_ex,self.domain.getNormal()))
561 # -------- get the solution ---------------------------
562 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
563 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
564 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
565 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
566 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
567 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
568 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
569 u=pde.getSolution()
570 # -------- test the solution ---------------------------
571 error=Lsup(u-u_ex)/Lsup(u_ex)
572 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
573
574 def test_StrongCoupled2(self):
575 x=Solution(self.domain).getX()
576 # --- set exact solution ----
577 u_ex=Data(1.,(2,),Solution(self.domain))
578 g_ex=Data(0.,(2,self.domain.getDim()), Solution(self.domain))
579 A=Data(0.,(2,self.domain.getDim(),2,self.domain.getDim()), Function(self.domain))
580 for i in xrange(self.domain.getDim()):
581 u_ex[0]+= 1*(i+1)*x[i]
582 g_ex[0,i]=1*(i+1)
583 u_ex[1]+= 2*(i+1)*x[i]
584 g_ex[1,i]=2*(i+1)
585 A[0,i,0,i]=1
586 A[1,i,1,i]=1
587
588 Y=Data(0.,(2,),Function(self.domain))
589 Y[0]= u_ex[1]
590 Y[1]=u_ex[0]
591 # create PDE:
592 pde=LinearPDE(self.domain,numEquations=2)
593 pde.setSymmetryOn()
594 mask=whereZero(x[0])*[1,1]
595 pde.setValue(r=u_ex,q=mask)
596 pde.setValue(A=A,
597 D=-20.*(numpy.ones((2,2))-kronecker(2)),
598 Y=-20.*Y,
599 y=matrixmult(g_ex,self.domain.getNormal()))
600 # -------- get the solution ---------------------------
601 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
602 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
603 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
604 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
605 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
606 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
607 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
608 u=pde.getSolution()
609 # -------- test the solution ---------------------------
610 error=Lsup(u-u_ex)/Lsup(u_ex)
611 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
612
613 def test_StrongCoupled3(self):
614 x=Solution(self.domain).getX()
615 # --- set exact solution ----
616 u_ex=Data(1.,(3,),Solution(self.domain))
617 g_ex=Data(0.,(3,self.domain.getDim()), Solution(self.domain))
618 A=Data(0.,(3,self.domain.getDim(),3,self.domain.getDim()), Function(self.domain))
619 for i in xrange(self.domain.getDim()):
620 u_ex[0]+= 1*(i+1)*x[i]
621 g_ex[0,i]=1*(i+1)
622 u_ex[1]+= 2*(i+1)*x[i]
623 g_ex[1,i]=2*(i+1)
624 u_ex[2]+= 3*(i+1)*x[i]
625 g_ex[2,i]=3*(i+1)
626 A[0,i,0,i]=1
627 A[1,i,1,i]=1
628 A[2,i,2,i]=1
629
630 Y=Data(0.,(3,),Function(self.domain))
631 Y[0]= u_ex[1]+u_ex[2]
632 Y[1]=u_ex[0]+ u_ex[2]
633 Y[2]=u_ex[0]+u_ex[1]
634 # create PDE:
635 pde=LinearPDE(self.domain,numEquations=3)
636 pde.setSymmetryOn()
637 mask=whereZero(x[0])*[1,1,1]
638 pde.setValue(r=u_ex,q=mask)
639 pde.setValue(A=A,
640 D=-20.*(numpy.ones((3,3))-kronecker(3)),
641 Y=-20.*Y,
642 y=matrixmult(g_ex,self.domain.getNormal()))
643 # -------- get the solution ---------------------------
644 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
645 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
646 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
647 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
648 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
649 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
650 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
651
652 u=pde.getSolution()
653 # -------- test the solution ---------------------------
654 error=Lsup(u-u_ex)/Lsup(u_ex)
655 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
656 def test_StrongCoupled4(self):
657 x=Solution(self.domain).getX()
658 # --- set exact solution ----
659 u_ex=Data(1.,(4,),Solution(self.domain))
660 g_ex=Data(0.,(4,self.domain.getDim()), Solution(self.domain))
661 A=Data(0.,(4,self.domain.getDim(),4,self.domain.getDim()), Function(self.domain))
662 for i in xrange(self.domain.getDim()):
663 u_ex[0]+= 1*(i+1)*x[i]
664 g_ex[0,i]=1*(i+1)
665 u_ex[1]+= 2*(i+1)*x[i]
666 g_ex[1,i]=2*(i+1)
667 u_ex[2]+= 3*(i+1)*x[i]
668 g_ex[2,i]=3*(i+1)
669 u_ex[3]+= 4*(i+1)*x[i]
670 g_ex[3,i]=4*(i+1)
671 A[0,i,0,i]=1
672 A[1,i,1,i]=1
673 A[2,i,2,i]=1
674 A[3,i,3,i]=1
675
676 Y=Data(0.,(4,),Function(self.domain))
677 Y[0]= u_ex[1]+u_ex[2]+u_ex[3]
678 Y[1]=u_ex[0]+ u_ex[2]+u_ex[3]
679 Y[2]=u_ex[0]+u_ex[1]+ u_ex[3]
680 Y[3]=u_ex[0]+u_ex[1]+u_ex[2]
681 # create PDE:
682 pde=LinearPDE(self.domain,numEquations=4)
683 pde.setSymmetryOn()
684 mask=whereZero(x[0])*[1,1,1,1]
685 pde.setValue(r=u_ex,q=mask)
686 pde.setValue(A=A,
687 D=-20.*(numpy.ones((4,4))-kronecker(4)),
688 Y=-20.*Y,
689 y=matrixmult(g_ex,self.domain.getNormal()))
690 # -------- get the solution ---------------------------
691 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
692 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
693 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
694 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
695 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
696 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
697 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
698 u=pde.getSolution()
699 # -------- test the solution ---------------------------
700 error=Lsup(u-u_ex)/Lsup(u_ex)
701 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
702
703 def test_Square(self):
704 # PDE constants
705 h1 = 0.5
706 h2 = 0.5
707 beta = 1.
708 alpha1 = 2.
709 alpha2 = 1.
710
711 # domain masks and domain-specific constants
712 x = Solution(self.domain).getX(); x0 = x[0]; x1 = x[1]
713 omega2 = wherePositive(x0-h1)*wherePositive(x1-h2)
714 omega1 = 1-omega2
715 ratio = alpha1/alpha2
716 alpha = alpha1*omega1 + alpha2*omega2
717
718 # --- set exact solution ----
719 a1 = 1.
720 d1 = 1.
721 b1 = -d1*h2
722 c1 = -d1*h1
723
724 a2 = a1 - (1.-ratio)*d1*h1*h2
725 b2 = ratio*b1
726 c2 = ratio*c1
727 d2 = ratio*d1
728
729 u_ex = omega1*(a1 + b1*x0 + c1*x1 + d1*x0*x1) + \
730 omega2*(a2 + b2*x0 + c2*x1 + d2*x0*x1)
731
732 # create PDE:
733 pde = LinearPDE(self.domain,numEquations=1)
734
735 # set the value to that of the solution on the boundary
736 q = whereZero(x0) + whereZero(x1) + \
737 whereZero(sup(x0)-x0) + whereZero(sup(x1)-x1)
738 pde.setValue(q=q,r=u_ex)
739
740 # create X points in the centre of the grid elements
741 xe = Function(self.domain).getX()
742 x0 = xe[0]
743 x1 = xe[1]
744
745 # redefine omega so that apha is more precise on the diagonal (?)
746 omega2 = wherePositive(x0-h1)*wherePositive(x1-h2)
747 omega1 = 1-omega2
748 ratio = alpha1/alpha2
749 alpha = alpha1*omega1 + alpha2*omega2
750
751 # set up PDE coefficients
752 pde.setValue(A=alpha*kronecker(self.domain), D=beta, Y=beta*u_ex)
753 pde.setSymmetryOn()
754
755 # -------- get the solution ---------------------------
756 pde.getSolverOptions().setTolerance(self.SOLVER_TOL)
757 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
758 if (USE_AMG): pde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
759 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
760 if MIN_MATRIX_SIZE!= None: pde.getSolverOptions().setMinCoarseMatrixSize(MIN_MATRIX_SIZE)
761 if MIN_SPARSITY!=None: pde.getSolverOptions().setMinCoarseMatrixSparsity(MIN_SPARSITY)
762 if MAX_LEVEL!=None: pde.getSolverOptions().setLevelMax(MAX_LEVEL)
763 u = pde.getSolution()
764
765 # -------- test the solution ---------------------------
766 error=Lsup(u-u_ex)/Lsup(u_ex)
767 self.assertTrue(error<self.RES_TOL, "solution error %s is too big."%error)
768
769
770 class Test_AMGOnFinleyHex2DOrder1(Test_AMG):
771 RES_TOL=5.e-7
772 SOLVER_TOL=1.e-8
773 def setUp(self):
774 NE=int(float(NE_TOTAL)**(1./2.)+0.5)
775 self.domain = Rectangle(NE,NE,1, optimize=OPTIMIZE)
776 def tearDown(self):
777 del self.domain
778
779 class Test_AMGOnFinleyHex3DOrder1(Test_AMG):
780 RES_TOL=5.e-7
781 SOLVER_TOL=1.e-8
782 def setUp(self):
783 NE=int(float(NE_TOTAL)**(1./3.)+0.5)
784 self.domain = Brick(NE,NE,NE,1, optimize=OPTIMIZE)
785 def tearDown(self):
786 del self.domain
787 if __name__ == '__main__':
788 suite = unittest.TestSuite()
789 suite.addTest(unittest.makeSuite(Test_AMGOnFinleyHex2DOrder1))
790 suite.addTest(unittest.makeSuite(Test_AMGOnFinleyHex3DOrder1))
791
792 s=unittest.TextTestRunner(verbosity=2).run(suite)
793 if not s.wasSuccessful(): sys.exit(1)

  ViewVC Help
Powered by ViewVC 1.1.26