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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1374 - (show annotations)
Tue Jan 8 09:37:55 2008 UTC (11 years, 5 months ago) by gross
File MIME type: text/x-python
File size: 15696 byte(s)
some changes to get things going on the cognac.ivec.org.
1 #
2 # $Id$
3 #
4 #######################################################
5 #
6 # Copyright 2003-2007 by ACceSS MNRF
7 # Copyright 2007 by University of Queensland
8 #
9 # http://esscc.uq.edu.au
10 # Primary Business: Queensland, Australia
11 # Licensed under the Open Software License version 3.0
12 # http://www.opensource.org/licenses/osl-3.0.php
13 #
14 #######################################################
15 #
16
17 """
18 Test suite for the linearPDE and pdetools test on finley
19
20 @remark:
21
22 @var __author__: name of author
23 @var __licence__: licence agreement
24 @var __url__: url entry point on documentation
25 @var __version__: version
26 @var __date__: date of the version
27 """
28
29 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
30 http://www.access.edu.au
31 Primary Business: Queensland, Australia"""
32 __license__="""Licensed under the Open Software License version 3.0
33 http://www.opensource.org/licenses/osl-3.0.php"""
34 __author__="Lutz Gross, l.gross@uq.edu.au"
35 __url__="http://www.iservo.edu.au/esys/escript"
36 __version__="$Revision: 859 $"
37 __date__="$Date: 2006-09-26 12:19:18 +1000 (Tue, 26 Sep 2006) $"
38
39 import unittest
40
41 from esys.escript import *
42 from esys.finley import Rectangle,Brick
43 from esys.escript.linearPDEs import LinearPDE
44 OPTIMIZE=False
45 SOLVER_VERBOSE=False
46
47 try:
48 FINLEY_TEST_DATA=os.environ['FINLEY_TEST_DATA']
49 except KeyError:
50 FINLEY_TEST_DATA='.'
51
52 FINLEY_TEST_MESH_PATH=FINLEY_TEST_DATA+"/data_meshes/"
53
54 # number of elements in the spatial directions
55 NE0=8
56 NE1=10
57 NE2=12
58
59 NE0=12
60 NE1=12
61 NE2=8
62
63 SOLVER_TOL=1.e-8
64 REL_TOL=1.e-6
65
66 FAC_DIAG=1.
67 FAC_OFFDIAG=-0.4
68
69 class SimpleSolve_Rectangle_Order1_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
70 def test_solve(self):
71 domain=Rectangle(NE0,NE1,1, optimize=OPTIMIZE)
72 x=Solution(domain).getX()
73 # --- set exact solution ----
74 u_ex=Scalar(0,Solution(domain))
75 u_ex=1.+2.*x[0]+3.*x[1]
76 # --- set exact gradient -----------
77 g_ex=Data(0.,(2,),Solution(domain))
78 g_ex[0]=2.
79 g_ex[1]=3.
80 # -------- test gradient --------------------------------
81 g=grad(u_ex)
82 self.failUnless(Lsup(g_ex-g)<REL_TOL*Lsup(g_ex))
83 # -------- set-up PDE -----------------------------------
84 pde=LinearPDE(domain,numEquations=1)
85 mask=whereZero(x[0])
86 pde.setValue(r=u_ex,q=mask)
87 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()))
88 # -------- get the solution ---------------------------
89 pde.setTolerance(SOLVER_TOL)
90 pde.setSolverMethod(pde.PCG,pde.JACOBI)
91 pde.setSolverPackage(pde.PASO)
92 u=pde.getSolution(verbose=SOLVER_VERBOSE)
93 # -------- test the solution ---------------------------
94 error=Lsup(u-u_ex)/Lsup(u_ex)
95 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
96 class SimpleSolve_Rectangle_Order1_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
97 def test_solve(self):
98 domain=Rectangle(NE0,NE1,1,optimize=OPTIMIZE)
99 x=Solution(domain).getX()
100 # --- set exact solution ----
101 u_ex=Vector(0,Solution(domain))
102 u_ex[0]=1.+2.*x[0]+3.*x[1]
103 u_ex[1]=-1.+3.*x[0]+2.*x[1]
104 # --- set exact gradient -----------
105 g_ex=Data(0.,(2,2),Solution(domain))
106 g_ex[0,0]=2.
107 g_ex[0,1]=3.
108 g_ex[1,0]=3.
109 g_ex[1,1]=2.
110 # -------- test gradient --------------------------------
111 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
112 # -------- set-up PDE -----------------------------------
113 pde=LinearPDE(domain,numEquations=2)
114 mask=whereZero(x[0])
115 pde.setValue(r=u_ex,q=mask*numarray.ones(2,))
116 A=Tensor4(0,Function(domain))
117 A[0,:,0,:]=kronecker(2)
118 A[1,:,1,:]=kronecker(2)
119 Y=Vector(0.,Function(domain))
120 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
121 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
122 pde.setValue(A=A,
123 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numarray.ones((2,2))*FAC_OFFDIAG,
124 Y=Y,
125 y=matrixmult(g_ex,domain.getNormal()))
126 # -------- get the solution ---------------------------
127 pde.setTolerance(SOLVER_TOL)
128 pde.setSolverMethod(pde.PCG,pde.JACOBI)
129 pde.setSolverPackage(pde.PASO)
130 u=pde.getSolution(verbose=SOLVER_VERBOSE)
131 # -------- test the solution ---------------------------
132 error=Lsup(u-u_ex)/Lsup(u_ex)
133 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
134 class SimpleSolve_Rectangle_Order2_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
135 def test_solve(self):
136 domain=Rectangle(NE0,NE1,2,l0=1.,l1=1,optimize=OPTIMIZE)
137 x=Solution(domain).getX()
138 # --- set exact solution ----
139 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
140 # --- set exact gradient -----------
141 g_ex=Data(0.,(2,),Solution(domain))
142 g_ex[0]=2.+8.*x[0]+5.*x[1]
143 g_ex[1]=3.+5.*x[0]+12.*x[1]
144 # -------- test gradient --------------------------------
145 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
146 # -------- set-up PDE -----------------------------------
147 pde=LinearPDE(domain,numEquations=1)
148 mask=whereZero(x[0])
149 pde.setValue(r=u_ex,q=mask)
150 pde.setValue(A=kronecker(2),y=inner(g_ex,domain.getNormal()),Y=-20.)
151 # -------- get the solution ---------------------------
152 pde.setTolerance(SOLVER_TOL)
153 pde.setSolverMethod(pde.PCG,pde.JACOBI)
154 pde.setSolverPackage(pde.PASO)
155 u=pde.getSolution(verbose=SOLVER_VERBOSE)
156 # -------- test the solution ---------------------------
157 error=Lsup(u-u_ex)/Lsup(u_ex)
158 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
159 class SimpleSolve_Rectangle_Order2_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
160 def test_solve(self):
161 domain=Rectangle(NE0,NE1,2,optimize=OPTIMIZE)
162 x=Solution(domain).getX()
163 # --- set exact solution ----
164 u_ex=Vector(0,Solution(domain))
165 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[0]**2+5.*x[1]*x[0]+6.*x[1]**2
166 u_ex[1]=-1.+4.*x[0]+2.*x[1]+1.*x[0]**2+6.*x[1]*x[0]+4.*x[1]**2
167 # --- set exact gradient -----------
168 g_ex=Data(0.,(2,2),Solution(domain))
169 g_ex[0,0]=2.+8.*x[0]+5.*x[1]
170 g_ex[0,1]=3.+5.*x[0]+12.*x[1]
171 g_ex[1,0]=4.+2.*x[0]+6.*x[1]
172 g_ex[1,1]=2.+6.*x[0]+8.*x[1]
173 # -------- test gradient --------------------------------
174 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
175 # -------- set-up PDE -----------------------------------
176 pde=LinearPDE(domain,numEquations=2)
177 mask=whereZero(x[0])
178 pde.setValue(r=u_ex,q=mask*numarray.ones(2,))
179 A=Tensor4(0,Function(domain))
180 A[0,:,0,:]=kronecker(2)
181 A[1,:,1,:]=kronecker(2)
182 Y=Vector(0.,Function(domain))
183 Y[0]=u_ex[0]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG
184 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG
185 pde.setValue(A=A,
186 D=kronecker(2)*(FAC_DIAG-FAC_OFFDIAG)+numarray.ones((2,2))*FAC_OFFDIAG,
187 Y=Y-[20.,10.],
188 y=matrixmult(g_ex,domain.getNormal()))
189 # -------- get the solution ---------------------------
190 pde.setTolerance(SOLVER_TOL)
191 pde.setSolverMethod(pde.PCG,pde.JACOBI)
192 pde.setSolverPackage(pde.PASO)
193 u=pde.getSolution(verbose=SOLVER_VERBOSE)
194 # -------- test the solution ---------------------------
195 error=Lsup(u-u_ex)/Lsup(u_ex)
196 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
197 class SimpleSolve_Brick_Order1_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
198 def test_solve(self):
199 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
200 x=Solution(domain).getX()
201 u_ex=1.+2.*x[0]+3.*x[1]+4.*x[2]
202 # --- set exact gradient -----------
203 g_ex=Data(0.,(3,),Solution(domain))
204 g_ex[0]=2.
205 g_ex[1]=3.
206 g_ex[2]=4.
207 # -------- test gradient --------------------------------
208 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
209 # -------- set-up PDE -----------------------------------
210 pde=LinearPDE(domain,numEquations=1)
211 mask=whereZero(x[0])
212 pde.setValue(r=u_ex,q=mask)
213 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()))
214 # -------- get the solution ---------------------------
215 pde.setTolerance(SOLVER_TOL)
216 pde.setSolverMethod(pde.PCG,pde.JACOBI)
217 pde.setSolverPackage(pde.PASO)
218 u=pde.getSolution(verbose=SOLVER_VERBOSE)
219 # -------- test the solution ---------------------------
220 error=Lsup(u-u_ex)/Lsup(u_ex)
221 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
222 class SimpleSolve_Brick_Order1_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
223 def test_solve(self):
224 domain=Brick(NE0,NE1,NE2,1,optimize=OPTIMIZE)
225 x=Solution(domain).getX()
226 # --- set exact solution ----
227 u_ex=Vector(0,Solution(domain))
228 u_ex[0]=1.+2.*x[0]+3.*x[1]+4.*x[2]
229 u_ex[1]=-1.+4.*x[0]+1.*x[1]-2.*x[2]
230 u_ex[2]=5.+8.*x[0]+4.*x[1]+5.*x[2]
231 # --- set exact gradient -----------
232 g_ex=Data(0.,(3,3),Solution(domain))
233 g_ex[0,0]=2.
234 g_ex[0,1]=3.
235 g_ex[0,2]=4.
236 g_ex[1,0]=4.
237 g_ex[1,1]=1.
238 g_ex[1,2]=-2.
239 g_ex[2,0]=8.
240 g_ex[2,1]=4.
241 g_ex[2,2]=5.
242 # -------- test gradient --------------------------------
243 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
244 # -------- set-up PDE -----------------------------------
245 pde=LinearPDE(domain,numEquations=3)
246 mask=whereZero(x[0])
247 pde.setValue(r=u_ex,q=mask*numarray.ones(3,))
248 A=Tensor4(0,Function(domain))
249 A[0,:,0,:]=kronecker(3)
250 A[1,:,1,:]=kronecker(3)
251 A[2,:,2,:]=kronecker(3)
252 Y=Vector(0.,Function(domain))
253 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
254 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
255 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
256 pde.setValue(A=A,
257 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numarray.ones((3,3))*FAC_OFFDIAG,
258 Y=Y,
259 y=matrixmult(g_ex,domain.getNormal()))
260 # -------- get the solution ---------------------------
261 pde.setTolerance(SOLVER_TOL)
262 pde.setSolverMethod(pde.PCG,pde.JACOBI)
263 pde.setSolverPackage(pde.PASO)
264 u=pde.getSolution(verbose=SOLVER_VERBOSE)
265 # -------- test the solution ---------------------------
266 error=Lsup(u-u_ex)/Lsup(u_ex)
267 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
268 class SimpleSolve_Brick_Order2_SinglePDE_Paso_PCG_Jacobi(unittest.TestCase):
269 def test_solve(self):
270 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
271 x=Solution(domain).getX()
272 # --- set exact solution ----
273 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
274 # --- set exact gradient -----------
275 g_ex=Data(0.,(3,),Solution(domain))
276 g_ex[0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
277 g_ex[1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
278 g_ex[2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
279 # -------- test gradient --------------------------------
280 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
281 # -------- set-up PDE -----------------------------------
282 pde=LinearPDE(domain,numEquations=1)
283 mask=whereZero(x[0])
284 pde.setValue(r=u_ex,q=mask)
285 pde.setValue(A=kronecker(3),y=inner(g_ex,domain.getNormal()),Y=-60.)
286 # -------- get the solution ---------------------------
287 pde.setTolerance(SOLVER_TOL)
288 pde.setSolverMethod(pde.PCG,pde.JACOBI)
289 pde.setSolverPackage(pde.PASO)
290 u=pde.getSolution(verbose=SOLVER_VERBOSE)
291 # -------- test the solution ---------------------------
292 error=Lsup(u-u_ex)/Lsup(u_ex)
293 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
294 class SimpleSolve_Brick_Order2_SystemPDE_Paso_PCG_Jacobi(unittest.TestCase):
295 def test_solve(self):
296 domain=Brick(NE0,NE1,NE2,2,optimize=OPTIMIZE)
297 x=Solution(domain).getX()
298 # --- set exact solution ----
299 u_ex=Vector(0,Solution(domain))
300 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
301 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
302 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
303 # --- set exact gradient -----------
304 g_ex=Data(0.,(3,3),Solution(domain))
305 g_ex[0,0]=2.+6.*x[1]+8.*x[2]+18.*x[0]
306 g_ex[0,1]=3.+6.*x[0]+7.*x[2]+20.*x[1]
307 g_ex[0,2]=4.+7.*x[1]+8.*x[0]+22.*x[2]
308 g_ex[1,0]=4.+3.*x[1]-8.*x[2]-4.*x[0]
309 g_ex[1,1]=1+3.*x[0]+2.*x[2]+14.*x[1]
310 g_ex[1,2]=-6.+2.*x[1]-8.*x[0]+10.*x[2]
311 g_ex[2,0]=7.-6.*x[1]+2.*x[2]+4.*x[0]
312 g_ex[2,1]=9.-6.*x[0]+8.*x[2]+16.*x[1]
313 g_ex[2,2]=2+8.*x[1]+2.*x[0]+2.*x[2]
314 # -------- test gradient --------------------------------
315 self.failUnless(Lsup(g_ex-grad(u_ex))<REL_TOL*Lsup(g_ex))
316 # -------- set-up PDE -----------------------------------
317 pde=LinearPDE(domain,numEquations=3)
318 mask=whereZero(x[0])
319 pde.setValue(r=u_ex,q=mask*numarray.ones(3,))
320 Y=Vector(0.,Function(domain))
321 Y[0]=u_ex[0]*FAC_DIAG+u_ex[2]*FAC_OFFDIAG+u_ex[1]*FAC_OFFDIAG
322 Y[1]=u_ex[1]*FAC_DIAG+u_ex[0]*FAC_OFFDIAG+u_ex[2]*FAC_OFFDIAG
323 Y[2]=u_ex[2]*FAC_DIAG+u_ex[1]*FAC_OFFDIAG+u_ex[0]*FAC_OFFDIAG
324 A=Tensor4(0,Function(domain))
325 A[0,:,0,:]=kronecker(3)
326 A[1,:,1,:]=kronecker(3)
327 A[2,:,2,:]=kronecker(3)
328 pde.setValue(A=A,
329 D=kronecker(3)*(FAC_DIAG-FAC_OFFDIAG)+numarray.ones((3,3))*FAC_OFFDIAG,
330 Y=Y-numarray.array([60.,20.,22.]),
331 y=matrixmult(g_ex,domain.getNormal()))
332 # -------- get the solution ---------------------------
333 pde.setTolerance(SOLVER_TOL)
334 pde.setSolverMethod(pde.PCG,pde.JACOBI)
335 pde.setSolverPackage(pde.PASO)
336 u=pde.getSolution(verbose=SOLVER_VERBOSE)
337 # -------- test the solution ---------------------------
338 error=Lsup(u-u_ex)/Lsup(u_ex)
339 self.failUnless(error<REL_TOL*Lsup(u_ex), "solution error %s is too big."%error)
340
341 if __name__ == '__main__':
342 suite = unittest.TestSuite()
343 suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_Order1_SinglePDE_Paso_PCG_Jacobi))
344 suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_Order1_SystemPDE_Paso_PCG_Jacobi))
345 suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_Order2_SinglePDE_Paso_PCG_Jacobi))
346 suite.addTest(unittest.makeSuite(SimpleSolve_Rectangle_Order2_SystemPDE_Paso_PCG_Jacobi))
347 suite.addTest(unittest.makeSuite(SimpleSolve_Brick_Order1_SinglePDE_Paso_PCG_Jacobi))
348 suite.addTest(unittest.makeSuite(SimpleSolve_Brick_Order1_SystemPDE_Paso_PCG_Jacobi))
349 suite.addTest(unittest.makeSuite(SimpleSolve_Brick_Order2_SinglePDE_Paso_PCG_Jacobi))
350 suite.addTest(unittest.makeSuite(SimpleSolve_Brick_Order2_SystemPDE_Paso_PCG_Jacobi))
351 s=unittest.TextTestRunner(verbosity=2).run(suite)

  ViewVC Help
Powered by ViewVC 1.1.26