/[escript]/trunk/escriptcore/test/python/test_linearPDEs.py
ViewVC logotype

Contents of /trunk/escriptcore/test/python/test_linearPDEs.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2455 - (show annotations)
Wed Jun 3 03:29:07 2009 UTC (10 years, 3 months ago) by jfenwick
Original Path: trunk/escript/test/python/test_linearPDEs.py
File MIME type: text/x-python
File size: 144489 byte(s)
Merging changes from numpy branch.

1
2 ########################################################
3 #
4 # Copyright (c) 2003-2008 by University of Queensland
5 # 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 __copyright__="""Copyright (c) 2003-2008 by University of Queensland
15 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 linearPDEs class
24
25 The tests must be linked with a Domain class object in the setUp method:
26
27 from esys.finley import Rectangle
28 class Test_LinearPDEOnFinley(Test_LinearPDE):
29 def setUp(self):
30 self.domain = Rectangle(10,10,2)
31 def tearDown(self):
32 del self.domain
33 suite = unittest.TestSuite()
34 suite.addTest(unittest.makeSuite(Test_LinearPDEOnFinley))
35 unittest.TextTestRunner(verbosity=2).run(suite)
36
37 @var __author__: name of author
38 @var __copyright__: copyrights
39 @var __license__: licence agreement
40 @var __url__: url entry point on documentation
41 @var __version__: version
42 @var __date__: date of the version
43 """
44
45 __author__="Lutz Gross, l.gross@uq.edu.au"
46
47 from esys.escript.util import Lsup,kronecker,interpolate,whereZero, outer, swap_axes
48 from esys.escript import Function,FunctionOnBoundary,FunctionOnContactZero,Solution,ReducedSolution,Vector,ContinuousFunction,Scalar, ReducedFunction,ReducedFunctionOnBoundary,ReducedFunctionOnContactZero,Data, Tensor4, Tensor
49 from esys.escript.linearPDEs import LinearPDE,IllegalCoefficientValue,Poisson, IllegalCoefficientFunctionSpace, TransportPDE, IllegalCoefficient, Helmholtz, LameEquation
50 import numpy
51 import unittest
52
53 class Test_linearPDEs(unittest.TestCase):
54 TOL=1.e-6
55 SOLVER_TOL=1.e-10
56 DEBUG=False
57 VERBOSE=False
58 def check(self,arg,ref_arg,tol=None):
59 """
60 checks if arg and ref_arg are nearly identical using the L{Lsup<esys.escript.util.Lsup>}
61 """
62 if tol==None: tol=self.TOL
63 return Lsup(arg-ref_arg)<=tol*Lsup(ref_arg)
64
65 class Test_LameEquation(Test_linearPDEs):
66
67 def test_config(self):
68 mypde=LameEquation(self.domain,debug=self.DEBUG)
69 d=self.domain.getDim()
70 self.failUnlessEqual((mypde.getNumEquations(),mypde.getNumSolutions(),mypde.isSymmetric()),(d,d,True),"set up incorrect")
71
72 def test_setCoefficient_q(self):
73 mypde=LameEquation(self.domain,debug=self.DEBUG)
74 x=self.domain.getX()
75 mypde.setValue(q=x)
76
77 q_ref=interpolate(x,Solution(self.domain))
78 self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
79 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
80 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
81 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
82 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
83 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
84 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
85 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
86 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
87 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
88 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
89 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
90 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
91 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
92 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
93 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
94 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
95 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
96 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
97 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
98 self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
99 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
100
101 def test_setCoefficient_r(self):
102 mypde=LameEquation(self.domain,debug=self.DEBUG)
103 x=self.domain.getX()
104 mypde.setValue(r=x)
105
106 r_ref=interpolate(x,Solution(self.domain))
107 self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
108 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
109 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
110 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
111 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
112 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
113 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
114 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
115 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
116 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
117 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
118 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
119 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
120 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
121 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
122 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
123 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
124 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
125 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
126 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
127 self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x")
128 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
129
130
131 def test_setCoefficient_F(self):
132 mypde=LameEquation(self.domain,debug=self.DEBUG)
133 x=self.domain.getX()
134 mypde.setValue(F=x)
135
136 Y_ref=interpolate(x,Function(self.domain))
137 self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
138 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
139 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
140 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
141 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
142 self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x")
143 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
144 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
145 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
146 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
147 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
148 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
149 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
150 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
151 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
152 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
153 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
154 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
155 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
156 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
157 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
158 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
159
160 def test_setCoefficient_f(self):
161 mypde=LameEquation(self.domain,debug=self.DEBUG)
162 x=self.domain.getX()
163 mypde.setValue(f=x)
164
165 y_ref=interpolate(x,FunctionOnBoundary(self.domain))
166 self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
167 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
168 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
169 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
170 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
171 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
172 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
173 self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"d is not x[0]")
174 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
175 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
176 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
177 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
178 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
179 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
180 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
181 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
182 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
183 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
184 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
185 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
186 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
187 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
188
189 def test_setCoefficient_sigma(self):
190 mypde=LameEquation(self.domain,debug=self.DEBUG)
191 x=self.domain.getX()
192 mypde.setValue(sigma=outer(x,x))
193
194 X_ref=interpolate(outer(x,x),Function(self.domain))
195 self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
196 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
197 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
198 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
199 self.failUnless(self.check(mypde.getCoefficient("X"),X_ref),"X is not x X x")
200 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
201 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
202 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
203 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
204 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
205 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
206 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
207 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
208 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
209 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
210 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
211 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
212 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
213 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
214 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
215 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
216 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
217
218 def test_setCoefficient_lambda(self):
219 mypde=LameEquation(self.domain,debug=self.DEBUG)
220 x=self.domain.getX()
221 mypde.setValue(lame_lambda=x[0])
222
223
224 k3=kronecker(Function(self.domain))
225 k3Xk3=outer(k3,k3)
226 A_ref=x[0]*k3Xk3
227
228 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
229 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
230 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
231 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
232 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
233 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
234 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
235 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
236 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
237 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
238 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
239 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
240 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
241 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
242 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
243 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
244 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
245 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
246 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
247 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
248 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
249 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
250
251 def test_setCoefficient_mu(self):
252 mypde=LameEquation(self.domain,debug=self.DEBUG)
253 x=self.domain.getX()
254 mypde.setValue(lame_mu=x[0])
255
256
257 k3=kronecker(Function(self.domain))
258 k3Xk3=outer(k3,k3)
259 A_ref=x[0]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
260
261 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
262 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
263 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
264 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
265 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
266 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
267 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
268 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
269 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
270 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
271 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
272 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
273 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
274 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
275 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
276 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
277 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
278 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
279 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
280 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
281 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
282 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
283
284 def test_setCoefficient_lambdamu(self):
285 mypde=LameEquation(self.domain,debug=self.DEBUG)
286 x=self.domain.getX()
287 mypde.setValue(lame_lambda=x[0], lame_mu=x[1])
288
289 k3=kronecker(Function(self.domain))
290 k3Xk3=outer(k3,k3)
291 A_ref=x[0]*k3Xk3+x[1]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
292
293 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
294 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
295 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
296 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
297 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
298 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
299 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
300 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
301 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
302 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
303 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
304 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
305 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
306 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
307 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
308 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
309 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
310 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
311 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
312 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
313 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
314 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
315
316 def test_solve(self):
317 d=self.domain.getDim()
318 mypde=LameEquation(self.domain)
319 cf=ContinuousFunction(self.domain)
320 x=cf.getX()
321 u_ex=x
322 msk=Vector(0.,cf)
323 for i in range(d): msk[i]=whereZero(x[i])
324 mypde.setValue(q=msk,r=u_ex,lame_mu=3,lame_lambda=50,f=(2*3+50*d)*FunctionOnBoundary(self.domain).getNormal())
325
326 u=mypde.getSolution()
327 self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
328
329 class Test_Helmholtz(Test_linearPDEs):
330
331 def test_config(self):
332 mypde=Helmholtz(self.domain,debug=self.DEBUG)
333 self.failUnlessEqual((mypde.getNumEquations(),mypde.getNumSolutions(),mypde.isSymmetric()),(1,1,True),"set up incorrect")
334 def test_setCoefficient_q(self):
335 mypde=Helmholtz(self.domain,debug=self.DEBUG)
336 x=self.domain.getX()
337 mypde.setValue(q=whereZero(x[0]))
338
339 q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
340 A_ref=kronecker(self.domain)
341
342 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
343 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
344 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
345 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
346 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
347 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
348 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
349 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
350 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
351 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
352 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
353 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
354 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
355 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
356 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
357 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
358 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
359 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
360 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
361 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
362 self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
363 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
364
365 def test_setCoefficient_r(self):
366 mypde=Helmholtz(self.domain,debug=self.DEBUG)
367 x=self.domain.getX()
368 mypde.setValue(r=x[0])
369
370 r_ref=interpolate(x[0],Solution(self.domain))
371 A_ref=kronecker(self.domain)
372 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
373 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
374 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
375 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
376 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
377 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
378 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
379 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
380 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
381 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
382 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
383 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
384 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
385 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
386 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
387 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
388 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
389 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
390 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
391 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
392 self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x[0]")
393 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
394
395
396 def test_setCoefficient_f(self):
397 mypde=Helmholtz(self.domain,debug=self.DEBUG)
398 x=self.domain.getX()
399 mypde.setValue(f=x[0])
400
401 Y_ref=interpolate(x[0],Function(self.domain))
402 A_ref=kronecker(self.domain)
403 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
404 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
405 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
406 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
407 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
408 self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
409 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
410 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
411 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
412 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
413 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
414 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
415 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
416 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
417 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
418 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
419 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
420 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
421 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
422 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
423 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
424 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
425
426 def test_setCoefficient_alpha(self):
427 mypde=Helmholtz(self.domain,debug=self.DEBUG)
428 x=self.domain.getX()
429 mypde.setValue(alpha=x[0])
430
431 d_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
432 A_ref=kronecker(self.domain)
433 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
434 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
435 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
436 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
437 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
438 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
439 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
440 self.failUnless(self.check(mypde.getCoefficient("d"),d_ref),"d is not x[0]")
441 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
442 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
443 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
444 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
445 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
446 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
447 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
448 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
449 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
450 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
451 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
452 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
453 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
454 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
455
456 def test_setCoefficient_g(self):
457 mypde=Helmholtz(self.domain,debug=self.DEBUG)
458 x=self.domain.getX()
459 mypde.setValue(g=x[0])
460
461 y_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
462 A_ref=kronecker(self.domain)
463 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
464 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
465 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
466 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
467 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
468 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
469 self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"y is not x[0]")
470 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
471 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
472 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
473 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
474 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
475 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
476 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
477 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
478 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
479 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
480 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
481 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
482 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
483 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
484 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
485
486 def test_setCoefficient_omega(self):
487 mypde=Helmholtz(self.domain,debug=self.DEBUG)
488 x=self.domain.getX()
489 mypde.setValue(omega=x[0])
490
491 D_ref=interpolate(x[0],Function(self.domain))
492 A_ref=kronecker(self.domain)
493 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
494 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
495 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
496 self.failUnless(self.check(mypde.getCoefficient("D"),D_ref),"D is not x[0]")
497 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
498 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
499 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
500 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
501 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
502 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
503 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
504 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
505 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
506 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
507 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
508 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
509 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
510 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
511 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
512 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
513 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
514 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
515
516 def test_solve(self):
517 d=self.domain.getDim()
518 cf=ContinuousFunction(self.domain)
519 u_ex=Scalar(1.,cf)
520 mypde=Helmholtz(self.domain)
521 mypde.setValue(f=3,omega=3,alpha=2,g=2)
522 u=mypde.getSolution()
523 self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
524
525 class Test_Poisson(Test_linearPDEs):
526
527 def test_config(self):
528 mypde=Poisson(self.domain,debug=self.DEBUG)
529 self.failUnlessEqual((mypde.getNumEquations(),mypde.getNumSolutions(),mypde.isSymmetric()),(1,1,True),"set up incorrect")
530 def test_setCoefficient_q(self):
531 mypde=Poisson(self.domain,debug=self.DEBUG)
532 x=self.domain.getX()
533 q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
534 A_ref=kronecker(self.domain)
535 mypde.setValue(q=whereZero(x[0]))
536 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
537 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
538 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
539 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
540 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
541 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
542 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
543 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
544 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
545 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
546 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
547 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
548 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
549 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
550 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
551 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
552 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
553 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
554 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
555 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
556 self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
557 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
558 def test_setCoefficient_f(self):
559 mypde=Poisson(self.domain,debug=self.DEBUG)
560 x=self.domain.getX()
561 Y_ref=interpolate(x[0],Function(self.domain))
562 A_ref=kronecker(self.domain)
563 mypde.setValue(f=x[0])
564 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
565 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
566 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
567 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
568 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
569 self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
570 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
571 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
572 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
573 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
574 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
575 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
576 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
577 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
578 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
579 self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
580 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
581 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
582 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
583 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
584 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
585 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
586 def test_setCoefficient_f_reduced(self):
587 mypde=Poisson(self.domain,debug=self.DEBUG)
588 x=self.domain.getX()
589 Y_ref=interpolate(x[0],ReducedFunction(self.domain))
590 A_ref=kronecker(self.domain)
591 mypde.setValue(f_reduced=x[0])
592 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
593 self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
594 self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
595 self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
596 self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
597 self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
598 self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
599 self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
600 self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
601 self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
602 self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
603 self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
604 self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
605 self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
606 self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
607 self.failUnless(self.check(mypde.getCoefficient("Y_reduced"),Y_ref),"Y_reduced is not x[0]")
608 self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
609 self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
610 self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
611 self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
612 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
613 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
614 def test_solve(self):
615 d=self.domain.getDim()
616 cf=ContinuousFunction(self.domain)
617 x=cf.getX()
618 #construct exact solution:
619 u_ex=Scalar(1.,cf)
620 for i in range(d):
621 u_ex*=x[i]*(2.-x[i])
622 #construct mask:
623 msk=Scalar(0.,cf)
624 for i in range(d):
625 msk+=whereZero(x[i])
626 #construct right hand side
627 f=Scalar(0,cf)
628 for i in range(d):
629 f_p=Scalar(1,cf)
630 for j in range(d):
631 if i==j:
632 f_p*=2.
633 else:
634 f_p*=x[j]*(2-x[j])
635 f+=f_p
636 mypde=Poisson(self.domain)
637 mypde.setValue(f=f,q=msk)
638 u=mypde.getSolution()
639 self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
640
641 class Test_LinearPDE_noLumping(Test_linearPDEs):
642 N=4
643 def test_setCoefficient_WithIllegalFunctionSpace(self):
644 mypde=LinearPDE(self.domain,debug=self.DEBUG)
645 self.failUnlessRaises(IllegalCoefficientFunctionSpace,mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))
646
647 def test_setCoefficient_WithWrongName(self):
648 mypde=LinearPDE(self.domain,debug=self.DEBUG)
649 self.failUnlessRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)
650
651 def test_resetCoefficient_WithWrongShape(self):
652 mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)
653 self.failUnlessRaises(IllegalCoefficientValue, mypde.setValue, C=0.)
654
655 def test_reducedOn(self):
656 mypde=LinearPDE(self.domain,debug=self.DEBUG)
657 x=self.domain.getX()
658 mypde.setReducedOrderOn()
659 mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])
660 u=mypde.getSolution()
661 self.failUnless(self.check(u,1.),'solution is wrong.')
662
663 def test_attemptToChangeOrderAfterDefinedCoefficient(self):
664 mypde=LinearPDE(self.domain,debug=self.DEBUG)
665 mypde.setValue(D=1.)
666 self.failUnlessRaises(RuntimeError,mypde.setReducedOrderOn)
667
668 def test_reducedOnConfig(self):
669 mypde=LinearPDE(self.domain,debug=self.DEBUG)
670 mypde.setReducedOrderOn()
671 self.failUnlessEqual((mypde.getFunctionSpaceForSolution(),mypde.getFunctionSpaceForEquation()),(ReducedSolution(self.domain),ReducedSolution(self.domain)),"reduced function spaces expected.")
672 #
673 # set coefficients for scalars:
674 #
675 def test_setCoefficient_A_Scalar(self):
676 d=self.domain.getDim()
677 mypde=LinearPDE(self.domain,debug=self.DEBUG)
678 mypde.setValue(A=numpy.ones((d,d)))
679 coeff=mypde.getCoefficient("A")
680 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,d),Function(self.domain),1,1))
681 def test_setCoefficient_B_Scalar(self):
682 d=self.domain.getDim()
683 mypde=LinearPDE(self.domain,debug=self.DEBUG)
684 mypde.setValue(B=numpy.ones((d,)))
685 coeff=mypde.getCoefficient("B")
686 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),Function(self.domain),1,1))
687 def test_setCoefficient_C_Scalar(self):
688 d=self.domain.getDim()
689 mypde=LinearPDE(self.domain,debug=self.DEBUG)
690 mypde.setValue(C=numpy.ones((d,)))
691 coeff=mypde.getCoefficient("C")
692 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),Function(self.domain),1,1))
693 def test_setCoefficient_D_Scalar(self):
694 d=self.domain.getDim()
695 mypde=LinearPDE(self.domain,debug=self.DEBUG)
696 mypde.setValue(D=1.)
697 coeff=mypde.getCoefficient("D")
698 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),Function(self.domain),1,1))
699 def test_setCoefficient_X_Scalar(self):
700 d=self.domain.getDim()
701 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
702 mypde.setValue(X=numpy.ones((d,)))
703 coeff=mypde.getCoefficient("X")
704 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((d,),Function(self.domain),1))
705 def test_setCoefficient_Y_Scalar(self):
706 d=self.domain.getDim()
707 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
708 mypde.setValue(Y=1.)
709 coeff=mypde.getCoefficient("Y")
710 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),Function(self.domain),1))
711 def test_setCoefficient_y_Scalar(self):
712 d=self.domain.getDim()
713 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
714 mypde.setValue(y=1.)
715 coeff=mypde.getCoefficient("y")
716 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))
717 def test_setCoefficient_d_Scalar(self):
718 d=self.domain.getDim()
719 mypde=LinearPDE(self.domain,debug=self.DEBUG)
720 mypde.setValue(d=1.)
721 coeff=mypde.getCoefficient("d")
722 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1,1))
723 def test_setCoefficient_d_contact_Scalar(self):
724 d=self.domain.getDim()
725 mypde=LinearPDE(self.domain,debug=self.DEBUG)
726 mypde.setValue(d_contact=1.)
727 coeff=mypde.getCoefficient("d_contact")
728 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))
729 def test_setCoefficient_y_contact_Scalar(self):
730 d=self.domain.getDim()
731 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
732 mypde.setValue(y_contact=1.)
733 coeff=mypde.getCoefficient("y_contact")
734 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))
735 def test_setCoefficient_A_reduced_Scalar(self):
736 d=self.domain.getDim()
737 mypde=LinearPDE(self.domain,debug=self.DEBUG)
738 mypde.setValue(A_reduced=numpy.ones((d,d)))
739 coeff=mypde.getCoefficient("A_reduced")
740 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))
741 def test_setCoefficient_B_reduced_Scalar(self):
742 d=self.domain.getDim()
743 mypde=LinearPDE(self.domain,debug=self.DEBUG)
744 mypde.setValue(B_reduced=numpy.ones((d,)))
745 coeff=mypde.getCoefficient("B_reduced")
746 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
747 def test_setCoefficient_C_reduced_Scalar(self):
748 d=self.domain.getDim()
749 mypde=LinearPDE(self.domain,debug=self.DEBUG)
750 mypde.setValue(C_reduced=numpy.ones((d,)))
751 coeff=mypde.getCoefficient("C_reduced")
752 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
753 def test_setCoefficient_D_reduced_Scalar(self):
754 d=self.domain.getDim()
755 mypde=LinearPDE(self.domain,debug=self.DEBUG)
756 mypde.setValue(D_reduced=1.)
757 coeff=mypde.getCoefficient("D_reduced")
758 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))
759 def test_setCoefficient_X_reduced_Scalar(self):
760 d=self.domain.getDim()
761 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
762 mypde.setValue(X_reduced=numpy.ones((d,)))
763 coeff=mypde.getCoefficient("X_reduced")
764 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
765 def test_setCoefficient_Y_reduced_Scalar(self):
766 d=self.domain.getDim()
767 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
768 mypde.setValue(Y_reduced=1.)
769 coeff=mypde.getCoefficient("Y_reduced")
770 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
771 def test_setCoefficient_y_reduced_Scalar(self):
772 d=self.domain.getDim()
773 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
774 mypde.setValue(y_reduced=1.)
775 coeff=mypde.getCoefficient("y_reduced")
776 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
777 def test_setCoefficient_d_reduced_Scalar(self):
778 d=self.domain.getDim()
779 mypde=LinearPDE(self.domain,debug=self.DEBUG)
780 mypde.setValue(d_reduced=1.)
781 coeff=mypde.getCoefficient("d_reduced")
782 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))
783 def test_setCoefficient_d_contact_reduced_Scalar(self):
784 d=self.domain.getDim()
785 mypde=LinearPDE(self.domain,debug=self.DEBUG)
786 mypde.setValue(d_contact_reduced=1.)
787 coeff=mypde.getCoefficient("d_contact_reduced")
788 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
789 def test_setCoefficient_y_contact_reduced_Scalar(self):
790 d=self.domain.getDim()
791 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
792 mypde.setValue(y_contact_reduced=1.)
793 coeff=mypde.getCoefficient("y_contact_reduced")
794 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
795 def test_setCoefficient_r_Scalar(self):
796 d=self.domain.getDim()
797 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
798 mypde.setValue(r=1.)
799 coeff=mypde.getCoefficient("r")
800 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),Solution(self.domain),1))
801 def test_setCoefficient_q_Scalar(self):
802 d=self.domain.getDim()
803 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
804 mypde.setValue(q=1.)
805 coeff=mypde.getCoefficient("q")
806 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),Solution(self.domain),1))
807 def test_setCoefficient_r_Scalar_reducedOn(self):
808 d=self.domain.getDim()
809 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
810 mypde.setReducedOrderOn()
811 mypde.setValue(r=1.)
812 coeff=mypde.getCoefficient("r")
813 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
814 def test_setCoefficient_q_Scalar_reducedOn(self):
815 d=self.domain.getDim()
816 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
817 mypde.setReducedOrderOn()
818 mypde.setValue(q=1.)
819 coeff=mypde.getCoefficient("q")
820 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
821
822 def test_setCoefficient_A_reduced_Scalar_usingA(self):
823 d=self.domain.getDim()
824 mypde=LinearPDE(self.domain,debug=self.DEBUG)
825 mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))
826 coeff=mypde.getCoefficient("A_reduced")
827 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))
828 def test_setCoefficient_B_reduced_Scalar_usingB(self):
829 d=self.domain.getDim()
830 mypde=LinearPDE(self.domain,debug=self.DEBUG)
831 mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
832 coeff=mypde.getCoefficient("B_reduced")
833 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
834 def test_setCoefficient_C_reduced_Scalar_usingC(self):
835 d=self.domain.getDim()
836 mypde=LinearPDE(self.domain,debug=self.DEBUG)
837 mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
838 coeff=mypde.getCoefficient("C_reduced")
839 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
840 def test_setCoefficient_D_reduced_Scalar_usingD(self):
841 d=self.domain.getDim()
842 mypde=LinearPDE(self.domain,debug=self.DEBUG)
843 mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))
844 coeff=mypde.getCoefficient("D_reduced")
845 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))
846 def test_setCoefficient_X_reduced_Scalar_usingX(self):
847 d=self.domain.getDim()
848 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
849 mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
850 coeff=mypde.getCoefficient("X_reduced")
851 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
852 def test_setCoefficient_Y_reduced_Scalar_usingY(self):
853 d=self.domain.getDim()
854 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
855 mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))
856 coeff=mypde.getCoefficient("Y_reduced")
857 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
858 def test_setCoefficient_y_reduced_Scalar_using_y(self):
859 d=self.domain.getDim()
860 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
861 mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
862 coeff=mypde.getCoefficient("y_reduced")
863 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
864 def test_setCoefficient_d_reduced_Scalar_using_d(self):
865 d=self.domain.getDim()
866 mypde=LinearPDE(self.domain,debug=self.DEBUG)
867 mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
868 coeff=mypde.getCoefficient("d_reduced")
869 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))
870 def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):
871 d=self.domain.getDim()
872 mypde=LinearPDE(self.domain,debug=self.DEBUG)
873 mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
874 coeff=mypde.getCoefficient("d_contact_reduced")
875 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
876 def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):
877 d=self.domain.getDim()
878 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
879 mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
880 coeff=mypde.getCoefficient("y_contact_reduced")
881 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
882 #
883 # set coefficients for systems:
884 #
885 def test_setCoefficient_A_System(self):
886 d=self.domain.getDim()
887 mypde=LinearPDE(self.domain,debug=self.DEBUG)
888 mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))
889 coeff=mypde.getCoefficient("A")
890 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N,d),Function(self.domain),self.N,self.N))
891 def test_setCoefficient_B_System(self):
892 d=self.domain.getDim()
893 mypde=LinearPDE(self.domain,debug=self.DEBUG)
894 mypde.setValue(B=numpy.ones((self.N,d,self.N)))
895 coeff=mypde.getCoefficient("B")
896 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N),Function(self.domain),self.N,self.N))
897 def test_setCoefficient_C_System(self):
898 d=self.domain.getDim()
899 mypde=LinearPDE(self.domain,debug=self.DEBUG)
900 mypde.setValue(C=numpy.ones((self.N,self.N,d)))
901 coeff=mypde.getCoefficient("C")
902 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N,d),Function(self.domain),self.N,self.N))
903 def test_setCoefficient_D_System(self):
904 d=self.domain.getDim()
905 mypde=LinearPDE(self.domain,debug=self.DEBUG)
906 mypde.setValue(D=numpy.ones((self.N,self.N)))
907 coeff=mypde.getCoefficient("D")
908 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),Function(self.domain),self.N,self.N))
909 def test_setCoefficient_X_System(self):
910 d=self.domain.getDim()
911 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
912 mypde.setValue(X=numpy.ones((self.N,d)))
913 coeff=mypde.getCoefficient("X")
914 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,d),Function(self.domain),self.N))
915 def test_setCoefficient_Y_System(self):
916 d=self.domain.getDim()
917 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
918 mypde.setValue(Y=numpy.ones((self.N,)))
919 coeff=mypde.getCoefficient("Y")
920 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),Function(self.domain),self.N))
921 def test_setCoefficient_y_System(self):
922 d=self.domain.getDim()
923 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
924 mypde.setValue(y=numpy.ones((self.N,)))
925 coeff=mypde.getCoefficient("y")
926 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),FunctionOnBoundary(self.domain),self.N))
927 def test_setCoefficient_d_System(self):
928 d=self.domain.getDim()
929 mypde=LinearPDE(self.domain,debug=self.DEBUG)
930 mypde.setValue(d=numpy.ones((self.N,self.N)))
931 coeff=mypde.getCoefficient("d")
932 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FunctionOnBoundary(self.domain),self.N,self.N))
933 def test_setCoefficient_d_contact_System(self):
934 d=self.domain.getDim()
935 mypde=LinearPDE(self.domain,debug=self.DEBUG)
936 mypde.setValue(d_contact=numpy.ones((self.N,self.N)))
937 coeff=mypde.getCoefficient("d_contact")
938 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FunctionOnContactZero(self.domain),self.N,self.N))
939 def test_setCoefficient_y_contact_System(self):
940 d=self.domain.getDim()
941 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
942 mypde.setValue(y_contact=numpy.ones((self.N,)))
943 coeff=mypde.getCoefficient("y_contact")
944 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))
945 def test_setCoefficient_A_reduced_System(self):
946 d=self.domain.getDim()
947 mypde=LinearPDE(self.domain,debug=self.DEBUG)
948 mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))
949 coeff=mypde.getCoefficient("A_reduced")
950 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))
951 def test_setCoefficient_B_reduced_System(self):
952 d=self.domain.getDim()
953 mypde=LinearPDE(self.domain,debug=self.DEBUG)
954 mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))
955 coeff=mypde.getCoefficient("B_reduced")
956 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))
957 def test_setCoefficient_C_reduced_System(self):
958 d=self.domain.getDim()
959 mypde=LinearPDE(self.domain,debug=self.DEBUG)
960 mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))
961 coeff=mypde.getCoefficient("C_reduced")
962 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))
963 def test_setCoefficient_D_System_reduced(self):
964 d=self.domain.getDim()
965 mypde=LinearPDE(self.domain,debug=self.DEBUG)
966 mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))
967 coeff=mypde.getCoefficient("D_reduced")
968 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))
969 def test_setCoefficient_X_System_reduced(self):
970 d=self.domain.getDim()
971 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
972 mypde.setValue(X_reduced=numpy.ones((self.N,d)))
973 coeff=mypde.getCoefficient("X_reduced")
974 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))
975 def test_setCoefficient_Y_System_reduced(self):
976 d=self.domain.getDim()
977 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
978 mypde.setValue(Y_reduced=numpy.ones((self.N,)))
979 coeff=mypde.getCoefficient("Y_reduced")
980 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))
981 def test_setCoefficient_y_System_reduced(self):
982 d=self.domain.getDim()
983 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
984 mypde.setValue(y_reduced=numpy.ones((self.N,)))
985 coeff=mypde.getCoefficient("y_reduced")
986 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))
987 def test_setCoefficient_d_reduced_System(self):
988 d=self.domain.getDim()
989 mypde=LinearPDE(self.domain,debug=self.DEBUG)
990 mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))
991 coeff=mypde.getCoefficient("d_reduced")
992 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))
993 def test_setCoefficient_d_contact_reduced_System(self):
994 d=self.domain.getDim()
995 mypde=LinearPDE(self.domain,debug=self.DEBUG)
996 mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))
997 coeff=mypde.getCoefficient("d_contact_reduced")
998 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
999 def test_setCoefficient_y_contact_reduced_System(self):
1000 d=self.domain.getDim()
1001 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1002 mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))
1003 coeff=mypde.getCoefficient("y_contact_reduced")
1004 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1005 def test_setCoefficient_r_System(self):
1006 d=self.domain.getDim()
1007 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1008 mypde.setValue(r=numpy.ones((self.N,)))
1009 coeff=mypde.getCoefficient("r")
1010 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1011 def test_setCoefficient_q_System(self):
1012 d=self.domain.getDim()
1013 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1014 mypde.setValue(q=numpy.ones((self.N,)))
1015 coeff=mypde.getCoefficient("q")
1016 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1017 def test_setCoefficient_r_System_reducedOn(self):
1018 d=self.domain.getDim()
1019 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1020 mypde.setReducedOrderOn()
1021 mypde.setValue(r=numpy.ones((self.N,)))
1022 coeff=mypde.getCoefficient("r")
1023 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1024 def test_setCoefficient_q_System_reducedOn(self):
1025 d=self.domain.getDim()
1026 mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1027 mypde.setReducedOrderOn()
1028 mypde.setValue(q=numpy.ones((self.N,)))
1029 coeff=mypde.getCoefficient("q")
1030 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1031
1032 def test_setCoefficient_A_reduced_System_using_A(self):
1033 d=self.domain.getDim()
1034 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1035 mypde.setValue(A=Data(numpy.ones((self.N,d,self.N,d)),ReducedFunction(self.domain)))
1036 coeff=mypde.getCoefficient("A_reduced")
1037 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1038 def test_setCoefficient_B_reduced_System_using_B(self):
1039 d=self.domain.getDim()
1040 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1041 mypde.setValue(B=Data(numpy.ones((self.N,d,self.N)),ReducedFunction(self.domain)))
1042 coeff=mypde.getCoefficient("B_reduced")
1043 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))
1044 def test_setCoefficient_C_reduced_System_using_C(self):
1045 d=self.domain.getDim()
1046 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1047 mypde.setValue(C=Data(numpy.ones((self.N,self.N,d)),ReducedFunction(self.domain)))
1048 coeff=mypde.getCoefficient("C_reduced")
1049 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1050 def test_setCoefficient_D_System_reduced_using_D(self):
1051 d=self.domain.getDim()
1052 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1053 mypde.setValue(D=Data(numpy.ones((self.N,self.N)),ReducedFunction(self.domain)))
1054 coeff=mypde.getCoefficient("D_reduced")
1055 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))
1056 def test_setCoefficient_X_System_reduced_using_X(self):
1057 d=self.domain.getDim()
1058 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1059 mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))
1060 coeff=mypde.getCoefficient("X_reduced")
1061 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))
1062 def test_setCoefficient_Y_System_reduced_using_Y(self):
1063 d=self.domain.getDim()
1064 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1065 mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))
1066 coeff=mypde.getCoefficient("Y_reduced")
1067 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))
1068 def test_setCoefficient_y_reduced_System_using_y(self):
1069 d=self.domain.getDim()
1070 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1071 mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))
1072 coeff=mypde.getCoefficient("y_reduced")
1073 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))
1074 def test_setCoefficient_d_reduced_System_using_d(self):
1075 d=self.domain.getDim()
1076 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1077 mypde.setValue(d=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnBoundary(self.domain)))
1078 coeff=mypde.getCoefficient("d_reduced")
1079 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))
1080 def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):
1081 d=self.domain.getDim()
1082 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1083 mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))
1084 coeff=mypde.getCoefficient("d_contact_reduced")
1085 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
1086 def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):
1087 d=self.domain.getDim()
1088 mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1089 mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))
1090 coeff=mypde.getCoefficient("y_contact_reduced")
1091 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1092 def test_resetCoefficient_HomogeneousConstraint(self):
1093 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1094 x=self.domain.getX()
1095 mypde.setValue(A=kronecker(self.domain),Y=1.,q=whereZero(x[0]))
1096 u1=mypde.getSolution()
1097 mypde.setValue(Y=2.)
1098 u2=mypde.getSolution()
1099 self.failUnless(self.check(u2,2*u1),'solution is wrong.')
1100
1101 def test_resetCoefficient_InHomogeneousConstraint(self):
1102 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1103 mypde.setSymmetryOn()
1104 x=self.domain.getX()
1105 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.,r=1,q=whereZero(x[0]))
1106 u1=mypde.getSolution(verbose=self.VERBOSE)
1107 mypde.setValue(Y=2.,D=2)
1108 u2=mypde.getSolution(verbose=self.VERBOSE)
1109 self.failUnless(self.check(u2,u1),'first solution is wrong.')
1110 u2=mypde.getSolution(verbose=self.VERBOSE)
1111 self.failUnless(self.check(u2,u1),'first solution is wrong.')
1112 mypde.setValue(r=2,Y=4.)
1113 u2=mypde.getSolution(verbose=self.VERBOSE)
1114 self.failUnless(self.check(u2,2*u1),'second solution is wrong.')
1115
1116 def test_symmetryCheckTrue_System(self):
1117 d=self.domain.getDim()
1118 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1119 A=numpy.ones((self.N,d,self.N,d))
1120 C=2*numpy.ones((self.N,self.N,d))
1121 B=2*numpy.ones((self.N,d,self.N))
1122 D=3*numpy.ones((self.N,self.N))
1123 d=4*numpy.ones((self.N,self.N))
1124 d_contact=5*numpy.ones((self.N,self.N))
1125 mypde.setValue(A=A,B=B,C=C,D=D,d=d,d_contact=d_contact,A_reduced=-A,B_reduced=-B,C_reduced=-C,D_reduced=-D,d_reduced=-d,d_contact_reduced=-d_contact)
1126 self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")
1127
1128 def test_symmetryCheckFalse_A_System(self):
1129 d=self.domain.getDim()
1130 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1131 A=numpy.ones((self.N,d,self.N,d))
1132 A[1,1,1,0]=0.
1133 mypde.setValue(A=A)
1134 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1135 def test_symmetryCheckFalse_BC_System(self):
1136 d=self.domain.getDim()
1137 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1138 C=2*numpy.ones((self.N,self.N,d))
1139 B=2*numpy.ones((self.N,d,self.N))
1140 B[0,0,1]=1.
1141 mypde.setValue(B=B,C=C)
1142 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1143
1144 def test_symmetryCheckFalse_D_System(self):
1145 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1146 D=3*numpy.ones((self.N,self.N))
1147 D[0,1]=0.
1148 mypde.setValue(D=D)
1149 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1150
1151 def test_symmetryCheckFalse_d_System(self):
1152 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1153 d=4*numpy.ones((self.N,self.N))
1154 d[0,1]=0.
1155 mypde.setValue(d=d)
1156 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1157
1158 def test_symmetryCheckFalse_d_contact_System(self):
1159 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1160 d_contact=5*numpy.ones((self.N,self.N))
1161 d_contact[0,1]=0.
1162 mypde.setValue(d_contact=d_contact)
1163 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1164
1165 def test_symmetryCheckFalse_A_reduced_System(self):
1166 d=self.domain.getDim()
1167 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1168 A=numpy.ones((self.N,d,self.N,d))
1169 A[1,1,1,0]=0.
1170 mypde.setValue(A_reduced=A)
1171 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1172 def test_symmetryCheckFalse_BC_reduced_System(self):
1173 d=self.domain.getDim()
1174 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1175 C=2*numpy.ones((self.N,self.N,d))
1176 B=2*numpy.ones((self.N,d,self.N))
1177 B[0,0,1]=1.
1178 mypde.setValue(B_reduced=B,C_reduced=C)
1179 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1180
1181 def test_symmetryCheckFalse_D_reduced_System(self):
1182 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1183 D=3*numpy.ones((self.N,self.N))
1184 D[0,1]=0.
1185 mypde.setValue(D_reduced=D)
1186 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1187
1188 def test_symmetryCheckFalse_d_reduced_System(self):
1189 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1190 d=4*numpy.ones((self.N,self.N))
1191 d[0,1]=0.
1192 mypde.setValue(d_reduced=d)
1193 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1194
1195 def test_symmetryCheckFalse_d_contact_reduced_System(self):
1196 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1197 d_contact=5*numpy.ones((self.N,self.N))
1198 d_contact[0,1]=0.
1199 mypde.setValue(d_contact_reduced=d_contact)
1200 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1201
1202 def test_symmetryCheckTrue_Scalar(self):
1203 d=self.domain.getDim()
1204 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1205 A=numpy.ones((d,d))
1206 C=2*numpy.ones((d,))
1207 B=2*numpy.ones((d,))
1208 D=3
1209 d=4
1210 d_contact=5
1211 mypde.setValue(A=A,B=B,C=C,D=D,d=d,d_contact=d_contact,A_reduced=-A,B_reduced=-B,C_reduced=-C,D_reduced=-D,d_reduced=-d,d_contact_reduced=-d_contact)
1212 self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")
1213
1214 def test_symmetryCheckFalse_A_Scalar(self):
1215 d=self.domain.getDim()
1216 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1217 A=numpy.ones((d,d))
1218 A[1,0]=0.
1219 mypde.setValue(A=A)
1220 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1221 def test_symmetryCheckFalse_BC_Scalar(self):
1222 d=self.domain.getDim()
1223 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1224 C=2*numpy.ones((d,))
1225 B=2*numpy.ones((d,))
1226 B[0]=1.
1227 mypde.setValue(B=B,C=C)
1228 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1229 def test_symmetryCheckFalse_A_reduced_Scalar(self):
1230 d=self.domain.getDim()
1231 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1232 A=numpy.ones((d,d))
1233 A[1,0]=0.
1234 mypde.setValue(A_reduced=A)
1235 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1236 def test_symmetryCheckFalse_BC_reduced_Scalar(self):
1237 d=self.domain.getDim()
1238 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1239 C=2*numpy.ones((d,))
1240 B=2*numpy.ones((d,))
1241 B[0]=1.
1242 mypde.setValue(B_reduced=B,C_reduced=C)
1243 self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1244 #
1245 # solver checks (single PDE)
1246 #
1247 def test_symmetryOnIterative(self):
1248 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1249 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1250 u=mypde.getSolution(verbose=self.VERBOSE)
1251 self.failUnless(self.check(u,1.),'solution is wrong.')
1252 def test_symmetryOnDirect(self):
1253 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1254 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1255 mypde.setSolverMethod(mypde.DIRECT)
1256 u=mypde.getSolution(verbose=self.VERBOSE)
1257 self.failUnless(self.check(u,1.),'solution is wrong.')
1258 def test_PCG_JACOBI(self):
1259 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1260 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1261 mypde.setSolverMethod(mypde.PCG,mypde.JACOBI)
1262 u=mypde.getSolution(verbose=self.VERBOSE)
1263 self.failUnless(self.check(u,1.),'solution is wrong.')
1264 def test_PCG_ILU0(self):
1265 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1266 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1267 mypde.setSolverMethod(mypde.PCG,mypde.ILU0)
1268 u=mypde.getSolution(verbose=self.VERBOSE)
1269 self.failUnless(self.check(u,1.),'solution is wrong.')
1270 def test_PCG_RILU(self):
1271 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1272 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1273 mypde.setSolverMethod(mypde.PCG,mypde.RILU)
1274 u=mypde.getSolution(verbose=self.VERBOSE)
1275 self.failUnless(self.check(u,1.),'solution is wrong.')
1276 def test_DIRECT(self):
1277 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1278 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1279 mypde.setSolverMethod(mypde.DIRECT)
1280 u=mypde.getSolution(verbose=self.VERBOSE)
1281 self.failUnless(self.check(u,1.),'solution is wrong.')
1282 def test_BICGSTAB_JACOBI(self):
1283 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1284 mypde.setSolverMethod(mypde.BICGSTAB,mypde.JACOBI)
1285 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1286 u=mypde.getSolution(verbose=self.VERBOSE)
1287 self.failUnless(self.check(u,1.),'solution is wrong.')
1288 def test_BICGSTAB_ILU0(self):
1289 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1290 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1291 mypde.setSolverMethod(mypde.BICGSTAB,mypde.ILU0)
1292 u=mypde.getSolution(verbose=self.VERBOSE)
1293 self.failUnless(self.check(u,1.),'solution is wrong.')
1294 def test_BICGSTAB_RILU(self):
1295 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1296 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1297 mypde.setSolverMethod(mypde.BICGSTAB,mypde.RILU)
1298 u=mypde.getSolution(verbose=self.VERBOSE)
1299 self.failUnless(self.check(u,1.),'solution is wrong.')
1300 def test_MINRES_JACOBI(self):
1301 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1302 mypde.setSolverMethod(mypde.MINRES,mypde.JACOBI)
1303 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1304 u=mypde.getSolution(verbose=self.VERBOSE)
1305 self.failUnless(self.check(u,1.),'solution is wrong.')
1306 def test_MINRES_ILU0(self):
1307 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1308 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1309 mypde.setSolverMethod(mypde.MINRES,mypde.ILU0)
1310 u=mypde.getSolution(verbose=self.VERBOSE)
1311 self.failUnless(self.check(u,1.),'solution is wrong.')
1312 def test_MINRES_RILU(self):
1313 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1314 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1315 mypde.setSolverMethod(mypde.MINRES,mypde.RILU)
1316 u=mypde.getSolution(verbose=self.VERBOSE)
1317 self.failUnless(self.check(u,1.),'solution is wrong.')
1318 def test_TFQMR_JACOBI(self):
1319 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1320 mypde.setSolverMethod(mypde.TFQMR,mypde.JACOBI)
1321 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1322 u=mypde.getSolution(verbose=self.VERBOSE)
1323 self.failUnless(self.check(u,1.),'solution is wrong.')
1324 def test_TFQMR_ILU0(self):
1325 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1326 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1327 mypde.setSolverMethod(mypde.TFQMR,mypde.ILU0)
1328 u=mypde.getSolution(verbose=self.VERBOSE)
1329 self.failUnless(self.check(u,1.),'solution is wrong.')
1330 def test_TFQMR_RILU(self):
1331 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1332 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1333 mypde.setSolverMethod(mypde.TFQMR,mypde.RILU)
1334 u=mypde.getSolution(verbose=self.VERBOSE)
1335 self.failUnless(self.check(u,1.),'solution is wrong.')
1336 def test_PRES20_JACOBI(self):
1337 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1338 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1339 mypde.setSolverMethod(mypde.PRES20,mypde.JACOBI)
1340 u=mypde.getSolution(verbose=self.VERBOSE)
1341 self.failUnless(self.check(u,1.),'solution is wrong.')
1342 def test_PRES20_ILU0(self):
1343 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1344 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1345 mypde.setSolverMethod(mypde.PRES20,mypde.ILU0)
1346 u=mypde.getSolution(verbose=self.VERBOSE)
1347 self.failUnless(self.check(u,1.),'solution is wrong.')
1348 def test_PRES20_RILU(self):
1349 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1350 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1351 mypde.setSolverMethod(mypde.PRES20,mypde.RILU)
1352 u=mypde.getSolution(verbose=self.VERBOSE)
1353 self.failUnless(self.check(u,1.),'solution is wrong.')
1354 def test_GMRESnoRestart_JACOBI(self):
1355 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1356 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1357 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1358 u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)
1359 self.failUnless(self.check(u,1.),'solution is wrong.')
1360 def test_GMRESnoRestart_ILU0(self):
1361 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1362 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1363 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1364 u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)
1365 u=mypde.getSolution(verbose=self.VERBOSE)
1366 self.failUnless(self.check(u,1.),'solution is wrong.')
1367 def test_GMRESnoRestart_RILU(self):
1368 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1369 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1370 mypde.setSolverMethod(mypde.GMRES,mypde.RILU)
1371 u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)
1372 u=mypde.getSolution(verbose=self.VERBOSE)
1373 self.failUnless(self.check(u,1.),'solution is wrong.')
1374 def test_GMRES_JACOBI(self):
1375 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1376 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1377 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1378 u=mypde.getSolution(verbose=self.VERBOSE)
1379 self.failUnless(self.check(u,1.),'solution is wrong.')
1380 def test_GMRES_ILU0(self):
1381 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1382 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1383 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1384 u=mypde.getSolution(verbose=self.VERBOSE)
1385 self.failUnless(self.check(u,1.),'solution is wrong.')
1386 def test_GMRES_RILU(self):
1387 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1388 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1389 mypde.setSolverMethod(mypde.GMRES,mypde.RILU)
1390 u=mypde.getSolution(verbose=self.VERBOSE)
1391 self.failUnless(self.check(u,1.),'solution is wrong.')
1392 def test_GMRES_truncation_restart_JACOBI(self):
1393 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1394 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1395 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1396 u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)
1397 self.failUnless(self.check(u,1.),'solution is wrong.')
1398 def test_GMRES_truncation_restart_ILU0(self):
1399 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1400 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1401 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1402 u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)
1403 self.failUnless(self.check(u,1.),'solution is wrong.')
1404 def test_GMRES_truncation_restart_RILU(self):
1405 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1406 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1407 mypde.setSolverMethod(mypde.GMRES,mypde.RILU)
1408 u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)
1409 self.failUnless(self.check(u,1.),'solution is wrong.')
1410 #
1411 # solver checks (PDE system)
1412 #
1413 def test_symmetryOnIterative_System(self):
1414 A=Tensor4(0.,Function(self.domain))
1415 D=Tensor(1.,Function(self.domain))
1416 Y=Vector(self.domain.getDim(),Function(self.domain))
1417 for i in range(self.domain.getDim()):
1418 A[i,:,i,:]=kronecker(self.domain)
1419 D[i,i]+=i
1420 Y[i]+=i
1421 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1422 mypde.setValue(A=A,D=D,Y=Y)
1423 u=mypde.getSolution(verbose=self.VERBOSE)
1424 self.failUnless(self.check(u,1.),'solution is wrong.')
1425 def test_symmetryOnDirect_System(self):
1426 A=Tensor4(0.,Function(self.domain))
1427 D=Tensor(1.,Function(self.domain))
1428 Y=Vector(self.domain.getDim(),Function(self.domain))
1429 for i in range(self.domain.getDim()):
1430 A[i,:,i,:]=kronecker(self.domain)
1431 D[i,i]+=i
1432 Y[i]+=i
1433 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1434 mypde.setValue(A=A,D=D,Y=Y)
1435 mypde.setSolverMethod(mypde.DIRECT)
1436 u=mypde.getSolution(verbose=self.VERBOSE)
1437 self.failUnless(self.check(u,1.),'solution is wrong.')
1438 def test_PCG_JACOBI_System(self):
1439 A=Tensor4(0.,Function(self.domain))
1440 D=Tensor(1.,Function(self.domain))
1441 Y=Vector(self.domain.getDim(),Function(self.domain))
1442 for i in range(self.domain.getDim()):
1443 A[i,:,i,:]=kronecker(self.domain)
1444 D[i,i]+=i
1445 Y[i]+=i
1446 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1447 mypde.setValue(A=A,D=D,Y=Y)
1448 mypde.setSolverMethod(mypde.PCG,mypde.JACOBI)
1449 u=mypde.getSolution(verbose=self.VERBOSE)
1450 self.failUnless(self.check(u,1.),'solution is wrong.')
1451 def test_PCG_ILU0_System(self):
1452 A=Tensor4(0.,Function(self.domain))
1453 D=Tensor(1.,Function(self.domain))
1454 Y=Vector(self.domain.getDim(),Function(self.domain))
1455 for i in range(self.domain.getDim()):
1456 A[i,:,i,:]=kronecker(self.domain)
1457 D[i,i]+=i
1458 Y[i]+=i
1459 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1460 mypde.setValue(A=A,D=D,Y=Y)
1461 mypde.setSolverMethod(mypde.PCG,mypde.ILU0)
1462 u=mypde.getSolution(verbose=self.VERBOSE)
1463 self.failUnless(self.check(u,1.),'solution is wrong.')
1464 def test_DIRECT_System(self):
1465 A=Tensor4(0.,Function(self.domain))
1466 D=Tensor(1.,Function(self.domain))
1467 Y=Vector(self.domain.getDim(),Function(self.domain))
1468 for i in range(self.domain.getDim()):
1469 A[i,:,i,:]=kronecker(self.domain)
1470 D[i,i]+=i
1471 Y[i]+=i
1472 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1473 mypde.setValue(A=A,D=D,Y=Y)
1474 mypde.setSolverMethod(mypde.DIRECT)
1475 u=mypde.getSolution(verbose=self.VERBOSE)
1476 self.failUnless(self.check(u,1.),'solution is wrong.')
1477 def test_BICGSTAB_JACOBI_System(self):
1478 A=Tensor4(0.,Function(self.domain))
1479 D=Tensor(1.,Function(self.domain))
1480 Y=Vector(self.domain.getDim(),Function(self.domain))
1481 for i in range(self.domain.getDim()):
1482 A[i,:,i,:]=kronecker(self.domain)
1483 D[i,i]+=i
1484 Y[i]+=i
1485 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1486 mypde.setValue(A=A,D=D,Y=Y)
1487 mypde.setSolverMethod(mypde.BICGSTAB,mypde.JACOBI)
1488 u=mypde.getSolution(verbose=self.VERBOSE)
1489 self.failUnless(self.check(u,1.),'solution is wrong.')
1490 def test_BICGSTAB_ILU0_System(self):
1491 A=Tensor4(0.,Function(self.domain))
1492 D=Tensor(1.,Function(self.domain))
1493 Y=Vector(self.domain.getDim(),Function(self.domain))
1494 for i in range(self.domain.getDim()):
1495 A[i,:,i,:]=kronecker(self.domain)
1496 D[i,i]+=i
1497 Y[i]+=i
1498 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1499 mypde.setValue(A=A,D=D,Y=Y)
1500 mypde.setSolverMethod(mypde.BICGSTAB,mypde.ILU0)
1501 u=mypde.getSolution(verbose=self.VERBOSE)
1502 self.failUnless(self.check(u,1.),'solution is wrong.')
1503 def test_PRES20_JACOBI_System(self):
1504 A=Tensor4(0.,Function(self.domain))
1505 D=Tensor(1.,Function(self.domain))
1506 Y=Vector(self.domain.getDim(),Function(self.domain))
1507 for i in range(self.domain.getDim()):
1508 A[i,:,i,:]=kronecker(self.domain)
1509 D[i,i]+=i
1510 Y[i]+=i
1511 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1512 mypde.setValue(A=A,D=D,Y=Y)
1513 mypde.setSolverMethod(mypde.PRES20,mypde.JACOBI)
1514 u=mypde.getSolution(verbose=self.VERBOSE)
1515 self.failUnless(self.check(u,1.),'solution is wrong.')
1516 def test_PRES20_ILU0_System(self):
1517 A=Tensor4(0.,Function(self.domain))
1518 D=Tensor(1.,Function(self.domain))
1519 Y=Vector(self.domain.getDim(),Function(self.domain))
1520 for i in range(self.domain.getDim()):
1521 A[i,:,i,:]=kronecker(self.domain)
1522 D[i,i]+=i
1523 Y[i]+=i
1524 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1525 mypde.setValue(A=A,D=D,Y=Y)
1526 mypde.setSolverMethod(mypde.PRES20,mypde.ILU0)
1527 u=mypde.getSolution(verbose=self.VERBOSE)
1528 self.failUnless(self.check(u,1.),'solution is wrong.')
1529 def test_GMRESnoRestart_JACOBI_System(self):
1530 A=Tensor4(0.,Function(self.domain))
1531 D=Tensor(1.,Function(self.domain))
1532 Y=Vector(self.domain.getDim(),Function(self.domain))
1533 for i in range(self.domain.getDim()):
1534 A[i,:,i,:]=kronecker(self.domain)
1535 D[i,i]+=i
1536 Y[i]+=i
1537 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1538 mypde.setValue(A=A,D=D,Y=Y)
1539 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1540 # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
1541 u=mypde.getSolution(verbose=self.VERBOSE)
1542 self.failUnless(self.check(u,1.),'solution is wrong.')
1543 def test_GMRESnoRestart_ILU0_System(self):
1544 A=Tensor4(0.,Function(self.domain))
1545 D=Tensor(1.,Function(self.domain))
1546 Y=Vector(self.domain.getDim(),Function(self.domain))
1547 for i in range(self.domain.getDim()):
1548 A[i,:,i,:]=kronecker(self.domain)
1549 D[i,i]+=i
1550 Y[i]+=i
1551 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1552 mypde.setValue(A=A,D=D,Y=Y)
1553 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1554 # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
1555 u=mypde.getSolution(verbose=self.VERBOSE)
1556 self.failUnless(self.check(u,1.),'solution is wrong.')
1557 def test_GMRES_JACOBI_System(self):
1558 A=Tensor4(0.,Function(self.domain))
1559 D=Tensor(1.,Function(self.domain))
1560 Y=Vector(self.domain.getDim(),Function(self.domain))
1561 for i in range(self.domain.getDim()):
1562 A[i,:,i,:]=kronecker(self.domain)
1563 D[i,i]+=i
1564 Y[i]+=i
1565 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1566 mypde.setValue(A=A,D=D,Y=Y)
1567 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1568 u=mypde.getSolution(verbose=self.VERBOSE)
1569 self.failUnless(self.check(u,1.),'solution is wrong.')
1570 def test_GMRES_ILU0_System(self):
1571 A=Tensor4(0.,Function(self.domain))
1572 D=Tensor(1.,Function(self.domain))
1573 Y=Vector(self.domain.getDim(),Function(self.domain))
1574 for i in range(self.domain.getDim()):
1575 A[i,:,i,:]=kronecker(self.domain)
1576 D[i,i]+=i
1577 Y[i]+=i
1578 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1579 mypde.setValue(A=A,D=D,Y=Y)
1580 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1581 u=mypde.getSolution(verbose=self.VERBOSE)
1582 self.failUnless(self.check(u,1.),'solution is wrong.')
1583 def test_GMRES_truncation_restart_JACOBI_System(self):
1584 A=Tensor4(0.,Function(self.domain))
1585 D=Tensor(1.,Function(self.domain))
1586 Y=Vector(self.domain.getDim(),Function(self.domain))
1587 for i in range(self.domain.getDim()):
1588 A[i,:,i,:]=kronecker(self.domain)
1589 D[i,i]+=i
1590 Y[i]+=i
1591 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1592 mypde.setValue(A=A,D=D,Y=Y)
1593 mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)
1594 u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)
1595 self.failUnless(self.check(u,1.),'solution is wrong.')
1596 def test_GMRES_truncation_restart_ILU0_System(self):
1597 A=Tensor4(0.,Function(self.domain))
1598 D=Tensor(1.,Function(self.domain))
1599 Y=Vector(self.domain.getDim(),Function(self.domain))
1600 for i in range(self.domain.getDim()):
1601 A[i,:,i,:]=kronecker(self.domain)
1602 D[i,i]+=i
1603 Y[i]+=i
1604 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1605 mypde.setValue(A=A,D=D,Y=Y)
1606 mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)
1607 u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)
1608 self.failUnless(self.check(u,1.),'solution is wrong.')
1609
1610 class Test_LinearPDE(Test_LinearPDE_noLumping):
1611 def test_Lumping_attemptToSetA(self):
1612 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1613 try:
1614 success=True
1615 mypde.setSolverMethod(mypde.LUMPING)
1616 mypde.setValue(A=kronecker(self.domain))
1617 u=mypde.getSolution(verbose=self.VERBOSE)
1618 except ValueError:
1619 success=False
1620 self.failUnless(not success,'error should be issued')
1621 def test_Lumping_attemptToSetB(self):
1622 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1623 try:
1624 success=True
1625 mypde.setSolverMethod(mypde.LUMPING)
1626 mypde.setValue(B=kronecker(self.domain)[0])
1627 u=mypde.getSolution(verbose=self.VERBOSE)
1628 except ValueError:
1629 success=False
1630 self.failUnless(not success,'error should be issued')
1631 def test_Lumping_attemptToSetC(self):
1632 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1633 try:
1634 success=True
1635 mypde.setSolverMethod(mypde.LUMPING)
1636 mypde.setValue(C=kronecker(self.domain)[0])
1637 u=mypde.getSolution(verbose=self.VERBOSE)
1638 except ValueError:
1639 success=False
1640 self.failUnless(not success,'error should be issued')
1641
1642 def test_Lumping_attemptToSetA_reduced(self):
1643 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1644 try:
1645 success=True
1646 mypde.setSolverMethod(mypde.LUMPING)
1647 mypde.setValue(A_reduced=kronecker(self.domain))
1648 u=mypde.getSolution(verbose=self.VERBOSE)
1649 except ValueError:
1650 success=False
1651 self.failUnless(not success,'error should be issued')
1652 def test_Lumping_attemptToSetB_reduced(self):
1653 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1654 try:
1655 success=True
1656 mypde.setSolverMethod(mypde.LUMPING)
1657 mypde.setValue(B_reduced=kronecker(self.domain)[0])
1658 u=mypde.getSolution(verbose=self.VERBOSE)
1659 except ValueError:
1660 success=False
1661 self.failUnless(not success,'error should be issued')
1662 def test_Lumping_attemptToSetC_reduced(self):
1663 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1664 try:
1665 success=True
1666 mypde.setSolverMethod(mypde.LUMPING)
1667 mypde.setValue(C_reduced=kronecker(self.domain)[0])
1668 u=mypde.getSolution(verbose=self.VERBOSE)
1669 except ValueError:
1670 success=False
1671 self.failUnless(not success,'error should be issued')
1672
1673 def test_Lumping(self):
1674 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1675 mypde.setSolverMethod(mypde.LUMPING)
1676 mypde.setValue(D=1.,Y=1.)
1677 u=mypde.getSolution(verbose=self.VERBOSE)
1678 self.failUnless(self.check(u,1.),'solution is wrong.')
1679 def test_Constrained_Lumping(self):
1680 x=self.domain.getX()
1681 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1682 mypde.setSolverMethod(mypde.LUMPING)
1683 mypde.setValue(D=1.,Y=1.,q=whereZero(x[0]),r=1.)
1684 u=mypde.getSolution(verbose=self.VERBOSE)
1685 self.failUnless(self.check(u,1.),'solution is wrong.')
1686
1687 def test_Lumping_System(self):
1688 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1689 mypde.setSolverMethod(mypde.LUMPING)
1690 mypde.setValue(D=numpy.array([[1.,0.],[0.,2.]]),Y=numpy.array([1.,2.]))
1691 u=mypde.getSolution(verbose=self.VERBOSE)
1692 self.failUnless(self.check(u,numpy.ones((2,))),'solution is wrong.')
1693 def test_Constrained_Lumping_System(self):
1694 x=self.domain.getX()
1695 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1696 mypde.setSolverMethod(mypde.LUMPING)
1697 mypde.setValue(D=numpy.array([[1.,0.],[0.,2.]]),Y=numpy.array([1.,2.]), \
1698 q=whereZero(x[0])*[0.,1],r=[0.,1.])
1699 u=mypde.getSolution(verbose=self.VERBOSE)
1700 self.failUnless(self.check(u,numpy.ones((2,))),'solution is wrong.')
1701
1702 def test_Lumping_updateRHS(self):
1703 x=self.domain.getX()
1704 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1705 mypde.setSolverMethod(mypde.LUMPING)
1706 mypde.setValue(D=1.,Y=1.)
1707 u=mypde.getSolution(verbose=self.VERBOSE)
1708 self.failUnless(self.check(u,1.),'first solution is wrong.')
1709 mypde.setValue(Y=2.,q=whereZero(x[0]),r=2.)
1710 u=mypde.getSolution(verbose=self.VERBOSE)
1711 self.failUnless(self.check(u,2.),'second solution is wrong.')
1712 def test_Lumping_updateOperator(self):
1713 x=self.domain.getX()
1714 mypde=LinearPDE(self.domain,debug=self.DEBUG)
1715 mypde.setSolverMethod(mypde.LUMPING)
1716 mypde.setValue(D=1.,Y=1.)
1717 u=mypde.getSolution(verbose=self.VERBOSE)
1718 mypde.setValue(D=2.)
1719 u=mypde.getSolution(verbose=self.VERBOSE)
1720 self.failUnless(self.check(u,0.5),'second solution is wrong.')
1721
1722
1723 class Test_TransportPDE(Test_linearPDEs):
1724 N=4
1725 def test_init_useBackwardEuler(self):
1726 mypde=TransportPDE(self.domain,debug=self.DEBUG, useBackwardEuler=True)
1727 self.failUnless(mypde.useBackwardEuler()==True,'backward Euler should be used')
1728 def test_init_donntUseBackwardEuler(self):
1729 mypde=TransportPDE(self.domain,debug=self.DEBUG, useBackwardEuler=False)
1730 self.failUnless(mypde.useBackwardEuler()==False,'backward Euler should not be used')
1731 def test_setCoefficient_WithWrongName(self):
1732 mypde=TransportPDE(self.domain,debug=self.DEBUG)
1733 self.failUnlessRaises(IllegalCoefficient,mypde.setValue, ROMA=Vector(0.,FunctionOnBoundary(self.domain)))
1734
1735 def test_setCoefficient_WithIllegalFunctionSpace(self):
1736 mypde=TransportPDE(self.domain,debug=self.DEBUG)
1737 self.failUnlessRaises(IllegalCoefficientFunctionSpace,mypde.setValue,C=Vector(0.,FunctionOnBoundary(self.domain)))
1738
1739 def test_resetCoefficient_WithWrongShape(self):
1740 mypde=TransportPDE(self.domain,numEquations=2,debug=self.DEBUG)
1741 self.failUnlessRaises(IllegalCoefficientValue, mypde.setValue, C=0.)
1742
1743 def test_setInitialSolution_scalar(self):
1744 mypde=TransportPDE(self.domain,numSolutions=1,debug=self.DEBUG)
1745 mypde.setInitialSolution(1.)
1746
1747 def test_setInitialSolution_scalar_negative(self):
1748 mypde=TransportPDE(self.domain,numSolutions=1,debug=self.DEBUG)
1749 self.failUnlessRaises(RuntimeError,mypde.setInitialSolution,-1.)
1750
1751 def test_setInitialSolution_scalar_WithWrongShape(self):
1752 mypde=TransportPDE(self.domain,numSolutions=1,debug=self.DEBUG)
1753 self.failUnlessRaises(ValueError,mypde.setInitialSolution,[1.,2.])
1754
1755 def test_setInitialSolution_system(self):
1756 mypde=TransportPDE(self.domain,numSolutions=2,debug=self.DEBUG)
1757 mypde.setInitialSolution([1.,2.])
1758
1759 def test_setInitialSolution_system(self):
1760 mypde=TransportPDE(self.domain,numSolutions=2,debug=self.DEBUG)
1761 self.failUnlessRaises(RuntimeError,mypde.setInitialSolution,[-1,2.])
1762
1763 def test_setInitialSolution_system_WithWrongShape(self):
1764 mypde=TransportPDE(self.domain,numSolutions=2,debug=self.DEBUG)
1765 self.failUnlessRaises(ValueError,mypde.setInitialSolution,1.)
1766
1767
1768 def test_attemptToChangeOrderAfterDefinedCoefficient(self):
1769 mypde=TransportPDE(self.domain,debug=self.DEBUG)
1770 mypde.setValue(D=1.)
1771 self.failUnlessRaises(RuntimeError, mypde.setReducedOrderOn)
1772
1773 def test_reducedOnConfig(self):
1774 mypde=TransportPDE(self.domain,debug=self.DEBUG)
1775 mypde.setReducedOrderOn()
1776 self.failUnlessEqual((mypde.getFunctionSpaceForSolution(),mypde.getFunctionSpaceForEquation()),(ReducedSolution(self.domain),ReducedSolution(self.domain)),"reduced function spaces expected.")
1777 #
1778 # set coefficients for scalars:
1779 #
1780 def test_setCoefficient_M_Scalar(self):
1781 d=self.domain.getDim()
1782 mypde=TransportPDE(self.domain,debug=self.DEBUG)
1783 mypde.setValue(M=1.)
1784 coeff=mypde.getCoefficient("M")
1785