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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3283 - (hide annotations)
Mon Oct 18 22:39:28 2010 UTC (9 years ago) by gross
Original Path: trunk/escript/test/python/test_linearPDEs.py
File MIME type: text/x-python
File size: 191774 byte(s)
AMG reengineered, BUG is Smoother fixed.


1 gross 3103 # -*- coding: utf-8 -*-
2 ksteube 1809
3     ########################################################
4 ksteube 1312 #
5 jfenwick 2881 # Copyright (c) 2003-2010 by University of Queensland
6 ksteube 1809 # Earth Systems Science Computational Center (ESSCC)
7     # http://www.uq.edu.au/esscc
8 ksteube 1312 #
9 ksteube 1809 # Primary Business: Queensland, Australia
10     # Licensed under the Open Software License version 3.0
11     # http://www.opensource.org/licenses/osl-3.0.php
12 ksteube 1312 #
13 ksteube 1809 ########################################################
14 jgs 149
15 jfenwick 2881 __copyright__="""Copyright (c) 2003-2010 by University of Queensland
16 ksteube 1809 Earth Systems Science Computational Center (ESSCC)
17     http://www.uq.edu.au/esscc
18     Primary Business: Queensland, Australia"""
19     __license__="""Licensed under the Open Software License version 3.0
20     http://www.opensource.org/licenses/osl-3.0.php"""
21 jfenwick 2344 __url__="https://launchpad.net/escript-finley"
22 ksteube 1809
23 jgs 149 """
24     Test suite for linearPDEs class
25    
26     """
27    
28     __author__="Lutz Gross, l.gross@uq.edu.au"
29    
30 gross 2325 from esys.escript.util import Lsup,kronecker,interpolate,whereZero, outer, swap_axes
31 ksteube 1312 from esys.escript import Function,FunctionOnBoundary,FunctionOnContactZero,Solution,ReducedSolution,Vector,ContinuousFunction,Scalar, ReducedFunction,ReducedFunctionOnBoundary,ReducedFunctionOnContactZero,Data, Tensor4, Tensor
32 gross 2470 from esys.escript.linearPDEs import LinearPDE,IllegalCoefficientValue,Poisson, IllegalCoefficientFunctionSpace, TransportPDE, IllegalCoefficient, Helmholtz, LameEquation, SolverOptions
33 jfenwick 2455 import numpy
34 jgs 149 import unittest
35    
36     class Test_linearPDEs(unittest.TestCase):
37     TOL=1.e-6
38 jgs 153 SOLVER_TOL=1.e-10
39 jgs 149 DEBUG=False
40     VERBOSE=False
41     def check(self,arg,ref_arg,tol=None):
42     """
43 jfenwick 2625 checks if arg and ref_arg are nearly identical using the `Lsup`
44 jgs 149 """
45     if tol==None: tol=self.TOL
46     return Lsup(arg-ref_arg)<=tol*Lsup(ref_arg)
47    
48 gross 2325 class Test_LameEquation(Test_linearPDEs):
49    
50     def test_config(self):
51     mypde=LameEquation(self.domain,debug=self.DEBUG)
52     d=self.domain.getDim()
53 gross 2474 self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(d,d,True),"set up incorrect")
54 gross 2325
55     def test_setCoefficient_q(self):
56     mypde=LameEquation(self.domain,debug=self.DEBUG)
57     x=self.domain.getX()
58     mypde.setValue(q=x)
59    
60     q_ref=interpolate(x,Solution(self.domain))
61     self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
62     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
63     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
64     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
65     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
66     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
67     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
68     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
69     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
70     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
71     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
72     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
73     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
74     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
75     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
76     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
77 jfenwick 3259 if self.domain.supportsContactElements():
78     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
79     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
80     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
81     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
82 gross 2325 self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
83     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
84    
85     def test_setCoefficient_r(self):
86     mypde=LameEquation(self.domain,debug=self.DEBUG)
87     x=self.domain.getX()
88     mypde.setValue(r=x)
89    
90     r_ref=interpolate(x,Solution(self.domain))
91     self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
92     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
93     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
94     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
95     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
96     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
97     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
98     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
99     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
100     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
101     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
102     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
103     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
104     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
105     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
106     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
107 jfenwick 3259 if self.domain.supportsContactElements():
108     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
109     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
110     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
111     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
112 gross 2325 self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x")
113     self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
114    
115    
116     def test_setCoefficient_F(self):
117     mypde=LameEquation(self.domain,debug=self.DEBUG)
118     x=self.domain.getX()
119     mypde.setValue(F=x)
120    
121     Y_ref=interpolate(x,Function(self.domain))
122     self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
123     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
124     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
125     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
126     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
127     self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x")
128     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
129     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
130     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
131     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
132     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
133     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
134     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
135     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
136     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
137     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
138 jfenwick 3259 if self.domain.supportsContactElements():
139     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
140     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
141     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
142     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
143     self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
144 gross 2325 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
145    
146     def test_setCoefficient_f(self):
147     mypde=LameEquation(self.domain,debug=self.DEBUG)
148     x=self.domain.getX()
149     mypde.setValue(f=x)
150    
151     y_ref=interpolate(x,FunctionOnBoundary(self.domain))
152     self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
153     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
154     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
155     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
156     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
157     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
158     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
159     self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"d is not x[0]")
160     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
161     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
162     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
163     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
164     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
165     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
166     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
167     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
168 jfenwick 3259 if self.domain.supportsContactElements():
169     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
170     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
171     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
172     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
173 gross 2325 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
174     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
175    
176     def test_setCoefficient_sigma(self):
177     mypde=LameEquation(self.domain,debug=self.DEBUG)
178     x=self.domain.getX()
179     mypde.setValue(sigma=outer(x,x))
180    
181     X_ref=interpolate(outer(x,x),Function(self.domain))
182     self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")
183     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
184     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
185     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
186     self.failUnless(self.check(mypde.getCoefficient("X"),X_ref),"X is not x X x")
187     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
188     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
189     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
190     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
191     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
192     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
193     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
194     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
195     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
196     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
197     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
198 jfenwick 3259 if self.domain.supportsContactElements():
199     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
200     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
201     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
202     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
203 gross 2325 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
204     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
205    
206     def test_setCoefficient_lambda(self):
207     mypde=LameEquation(self.domain,debug=self.DEBUG)
208     x=self.domain.getX()
209     mypde.setValue(lame_lambda=x[0])
210    
211    
212     k3=kronecker(Function(self.domain))
213     k3Xk3=outer(k3,k3)
214     A_ref=x[0]*k3Xk3
215    
216     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
217     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
218     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
219     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
220     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
221     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
222     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
223     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
224     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
225     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
226     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
227     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
228     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
229     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
230     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
231     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
232 jfenwick 3259 if self.domain.supportsContactElements():
233     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
234     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
235     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
236     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
237 gross 2325 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
238     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
239    
240     def test_setCoefficient_mu(self):
241     mypde=LameEquation(self.domain,debug=self.DEBUG)
242     x=self.domain.getX()
243     mypde.setValue(lame_mu=x[0])
244    
245    
246     k3=kronecker(Function(self.domain))
247     k3Xk3=outer(k3,k3)
248     A_ref=x[0]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
249    
250     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
251     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
252     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
253     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
254     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
255     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
256     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
257     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
258     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
259     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
260     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
261     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
262     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
263     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
264     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
265     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
266 jfenwick 3259 if self.domain.supportsContactElements():
267     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
268     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced 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("q").isEmpty(),"q is not empty")
272 gross 2325 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
273    
274     def test_setCoefficient_lambdamu(self):
275     mypde=LameEquation(self.domain,debug=self.DEBUG)
276     x=self.domain.getX()
277     mypde.setValue(lame_lambda=x[0], lame_mu=x[1])
278    
279     k3=kronecker(Function(self.domain))
280     k3Xk3=outer(k3,k3)
281     A_ref=x[0]*k3Xk3+x[1]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
282    
283     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
284     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
285     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
286     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
287     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
288     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
289     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
290     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
291     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
292     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
293     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
294     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
295     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
296     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
297     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
298     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
299 jfenwick 3259 if self.domain.supportsContactElements():
300     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
301     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
302     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
303     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
304     self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
305 gross 2325 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
306    
307     def test_solve(self):
308     d=self.domain.getDim()
309     mypde=LameEquation(self.domain)
310     cf=ContinuousFunction(self.domain)
311     x=cf.getX()
312     u_ex=x
313     msk=Vector(0.,cf)
314     for i in range(d): msk[i]=whereZero(x[i])
315     mypde.setValue(q=msk,r=u_ex,lame_mu=3,lame_lambda=50,f=(2*3+50*d)*FunctionOnBoundary(self.domain).getNormal())
316    
317     u=mypde.getSolution()
318     self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
319    
320 gross 2323 class Test_Helmholtz(Test_linearPDEs):
321    
322     def test_config(self):
323     mypde=Helmholtz(self.domain,debug=self.DEBUG)
324 gross 2474 self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")
325 gross 2323 def test_setCoefficient_q(self):
326     mypde=Helmholtz(self.domain,debug=self.DEBUG)
327     x=self.domain.getX()
328     mypde.setValue(q=whereZero(x[0]))
329    
330     q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
331     A_ref=kronecker(self.domain)
332    
333     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
334     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
335     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
336     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
337     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
338     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
339     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
340     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
341     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
342     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
343     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
344     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
345     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
346     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
347     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
348     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
349 jfenwick 3259 if self.domain.supportsContactElements():
350     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
351     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
352     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
353     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
354 gross 2323 self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
355     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
356    
357     def test_setCoefficient_r(self):
358     mypde=Helmholtz(self.domain,debug=self.DEBUG)
359     x=self.domain.getX()
360     mypde.setValue(r=x[0])
361    
362     r_ref=interpolate(x[0],Solution(self.domain))
363     A_ref=kronecker(self.domain)
364     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
365     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
366     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
367     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
368     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
369     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
370     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
371     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
372     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
373     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
374     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
375     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
376     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
377     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
378     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
379     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
380 jfenwick 3259 if self.domain.supportsContactElements():
381     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
382     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
383     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
384     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
385     self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x[0]")
386 gross 2323 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
387    
388    
389     def test_setCoefficient_f(self):
390     mypde=Helmholtz(self.domain,debug=self.DEBUG)
391     x=self.domain.getX()
392     mypde.setValue(f=x[0])
393    
394     Y_ref=interpolate(x[0],Function(self.domain))
395     A_ref=kronecker(self.domain)
396     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
397     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
398     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
399     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
400     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
401     self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
402     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
403     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
404     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
405     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
406     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
407     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
408     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
409     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
410     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
411     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
412 jfenwick 3259 if self.domain.supportsContactElements():
413     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
414     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
415     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
416     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
417 gross 2323 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
418     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
419    
420     def test_setCoefficient_alpha(self):
421     mypde=Helmholtz(self.domain,debug=self.DEBUG)
422     x=self.domain.getX()
423     mypde.setValue(alpha=x[0])
424    
425     d_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
426     A_ref=kronecker(self.domain)
427     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
428     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
429     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
430     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
431     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
432     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
433     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
434     self.failUnless(self.check(mypde.getCoefficient("d"),d_ref),"d is not x[0]")
435     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
436     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
437     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
438     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
439     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
440     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")
441     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
442     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
443 jfenwick 3259 if self.domain.supportsContactElements():
444     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
445     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
446     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
447     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
448 gross 2323 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
449     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
450    
451     def test_setCoefficient_g(self):
452     mypde=Helmholtz(self.domain,debug=self.DEBUG)
453     x=self.domain.getX()
454     mypde.setValue(g=x[0])
455    
456     y_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
457     A_ref=kronecker(self.domain)
458     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
459     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
460     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
461     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
462     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
463     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
464     self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"y is not x[0]")
465     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
466     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
467     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
468     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
469     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
470     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
471     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
472     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
473     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
474 jfenwick 3259 if self.domain.supportsContactElements():
475     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
476     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
477     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
478     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
479 gross 2323 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
480     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
481    
482     def test_setCoefficient_omega(self):
483     mypde=Helmholtz(self.domain,debug=self.DEBUG)
484     x=self.domain.getX()
485     mypde.setValue(omega=x[0])
486    
487     D_ref=interpolate(x[0],Function(self.domain))
488     A_ref=kronecker(self.domain)
489     self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
490     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
491     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
492     self.failUnless(self.check(mypde.getCoefficient("D"),D_ref),"D is not x[0]")
493     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
494     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
495     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
496     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
497     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
498     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
499     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
500     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
501     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
502     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
503     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
504     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
505 jfenwick 3259 if self.domain.supportsContactElements():
506     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
507     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
508     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
509     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
510     self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
511 gross 2323 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
512    
513     def test_solve(self):
514     d=self.domain.getDim()
515     cf=ContinuousFunction(self.domain)
516     u_ex=Scalar(1.,cf)
517     mypde=Helmholtz(self.domain)
518     mypde.setValue(f=3,omega=3,alpha=2,g=2)
519     u=mypde.getSolution()
520     self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
521    
522 jgs 149 class Test_Poisson(Test_linearPDEs):
523    
524     def test_config(self):
525     mypde=Poisson(self.domain,debug=self.DEBUG)
526 gross 2474 self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")
527 jgs 149 def test_setCoefficient_q(self):
528     mypde=Poisson(self.domain,debug=self.DEBUG)
529     x=self.domain.getX()
530 gross 304 q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
531 jgs 149 A_ref=kronecker(self.domain)
532 gross 304 mypde.setValue(q=whereZero(x[0]))
533 gross 1841 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
534     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
535     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
536     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
537     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
538     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
539     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
540     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
541     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
542     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
543     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
544     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
545     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
546     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
547     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
548     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
549 jfenwick 3259 if self.domain.supportsContactElements():
550     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
551     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
552     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
553     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
554     self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
555 gross 1841 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
556 jgs 149 def test_setCoefficient_f(self):
557     mypde=Poisson(self.domain,debug=self.DEBUG)
558     x=self.domain.getX()
559     Y_ref=interpolate(x[0],Function(self.domain))
560     A_ref=kronecker(self.domain)
561     mypde.setValue(f=x[0])
562 gross 1841 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
563     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
564     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
565     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
566     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
567     self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
568     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
569     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
570     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
571     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
572     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
573     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
574     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
575     self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")
576     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
577     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
578 jfenwick 3259 if self.domain.supportsContactElements():
579     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
580     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
581     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
582     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
583 gross 1841 self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
584     self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
585 gross 1072 def test_setCoefficient_f_reduced(self):
586     mypde=Poisson(self.domain,debug=self.DEBUG)
587     x=self.domain.getX()
588     Y_ref=interpolate(x[0],ReducedFunction(self.domain))
589     A_ref=kronecker(self.domain)
590     mypde.setValue(f_reduced=x[0])
591 gross 1841 self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
592     self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")
593     self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")
594     self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")
595     self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")
596     self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")
597     self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")
598     self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")
599     self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")
600     self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")
601     self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")
602     self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")
603     self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")
604     self.failUnless(self.check(mypde.getCoefficient("Y_reduced"),Y_ref),"Y_reduced is not x[0]")
605     self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")
606     self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")
607 jfenwick 3259 if self.domain.supportsContactElements():
608     self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")
609     self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")
610     self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")
611     self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")
612     self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")
613 gross 1841 self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")
614 jgs 149 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 gross 304 msk+=whereZero(x[i])
626 jgs 149 #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 jfenwick 3259 print u,u_ex, self.TOL
640 jgs 149 self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")
641    
642 gross 855 class Test_LinearPDE_noLumping(Test_linearPDEs):
643 jgs 149 N=4
644 gross 2470 def test_SolverOptions(self):
645     so=SolverOptions()
646 artak 2828
647 artak 2835 self.failUnless(so.getSmoother() == 28, "initial Smoother is wrong.")
648 artak 2828 self.failUnlessRaises(ValueError,so.setSmoother,-1)
649     so.setSmoother(so.GAUSS_SEIDEL)
650     self.failUnless(so.getSmoother() == 28, "Gauss-Seidel smoother is not set.")
651 artak 2835 so.setSmoother(so.JACOBI)
652     self.failUnless(so.getSmoother() == 10, "Jacobi smoother is not set.")
653 gross 2470
654 gross 3283 self.failUnless(so.getLevelMax() == 100, "initial LevelMax is wrong.")
655 gross 2470 self.failUnlessRaises(ValueError,so.setLevelMax,-1)
656 artak 2688 so.setLevelMax(20)
657     self.failUnless(so.getLevelMax() == 20, "LevelMax is wrong.")
658 gross 2470
659 artak 2688 self.failUnless(so.getCoarseningThreshold() == 0.25, "initial CoarseningThreshold is wrong.")
660 gross 2470 self.failUnlessRaises(ValueError,so.setCoarseningThreshold,-1)
661     so.setCoarseningThreshold(0.1)
662     self.failUnless(so.getCoarseningThreshold() == 0.1, "CoarseningThreshold is wrong.")
663 artak 2524
664     self.failUnless(so.getMinCoarseMatrixSize() == 500, "initial Minimum Coarse Matrix Size is wrong.")
665     self.failUnlessRaises(ValueError,so.setMinCoarseMatrixSize,-1)
666     so.setMinCoarseMatrixSize(1000)
667     self.failUnless(so.getMinCoarseMatrixSize() == 1000, "Minimum Coarse Matrix Size is wrong.")
668 gross 2470
669 gross 3160 self.failUnless(so.getNumSweeps() == 1, "initial Sweeps is wrong.")
670 gross 2470 self.failUnlessRaises(ValueError,so.setNumSweeps,-1)
671     so.setNumSweeps(3)
672     self.failUnless(so.getNumSweeps() == 3, "Sweeps is wrong.")
673    
674     self.failUnless(so.getNumPreSweeps() == 2, "initial PreSweeps is wrong.")
675     self.failUnlessRaises(ValueError,so.setNumPreSweeps,-1)
676     so.setNumPreSweeps(4)
677     self.failUnless(so.getNumPreSweeps() == 4, "PreSweeps is wrong.")
678    
679     self.failUnless(so.getNumPostSweeps() == 2, "initial PreSweeps is wrong.")
680     self.failUnlessRaises(ValueError,so.setNumPostSweeps,-1)
681     so.setNumPostSweeps(5)
682     self.failUnless(so.getNumPostSweeps() == 5, "PreSweeps is wrong.")
683    
684     self.failUnless(so.getTolerance() == 1.e-8, "initial Tolerance is wrong.")
685     self.failUnlessRaises(ValueError,so.setTolerance,-1)
686     so.setTolerance(0.2)
687     self.failUnless(so.getTolerance() == 0.2, "Tolerance is wrong.")
688    
689     self.failUnless(so.getAbsoluteTolerance() == 0., "initial AbsoluteTolerance is wrong.")
690     self.failUnlessRaises(ValueError,so.setAbsoluteTolerance,-1)
691     so.setAbsoluteTolerance(0.3)
692     self.failUnless(so.getAbsoluteTolerance() == 0.3, "AbsoluteTolerance is wrong.")
693    
694     self.failUnless(so.getInnerTolerance() == 0.9, "initial InnerTolerance is wrong.")
695     self.failUnlessRaises(ValueError,so.setInnerTolerance,-1)
696     so.setInnerTolerance(0.4)
697     self.failUnless(so.getInnerTolerance() == 0.4, "InnerTolerance is wrong.")
698    
699     self.failUnless(so.getDropTolerance() == 0.01, "initial DropTolerance is wrong.")
700     self.failUnlessRaises(ValueError,so.setDropTolerance,-1)
701     so.setDropTolerance(0.5)
702     self.failUnless(so.getDropTolerance() == 0.5, "DropDropTolerance is wrong.")
703    
704     self.failUnless(so.getDropStorage() == 2., "initial DropStorage is wrong.")
705     self.failUnlessRaises(ValueError,so.setDropStorage,-1)
706     so.setDropStorage(10)
707     self.failUnless(so.getDropStorage() == 10, "DropStorage is wrong.")
708    
709     self.failUnless(so.getRelaxationFactor() == 0.3, "initial RelaxationFactor is wrong.")
710     self.failUnlessRaises(ValueError,so.setRelaxationFactor,-1)
711     so.setRelaxationFactor(0.1)
712     self.failUnless(so.getRelaxationFactor() == 0.1, "Relaxation is wrong.")
713    
714    
715 gross 2480 self.failUnless(so.getIterMax() == 100000, "initial IterMax is wrong.")
716 gross 2470 self.failUnlessRaises(ValueError,so.setIterMax,0)
717     so.setIterMax(11)
718     self.failUnless(so.getIterMax() == 11, "IterMax is wrong.")
719    
720     self.failUnless(so.getInnerIterMax() == 10, "initial InnerIterMax is wrong.")
721     self.failUnlessRaises(ValueError,so.setInnerIterMax,0)
722     so.setInnerIterMax(12)
723     self.failUnless(so.getInnerIterMax() == 12, "InnerIterMax is wrong.")
724    
725     self.failUnless(so.getTruncation() == 20, "initial Truncation is wrong.")
726     self.failUnlessRaises(ValueError,so.setTruncation,0)
727     so.setTruncation(13)
728     self.failUnless(so.getTruncation() == 13, "Truncation is wrong.")
729    
730     self.failUnless(so.getRestart() == None, "initial Truncation is wrong.")
731     self.failUnlessRaises(ValueError,so.setTruncation,0)
732     so.setRestart(14)
733     self.failUnless(so.getRestart() == 14, "Truncation is wrong.")
734     so.setRestart(None)
735     self.failUnless(so.getRestart() == None, "Truncation is wrong.")
736 gross 2474
737     self.failUnless(not so.isVerbose(), "initial verbosity flag is wrong.")
738     so.setVerbosityOn()
739     self.failUnless(so.isVerbose(), "verbosity (1) flag is wrong.")
740     so.setVerbosityOff()
741     self.failUnless(not so.isVerbose(), "verbosity (2) flag is wrong.")
742     so.setVerbosity(verbose=True)
743     self.failUnless(so.isVerbose(), "verbosity (3) flag is wrong.")
744     so.setVerbosity(verbose=False)
745     self.failUnless(not so.isVerbose(), "verbosity (4) flag is wrong.")
746    
747 gross 2470 self.failUnless(not so.isSymmetric(), "initial symmetry flag is wrong.")
748     so.setSymmetryOn()
749     self.failUnless(so.isSymmetric(), "symmetry (1) flag is wrong.")
750     so.setSymmetryOff()
751     self.failUnless(not so.isSymmetric(), "symmetry (2) flag is wrong.")
752     so.setSymmetry(flag=True)
753     self.failUnless(so.isSymmetric(), "symmetry (3) flag is wrong.")
754     so.setSymmetry(flag=False)
755     self.failUnless(not so.isSymmetric(), "symmetry (4) flag is wrong.")
756    
757     self.failUnless(so.adaptInnerTolerance(), "initial InnerToleranceAdaption flag is wrong.")
758     so.setInnerToleranceAdaptionOn()
759     self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (1) flag is wrong.")
760     so.setInnerToleranceAdaptionOff()
761     self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (2) flag is wrong.")
762     so.setInnerToleranceAdaption(adapt=True)
763     self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (3) flag is wrong.")
764     so.setInnerToleranceAdaption(adapt=False)
765     self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (4) flag is wrong.")
766    
767     self.failUnless(not so.acceptConvergenceFailure(), "initial acceptConvergenceFailure flag is wrong.")
768     so.setAcceptanceConvergenceFailureOn()
769     self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (1) flag is wrong.")
770     so.setAcceptanceConvergenceFailureOff()
771     self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (2) flag is wrong.")
772     so.setAcceptanceConvergenceFailure(accept=True)
773     self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (3) flag is wrong.")
774     so.setAcceptanceConvergenceFailure(accept=False)
775     self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (4) flag is wrong.")
776    
777     self.failUnless(so.getReordering() == 30, "initial Reordering is wrong.")
778     self.failUnlessRaises(ValueError,so.setReordering,-1)
779     so.setReordering(so.NO_REORDERING)
780     self.failUnless(so.getReordering() == 17, "NO_REORDERING is not set.")
781     so.setReordering(so.MINIMUM_FILL_IN)
782     self.failUnless(so.getReordering() == 18, "MINIMUM_FILL_IN is not set.")
783     so.setReordering(so.NESTED_DISSECTION)
784     self.failUnless(so.getReordering() == 19, "NESTED_DISSECTION is not set.")
785     so.setReordering(so.DEFAULT_REORDERING)
786     self.failUnless(so.getReordering() == 30, "DEFAULT_REORDERING is not set.")
787    
788     self.failUnless(so.getPackage() == 0, "initial solver package is wrong.")
789     self.failUnlessRaises(ValueError,so.setPackage,-1)
790     so.setPackage(so.PASO)
791     self.failUnless(so.getPackage() == 21, "PASO is not set.")
792     so.setPackage(so.SUPER_LU)
793     self.failUnless(so.getPackage() == 31, "SUPER_LU is not set.")
794     so.setPackage(so.PASTIX)
795     self.failUnless(so.getPackage() == 32, "PASTIX is not set.")
796     so.setPackage(so.MKL)
797     self.failUnless(so.getPackage() == 15, "MKL is not set.")
798     so.setPackage(so.UMFPACK)
799     self.failUnless(so.getPackage() == 16, "UMFPACK is not set.")
800     so.setPackage(so.TRILINOS)
801     self.failUnless(so.getPackage() == 24, "TRILINOS is not set.")
802    
803     self.failUnless(so.getSolverMethod() == 0, "initial SolverMethod is wrong.")
804     self.failUnlessRaises(ValueError,so.setSolverMethod,-1)
805     so.setSolverMethod(so.DIRECT)
806     self.failUnless(so.getSolverMethod() == 1, "DIRECT is not set.")
807     so.setSolverMethod(so.CHOLEVSKY)
808     self.failUnless(so.getSolverMethod() == 2, "CHOLEVSKY is not set.")
809     so.setSolverMethod(so.PCG)
810     self.failUnless(so.getSolverMethod() == 3, "PCG is not set.")
811     so.setSolverMethod(so.CR)
812     self.failUnless(so.getSolverMethod() == 4, "CR is not set.")
813     so.setSolverMethod(so.CGS)
814     self.failUnless(so.getSolverMethod() == 5, "CGS is not set.")
815     so.setSolverMethod(so.BICGSTAB)
816     self.failUnless(so.getSolverMethod() == 6, "BICGSTAB is not set.")
817     so.setSolverMethod(so.GMRES)
818     self.failUnless(so.getSolverMethod() == 11, "GMRES is not set.")
819     so.setSolverMethod(so.PRES20)
820     self.failUnless(so.getSolverMethod() == 12, "PRES20 is not set.")
821     so.setSolverMethod(so.LUMPING)
822     self.failUnless(so.getSolverMethod() == 13, "LUMPING is not set.")
823     so.setSolverMethod(so.ITERATIVE)
824     self.failUnless(so.getSolverMethod() == 20, "ITERATIVE is not set.")
825     so.setSolverMethod(so.NONLINEAR_GMRES)
826     self.failUnless(so.getSolverMethod() == 25, "NONLINEAR_GMRES is not set.")
827     so.setSolverMethod(so.TFQMR)
828     self.failUnless(so.getSolverMethod() == 26, "TFQMR is not set.")
829     so.setSolverMethod(so.MINRES)
830     self.failUnless(so.getSolverMethod() == 27, "MINRES is not set.")
831     so.setSolverMethod(so.DEFAULT)
832     self.failUnless(so.getSolverMethod() == 0, "DEFAULT is not set.")
833    
834     self.failUnless(so.getPreconditioner() == 10, "initial Preconditioner is wrong.")
835     self.failUnlessRaises(ValueError,so.setPreconditioner,-1)
836     so.setPreconditioner(so.ILU0)
837     self.failUnless(so.getPreconditioner() == 8, "ILU0 is not set.")
838     so.setPreconditioner(so.ILUT)
839     self.failUnless(so.getPreconditioner() == 9, "ILUT is not set.")
840     so.setPreconditioner(so.JACOBI)
841     self.failUnless(so.getPreconditioner() == 10, "JACOBI is not set.")
842     so.setPreconditioner(so.AMG)
843     self.failUnless(so.getPreconditioner() == 22, "AMG is not set.")
844     so.setPreconditioner(so.REC_ILU)
845     self.failUnless(so.getPreconditioner() == 23, "REC_ILU is not set.")
846     so.setPreconditioner(so.GAUSS_SEIDEL)
847     self.failUnless(so.getPreconditioner() == 28, "GAUSS_SEIDEL is not set.")
848     so.setPreconditioner(so.RILU)
849     self.failUnless(so.getPreconditioner() == 29, "RILU is not set.")
850 gross 3103 so.setPreconditioner(so.AMLI)
851     self.failUnless(so.getPreconditioner() == 38, "AMLI is not set.")
852 gross 2470 so.setPreconditioner(so.NO_PRECONDITIONER)
853     self.failUnless(so.getPreconditioner() == 36, "NO_PRECONDITIONER is not set.")
854    
855     self.failUnless(so.getCoarsening() == 0, "initial Coarseningr is wrong.")
856     self.failUnlessRaises(ValueError,so.setCoarsening,-1)
857     so.setCoarsening(so.YAIR_SHAPIRA_COARSENING)
858     self.failUnless(so.getCoarsening() == 33, "YAIR_SHAPIRA_COARSENING is not set.")
859     so.setCoarsening(so.RUGE_STUEBEN_COARSENING)
860     self.failUnless(so.getCoarsening() == 34, "RUGE_STUEBEN_COARSENING is not set.")
861     so.setCoarsening(so.AGGREGATION_COARSENING)
862     self.failUnless(so.getCoarsening() == 35, "AGREGATION_COARSENING is not set.")
863 artak 2817 so.setCoarsening(so.STANDARD_COARSENING)
864     self.failUnless(so.getCoarsening() == 39, "STANDARD_COARSENING is not set.")
865 gross 2470 so.setCoarsening(so.DEFAULT)
866     self.failUnless(so.getCoarsening() == 0, "DEFAULT is not set.")
867    
868     self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
869     self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
870     self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")
871     self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
872     self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
873     self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")
874     self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")
875     self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
876     self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
877     self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
878     self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
879    
880     so._updateDiagnostics("num_iter",1)
881     so._updateDiagnostics("num_inner_iter",2)
882     so._updateDiagnostics("time",3)
883     so._updateDiagnostics("set_up_time",4)
884     so._updateDiagnostics("residual_norm",5)
885     so._updateDiagnostics("converged",True)
886    
887     self.failUnless(so.getDiagnostics("num_iter") == 1, "num_iter is wrong.")
888     self.failUnless(so.getDiagnostics("num_inner_iter") == 2, "num_inner_iter is wrong.")
889     self.failUnless(so.getDiagnostics("time") == 3, "time is wrong.")
890     self.failUnless(so.getDiagnostics("set_up_time") == 4, "set_up_time is wrong.")
891     self.failUnless(so.getDiagnostics("residual_norm") == 5, "residual_norm is wrong.")
892     self.failUnless(so.getDiagnostics("converged"), "converged is wrong.")
893     self.failUnless(so.hasConverged(), "convergence flag is wrong.")
894     self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
895     self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
896     self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")
897     self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")
898    
899     so.resetDiagnostics()
900     self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
901     self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
902     self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")
903     self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
904     self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
905     self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")
906     self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong")
907     self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
908     self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
909     self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")
910     self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")
911    
912     so._updateDiagnostics("num_iter",10)
913     so._updateDiagnostics("num_inner_iter",20)
914     so._updateDiagnostics("time",30)
915     so._updateDiagnostics("set_up_time",40)
916     so._updateDiagnostics("residual_norm",50)
917     so._updateDiagnostics("converged",False)
918    
919     self.failUnless(so.getDiagnostics("num_iter") == 10, "num_iter is wrong.")
920     self.failUnless(so.getDiagnostics("num_inner_iter") == 20, "num_inner_iter is wrong.")
921     self.failUnless(so.getDiagnostics("time") == 30, "time is wrong.")
922     self.failUnless(so.getDiagnostics("set_up_time") == 40, "set_up_time is wrong.")
923     self.failUnless(so.getDiagnostics("residual_norm") == 50, "residual_norm is wrong.")
924     self.failUnless(not so.getDiagnostics("converged"), "converged is wrong.")
925     self.failUnless(not so.hasConverged(), "convergence flag is wrong.")
926     self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 22, "cum_num_inner_iter is wrong.")
927     self.failUnless(so.getDiagnostics("cum_num_iter") == 11, "cum_num_iter is wrong.")
928     self.failUnless(so.getDiagnostics("cum_time") ==33, "cum_time is wrong.")
929     self.failUnless(so.getDiagnostics("cum_set_up_time") == 44, "cum_set_up_time is wrong.")
930    
931     so.resetDiagnostics(all=True)
932     self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
933     self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
934     self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")
935     self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
936     self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
937     self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")
938     self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")
939     self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
940     self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
941     self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
942     self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
943    
944 jgs 149 def test_setCoefficient_WithIllegalFunctionSpace(self):
945     mypde=LinearPDE(self.domain,debug=self.DEBUG)
946 gross 2474 self.failUnlessRaises(IllegalCoefficientFunctionSpace, mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))
947 gross 1859
948     def test_setCoefficient_WithWrongName(self):
949     mypde=LinearPDE(self.domain,debug=self.DEBUG)
950     self.failUnlessRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)
951    
952 jgs 149 def test_resetCoefficient_WithWrongShape(self):
953     mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)
954 gross 1859 self.failUnlessRaises(IllegalCoefficientValue, mypde.setValue, C=0.)
955    
956 jgs 149 def test_reducedOn(self):
957     mypde=LinearPDE(self.domain,debug=self.DEBUG)
958     x=self.domain.getX()
959     mypde.setReducedOrderOn()
960     mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])
961     u=mypde.getSolution()
962     self.failUnless(self.check(u,1.),'solution is wrong.')
963    
964     def test_attemptToChangeOrderAfterDefinedCoefficient(self):
965     mypde=LinearPDE(self.domain,debug=self.DEBUG)
966     mypde.setValue(D=1.)
967 gross 2474 self.failUnlessRaises(RuntimeError, mypde.setReducedOrderOn)
968 jgs 149
969     def test_reducedOnConfig(self):
970     mypde=LinearPDE(self.domain,debug=self.DEBUG)
971     mypde.setReducedOrderOn()
972 gross 2474 self.failUnlessEqual((mypde.getFunctionSpaceForSolution(), mypde.getFunctionSpaceForEquation()),(ReducedSolution(self.domain),ReducedSolution(self.domain)),"reduced function spaces expected.")
973 jgs 149 #
974     # set coefficients for scalars:
975     #
976     def test_setCoefficient_A_Scalar(self):
977     d=self.domain.getDim()
978     mypde=LinearPDE(self.domain,debug=self.DEBUG)
979 jfenwick 2455 mypde.setValue(A=numpy.ones((d,d)))
980 gross 1841 coeff=mypde.getCoefficient("A")
981 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),Function(self.domain),1,1))
982 jgs 149 def test_setCoefficient_B_Scalar(self):
983     d=self.domain.getDim()
984     mypde=LinearPDE(self.domain,debug=self.DEBUG)
985 jfenwick 2455 mypde.setValue(B=numpy.ones((d,)))
986 gross 1841 coeff=mypde.getCoefficient("B")
987 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))
988 jgs 149 def test_setCoefficient_C_Scalar(self):
989     d=self.domain.getDim()
990     mypde=LinearPDE(self.domain,debug=self.DEBUG)
991 jfenwick 2455 mypde.setValue(C=numpy.ones((d,)))
992 gross 1841 coeff=mypde.getCoefficient("C")
993 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))
994 jgs 149 def test_setCoefficient_D_Scalar(self):
995     d=self.domain.getDim()
996     mypde=LinearPDE(self.domain,debug=self.DEBUG)
997     mypde.setValue(D=1.)
998 gross 1841 coeff=mypde.getCoefficient("D")
999 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),Function(self.domain),1,1))
1000 jgs 149 def test_setCoefficient_X_Scalar(self):
1001     d=self.domain.getDim()
1002     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1003 jfenwick 2455 mypde.setValue(X=numpy.ones((d,)))
1004 gross 1841 coeff=mypde.getCoefficient("X")
1005 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),Function(self.domain),1))
1006 jgs 149 def test_setCoefficient_Y_Scalar(self):
1007     d=self.domain.getDim()
1008     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1009     mypde.setValue(Y=1.)
1010 gross 1841 coeff=mypde.getCoefficient("Y")
1011 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),Function(self.domain),1))
1012 jgs 149 def test_setCoefficient_y_Scalar(self):
1013     d=self.domain.getDim()
1014     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1015     mypde.setValue(y=1.)
1016 gross 1841 coeff=mypde.getCoefficient("y")
1017 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))
1018 jgs 149 def test_setCoefficient_d_Scalar(self):
1019     d=self.domain.getDim()
1020     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1021     mypde.setValue(d=1.)
1022 gross 1841 coeff=mypde.getCoefficient("d")
1023 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1,1))
1024 jgs 149 def test_setCoefficient_d_contact_Scalar(self):
1025 jfenwick 3259 if self.domain.supportsContactElements():
1026     d=self.domain.getDim()
1027     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1028     mypde.setValue(d_contact=1.)
1029     coeff=mypde.getCoefficient("d_contact")
1030     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))
1031 jgs 149 def test_setCoefficient_y_contact_Scalar(self):
1032     d=self.domain.getDim()
1033 jfenwick 3259 if self.domain.supportsContactElements():
1034     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1035     mypde.setValue(y_contact=1.)
1036     coeff=mypde.getCoefficient("y_contact")
1037     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))
1038 gross 1072 def test_setCoefficient_A_reduced_Scalar(self):
1039     d=self.domain.getDim()
1040     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1041 jfenwick 2455 mypde.setValue(A_reduced=numpy.ones((d,d)))
1042 gross 1841 coeff=mypde.getCoefficient("A_reduced")
1043 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))
1044 gross 1072 def test_setCoefficient_B_reduced_Scalar(self):
1045     d=self.domain.getDim()
1046     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1047 jfenwick 2455 mypde.setValue(B_reduced=numpy.ones((d,)))
1048 gross 1841 coeff=mypde.getCoefficient("B_reduced")
1049 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
1050 gross 1072 def test_setCoefficient_C_reduced_Scalar(self):
1051     d=self.domain.getDim()
1052     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1053 jfenwick 2455 mypde.setValue(C_reduced=numpy.ones((d,)))
1054 gross 1841 coeff=mypde.getCoefficient("C_reduced")
1055 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
1056 gross 1072 def test_setCoefficient_D_reduced_Scalar(self):
1057     d=self.domain.getDim()
1058     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1059     mypde.setValue(D_reduced=1.)
1060 gross 1841 coeff=mypde.getCoefficient("D_reduced")
1061 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))
1062 gross 1072 def test_setCoefficient_X_reduced_Scalar(self):
1063     d=self.domain.getDim()
1064     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1065 jfenwick 2455 mypde.setValue(X_reduced=numpy.ones((d,)))
1066 gross 1841 coeff=mypde.getCoefficient("X_reduced")
1067 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
1068 gross 1072 def test_setCoefficient_Y_reduced_Scalar(self):
1069     d=self.domain.getDim()
1070     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1071     mypde.setValue(Y_reduced=1.)
1072 gross 1841 coeff=mypde.getCoefficient("Y_reduced")
1073 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
1074 gross 1072 def test_setCoefficient_y_reduced_Scalar(self):
1075     d=self.domain.getDim()
1076     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1077     mypde.setValue(y_reduced=1.)
1078 gross 1841 coeff=mypde.getCoefficient("y_reduced")
1079 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
1080 gross 1072 def test_setCoefficient_d_reduced_Scalar(self):
1081     d=self.domain.getDim()
1082     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1083     mypde.setValue(d_reduced=1.)
1084 gross 1841 coeff=mypde.getCoefficient("d_reduced")
1085 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))
1086 gross 1072 def test_setCoefficient_d_contact_reduced_Scalar(self):
1087 jfenwick 3259 if self.domain.supportsContactElements():
1088     d=self.domain.getDim()
1089     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1090     mypde.setValue(d_contact_reduced=1.)
1091     coeff=mypde.getCoefficient("d_contact_reduced")
1092     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
1093 gross 1072 def test_setCoefficient_y_contact_reduced_Scalar(self):
1094 jfenwick 3259 if self.domain.supportsContactElements():
1095     d=self.domain.getDim()
1096     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1097     mypde.setValue(y_contact_reduced=1.)
1098     coeff=mypde.getCoefficient("y_contact_reduced")
1099     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
1100 jgs 149 def test_setCoefficient_r_Scalar(self):
1101     d=self.domain.getDim()
1102     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1103     mypde.setValue(r=1.)
1104 gross 1841 coeff=mypde.getCoefficient("r")
1105 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
1106 jgs 149 def test_setCoefficient_q_Scalar(self):
1107     d=self.domain.getDim()
1108     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1109     mypde.setValue(q=1.)
1110 gross 1841 coeff=mypde.getCoefficient("q")
1111 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
1112 jgs 149 def test_setCoefficient_r_Scalar_reducedOn(self):
1113     d=self.domain.getDim()
1114     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1115     mypde.setReducedOrderOn()
1116     mypde.setValue(r=1.)
1117 gross 1841 coeff=mypde.getCoefficient("r")
1118 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
1119 jgs 149 def test_setCoefficient_q_Scalar_reducedOn(self):
1120     d=self.domain.getDim()
1121     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1122     mypde.setReducedOrderOn()
1123     mypde.setValue(q=1.)
1124 gross 1841 coeff=mypde.getCoefficient("q")
1125 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
1126 jgs 149
1127 gross 1072 def test_setCoefficient_A_reduced_Scalar_usingA(self):
1128     d=self.domain.getDim()
1129     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1130 jfenwick 2455 mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))
1131 gross 1841 coeff=mypde.getCoefficient("A_reduced")
1132 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))
1133 gross 1072 def test_setCoefficient_B_reduced_Scalar_usingB(self):
1134     d=self.domain.getDim()
1135     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1136 jfenwick 2455 mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1137 gross 1841 coeff=mypde.getCoefficient("B_reduced")
1138 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
1139 gross 1072 def test_setCoefficient_C_reduced_Scalar_usingC(self):
1140     d=self.domain.getDim()
1141     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1142 jfenwick 2455 mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1143 gross 1841 coeff=mypde.getCoefficient("C_reduced")
1144 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
1145 gross 1072 def test_setCoefficient_D_reduced_Scalar_usingD(self):
1146     d=self.domain.getDim()
1147     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1148     mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))
1149 gross 1841 coeff=mypde.getCoefficient("D_reduced")
1150 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))
1151 gross 1072 def test_setCoefficient_X_reduced_Scalar_usingX(self):
1152     d=self.domain.getDim()
1153     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1154 jfenwick 2455 mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1155 gross 1841 coeff=mypde.getCoefficient("X_reduced")
1156 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
1157 gross 1072 def test_setCoefficient_Y_reduced_Scalar_usingY(self):
1158     d=self.domain.getDim()
1159     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1160     mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))
1161 gross 1841 coeff=mypde.getCoefficient("Y_reduced")
1162 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
1163 gross 1072 def test_setCoefficient_y_reduced_Scalar_using_y(self):
1164     d=self.domain.getDim()
1165     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1166     mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1167 gross 1841 coeff=mypde.getCoefficient("y_reduced")
1168 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
1169 gross 1072 def test_setCoefficient_d_reduced_Scalar_using_d(self):
1170     d=self.domain.getDim()
1171     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1172     mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1173 gross 1841 coeff=mypde.getCoefficient("d_reduced")
1174 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))
1175 gross 1072 def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):
1176 jfenwick 3259 if self.domain.supportsContactElements():
1177     d=self.domain.getDim()
1178     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1179     mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1180     coeff=mypde.getCoefficient("d_contact_reduced")
1181     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
1182 gross 1072 def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):
1183 jfenwick 3259 if self.domain.supportsContactElements():
1184     d=self.domain.getDim()
1185     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1186     mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1187     coeff=mypde.getCoefficient("y_contact_reduced")
1188     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
1189 jgs 149 #
1190     # set coefficients for systems:
1191     #
1192     def test_setCoefficient_A_System(self):
1193     d=self.domain.getDim()
1194     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1195 jfenwick 2455 mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))
1196 gross 1841 coeff=mypde.getCoefficient("A")
1197 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),Function(self.domain),self.N,self.N))
1198 jgs 149 def test_setCoefficient_B_System(self):
1199     d=self.domain.getDim()
1200     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1201 jfenwick 2455 mypde.setValue(B=numpy.ones((self.N,d,self.N)))
1202 gross 1841 coeff=mypde.getCoefficient("B")
1203 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),Function(self.domain),self.N,self.N))
1204 jgs 149 def test_setCoefficient_C_System(self):
1205     d=self.domain.getDim()
1206     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1207 jfenwick 2455 mypde.setValue(C=numpy.ones((self.N,self.N,d)))
1208 gross 1841 coeff=mypde.getCoefficient("C")
1209 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),Function(self.domain),self.N,self.N))
1210 jgs 149 def test_setCoefficient_D_System(self):
1211     d=self.domain.getDim()
1212     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1213 jfenwick 2455 mypde.setValue(D=numpy.ones((self.N,self.N)))
1214 gross 1841 coeff=mypde.getCoefficient("D")
1215 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),Function(self.domain),self.N,self.N))
1216 jgs 149 def test_setCoefficient_X_System(self):
1217     d=self.domain.getDim()
1218     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1219 jfenwick 2455 mypde.setValue(X=numpy.ones((self.N,d)))
1220 gross 1841 coeff=mypde.getCoefficient("X")
1221 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),Function(self.domain),self.N))
1222 jgs 149 def test_setCoefficient_Y_System(self):
1223     d=self.domain.getDim()
1224     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1225 jfenwick 2455 mypde.setValue(Y=numpy.ones((self.N,)))
1226 gross 1841 coeff=mypde.getCoefficient("Y")
1227 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),Function(self.domain),self.N))
1228 jgs 149 def test_setCoefficient_y_System(self):
1229     d=self.domain.getDim()
1230     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1231 jfenwick 2455 mypde.setValue(y=numpy.ones((self.N,)))
1232 gross 1841 coeff=mypde.getCoefficient("y")
1233 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnBoundary(self.domain),self.N))
1234 jgs 149 def test_setCoefficient_d_System(self):
1235     d=self.domain.getDim()
1236     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1237 jfenwick 2455 mypde.setValue(d=numpy.ones((self.N,self.N)))
1238 gross 1841 coeff=mypde.getCoefficient("d")
1239 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnBoundary(self.domain),self.N,self.N))
1240 jgs 149 def test_setCoefficient_d_contact_System(self):
1241 jfenwick 3259 if self.domain.supportsContactElements():
1242     d=self.domain.getDim()
1243     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1244     mypde.setValue(d_contact=numpy.ones((self.N,self.N)))
1245     coeff=mypde.getCoefficient("d_contact")
1246     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnContactZero(self.domain),self.N,self.N))
1247 jgs 149 def test_setCoefficient_y_contact_System(self):
1248 jfenwick 3259 if self.domain.supportsContactElements():
1249     d=self.domain.getDim()
1250     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1251     mypde.setValue(y_contact=numpy.ones((self.N,)))
1252     coeff=mypde.getCoefficient("y_contact")
1253     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))
1254 gross 1072 def test_setCoefficient_A_reduced_System(self):
1255     d=self.domain.getDim()
1256     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1257 jfenwick 2455 mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))
1258 gross 1841 coeff=mypde.getCoefficient("A_reduced")
1259 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1260 gross 1072 def test_setCoefficient_B_reduced_System(self):
1261     d=self.domain.getDim()
1262     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1263 jfenwick 2455 mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))
1264 gross 1841 coeff=mypde.getCoefficient("B_reduced")
1265 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))
1266 gross 1072 def test_setCoefficient_C_reduced_System(self):
1267     d=self.domain.getDim()
1268     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1269 jfenwick 2455 mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))
1270 gross 1841 coeff=mypde.getCoefficient("C_reduced")
1271 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1272 gross 1072 def test_setCoefficient_D_System_reduced(self):
1273     d=self.domain.getDim()
1274     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1275 jfenwick 2455 mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))
1276 gross 1841 coeff=mypde.getCoefficient("D_reduced")
1277 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))
1278 gross 1072 def test_setCoefficient_X_System_reduced(self):
1279     d=self.domain.getDim()
1280     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1281 jfenwick 2455 mypde.setValue(X_reduced=numpy.ones((self.N,d)))
1282 gross 1841 coeff=mypde.getCoefficient("X_reduced")
1283 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))
1284 gross 1072 def test_setCoefficient_Y_System_reduced(self):
1285     d=self.domain.getDim()
1286     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1287 jfenwick 2455 mypde.setValue(Y_reduced=numpy.ones((self.N,)))
1288 gross 1841 coeff=mypde.getCoefficient("Y_reduced")
1289 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))
1290 gross 1072 def test_setCoefficient_y_System_reduced(self):
1291     d=self.domain.getDim()
1292     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1293 jfenwick 2455 mypde.setValue(y_reduced=numpy.ones((self.N,)))
1294 gross 1841 coeff=mypde.getCoefficient("y_reduced")
1295 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))
1296 gross 1072 def test_setCoefficient_d_reduced_System(self):
1297     d=self.domain.getDim()
1298     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1299 jfenwick 2455 mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))
1300 gross 1841 coeff=mypde.getCoefficient("d_reduced")
1301 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))
1302 gross 1072 def test_setCoefficient_d_contact_reduced_System(self):
1303 jfenwick 3259 if self.domain.supportsContactElements():
1304     d=self.domain.getDim()
1305     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1306     mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))
1307     coeff=mypde.getCoefficient("d_contact_reduced")
1308     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
1309 gross 1072 def test_setCoefficient_y_contact_reduced_System(self):
1310 jfenwick 3259 if self.domain.supportsContactElements():
1311     d=self.domain.getDim()
1312     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1313     mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))
1314     coeff=mypde.getCoefficient("y_contact_reduced")
1315     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1316 jgs 149 def test_setCoefficient_r_System(self):
1317     d=self.domain.getDim()
1318     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1319 jfenwick 2455 mypde.setValue(r=numpy.ones((self.N,)))
1320 gross 1841 coeff=mypde.getCoefficient("r")
1321 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1322 jgs 149 def test_setCoefficient_q_System(self):
1323     d=self.domain.getDim()
1324     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1325 jfenwick 2455 mypde.setValue(q=numpy.ones((self.N,)))
1326 gross 1841 coeff=mypde.getCoefficient("q")
1327 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1328 jgs 149 def test_setCoefficient_r_System_reducedOn(self):
1329     d=self.domain.getDim()
1330     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1331     mypde.setReducedOrderOn()
1332 jfenwick 2455 mypde.setValue(r=numpy.ones((self.N,)))
1333 gross 1841 coeff=mypde.getCoefficient("r")
1334 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1335 jgs 149 def test_setCoefficient_q_System_reducedOn(self):
1336     d=self.domain.getDim()
1337     mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1338     mypde.setReducedOrderOn()
1339 jfenwick 2455 mypde.setValue(q=numpy.ones((self.N,)))
1340 gross 1841 coeff=mypde.getCoefficient("q")
1341 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1342 jgs 149
1343 gross 1072 def test_setCoefficient_A_reduced_System_using_A(self):
1344     d=self.domain.getDim()
1345     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1346 jfenwick 2455 mypde.setValue(A=Data(numpy.ones((self.N,d,self.N,d)),ReducedFunction(self.domain)))
1347 gross 1841 coeff=mypde.getCoefficient("A_reduced")
1348 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1349 gross 1072 def test_setCoefficient_B_reduced_System_using_B(self):
1350     d=self.domain.getDim()
1351     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1352 jfenwick 2455 mypde.setValue(B=Data(numpy.ones((self.N,d,self.N)),ReducedFunction(self.domain)))
1353 gross 1841 coeff=mypde.getCoefficient("B_reduced")
1354 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))
1355 gross 1072 def test_setCoefficient_C_reduced_System_using_C(self):
1356     d=self.domain.getDim()
1357     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1358 jfenwick 2455 mypde.setValue(C=Data(numpy.ones((self.N,self.N,d)),ReducedFunction(self.domain)))
1359 gross 1841 coeff=mypde.getCoefficient("C_reduced")
1360 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1361 gross 1072 def test_setCoefficient_D_System_reduced_using_D(self):
1362     d=self.domain.getDim()
1363     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1364 jfenwick 2455 mypde.setValue(D=Data(numpy.ones((self.N,self.N)),ReducedFunction(self.domain)))
1365 gross 1841 coeff=mypde.getCoefficient("D_reduced")
1366 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))
1367 gross 1072 def test_setCoefficient_X_System_reduced_using_X(self):
1368     d=self.domain.getDim()
1369     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1370 jfenwick 2455 mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))
1371 gross 1841 coeff=mypde.getCoefficient("X_reduced")
1372 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))
1373 gross 1072 def test_setCoefficient_Y_System_reduced_using_Y(self):
1374     d=self.domain.getDim()
1375     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1376 jfenwick 2455 mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))
1377 gross 1841 coeff=mypde.getCoefficient("Y_reduced")
1378 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))
1379 gross 1072 def test_setCoefficient_y_reduced_System_using_y(self):
1380     d=self.domain.getDim()
1381     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1382 jfenwick 2455 mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))
1383 gross 1841 coeff=mypde.getCoefficient("y_reduced")
1384 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))
1385 gross 1072 def test_setCoefficient_d_reduced_System_using_d(self):
1386     d=self.domain.getDim()
1387     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1388 jfenwick 2455 mypde.setValue(d=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnBoundary(self.domain)))
1389 gross 1841 coeff=mypde.getCoefficient("d_reduced")
1390 gross 2474 self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))
1391 gross 1072 def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):
1392 jfenwick 3259 if self.domain.supportsContactElements():
1393     d=self.domain.getDim()
1394     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1395     mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))
1396     coeff=mypde.getCoefficient("d_contact_reduced")
1397     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
1398 gross 1072 def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):
1399 jfenwick 3259 if self.domain.supportsContactElements():
1400     d=self.domain.getDim()
1401     mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1402     mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))
1403     coeff=mypde.getCoefficient("y_contact_reduced")
1404     self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1405 jgs 149 def test_resetCoefficient_HomogeneousConstraint(self):
1406     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1407     x=self.domain.getX()
1408 gross 304 mypde.setValue(A=kronecker(self.domain),Y=1.,q=whereZero(x[0]))
1409 jgs 149 u1=mypde.getSolution()
1410     mypde.setValue(Y=2.)
1411     u2=mypde.getSolution()
1412     self.failUnless(self.check(u2,2*u1),'solution is wrong.')
1413    
1414     def test_resetCoefficient_InHomogeneousConstraint(self):
1415     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1416 jgs 154 mypde.setSymmetryOn()
1417 gross 2474 mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1418 jgs 149 x=self.domain.getX()
1419 gross 304 mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.,r=1,q=whereZero(x[0]))
1420 gross 2474 u1=mypde.getSolution()
1421 jgs 149 mypde.setValue(Y=2.,D=2)
1422 gross 2474 u2=mypde.getSolution()
1423 jgs 154 self.failUnless(self.check(u2,u1),'first solution is wrong.')
1424 gross 2474 u2=mypde.getSolution()
1425 jgs 154 self.failUnless(self.check(u2,u1),'first solution is wrong.')
1426 jgs 149 mypde.setValue(r=2,Y=4.)
1427 gross 2474 u2=mypde.getSolution()
1428 jgs 154 self.failUnless(self.check(u2,2*u1),'second solution is wrong.')
1429 jgs 149
1430 gross 2535 def test_Status(self):
1431     DIM=self.domain.getDim()
1432     x=self.domain.getX()
1433     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1434     mypde.getSolverOptions().setSymmetryOn()
1435     mypde.getSolverOptions().setTolerance(self.RES_TOL)
1436     mypde.setValue(A=kronecker(self.domain), q=whereZero(x[0])+whereZero(x[0]-1.), Y=2.)
1437     x1=self.domain.getX()
1438     u1_ref=x1[0]*(1.-x1[0])
1439     u1=mypde.getSolution()
1440     error1=Lsup(u1-u1_ref)/Lsup(u1_ref)
1441     self.failUnless(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of first pde does not match domain status.")
1442     self.domain.setX(x*5)
1443     self.failUnless(mypde.getDomainStatus() != mypde.getSystemStatus(), "status of first pde matches updated domain status.")
1444     x2=self.domain.getX()
1445     u2_ref=x2[0]*(5.-x2[0])
1446     u2=mypde.getSolution()
1447     error2=Lsup(u2-u2_ref)/Lsup(u2_ref)
1448     self.failUnless(error2 <= max(error1,self.RES_TOL)*10., "solution of second PDE wrong.")
1449     self.failUnless(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of second pde does not match domain status.")
1450    
1451 jgs 149 def test_symmetryCheckTrue_System(self):
1452     d=self.domain.getDim()
1453     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1454 jfenwick 2455 A=numpy.ones((self.N,d,self.N,d))
1455     C=2*numpy.ones((self.N,self.N,d))
1456     B=2*numpy.ones((self.N,d,self.N))
1457     D=3*numpy.ones((self.N,self.N))
1458     d=4*numpy.ones((self.N,self.N))
1459     d_contact=5*numpy.ones((self.N,self.N))
1460 jfenwick 3259 pars={"A":A, "B":B, "C":C, "D":D, "d":d, "A_reduced":-A, "B_reduced":-B, "C_reduced":-C, "D_reduced":-D, "d_reduced":-d}
1461     if self.domain.supportsContactElements():
1462     pars["d_contact"]=d_contact
1463     pars["d_contact_reduced"]=-d_contact
1464     mypde.setValue(**pars)
1465 jgs 149 self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")
1466    
1467     def test_symmetryCheckFalse_A_System(self):
1468     d=self.domain.getDim()
1469     mypde=LinearPDE(self.domain,debug=self.DEBUG)
1470 jfenwick 2455 A=numpy.ones((self.N,d,self.N,d))
1471 jgs 149 A[1,1,1,0]=0.
1472     mypde.setValue(A=A