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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4657 - (show annotations)
Thu Feb 6 06:12:20 2014 UTC (5 years, 7 months ago) by jfenwick
File MIME type: text/x-python
File size: 229416 byte(s)
I changed some files.
Updated copyright notices, added GeoComp.



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