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

Annotation of /temp/finley/test/python/run_simplesolve.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1374 - (hide annotations)
Tue Jan 8 09:37:55 2008 UTC (12 years, 10 months ago) by gross
Original Path: trunk/finley/test/python/run_simplesolve.py
File MIME type: text/x-python
File size: 15696 byte(s)
some changes to get things going on the cognac.ivec.org.
1 ksteube 1310 #
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