/[escript]/trunk/ripley/test/python/run_simplesolve.py
ViewVC logotype

Annotation of /trunk/ripley/test/python/run_simplesolve.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3781 - (hide annotations)
Mon Jan 23 22:18:49 2012 UTC (7 years, 7 months ago) by caltinay
Original Path: branches/ripleygmg_from_3668/ripley/test/python/run_simplesolve.py
File MIME type: text/x-python
File size: 25860 byte(s)
Since setX() is not supported the constructors now take location and extent.
Python wrapper allows construction the old way with l0,l1 or with
tuples (x0,x1),(y0,y1).

1 caltinay 3776
2     ########################################################
3     #
4 caltinay 3781 # Copyright (c) 2003-2012 by University of Queensland
5 caltinay 3776 # Earth Systems Science Computational Center (ESSCC)
6     # http://www.uq.edu.au/esscc
7     #
8     # Primary Business: Queensland, Australia
9     # Licensed under the Open Software License version 3.0
10     # http://www.opensource.org/licenses/osl-3.0.php
11     #
12     ########################################################
13    
14 caltinay 3781 __copyright__="""Copyright (c) 2003-2012 by University of Queensland
15 caltinay 3776 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     http://www.opensource.org/licenses/osl-3.0.php"""
20     __url__="https://launchpad.net/escript-finley"
21    
22     """
23     Test suite for the linearPDE and pdetools test on ripley
24    
25     :remark:
26    
27     :var __author__: name of author
28     :var __licence__: licence agreement
29     :var __url__: url entry point on documentation
30     :var __version__: version
31     :var __date__: date of the version
32     """
33    
34     __author__="Lutz Gross, l.gross@uq.edu.au"
35    
36     import unittest, sys
37    
38     from esys.escript import *
39     from esys.ripley import Rectangle,Brick
40     from esys.escript.linearPDEs import LinearPDE, SolverOptions
41     import numpy
42     SOLVER_VERBOSE=False
43     # setNumberOfThreads(2)
44    
45     try:
46     RIPLEY_TEST_DATA=os.environ['RIPLEY_TEST_DATA']
47     except KeyError:
48     RIPLEY_TEST_DATA='.'
49    
50     # number of elements in the spatial directions
51     NE0=8
52     NE1=10
53     NE2=12
54    
55     NE0=12
56     NE1=12
57     NE2=8
58    
59     SOLVER_TOL=1.e-8
60     REL_TOL=1.e-6
61    
62     FAC_DIAG=1.
63     FAC_OFFDIAG=-0.4
64    
65    
66     class SimpleSolve_Rectangle_SinglePDE_Paso_BICGSTAB_Jacobi(unittest.TestCase):
67     def test_solve(self):
68     # Tell about how many MPI CPUs and OpenMP threads
69     domain=Rectangle(NE0,NE1)
70     x=Solution(domain).getX()
71     # --- set exact solution ----
72     u_ex=Scalar(0,Solution(domain))
73     u_ex=1.+2.*x[0]+3.*x[1]
74     # --- set exact gradient -----------
75     g_ex=Data(0.,(2,),Solution(domain))
76     g_ex[0]=2.
77     g_ex[1]=3.
78     # -------- test gradient --------------------------------
79     g=grad(u_ex)
80     self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
81     # -------- set-up PDE -----------------------------------
82     pde=LinearPDE(domain,numEquations=1)
83     mask=whereZero(x[0])
84     pde.setValue(r=u_ex,q=mask)
85     pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
86     # -------- get the solution ---------------------------
87     pde.getSolverOptions().setTolerance(SOLVER_TOL)
88     pde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
89     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
90     pde.getSolverOptions().setPackage(SolverOptions.PASO)
91     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
92     u=pde.getSolution()
93     # -------- test the solution ---------------------------
94     error=Lsup(u-u_ex)
95     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
96     class SimpleSolve_Rectangle_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
97     def test_solve(self):
98     domain=Rectangle(NE0,NE1)
99     x=Solution(domain).getX()
100     # --- set exact solution ----
101     u_ex=Scalar(0,Solution(domain))
102     u_ex=1.+2.*x[0]+3.*x[1]
103     # --- set exact gradient -----------
104     g_ex=Data(0.,(2,),Solution(domain))
105     g_ex[0]=2.
106     g_ex[1]=3.
107     # -------- test gradient --------------------------------
108     g=grad(u_ex)
109     self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
110     # -------- set-up PDE -----------------------------------
111     pde=LinearPDE(domain,numEquations=1)
112     mask=whereZero(x[0])
113     pde.setValue(r=u_ex,q=mask)
114     pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
115     # -------- get the solution ---------------------------
116     pde.getSolverOptions().setTolerance(SOLVER_TOL)
117     pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
118     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
119     pde.getSolverOptions().setPackage(SolverOptions.PASO)
120     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
121     u=pde.getSolution()
122     # -------- test the solution ---------------------------
123     error=Lsup(u-u_ex)
124     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
125     class SimpleSolve_Rectangle_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
126     def test_solve(self):
127     domain=Rectangle(NE0,NE1)
128     x=Solution(domain).getX()
129     # --- set exact solution ----
130     u_ex=Vector(0,Solution(domain))
131     u_ex[0]=1.+2.*x[0]+3.*x[1]
132     u_ex[1]=-1.+3.*x[0]+2.*x[1]
133     # --- set exact gradient -----------
134     g_ex=Data(0.,(2,2),Solution(domain))
135     g_ex[0,0]=2.
136     g_ex[0,1]=3.
137     g_ex[1,0]=3.
138     g_ex[1,1]=2.
139     # -------- test gradient --------------------------------
140     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
141     # -------- set-up PDE -----------------------------------
142     pde=LinearPDE(domain,numEquations=2)
143     mask=whereZero(x[0])
144     pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
145     A=Tensor4(0,Function(domain))
146     A[0,:,0,:]=kronecker(2)
147     A[1,:,1,:]=kronecker(2)
148     Y=Vector(0.,Function(domain))
149     Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
150     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
151     pde.setValue(A=A,
152     D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
153     Y=Y,
154     y=matrixmult(g_ex,domain.getNormal()))
155     # -------- get the solution ---------------------------
156     pde.getSolverOptions().setTolerance(SOLVER_TOL)
157     pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
158     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
159     pde.getSolverOptions().setPackage(SolverOptions.PASO)
160     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
161     u=pde.getSolution()
162     # -------- test the solution ---------------------------
163     error=Lsup(u-u_ex)
164     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
165     class SimpleSolve_Rectangle_Order2_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
166     def test_solve(self):
167     domain=Rectangle(NE0,NE1,2,l0=1.,l1=1)
168     x=Solution(domain).getX()
169     # --- set exact solution ----
170     u_ex=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
171     # --- set exact gradient -----------
172     g_ex=Data(0.,(2,),Solution(domain))
173     g_ex[0]=2.+8.*x[0]+5.*x[1]
174     g_ex[1]=3.+5.*x[0]+12.*x[1]
175     # -------- test gradient --------------------------------
176     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
177     # -------- set-up PDE -----------------------------------
178     pde=LinearPDE(domain,numEquations=1)
179     mask=whereZero(x[0])
180     pde.setValue(r=u_ex,q=mask)
181     pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()),Y=-20.)
182     # -------- get the solution ---------------------------
183     pde.getSolverOptions().setTolerance(SOLVER_TOL)
184     pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
185     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
186     pde.getSolverOptions().setPackage(SolverOptions.PASO)
187     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
188     u=pde.getSolution()
189     # -------- test the solution ---------------------------
190     error=Lsup(u-u_ex)
191     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
192     class SimpleSolve_Brick_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
193     def test_solve(self):
194     domain=Brick(NE0,NE1,NE2)
195     x=Solution(domain).getX()
196     u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
197     # --- set exact gradient -----------
198     g_ex=Data(0.,(3,),Solution(domain))
199     g_ex[0]=2.
200     g_ex[1]=3.
201     g_ex[2]=4.
202     # -------- test gradient --------------------------------
203     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
204     # -------- set-up PDE -----------------------------------
205     pde=LinearPDE(domain,numEquations=1)
206     mask=whereZero(x[0])
207     pde.setValue(r=u_ex,q=mask)
208     pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
209     # -------- get the solution ---------------------------
210     pde.getSolverOptions().setTolerance(SOLVER_TOL)
211     pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
212     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
213     pde.getSolverOptions().setPackage(SolverOptions.PASO)
214     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
215     u=pde.getSolution()
216     # -------- test the solution ---------------------------
217     error=Lsup(u-u_ex)
218     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
219     class SimpleSolve_Brick_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
220     def test_solve(self):
221     domain=Brick(NE0,NE1,NE2)
222     x=Solution(domain).getX()
223     # --- set exact solution ----
224     u_ex=Vector(0,Solution(domain))
225     u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
226     u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
227     u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
228     # --- set exact gradient -----------
229     g_ex=Data(0.,(3,3),Solution(domain))
230     g_ex[0,0]=2.
231     g_ex[0,1]=3.
232     g_ex[0,2]=4.
233     g_ex[1,0]=4.
234     g_ex[1,1]=1.
235     g_ex[1,2]=-2.
236     g_ex[2,0]=8.
237     g_ex[2,1]=4.
238     g_ex[2,2]=5.
239     # -------- test gradient --------------------------------
240     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
241     # -------- set-up PDE -----------------------------------
242     pde=LinearPDE(domain,numEquations=3)
243     mask=whereZero(x[0])
244     pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
245     A=Tensor4(0,Function(domain))
246     A[0,:,0,:]=kronecker(3)
247     A[1,:,1,:]=kronecker(3)
248     A[2,:,2,:]=kronecker(3)
249     Y=Vector(0.,Function(domain))
250     Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
251     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
252     Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
253     pde.setValue(A=A,
254     D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
255     Y=Y,
256     y=matrixmult(g_ex,domain.getNormal()))
257     # -------- get the solution ---------------------------
258     pde.getSolverOptions().setTolerance(SOLVER_TOL)
259     pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
260     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
261     pde.getSolverOptions().setPackage(SolverOptions.PASO)
262     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
263     u=pde.getSolution()
264     # -------- test the solution ---------------------------
265     error=Lsup(u-u_ex)
266     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
267    
268     class SimpleSolve_Rectangle_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
269     def test_solve(self):
270     domain=Rectangle(NE0,NE1)
271     x=Solution(domain).getX()
272     # --- set exact solution ----
273     u_ex=Scalar(0,Solution(domain))
274     u_ex=1.+2.*x[0]+3.*x[1]
275     # --- set exact gradient -----------
276     g_ex=Data(0.,(2,),Solution(domain))
277     g_ex[0]=2.
278     g_ex[1]=3.
279     # -------- test gradient --------------------------------
280     g=grad(u_ex)
281     self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
282     # -------- set-up PDE -----------------------------------
283     pde=LinearPDE(domain,numEquations=1)
284     mask=whereZero(x[0])
285     pde.setValue(r=u_ex,q=mask)
286     pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
287     # -------- get the solution ---------------------------
288     pde.getSolverOptions().setTolerance(SOLVER_TOL)
289     pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
290     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
291     pde.getSolverOptions().setPackage(SolverOptions.PASO)
292     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
293     u=pde.getSolution()
294     # -------- test the solution ---------------------------
295     error=Lsup(u-u_ex)
296     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
297    
298     class SimpleSolve_Rectangle_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
299     def test_solve(self):
300     domain=Rectangle(NE0,NE1)
301     x=Solution(domain).getX()
302     # --- set exact solution ----
303     u_ex=Vector(0,Solution(domain))
304     u_ex[0]=1.+2.*x[0]+3.*x[1]
305     u_ex[1]=-1.+3.*x[0]+2.*x[1]
306     # --- set exact gradient -----------
307     g_ex=Data(0.,(2,2),Solution(domain))
308     g_ex[0,0]=2.
309     g_ex[0,1]=3.
310     g_ex[1,0]=3.
311     g_ex[1,1]=2.
312     # -------- test gradient --------------------------------
313     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
314     # -------- set-up PDE -----------------------------------
315     pde=LinearPDE(domain,numEquations=2)
316     mask=whereZero(x[0])
317     pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
318     A=Tensor4(0,Function(domain))
319     A[0,:,0,:]=kronecker(2)
320     A[1,:,1,:]=kronecker(2)
321     Y=Vector(0.,Function(domain))
322     Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
323     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
324     pde.setValue(A=A,
325     D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
326     Y=Y,
327     y=matrixmult(g_ex,domain.getNormal()))
328     # -------- get the solution ---------------------------
329     pde.getSolverOptions().setTolerance(SOLVER_TOL)
330     pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
331     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
332     pde.getSolverOptions().setPackage(SolverOptions.PASO)
333     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
334     u=pde.getSolution()
335     # -------- test the solution ---------------------------
336     error=Lsup(u-u_ex)
337     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
338     class SimpleSolve_Brick_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
339     def test_solve(self):
340     domain=Brick(NE0,NE1,NE2)
341     x=Solution(domain).getX()
342     u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
343     # --- set exact gradient -----------
344     g_ex=Data(0.,(3,),Solution(domain))
345     g_ex[0]=2.
346     g_ex[1]=3.
347     g_ex[2]=4.
348     # -------- test gradient --------------------------------
349     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
350     # -------- set-up PDE -----------------------------------
351     pde=LinearPDE(domain,numEquations=1)
352     mask=whereZero(x[0])
353     pde.setValue(r=u_ex,q=mask)
354     pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
355     # -------- get the solution ---------------------------
356     pde.getSolverOptions().setTolerance(SOLVER_TOL)
357     pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
358     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
359     pde.getSolverOptions().setPackage(SolverOptions.PASO)
360     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
361     u=pde.getSolution()
362     # -------- test the solution ---------------------------
363     error=Lsup(u-u_ex)
364     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
365    
366     class SimpleSolve_Brick_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
367     def test_solve(self):
368     domain=Brick(NE0,NE1,NE2)
369     x=Solution(domain).getX()
370     # --- set exact solution ----
371     u_ex=Vector(0,Solution(domain))
372     u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
373     u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
374     u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
375     # --- set exact gradient -----------
376     g_ex=Data(0.,(3,3),Solution(domain))
377     g_ex[0,0]=2.
378     g_ex[0,1]=3.
379     g_ex[0,2]=4.
380     g_ex[1,0]=4.
381     g_ex[1,1]=1.
382     g_ex[1,2]=-2.
383     g_ex[2,0]=8.
384     g_ex[2,1]=4.
385     g_ex[2,2]=5.
386     # -------- test gradient --------------------------------
387     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
388     # -------- set-up PDE -----------------------------------
389     pde=LinearPDE(domain,numEquations=3)
390     mask=whereZero(x[0])
391     pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
392     A=Tensor4(0,Function(domain))
393     A[0,:,0,:]=kronecker(3)
394     A[1,:,1,:]=kronecker(3)
395     A[2,:,2,:]=kronecker(3)
396     Y=Vector(0.,Function(domain))
397     Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
398     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
399     Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
400     pde.setValue(A=A,
401     D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
402     Y=Y,
403     y=matrixmult(g_ex,domain.getNormal()))
404     # -------- get the solution ---------------------------
405     pde.getSolverOptions().setTolerance(SOLVER_TOL)
406     pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
407     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
408     pde.getSolverOptions().setPackage(SolverOptions.PASO)
409     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
410     u=pde.getSolution()
411     # -------- test the solution ---------------------------
412     error=Lsup(u-u_ex)
413     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
414    
415     class SimpleSolve_Rectangle_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
416     def test_solve(self):
417     domain=Rectangle(NE0,NE1)
418     x=Solution(domain).getX()
419     # --- set exact solution ----
420     u_ex=Scalar(0,Solution(domain))
421     u_ex=1.+2.*x[0]+3.*x[1]
422     # --- set exact gradient -----------
423     g_ex=Data(0.,(2,),Solution(domain))
424     g_ex[0]=2.
425     g_ex[1]=3.
426     # -------- test gradient --------------------------------
427     g=grad(u_ex)
428     self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
429     # -------- set-up PDE -----------------------------------
430     pde=LinearPDE(domain,numEquations=1)
431     mask=whereZero(x[0])
432     pde.setValue(r=u_ex,q=mask)
433     pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
434     # -------- get the solution ---------------------------
435     pde.getSolverOptions().setTolerance(SOLVER_TOL)
436     pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
437     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
438     pde.getSolverOptions().setPackage(SolverOptions.PASO)
439    
440     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
441     u=pde.getSolution()
442     # -------- test the solution ---------------------------
443     error=Lsup(u-u_ex)
444     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
445    
446     class SimpleSolve_Rectangle_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
447     def test_solve(self):
448     domain=Rectangle(NE0,NE1)
449     x=Solution(domain).getX()
450     # --- set exact solution ----
451     u_ex=Vector(0,Solution(domain))
452     u_ex[0]=1.+2.*x[0]+3.*x[1]
453     u_ex[1]=-1.+3.*x[0]+2.*x[1]
454     # --- set exact gradient -----------
455     g_ex=Data(0.,(2,2),Solution(domain))
456     g_ex[0,0]=2.
457     g_ex[0,1]=3.
458     g_ex[1,0]=3.
459     g_ex[1,1]=2.
460     # -------- test gradient --------------------------------
461     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
462     # -------- set-up PDE -----------------------------------
463     pde=LinearPDE(domain,numEquations=2)
464     mask=whereZero(x[0])
465     pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
466     A=Tensor4(0,Function(domain))
467     A[0,:,0,:]=kronecker(2)
468     A[1,:,1,:]=kronecker(2)
469     Y=Vector(0.,Function(domain))
470     Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
471     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
472     pde.setValue(A=A,
473     D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
474     Y=Y,
475     y=matrixmult(g_ex,domain.getNormal()))
476     # -------- get the solution ---------------------------
477     pde.getSolverOptions().setTolerance(SOLVER_TOL)
478     pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
479     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
480     pde.getSolverOptions().setPackage(SolverOptions.PASO)
481     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
482     u=pde.getSolution()
483     # -------- test the solution ---------------------------
484     error=Lsup(u-u_ex)
485     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
486     class SimpleSolve_Brick_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
487     def test_solve(self):
488     domain=Brick(NE0,NE1,NE2)
489     x=Solution(domain).getX()
490     u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
491     # --- set exact gradient -----------
492     g_ex=Data(0.,(3,),Solution(domain))
493     g_ex[0]=2.
494     g_ex[1]=3.
495     g_ex[2]=4.
496     # -------- test gradient --------------------------------
497     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
498     # -------- set-up PDE -----------------------------------
499     pde=LinearPDE(domain,numEquations=1)
500     mask=whereZero(x[0])
501     pde.setValue(r=u_ex,q=mask)
502     pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
503     # -------- get the solution ---------------------------
504     pde.getSolverOptions().setTolerance(SOLVER_TOL)
505     pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
506     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
507     pde.getSolverOptions().setPackage(SolverOptions.PASO)
508     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
509     u=pde.getSolution()
510     # -------- test the solution ---------------------------
511     error=Lsup(u-u_ex)
512     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
513    
514     class SimpleSolve_Brick_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
515     def test_solve(self):
516     domain=Brick(NE0,NE1,NE2)
517     x=Solution(domain).getX()
518     # --- set exact solution ----
519     u_ex=Vector(0,Solution(domain))
520     u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
521     u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
522     u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
523     # --- set exact gradient -----------
524     g_ex=Data(0.,(3,3),Solution(domain))
525     g_ex[0,0]=2.
526     g_ex[0,1]=3.
527     g_ex[0,2]=4.
528     g_ex[1,0]=4.
529     g_ex[1,1]=1.
530     g_ex[1,2]=-2.
531     g_ex[2,0]=8.
532     g_ex[2,1]=4.
533     g_ex[2,2]=5.
534     # -------- test gradient --------------------------------
535     self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
536     # -------- set-up PDE -----------------------------------
537     pde=LinearPDE(domain,numEquations=3)
538     mask=whereZero(x[0])
539     pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
540     A=Tensor4(0,Function(domain))
541     A[0,:,0,:]=kronecker(3)
542     A[1,:,1,:]=kronecker(3)
543     A[2,:,2,:]=kronecker(3)
544     Y=Vector(0.,Function(domain))
545     Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
546     Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
547     Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
548     pde.setValue(A=A,
549     D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
550     Y=Y,
551     y=matrixmult(g_ex,domain.getNormal()))
552     # -------- get the solution ---------------------------
553     pde.getSolverOptions().setTolerance(SOLVER_TOL)
554     pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
555     pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
556     pde.getSolverOptions().setPackage(SolverOptions.PASO)
557     pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
558     u=pde.getSolution()
559     # -------- test the solution ---------------------------
560     error=Lsup(u-u_ex)
561     self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
562    
563    
564     if __name__ == '__main__':
565     suite = unittest.TestSuite()
566     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SinglePDE_Paso_BICGSTAB_Jacobi))
567     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SinglePDE_Paso_PCG_Jacobi))
568     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SystemPDE_Paso_PCG_Jacobi))
569     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SinglePDE_Paso_PCG_Jacobi))
570     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SystemPDE_Paso_PCG_Jacobi))
571     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SinglePDE_Paso_TFQMR_Jacobi))
572     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SystemPDE_Paso_TFQMR_Jacobi))
573     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SinglePDE_Paso_TFQMR_Jacobi))
574     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SystemPDE_Paso_TFQMR_Jacobi))
575     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SinglePDE_Paso_MINRES_Jacobi))
576     suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_SinglePDE_Paso_MINRES_Jacobi))
577     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SinglePDE_Paso_MINRES_Jacobi))
578     suite.addTest(unittest.makeSuite(SimpleSolve_Brick_SystemPDE_Paso_MINRES_Jacobi))
579    
580     s=unittest.TextTestRunner(verbosity=2).run(suite)
581     if not s.wasSuccessful(): sys.exit(1)

  ViewVC Help
Powered by ViewVC 1.1.26