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

Contents of /trunk/dudley/test/python/run_simplesolve.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5706 - (show annotations)
Mon Jun 29 03:41:36 2015 UTC (4 years, 2 months ago) by sshaw
File MIME type: text/x-python
File size: 47513 byte(s)
all python files now force use of python3 prints and division syntax to stop sneaky errors appearing in py3 environs
1
2 ##############################################################################
3 #
4 # Copyright (c) 2003-2015 by The University of Queensland
5 # http://www.uq.edu.au
6 #
7 # Primary Business: Queensland, Australia
8 # Licensed under the Open Software License version 3.0
9 # http://www.opensource.org/licenses/osl-3.0.php
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 print_function, division
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 http://www.opensource.org/licenses/osl-3.0.php"""
24 __url__="https://launchpad.net/escript-finley"
25
26 """
27 Test suite for the linearPDE and pdetools test on dudley
28
29 :remark:
30
31 :var __author__: name of author
32 :var __licence__: licence agreement
33 :var __url__: url entry point on documentation
34 :var __version__: version
35 :var __date__: date of the version
36 """
37
38 __author__="Lutz Gross, l.gross@uq.edu.au"
39
40 import esys.escriptcore.utestselect as unittest, sys
41 from esys.escriptcore.testing import *
42 from esys.escript import *
43 from esys.dudley import Rectangle,Brick
44 from esys.escript.linearPDEs import LinearPDE, SolverOptions
45 import numpy
46 OPTIMIZE=True
47 SOLVER_VERBOSE=False
48 # setNumberOfThreads(2)
49
50 try:
51 DUDLEY_TEST_DATA=os.environ['DUDLEY_TEST_DATA']
52 except KeyError:
53 DUDLEY_TEST_DATA='.'
54
55 DUDLEY_TEST_MESH_PATH=os.path.join(DUDLEY_TEST_DATA,"data_meshes")
56
57 # number of elements in the spatial directions
58 NE0=8
59 NE1=10
60 NE2=12
61
62 NE0=12
63 NE1=13
64 NE2=8
65
66 SOLVER_TOL=1.e-8
67 REL_TOL=1.e-6
68
69 FAC_DIAG=1.
70 FAC_OFFDIAG=-0.4
71
72
73 class Test_SimpleSolve_Rectangle_Order1_SinglePDE_Paso_BICGSTAB_Jacobi(unittest.TestCase):
74 def test_solve(self):
75 # Tell about how many MPI CPUs and OpenMP threads
76 domain=Rectangle(NE0,NE1,1, optimize=OPTIMIZE)
77 x=Solution(domain).getX()
78 # --- set exact solution ----
79 u_ex=Scalar(0,Solution(domain))
80 u_ex=1.+2.*x[0]+3.*x[1]
81 # --- set exact gradient -----------
82 g_ex=Data(0.,(2,),Solution(domain))
83 g_ex[0]=2.
84 g_ex[1]=3.
85 # -------- test gradient --------------------------------
86 g=grad(u_ex)
87 self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
88 # -------- set-up PDE -----------------------------------
89 pde=LinearPDE(domain,numEquations=1)
90 mask=whereZero(x[0])
91 pde.setValue(r=u_ex,q=mask)
92 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
93 # -------- get the solution ---------------------------
94 pde.getSolverOptions().setTolerance(SOLVER_TOL)
95 pde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
96 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
97 pde.getSolverOptions().setPackage(SolverOptions.PASO)
98 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
99 u=pde.getSolution()
100 # -------- test the solution ---------------------------
101 error=Lsup(u-u_ex)
102 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
103
104 class Test_SimpleSolve_Rectangle_Order1_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
105 def test_solve(self):
106 domain=Rectangle(NE0,NE1,1, optimize=OPTIMIZE)
107 x=Solution(domain).getX()
108 # --- set exact solution ----
109 u_ex=Scalar(0,Solution(domain))
110 u_ex=1.+2.*x[0]+3.*x[1]
111 # --- set exact gradient -----------
112 g_ex=Data(0.,(2,),Solution(domain))
113 g_ex[0]=2.
114 g_ex[1]=3.
115 # -------- test gradient --------------------------------
116 g=grad(u_ex)
117 self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
118 # -------- set-up PDE -----------------------------------
119 pde=LinearPDE(domain,numEquations=1)
120 mask=whereZero(x[0])
121 pde.setValue(r=u_ex,q=mask)
122 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
123 # -------- get the solution ---------------------------
124 pde.getSolverOptions().setTolerance(SOLVER_TOL)
125 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
126 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
127 pde.getSolverOptions().setPackage(SolverOptions.PASO)
128 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
129 u=pde.getSolution()
130 # -------- test the solution ---------------------------
131 error=Lsup(u-u_ex)
132 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
133
134 class Test_SimpleSolve_Rectangle_Order1_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
135 def test_solve(self):
136 domain=Rectangle(NE0,NE1,1,optimize=OPTIMIZE)
137 x=Solution(domain).getX()
138 # --- set exact solution ----
139 u_ex=Vector(0,Solution(domain))
140 u_ex[0]=1.+2.*x[0]+3.*x[1]
141 u_ex[1]=-1.+3.*x[0]+2.*x[1]
142 # --- set exact gradient -----------
143 g_ex=Data(0.,(2,2),Solution(domain))
144 g_ex[0,0]=2.
145 g_ex[0,1]=3.
146 g_ex[1,0]=3.
147 g_ex[1,1]=2.
148 # -------- test gradient --------------------------------
149 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
150 # -------- set-up PDE -----------------------------------
151 pde=LinearPDE(domain,numEquations=2)
152 mask=whereZero(x[0])
153 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
154 A=Tensor4(0,Function(domain))
155 A[0,:,0,:]=kronecker(2)
156 A[1,:,1,:]=kronecker(2)
157 Y=Vector(0.,Function(domain))
158 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
159 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
160 pde.setValue(A=A,
161 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
162 Y=Y,
163 y=matrixmult(g_ex,domain.getNormal()))
164 # -------- get the solution ---------------------------
165 pde.getSolverOptions().setTolerance(SOLVER_TOL)
166 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
167 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
168 pde.getSolverOptions().setPackage(SolverOptions.PASO)
169 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
170 u=pde.getSolution()
171 # -------- test the solution ---------------------------
172 error=Lsup(u-u_ex)
173 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
174
175 @unittest.skip("Test should be checked")
176 class Test_SimpleSolve_Rectangle_Order2_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
177 def test_solve(self):
178 domain=Rectangle(NE0,NE1,2,l0=1.,l1=1,optimize=OPTIMIZE)
179 x=Solution(domain).getX()
180 # --- set exact solution ----
181 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
182 # --- set exact gradient -----------
183 g_ex=Data(0.,(2,),Solution(domain))
184 g_ex[0]=2.+8.*x[0]+5.*x[1]
185 g_ex[1]=3.+5.*x[0]+12.*x[1]
186 # -------- test gradient --------------------------------
187 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
188 # -------- set-up PDE -----------------------------------
189 pde=LinearPDE(domain,numEquations=1)
190 mask=whereZero(x[0])
191 pde.setValue(r=u_ex,q=mask)
192 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()),Y=-20.)
193 # -------- get the solution ---------------------------
194 pde.getSolverOptions().setTolerance(SOLVER_TOL)
195 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
196 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
197 pde.getSolverOptions().setPackage(SolverOptions.PASO)
198 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
199 u=pde.getSolution()
200 # -------- test the solution ---------------------------
201 error=Lsup(u-u_ex)
202 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
203
204 @unittest.skip("Test should be checked")
205 class Test_SimpleSolve_Rectangle_Order2_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
206 def test_solve(self):
207 domain=Rectangle(NE0,NE1,2,optimize=OPTIMIZE)
208 x=Solution(domain).getX()
209 # --- set exact solution ----
210 u_ex=Vector(0,Solution(domain))
211 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
212 u_ex[1]=-1.+4.*x[0]+2.*x[1]+1.*x[0]**2+6.*x[1]*x[0]+4.*x[1]**2
213 # --- set exact gradient -----------
214 g_ex=Data(0.,(2,2),Solution(domain))
215 g_ex[0,0]=2.+8.*x[0]+5.*x[1]
216 g_ex[0,1]=3.+5.*x[0]+12.*x[1]
217 g_ex[1,0]=4.+2.*x[0]+6.*x[1]
218 g_ex[1,1]=2.+6.*x[0]+8.*x[1]
219 # -------- test gradient --------------------------------
220 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
221 # -------- set-up PDE -----------------------------------
222 pde=LinearPDE(domain,numEquations=2)
223 mask=whereZero(x[0])
224 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
225 A=Tensor4(0,Function(domain))
226 A[0,:,0,:]=kronecker(2)
227 A[1,:,1,:]=kronecker(2)
228 Y=Vector(0.,Function(domain))
229 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
230 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
231 pde.setValue(A=A,
232 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
233 Y=Y-[20.,10.],
234 y=matrixmult(g_ex,domain.getNormal()))
235 # -------- get the solution ---------------------------
236 pde.getSolverOptions().setTolerance(SOLVER_TOL)
237 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
238 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
239 pde.getSolverOptions().setPackage(SolverOptions.PASO)
240 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
241 u=pde.getSolution()
242 # -------- test the solution ---------------------------
243 error=Lsup(u-u_ex)
244 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
245
246 class Test_SimpleSolve_Brick_Order1_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
247 def test_solve(self):
248 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
249 x=Solution(domain).getX()
250 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
251 # --- set exact gradient -----------
252 g_ex=Data(0.,(3,),Solution(domain))
253 g_ex[0]=2.
254 g_ex[1]=3.
255 g_ex[2]=4.
256 # -------- test gradient --------------------------------
257 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
258 # -------- set-up PDE -----------------------------------
259 pde=LinearPDE(domain,numEquations=1)
260 mask=whereZero(x[0])
261 pde.setValue(r=u_ex,q=mask)
262 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
263 # -------- get the solution ---------------------------
264 pde.getSolverOptions().setTolerance(SOLVER_TOL)
265 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
266 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
267 pde.getSolverOptions().setPackage(SolverOptions.PASO)
268 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
269 u=pde.getSolution()
270 # -------- test the solution ---------------------------
271 error=Lsup(u-u_ex)
272 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
273
274 class Test_SimpleSolve_Brick_Order1_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
275 def test_solve(self):
276 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
277 x=Solution(domain).getX()
278 # --- set exact solution ----
279 u_ex=Vector(0,Solution(domain))
280 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
281 u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
282 u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
283 # --- set exact gradient -----------
284 g_ex=Data(0.,(3,3),Solution(domain))
285 g_ex[0,0]=2.
286 g_ex[0,1]=3.
287 g_ex[0,2]=4.
288 g_ex[1,0]=4.
289 g_ex[1,1]=1.
290 g_ex[1,2]=-2.
291 g_ex[2,0]=8.
292 g_ex[2,1]=4.
293 g_ex[2,2]=5.
294 # -------- test gradient --------------------------------
295 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
296 # -------- set-up PDE -----------------------------------
297 pde=LinearPDE(domain,numEquations=3)
298 mask=whereZero(x[0])
299 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
300 A=Tensor4(0,Function(domain))
301 A[0,:,0,:]=kronecker(3)
302 A[1,:,1,:]=kronecker(3)
303 A[2,:,2,:]=kronecker(3)
304 Y=Vector(0.,Function(domain))
305 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
306 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
307 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
308 pde.setValue(A=A,
309 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
310 Y=Y,
311 y=matrixmult(g_ex,domain.getNormal()))
312 # -------- get the solution ---------------------------
313 pde.getSolverOptions().setTolerance(SOLVER_TOL)
314 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
315 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
316 pde.getSolverOptions().setPackage(SolverOptions.PASO)
317 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
318 u=pde.getSolution()
319 # -------- test the solution ---------------------------
320 error=Lsup(u-u_ex)
321 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
322
323 @unittest.skip("Test should be checked")
324 class Test_SimpleSolve_Brick_Order2_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
325 def test_solve(self):
326 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
327 x=Solution(domain).getX()
328 # --- set exact solution ----
329 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
330 # --- set exact gradient -----------
331 g_ex=Data(0.,(3,),Solution(domain))
332 g_ex[0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
333 g_ex[1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
334 g_ex[2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
335 # -------- test gradient --------------------------------
336 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
337 # -------- set-up PDE -----------------------------------
338 pde=LinearPDE(domain,numEquations=1)
339 mask=whereZero(x[0])
340 pde.setValue(r=u_ex,q=mask)
341 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()),Y=-60.)
342 # -------- get the solution ---------------------------
343 pde.getSolverOptions().setTolerance(SOLVER_TOL)
344 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
345 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
346 pde.getSolverOptions().setPackage(SolverOptions.PASO)
347 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
348 u=pde.getSolution()
349 # -------- test the solution ---------------------------
350 error=Lsup(u-u_ex)
351 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
352
353 @unittest.skip("Test should be checked")
354 class Test_SimpleSolve_Brick_Order2_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
355 def test_solve(self):
356 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
357 x=Solution(domain).getX()
358 # --- set exact solution ----
359 u_ex=Vector(0,Solution(domain))
360 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
361 u_ex[1]=2.+4.*x[0]+1.*x[1]-6.*x[2]+3.*x[0]*x[1]+2.*x[1]*x[2]-8.*x[2]*x[0]-2.*x[0]**2+7.*x[1]**2+5.*x[2]**2
362 u_ex[2]=-2.+7.*x[0]+9.*x[1]+2*x[2]-6.*x[0]*x[1]+8.*x[1]*x[2]+2.*x[2]*x[0]+2.*x[0]**2+8.*x[1]**2+1.*x[2]**2
363 # --- set exact gradient -----------
364 g_ex=Data(0.,(3,3),Solution(domain))
365 g_ex[0,0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
366 g_ex[0,1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
367 g_ex[0,2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
368 g_ex[1,0]=4.+3.*x[1]-8.*x[2]-4.*x[0]
369 g_ex[1,1]=1+3.*x[0]+2.*x[2]+14.*x[1]
370 g_ex[1,2]=-6.+2.*x[1]-8.*x[0]+10.*x[2]
371 g_ex[2,0]=7.-6.*x[1]+2.*x[2]+4.*x[0]
372 g_ex[2,1]=9.-6.*x[0]+8.*x[2]+16.*x[1]
373 g_ex[2,2]=2+8.*x[1]+2.*x[0]+2.*x[2]
374 # -------- test gradient --------------------------------
375 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
376 # -------- set-up PDE -----------------------------------
377 pde=LinearPDE(domain,numEquations=3)
378 mask=whereZero(x[0])
379 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
380 Y=Vector(0.,Function(domain))
381 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
382 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
383 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
384 A=Tensor4(0,Function(domain))
385 A[0,:,0,:]=kronecker(3)
386 A[1,:,1,:]=kronecker(3)
387 A[2,:,2,:]=kronecker(3)
388 pde.setValue(A=A,
389 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
390 Y=Y-numpy.array([60.,20.,22.]),
391 y=matrixmult(g_ex,domain.getNormal()))
392 # -------- get the solution ---------------------------
393 pde.getSolverOptions().setTolerance(SOLVER_TOL)
394 pde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
395 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
396 pde.getSolverOptions().setPackage(SolverOptions.PASO)
397 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
398 u=pde.getSolution()
399 # -------- test the solution ---------------------------
400 error=Lsup(u-u_ex)
401 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
402
403 class Test_SimpleSolve_Rectangle_Order1_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
404 def test_solve(self):
405 domain=Rectangle(NE0,NE1,1, optimize=OPTIMIZE)
406 x=Solution(domain).getX()
407 # --- set exact solution ----
408 u_ex=Scalar(0,Solution(domain))
409 u_ex=1.+2.*x[0]+3.*x[1]
410 # --- set exact gradient -----------
411 g_ex=Data(0.,(2,),Solution(domain))
412 g_ex[0]=2.
413 g_ex[1]=3.
414 # -------- test gradient --------------------------------
415 g=grad(u_ex)
416 self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
417 # -------- set-up PDE -----------------------------------
418 pde=LinearPDE(domain,numEquations=1)
419 mask=whereZero(x[0])
420 pde.setValue(r=u_ex,q=mask)
421 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
422 # -------- get the solution ---------------------------
423 pde.getSolverOptions().setTolerance(SOLVER_TOL)
424 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
425 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
426 pde.getSolverOptions().setPackage(SolverOptions.PASO)
427 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
428 u=pde.getSolution()
429 # -------- test the solution ---------------------------
430 error=Lsup(u-u_ex)
431 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
432
433 @unittest.skip("Test should be checked")
434 class Test_SimpleSolve_Rectangle_Order2_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
435 def test_solve(self):
436 domain=Rectangle(NE0,NE1,2,l0=1.,l1=1,optimize=OPTIMIZE)
437 x=Solution(domain).getX()
438 # --- set exact solution ----
439 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
440 # --- set exact gradient -----------
441 g_ex=Data(0.,(2,),Solution(domain))
442 g_ex[0]=2.+8.*x[0]+5.*x[1]
443 g_ex[1]=3.+5.*x[0]+12.*x[1]
444 # -------- test gradient --------------------------------
445 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
446 # -------- set-up PDE -----------------------------------
447 pde=LinearPDE(domain,numEquations=1)
448 mask=whereZero(x[0])
449 pde.setValue(r=u_ex,q=mask)
450 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()),Y=-20.)
451 # -------- get the solution ---------------------------
452 pde.getSolverOptions().setTolerance(SOLVER_TOL)
453 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
454 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
455 pde.getSolverOptions().setPackage(SolverOptions.PASO)
456 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
457 u=pde.getSolution()
458 # -------- test the solution ---------------------------
459 error=Lsup(u-u_ex)
460 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
461
462 class Test_SimpleSolve_Rectangle_Order1_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
463 def test_solve(self):
464 domain=Rectangle(NE0,NE1,1,optimize=OPTIMIZE)
465 x=Solution(domain).getX()
466 # --- set exact solution ----
467 u_ex=Vector(0,Solution(domain))
468 u_ex[0]=1.+2.*x[0]+3.*x[1]
469 u_ex[1]=-1.+3.*x[0]+2.*x[1]
470 # --- set exact gradient -----------
471 g_ex=Data(0.,(2,2),Solution(domain))
472 g_ex[0,0]=2.
473 g_ex[0,1]=3.
474 g_ex[1,0]=3.
475 g_ex[1,1]=2.
476 # -------- test gradient --------------------------------
477 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
478 # -------- set-up PDE -----------------------------------
479 pde=LinearPDE(domain,numEquations=2)
480 mask=whereZero(x[0])
481 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
482 A=Tensor4(0,Function(domain))
483 A[0,:,0,:]=kronecker(2)
484 A[1,:,1,:]=kronecker(2)
485 Y=Vector(0.,Function(domain))
486 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
487 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
488 pde.setValue(A=A,
489 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
490 Y=Y,
491 y=matrixmult(g_ex,domain.getNormal()))
492 # -------- get the solution ---------------------------
493 pde.getSolverOptions().setTolerance(SOLVER_TOL)
494 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
495 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
496 pde.getSolverOptions().setPackage(SolverOptions.PASO)
497 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
498 u=pde.getSolution()
499 # -------- test the solution ---------------------------
500 error=Lsup(u-u_ex)
501 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
502
503 @unittest.skip("Test should be checked")
504 class Test_SimpleSolve_Rectangle_Order2_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
505 def test_solve(self):
506 domain=Rectangle(NE0,NE1,2,optimize=OPTIMIZE)
507 x=Solution(domain).getX()
508 # --- set exact solution ----
509 u_ex=Vector(0,Solution(domain))
510 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
511 u_ex[1]=-1.+4.*x[0]+2.*x[1]+1.*x[0]**2+6.*x[1]*x[0]+4.*x[1]**2
512 # --- set exact gradient -----------
513 g_ex=Data(0.,(2,2),Solution(domain))
514 g_ex[0,0]=2.+8.*x[0]+5.*x[1]
515 g_ex[0,1]=3.+5.*x[0]+12.*x[1]
516 g_ex[1,0]=4.+2.*x[0]+6.*x[1]
517 g_ex[1,1]=2.+6.*x[0]+8.*x[1]
518 # -------- test gradient --------------------------------
519 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
520 # -------- set-up PDE -----------------------------------
521 pde=LinearPDE(domain,numEquations=2)
522 mask=whereZero(x[0])
523 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
524 A=Tensor4(0,Function(domain))
525 A[0,:,0,:]=kronecker(2)
526 A[1,:,1,:]=kronecker(2)
527 Y=Vector(0.,Function(domain))
528 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
529 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
530 pde.setValue(A=A,
531 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
532 Y=Y-[20.,10.],
533 y=matrixmult(g_ex,domain.getNormal()))
534 # -------- get the solution ---------------------------
535 pde.getSolverOptions().setTolerance(SOLVER_TOL)
536 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
537 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
538 pde.getSolverOptions().setPackage(SolverOptions.PASO)
539 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
540 u=pde.getSolution()
541 # -------- test the solution ---------------------------
542 error=Lsup(u-u_ex)
543 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
544
545 class Test_SimpleSolve_Brick_Order1_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
546 def test_solve(self):
547 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
548 x=Solution(domain).getX()
549 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
550 # --- set exact gradient -----------
551 g_ex=Data(0.,(3,),Solution(domain))
552 g_ex[0]=2.
553 g_ex[1]=3.
554 g_ex[2]=4.
555 # -------- test gradient --------------------------------
556 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
557 # -------- set-up PDE -----------------------------------
558 pde=LinearPDE(domain,numEquations=1)
559 mask=whereZero(x[0])
560 pde.setValue(r=u_ex,q=mask)
561 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
562 # -------- get the solution ---------------------------
563 pde.getSolverOptions().setTolerance(SOLVER_TOL)
564 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
565 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
566 pde.getSolverOptions().setPackage(SolverOptions.PASO)
567 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
568 u=pde.getSolution()
569 # -------- test the solution ---------------------------
570 error=Lsup(u-u_ex)
571 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
572
573 class Test_SimpleSolve_Brick_Order1_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
574 def test_solve(self):
575 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
576 x=Solution(domain).getX()
577 # --- set exact solution ----
578 u_ex=Vector(0,Solution(domain))
579 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
580 u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
581 u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
582 # --- set exact gradient -----------
583 g_ex=Data(0.,(3,3),Solution(domain))
584 g_ex[0,0]=2.
585 g_ex[0,1]=3.
586 g_ex[0,2]=4.
587 g_ex[1,0]=4.
588 g_ex[1,1]=1.
589 g_ex[1,2]=-2.
590 g_ex[2,0]=8.
591 g_ex[2,1]=4.
592 g_ex[2,2]=5.
593 # -------- test gradient --------------------------------
594 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
595 # -------- set-up PDE -----------------------------------
596 pde=LinearPDE(domain,numEquations=3)
597 mask=whereZero(x[0])
598 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
599 A=Tensor4(0,Function(domain))
600 A[0,:,0,:]=kronecker(3)
601 A[1,:,1,:]=kronecker(3)
602 A[2,:,2,:]=kronecker(3)
603 Y=Vector(0.,Function(domain))
604 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
605 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
606 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
607 pde.setValue(A=A,
608 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
609 Y=Y,
610 y=matrixmult(g_ex,domain.getNormal()))
611 # -------- get the solution ---------------------------
612 pde.getSolverOptions().setTolerance(SOLVER_TOL)
613 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
614 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
615 pde.getSolverOptions().setPackage(SolverOptions.PASO)
616 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
617 u=pde.getSolution()
618 # -------- test the solution ---------------------------
619 error=Lsup(u-u_ex)
620 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
621
622 @unittest.skip("Test should be checked")
623 class Test_SimpleSolve_Brick_Order2_SinglePDE_Paso_TFQMR_Jacobi(unittest.TestCase):
624 def test_solve(self):
625 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
626 x=Solution(domain).getX()
627 # --- set exact solution ----
628 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
629 # --- set exact gradient -----------
630 g_ex=Data(0.,(3,),Solution(domain))
631 g_ex[0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
632 g_ex[1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
633 g_ex[2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
634 # -------- test gradient --------------------------------
635 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
636 # -------- set-up PDE -----------------------------------
637 pde=LinearPDE(domain,numEquations=1)
638 mask=whereZero(x[0])
639 pde.setValue(r=u_ex,q=mask)
640 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()),Y=-60.)
641 # -------- get the solution ---------------------------
642 pde.getSolverOptions().setTolerance(SOLVER_TOL)
643 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
644 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
645 pde.getSolverOptions().setPackage(SolverOptions.PASO)
646 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
647 u=pde.getSolution()
648 # -------- test the solution ---------------------------
649 error=Lsup(u-u_ex)
650 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
651
652 @unittest.skip("Test should be checked")
653 class Test_SimpleSolve_Brick_Order2_SystemPDE_Paso_TFQMR_Jacobi(unittest.TestCase):
654 def test_solve(self):
655 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
656 x=Solution(domain).getX()
657 # --- set exact solution ----
658 u_ex=Vector(0,Solution(domain))
659 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
660 u_ex[1]=2.+4.*x[0]+1.*x[1]-6.*x[2]+3.*x[0]*x[1]+2.*x[1]*x[2]-8.*x[2]*x[0]-2.*x[0]**2+7.*x[1]**2+5.*x[2]**2
661 u_ex[2]=-2.+7.*x[0]+9.*x[1]+2*x[2]-6.*x[0]*x[1]+8.*x[1]*x[2]+2.*x[2]*x[0]+2.*x[0]**2+8.*x[1]**2+1.*x[2]**2
662 # --- set exact gradient -----------
663 g_ex=Data(0.,(3,3),Solution(domain))
664 g_ex[0,0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
665 g_ex[0,1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
666 g_ex[0,2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
667 g_ex[1,0]=4.+3.*x[1]-8.*x[2]-4.*x[0]
668 g_ex[1,1]=1+3.*x[0]+2.*x[2]+14.*x[1]
669 g_ex[1,2]=-6.+2.*x[1]-8.*x[0]+10.*x[2]
670 g_ex[2,0]=7.-6.*x[1]+2.*x[2]+4.*x[0]
671 g_ex[2,1]=9.-6.*x[0]+8.*x[2]+16.*x[1]
672 g_ex[2,2]=2+8.*x[1]+2.*x[0]+2.*x[2]
673 # -------- test gradient --------------------------------
674 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
675 # -------- set-up PDE -----------------------------------
676 pde=LinearPDE(domain,numEquations=3)
677 mask=whereZero(x[0])
678 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
679 Y=Vector(0.,Function(domain))
680 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
681 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
682 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
683 A=Tensor4(0,Function(domain))
684 A[0,:,0,:]=kronecker(3)
685 A[1,:,1,:]=kronecker(3)
686 A[2,:,2,:]=kronecker(3)
687 pde.setValue(A=A,
688 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
689 Y=Y-numpy.array([60.,20.,22.]),
690 y=matrixmult(g_ex,domain.getNormal()))
691 # -------- get the solution ---------------------------
692 pde.getSolverOptions().setTolerance(SOLVER_TOL)
693 pde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
694 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
695 pde.getSolverOptions().setPackage(SolverOptions.PASO)
696 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
697 u=pde.getSolution()
698 # -------- test the solution ---------------------------
699 error=Lsup(u-u_ex)
700 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
701
702
703 class Test_SimpleSolve_Rectangle_Order1_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
704 def test_solve(self):
705 domain=Rectangle(NE0,NE1,1, optimize=OPTIMIZE)
706 x=Solution(domain).getX()
707 # --- set exact solution ----
708 u_ex=Scalar(0,Solution(domain))
709 u_ex=1.+2.*x[0]+3.*x[1]
710 # --- set exact gradient -----------
711 g_ex=Data(0.,(2,),Solution(domain))
712 g_ex[0]=2.
713 g_ex[1]=3.
714 # -------- test gradient --------------------------------
715 g=grad(u_ex)
716 self.assertTrue(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
717 # -------- set-up PDE -----------------------------------
718 pde=LinearPDE(domain,numEquations=1)
719 mask=whereZero(x[0])
720 pde.setValue(r=u_ex,q=mask)
721 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
722 # -------- get the solution ---------------------------
723 pde.getSolverOptions().setTolerance(SOLVER_TOL)
724 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
725 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
726 pde.getSolverOptions().setPackage(SolverOptions.PASO)
727
728 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
729 u=pde.getSolution()
730 # -------- test the solution ---------------------------
731 error=Lsup(u-u_ex)
732 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
733
734 @unittest.skip("Test should be checked")
735 class Test_SimpleSolve_Rectangle_Order2_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
736 def test_solve(self):
737 domain=Rectangle(NE0,NE1,2,l0=1.,l1=1,optimize=OPTIMIZE)
738 x=Solution(domain).getX()
739 # --- set exact solution ----
740 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
741 # --- set exact gradient -----------
742 g_ex=Data(0.,(2,),Solution(domain))
743 g_ex[0]=2.+8.*x[0]+5.*x[1]
744 g_ex[1]=3.+5.*x[0]+12.*x[1]
745 # -------- test gradient --------------------------------
746 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
747 # -------- set-up PDE -----------------------------------
748 pde=LinearPDE(domain,numEquations=1)
749 mask=whereZero(x[0])
750 pde.setValue(r=u_ex,q=mask)
751 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()),Y=-20.)
752 # -------- get the solution ---------------------------
753 pde.getSolverOptions().setTolerance(SOLVER_TOL)
754 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
755 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
756 pde.getSolverOptions().setPackage(SolverOptions.PASO)
757 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
758 u=pde.getSolution()
759 # -------- test the solution ---------------------------
760 error=Lsup(u-u_ex)
761 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
762
763 class Test_SimpleSolve_Rectangle_Order1_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
764 def test_solve(self):
765 domain=Rectangle(NE0,NE1,1,optimize=OPTIMIZE)
766 x=Solution(domain).getX()
767 # --- set exact solution ----
768 u_ex=Vector(0,Solution(domain))
769 u_ex[0]=1.+2.*x[0]+3.*x[1]
770 u_ex[1]=-1.+3.*x[0]+2.*x[1]
771 # --- set exact gradient -----------
772 g_ex=Data(0.,(2,2),Solution(domain))
773 g_ex[0,0]=2.
774 g_ex[0,1]=3.
775 g_ex[1,0]=3.
776 g_ex[1,1]=2.
777 # -------- test gradient --------------------------------
778 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
779 # -------- set-up PDE -----------------------------------
780 pde=LinearPDE(domain,numEquations=2)
781 mask=whereZero(x[0])
782 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
783 A=Tensor4(0,Function(domain))
784 A[0,:,0,:]=kronecker(2)
785 A[1,:,1,:]=kronecker(2)
786 Y=Vector(0.,Function(domain))
787 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
788 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
789 pde.setValue(A=A,
790 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
791 Y=Y,
792 y=matrixmult(g_ex,domain.getNormal()))
793 # -------- get the solution ---------------------------
794 pde.getSolverOptions().setTolerance(SOLVER_TOL)
795 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
796 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
797 pde.getSolverOptions().setPackage(SolverOptions.PASO)
798 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
799 u=pde.getSolution()
800 # -------- test the solution ---------------------------
801 error=Lsup(u-u_ex)
802 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
803
804 @unittest.skip("Test should be checked")
805 class Test_SimpleSolve_Rectangle_Order2_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
806 def test_solve(self):
807 domain=Rectangle(NE0,NE1,2,optimize=OPTIMIZE)
808 x=Solution(domain).getX()
809 # --- set exact solution ----
810 u_ex=Vector(0,Solution(domain))
811 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
812 u_ex[1]=-1.+4.*x[0]+2.*x[1]+1.*x[0]**2+6.*x[1]*x[0]+4.*x[1]**2
813 # --- set exact gradient -----------
814 g_ex=Data(0.,(2,2),Solution(domain))
815 g_ex[0,0]=2.+8.*x[0]+5.*x[1]
816 g_ex[0,1]=3.+5.*x[0]+12.*x[1]
817 g_ex[1,0]=4.+2.*x[0]+6.*x[1]
818 g_ex[1,1]=2.+6.*x[0]+8.*x[1]
819 # -------- test gradient --------------------------------
820 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
821 # -------- set-up PDE -----------------------------------
822 pde=LinearPDE(domain,numEquations=2)
823 mask=whereZero(x[0])
824 pde.setValue(r=u_ex,q=mask*numpy.ones(2,))
825 A=Tensor4(0,Function(domain))
826 A[0,:,0,:]=kronecker(2)
827 A[1,:,1,:]=kronecker(2)
828 Y=Vector(0.,Function(domain))
829 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
830 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
831 pde.setValue(A=A,
832 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((2,2))*FAC_OFFDIAG,
833 Y=Y-[20.,10.],
834 y=matrixmult(g_ex,domain.getNormal()))
835 # -------- get the solution ---------------------------
836 pde.getSolverOptions().setTolerance(SOLVER_TOL)
837 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
838 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
839 pde.getSolverOptions().setPackage(SolverOptions.PASO)
840 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
841 u=pde.getSolution()
842 # -------- test the solution ---------------------------
843 error=Lsup(u-u_ex)
844 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
845
846 class Test_SimpleSolve_Brick_Order1_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
847 def test_solve(self):
848 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
849 x=Solution(domain).getX()
850 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
851 # --- set exact gradient -----------
852 g_ex=Data(0.,(3,),Solution(domain))
853 g_ex[0]=2.
854 g_ex[1]=3.
855 g_ex[2]=4.
856 # -------- test gradient --------------------------------
857 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
858 # -------- set-up PDE -----------------------------------
859 pde=LinearPDE(domain,numEquations=1)
860 mask=whereZero(x[0])
861 pde.setValue(r=u_ex,q=mask)
862 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
863 # -------- get the solution ---------------------------
864 pde.getSolverOptions().setTolerance(SOLVER_TOL)
865 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
866 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
867 pde.getSolverOptions().setPackage(SolverOptions.PASO)
868 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
869 u=pde.getSolution()
870 # -------- test the solution ---------------------------
871 error=Lsup(u-u_ex)
872 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
873
874 class Test_SimpleSolve_Brick_Order1_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
875 def test_solve(self):
876 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
877 x=Solution(domain).getX()
878 # --- set exact solution ----
879 u_ex=Vector(0,Solution(domain))
880 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
881 u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
882 u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
883 # --- set exact gradient -----------
884 g_ex=Data(0.,(3,3),Solution(domain))
885 g_ex[0,0]=2.
886 g_ex[0,1]=3.
887 g_ex[0,2]=4.
888 g_ex[1,0]=4.
889 g_ex[1,1]=1.
890 g_ex[1,2]=-2.
891 g_ex[2,0]=8.
892 g_ex[2,1]=4.
893 g_ex[2,2]=5.
894 # -------- test gradient --------------------------------
895 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
896 # -------- set-up PDE -----------------------------------
897 pde=LinearPDE(domain,numEquations=3)
898 mask=whereZero(x[0])
899 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
900 A=Tensor4(0,Function(domain))
901 A[0,:,0,:]=kronecker(3)
902 A[1,:,1,:]=kronecker(3)
903 A[2,:,2,:]=kronecker(3)
904 Y=Vector(0.,Function(domain))
905 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
906 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
907 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
908 pde.setValue(A=A,
909 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
910 Y=Y,
911 y=matrixmult(g_ex,domain.getNormal()))
912 # -------- get the solution ---------------------------
913 pde.getSolverOptions().setTolerance(SOLVER_TOL)
914 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
915 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
916 pde.getSolverOptions().setPackage(SolverOptions.PASO)
917 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
918 u=pde.getSolution()
919 # -------- test the solution ---------------------------
920 error=Lsup(u-u_ex)
921 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
922
923 @unittest.skip("Test should be checked")
924 class Test_SimpleSolve_Brick_Order2_SinglePDE_Paso_MINRES_Jacobi(unittest.TestCase):
925 def test_solve(self):
926 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
927 x=Solution(domain).getX()
928 # --- set exact solution ----
929 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
930 # --- set exact gradient -----------
931 g_ex=Data(0.,(3,),Solution(domain))
932 g_ex[0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
933 g_ex[1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
934 g_ex[2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
935 # -------- test gradient --------------------------------
936 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
937 # -------- set-up PDE -----------------------------------
938 pde=LinearPDE(domain,numEquations=1)
939 mask=whereZero(x[0])
940 pde.setValue(r=u_ex,q=mask)
941 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()),Y=-60.)
942 # -------- get the solution ---------------------------
943 pde.getSolverOptions().setTolerance(SOLVER_TOL)
944 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
945 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
946 pde.getSolverOptions().setPackage(SolverOptions.PASO)
947 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
948 u=pde.getSolution()
949 # -------- test the solution ---------------------------
950 error=Lsup(u-u_ex)
951 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
952
953 @unittest.skip("Test should be checked")
954 class Test_SimpleSolve_Brick_Order2_SystemPDE_Paso_MINRES_Jacobi(unittest.TestCase):
955 def test_solve(self):
956 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
957 x=Solution(domain).getX()
958 # --- set exact solution ----
959 u_ex=Vector(0,Solution(domain))
960 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]+6.*x[0]*x[1]+7.*x[1]*x[2]+8.*x[2]*x[0]+9.*x[0]**2+10.*x[1]**2+11.*x[2]**2
961 u_ex[1]=2.+4.*x[0]+1.*x[1]-6.*x[2]+3.*x[0]*x[1]+2.*x[1]*x[2]-8.*x[2]*x[0]-2.*x[0]**2+7.*x[1]**2+5.*x[2]**2
962 u_ex[2]=-2.+7.*x[0]+9.*x[1]+2*x[2]-6.*x[0]*x[1]+8.*x[1]*x[2]+2.*x[2]*x[0]+2.*x[0]**2+8.*x[1]**2+1.*x[2]**2
963 # --- set exact gradient -----------
964 g_ex=Data(0.,(3,3),Solution(domain))
965 g_ex[0,0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
966 g_ex[0,1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
967 g_ex[0,2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
968 g_ex[1,0]=4.+3.*x[1]-8.*x[2]-4.*x[0]
969 g_ex[1,1]=1+3.*x[0]+2.*x[2]+14.*x[1]
970 g_ex[1,2]=-6.+2.*x[1]-8.*x[0]+10.*x[2]
971 g_ex[2,0]=7.-6.*x[1]+2.*x[2]+4.*x[0]
972 g_ex[2,1]=9.-6.*x[0]+8.*x[2]+16.*x[1]
973 g_ex[2,2]=2+8.*x[1]+2.*x[0]+2.*x[2]
974 # -------- test gradient --------------------------------
975 self.assertTrue(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
976 # -------- set-up PDE -----------------------------------
977 pde=LinearPDE(domain,numEquations=3)
978 mask=whereZero(x[0])
979 pde.setValue(r=u_ex,q=mask*numpy.ones(3,))
980 Y=Vector(0.,Function(domain))
981 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
982 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
983 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
984 A=Tensor4(0,Function(domain))
985 A[0,:,0,:]=kronecker(3)
986 A[1,:,1,:]=kronecker(3)
987 A[2,:,2,:]=kronecker(3)
988 pde.setValue(A=A,
989 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numpy.ones((3,3))*FAC_OFFDIAG,
990 Y=Y-numpy.array([60.,20.,22.]),
991 y=matrixmult(g_ex,domain.getNormal()))
992 # -------- get the solution ---------------------------
993 pde.getSolverOptions().setTolerance(SOLVER_TOL)
994 pde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
995 pde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
996 pde.getSolverOptions().setPackage(SolverOptions.PASO)
997 pde.getSolverOptions().setVerbosity(SOLVER_VERBOSE)
998 u=pde.getSolution()
999 # -------- test the solution ---------------------------
1000 error=Lsup(u-u_ex)
1001 self.assertTrue(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
1002
1003 if __name__ == '__main__':
1004 run_tests(__name__, exit_on_failure=True)

  ViewVC Help
Powered by ViewVC 1.1.26