/[escript]/trunk/escript/py_src/generateutil
ViewVC logotype

Contents of /trunk/escript/py_src/generateutil

Parent Directory Parent Directory | Revision Log Revision Log


Revision 441 - (show annotations)
Fri Jan 20 03:40:39 2006 UTC (13 years, 9 months ago) by gross
File size: 123349 byte(s)
and finally test for the divergence operator added. div() has been modified to take Symbol arguments
1 #!/usr/bin/python
2 # $Id$
3
4 """
5 program generates parts of the util.py and the test_util.py script
6 """
7 test_header=""
8 test_header+="import unittest\n"
9 test_header+="import numarray\n"
10 test_header+="from esys.escript import *\n"
11 test_header+="from esys.finley import Rectangle\n"
12 test_header+="class Test_util2(unittest.TestCase):\n"
13 test_header+=" RES_TOL=1.e-7\n"
14 test_header+=" def setUp(self):\n"
15 test_header+=" self.__dom =Rectangle(11,11,2)\n"
16 test_header+=" self.functionspace = FunctionOnBoundary(self.__dom)\n"
17 test_tail=""
18 test_tail+="suite = unittest.TestSuite()\n"
19 test_tail+="suite.addTest(unittest.makeSuite(Test_util2))\n"
20 test_tail+="unittest.TextTestRunner(verbosity=2).run(suite)\n"
21
22 case_set=["float","array","constData","taggedData","expandedData","Symbol"]
23 shape_set=[ (),(2,), (4,5), (6,2,2),(3,2,3,4)]
24
25 t_prog=""
26 t_prog_with_tags=""
27 t_prog_failing=""
28 u_prog=""
29
30 def wherepos(arg):
31 if arg>0.:
32 return 1.
33 else:
34 return 0.
35
36
37 class OPERATOR:
38 def __init__(self,nickname,rng=[-1000.,1000],test_expr="",math_expr=None,
39 numarray_expr="",symbol_expr=None,diff=None,name=""):
40 self.nickname=nickname
41 self.rng=rng
42 self.test_expr=test_expr
43 if math_expr==None:
44 self.math_expr=test_expr.replace("%a1%","arg")
45 else:
46 self.math_expr=math_expr
47 self.numarray_expr=numarray_expr
48 self.symbol_expr=symbol_expr
49 self.diff=diff
50 self.name=name
51
52 import random
53 import numarray
54 import math
55 finc=1.e-6
56
57 def getResultCaseForBin(case0,case1):
58 if case0=="float":
59 if case1=="float":
60 case="float"
61 elif case1=="array":
62 case="array"
63 elif case1=="constData":
64 case="constData"
65 elif case1=="taggedData":
66 case="taggedData"
67 elif case1=="expandedData":
68 case="expandedData"
69 elif case1=="Symbol":
70 case="Symbol"
71 else:
72 raise ValueError,"unknown case1=%s"%case1
73 elif case0=="array":
74 if case1=="float":
75 case="array"
76 elif case1=="array":
77 case="array"
78 elif case1=="constData":
79 case="constData"
80 elif case1=="taggedData":
81 case="taggedData"
82 elif case1=="expandedData":
83 case="expandedData"
84 elif case1=="Symbol":
85 case="Symbol"
86 else:
87 raise ValueError,"unknown case1=%s"%case1
88 elif case0=="constData":
89 if case1=="float":
90 case="constData"
91 elif case1=="array":
92 case="constData"
93 elif case1=="constData":
94 case="constData"
95 elif case1=="taggedData":
96 case="taggedData"
97 elif case1=="expandedData":
98 case="expandedData"
99 elif case1=="Symbol":
100 case="Symbol"
101 else:
102 raise ValueError,"unknown case1=%s"%case1
103 elif case0=="taggedData":
104 if case1=="float":
105 case="taggedData"
106 elif case1=="array":
107 case="taggedData"
108 elif case1=="constData":
109 case="taggedData"
110 elif case1=="taggedData":
111 case="taggedData"
112 elif case1=="expandedData":
113 case="expandedData"
114 elif case1=="Symbol":
115 case="Symbol"
116 else:
117 raise ValueError,"unknown case1=%s"%case1
118 elif case0=="expandedData":
119 if case1=="float":
120 case="expandedData"
121 elif case1=="array":
122 case="expandedData"
123 elif case1=="constData":
124 case="expandedData"
125 elif case1=="taggedData":
126 case="expandedData"
127 elif case1=="expandedData":
128 case="expandedData"
129 elif case1=="Symbol":
130 case="Symbol"
131 else:
132 raise ValueError,"unknown case1=%s"%case1
133 elif case0=="Symbol":
134 if case1=="float":
135 case="Symbol"
136 elif case1=="array":
137 case="Symbol"
138 elif case1=="constData":
139 case="Symbol"
140 elif case1=="taggedData":
141 case="Symbol"
142 elif case1=="expandedData":
143 case="Symbol"
144 elif case1=="Symbol":
145 case="Symbol"
146 else:
147 raise ValueError,"unknown case1=%s"%case1
148 else:
149 raise ValueError,"unknown case0=%s"%case0
150 return case
151
152
153 def makeArray(shape,rng):
154 l=rng[1]-rng[0]
155 out=numarray.zeros(shape,numarray.Float64)
156 if len(shape)==0:
157 out=l*random.random()+rng[0]
158 elif len(shape)==1:
159 for i0 in range(shape[0]):
160 out[i0]=l*random.random()+rng[0]
161 elif len(shape)==2:
162 for i0 in range(shape[0]):
163 for i1 in range(shape[1]):
164 out[i0,i1]=l*random.random()+rng[0]
165 elif len(shape)==3:
166 for i0 in range(shape[0]):
167 for i1 in range(shape[1]):
168 for i2 in range(shape[2]):
169 out[i0,i1,i2]=l*random.random()+rng[0]
170 elif len(shape)==4:
171 for i0 in range(shape[0]):
172 for i1 in range(shape[1]):
173 for i2 in range(shape[2]):
174 for i3 in range(shape[3]):
175 out[i0,i1,i2,i3]=l*random.random()+rng[0]
176 elif len(shape)==5:
177 for i0 in range(shape[0]):
178 for i1 in range(shape[1]):
179 for i2 in range(shape[2]):
180 for i3 in range(shape[3]):
181 for i4 in range(shape[4]):
182 out[i0,i1,i2,i3,i4]=l*random.random()+rng[0]
183 else:
184 raise SystemError,"rank is restricted to 5"
185 return out
186
187
188 def makeResult(val,test_expr):
189 if isinstance(val,float):
190 out=eval(test_expr.replace("%a1%","val"))
191 elif len(val.shape)==0:
192 out=eval(test_expr.replace("%a1%","val"))
193 elif len(val.shape)==1:
194 out=numarray.zeros(val.shape,numarray.Float64)
195 for i0 in range(val.shape[0]):
196 out[i0]=eval(test_expr.replace("%a1%","val[i0]"))
197 elif len(val.shape)==2:
198 out=numarray.zeros(val.shape,numarray.Float64)
199 for i0 in range(val.shape[0]):
200 for i1 in range(val.shape[1]):
201 out[i0,i1]=eval(test_expr.replace("%a1%","val[i0,i1]"))
202 elif len(val.shape)==3:
203 out=numarray.zeros(val.shape,numarray.Float64)
204 for i0 in range(val.shape[0]):
205 for i1 in range(val.shape[1]):
206 for i2 in range(val.shape[2]):
207 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val[i0,i1,i2]"))
208 elif len(val.shape)==4:
209 out=numarray.zeros(val.shape,numarray.Float64)
210 for i0 in range(val.shape[0]):
211 for i1 in range(val.shape[1]):
212 for i2 in range(val.shape[2]):
213 for i3 in range(val.shape[3]):
214 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val[i0,i1,i2,i3]"))
215 else:
216 raise SystemError,"rank is restricted to 4"
217 return out
218
219 def makeResult2(val0,val1,test_expr):
220 if isinstance(val0,float):
221 if isinstance(val1,float):
222 out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1"))
223 elif len(val1.shape)==0:
224 out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1"))
225 elif len(val1.shape)==1:
226 out=numarray.zeros(val1.shape,numarray.Float64)
227 for i0 in range(val1.shape[0]):
228 out[i0]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0]"))
229 elif len(val1.shape)==2:
230 out=numarray.zeros(val1.shape,numarray.Float64)
231 for i0 in range(val1.shape[0]):
232 for i1 in range(val1.shape[1]):
233 out[i0,i1]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1]"))
234 elif len(val1.shape)==3:
235 out=numarray.zeros(val1.shape,numarray.Float64)
236 for i0 in range(val1.shape[0]):
237 for i1 in range(val1.shape[1]):
238 for i2 in range(val1.shape[2]):
239 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2]"))
240 elif len(val1.shape)==4:
241 out=numarray.zeros(val1.shape,numarray.Float64)
242 for i0 in range(val1.shape[0]):
243 for i1 in range(val1.shape[1]):
244 for i2 in range(val1.shape[2]):
245 for i3 in range(val1.shape[3]):
246 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2,i3]"))
247 else:
248 raise SystemError,"rank of val1 is restricted to 4"
249 elif len(val0.shape)==0:
250 if isinstance(val1,float):
251 out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1"))
252 elif len(val1.shape)==0:
253 out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1"))
254 elif len(val1.shape)==1:
255 out=numarray.zeros(val1.shape,numarray.Float64)
256 for i0 in range(val1.shape[0]):
257 out[i0]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0]"))
258 elif len(val1.shape)==2:
259 out=numarray.zeros(val1.shape,numarray.Float64)
260 for i0 in range(val1.shape[0]):
261 for i1 in range(val1.shape[1]):
262 out[i0,i1]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1]"))
263 elif len(val1.shape)==3:
264 out=numarray.zeros(val1.shape,numarray.Float64)
265 for i0 in range(val1.shape[0]):
266 for i1 in range(val1.shape[1]):
267 for i2 in range(val1.shape[2]):
268 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2]"))
269 elif len(val1.shape)==4:
270 out=numarray.zeros(val1.shape,numarray.Float64)
271 for i0 in range(val1.shape[0]):
272 for i1 in range(val1.shape[1]):
273 for i2 in range(val1.shape[2]):
274 for i3 in range(val1.shape[3]):
275 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2,i3]"))
276 else:
277 raise SystemError,"rank of val1 is restricted to 4"
278 elif len(val0.shape)==1:
279 if isinstance(val1,float):
280 out=numarray.zeros(val0.shape,numarray.Float64)
281 for i0 in range(val0.shape[0]):
282 out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1"))
283 elif len(val1.shape)==0:
284 out=numarray.zeros(val0.shape,numarray.Float64)
285 for i0 in range(val0.shape[0]):
286 out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1"))
287 elif len(val1.shape)==1:
288 out=numarray.zeros(val0.shape,numarray.Float64)
289 for i0 in range(val0.shape[0]):
290 out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0]"))
291 elif len(val1.shape)==2:
292 out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64)
293 for i0 in range(val0.shape[0]):
294 for j1 in range(val1.shape[1]):
295 out[i0,j1]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1]"))
296 elif len(val1.shape)==3:
297 out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64)
298 for i0 in range(val0.shape[0]):
299 for j1 in range(val1.shape[1]):
300 for j2 in range(val1.shape[2]):
301 out[i0,j1,j2]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1,j2]"))
302 elif len(val1.shape)==4:
303 out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64)
304 for i0 in range(val0.shape[0]):
305 for j1 in range(val1.shape[1]):
306 for j2 in range(val1.shape[2]):
307 for j3 in range(val1.shape[3]):
308 out[i0,j1,j2,j3]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1,j2,j3]"))
309 else:
310 raise SystemError,"rank of val1 is restricted to 4"
311 elif len(val0.shape)==2:
312 if isinstance(val1,float):
313 out=numarray.zeros(val0.shape,numarray.Float64)
314 for i0 in range(val0.shape[0]):
315 for i1 in range(val0.shape[1]):
316 out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1"))
317 elif len(val1.shape)==0:
318 out=numarray.zeros(val0.shape,numarray.Float64)
319 for i0 in range(val0.shape[0]):
320 for i1 in range(val0.shape[1]):
321 out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1"))
322 elif len(val1.shape)==1:
323 out=numarray.zeros(val0.shape,numarray.Float64)
324 for i0 in range(val0.shape[0]):
325 for i1 in range(val0.shape[1]):
326 out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0]"))
327 elif len(val1.shape)==2:
328 out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64)
329 for i0 in range(val0.shape[0]):
330 for i1 in range(val0.shape[1]):
331 out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1]"))
332 elif len(val1.shape)==3:
333 out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64)
334 for i0 in range(val0.shape[0]):
335 for i1 in range(val0.shape[1]):
336 for j2 in range(val1.shape[2]):
337 out[i0,i1,j2]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1,j2]"))
338 elif len(val1.shape)==4:
339 out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64)
340 for i0 in range(val0.shape[0]):
341 for i1 in range(val0.shape[1]):
342 for j2 in range(val1.shape[2]):
343 for j3 in range(val1.shape[3]):
344 out[i0,i1,j2,j3]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1,j2,j3]"))
345 else:
346 raise SystemError,"rank of val1 is restricted to 4"
347 elif len(val0.shape)==3:
348 if isinstance(val1,float):
349 out=numarray.zeros(val0.shape,numarray.Float64)
350 for i0 in range(val0.shape[0]):
351 for i1 in range(val0.shape[1]):
352 for i2 in range(val0.shape[2]):
353 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1"))
354 elif len(val1.shape)==0:
355 out=numarray.zeros(val0.shape,numarray.Float64)
356 for i0 in range(val0.shape[0]):
357 for i1 in range(val0.shape[1]):
358 for i2 in range(val0.shape[2]):
359 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1"))
360 elif len(val1.shape)==1:
361 out=numarray.zeros(val0.shape,numarray.Float64)
362 for i0 in range(val0.shape[0]):
363 for i1 in range(val0.shape[1]):
364 for i2 in range(val0.shape[2]):
365 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0]"))
366 elif len(val1.shape)==2:
367 out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64)
368 for i0 in range(val0.shape[0]):
369 for i1 in range(val0.shape[1]):
370 for i2 in range(val0.shape[2]):
371 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1]"))
372 elif len(val1.shape)==3:
373 out=numarray.zeros(val0.shape,numarray.Float64)
374 for i0 in range(val0.shape[0]):
375 for i1 in range(val0.shape[1]):
376 for i2 in range(val0.shape[2]):
377 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1,i2]"))
378 elif len(val1.shape)==4:
379 out=numarray.zeros(val0.shape+val1.shape[3:],numarray.Float64)
380 for i0 in range(val0.shape[0]):
381 for i1 in range(val0.shape[1]):
382 for i2 in range(val0.shape[2]):
383 for j3 in range(val1.shape[3]):
384 out[i0,i1,i2,j3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1,i2,j3]"))
385 else:
386 raise SystemError,"rank of val1 is rargs[1]estricted to 4"
387 elif len(val0.shape)==4:
388 if isinstance(val1,float):
389 out=numarray.zeros(val0.shape,numarray.Float64)
390 for i0 in range(val0.shape[0]):
391 for i1 in range(val0.shape[1]):
392 for i2 in range(val0.shape[2]):
393 for i3 in range(val0.shape[3]):
394 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1"))
395 elif len(val1.shape)==0:
396 out=numarray.zeros(val0.shape,numarray.Float64)
397 for i0 in range(val0.shape[0]):
398 for i1 in range(val0.shape[1]):
399 for i2 in range(val0.shape[2]):
400 for i3 in range(val0.shape[3]):
401 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1"))
402 elif len(val1.shape)==1:
403 out=numarray.zeros(val0.shape,numarray.Float64)
404 for i0 in range(val0.shape[0]):
405 for i1 in range(val0.shape[1]):
406 for i2 in range(val0.shape[2]):
407 for i3 in range(val0.shape[3]):
408 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0]"))
409 elif len(val1.shape)==2:
410 out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64)
411 for i0 in range(val0.shape[0]):
412 for i1 in range(val0.shape[1]):
413 for i2 in range(val0.shape[2]):
414 for i3 in range(val0.shape[3]):
415 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1]"))
416 elif len(val1.shape)==3:
417 out=numarray.zeros(val0.shape,numarray.Float64)
418 for i0 in range(val0.shape[0]):
419 for i1 in range(val0.shape[1]):
420 for i2 in range(val0.shape[2]):
421 for i3 in range(val0.shape[3]):
422 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1,i2]"))
423 elif len(val1.shape)==4:
424 out=numarray.zeros(val0.shape,numarray.Float64)
425 for i0 in range(val0.shape[0]):
426 for i1 in range(val0.shape[1]):
427 for i2 in range(val0.shape[2]):
428 for i3 in range(val0.shape[3]):
429 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1,i2,i3]"))
430 else:
431 raise SystemError,"rank of val1 is restricted to 4"
432 else:
433 raise SystemError,"rank is restricted to 4"
434 return out
435
436
437 def mkText(case,name,a,a1=None,use_tagging_for_expanded_data=False):
438 t_out=""
439 if case=="float":
440 if isinstance(a,float):
441 t_out+=" %s=%s\n"%(name,a)
442 elif a.rank==0:
443 t_out+=" %s=%s\n"%(name,a)
444 else:
445 t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist())
446 elif case=="array":
447 if isinstance(a,float):
448 t_out+=" %s=numarray.array(%s)\n"%(name,a)
449 elif a.rank==0:
450 t_out+=" %s=numarray.array(%s)\n"%(name,a)
451 else:
452 t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist())
453 elif case=="constData":
454 if isinstance(a,float):
455 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
456 elif a.rank==0:
457 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
458 else:
459 t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist())
460 elif case=="taggedData":
461 if isinstance(a,float):
462 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
463 t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1)
464 elif a.rank==0:
465 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
466 t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1)
467 else:
468 t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist())
469 t_out+=" %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist())
470 elif case=="expandedData":
471 if use_tagging_for_expanded_data:
472 if isinstance(a,float):
473 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
474 t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1)
475 elif a.rank==0:
476 t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
477 t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1)
478 else:
479 t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist())
480 t_out+=" %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist())
481 t_out+=" %s.expand()\n"%name
482 else:
483 t_out+=" msk_%s=whereNegative(self.functionspace.getX()[0]-0.5)\n"%name
484 if isinstance(a,float):
485 t_out+=" %s=msk_%s*(%s)+(1.-msk_%s)*(%s)\n"%(name,name,a,name,a1)
486 elif a.rank==0:
487 t_out+=" %s=msk_%s*numarray.array(%s)+(1.-msk_%s)*numarray.array(%s)\n"%(name,name,a,name,a1)
488 else:
489 t_out+=" %s=msk_%s*numarray.array(%s)+(1.-msk_%s)*numarray.array(%s)\n"%(name,name,a.tolist(),name,a1.tolist())
490 elif case=="Symbol":
491 if isinstance(a,float):
492 t_out+=" %s=Symbol(shape=())\n"%(name)
493 elif a.rank==0:
494 t_out+=" %s=Symbol(shape=())\n"%(name)
495 else:
496 t_out+=" %s=Symbol(shape=%s)\n"%(name,str(a.shape))
497
498 return t_out
499
500 def mkTypeAndShapeTest(case,sh,argstr):
501 text=""
502 if case=="float":
503 text+=" self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%argstr
504 elif case=="array":
505 text+=" self.failUnless(isinstance(%s,numarray.NumArray),\"wrong type of result.\")\n"%argstr
506 text+=" self.failUnlessEqual(%s.shape,%s,\"wrong shape of result.\")\n"%(argstr,str(sh))
507 elif case in ["constData","taggedData","expandedData"]:
508 text+=" self.failUnless(isinstance(%s,Data),\"wrong type of result.\")\n"%argstr
509 text+=" self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result.\")\n"%(argstr,str(sh))
510 elif case=="Symbol":
511 text+=" self.failUnless(isinstance(%s,Symbol),\"wrong type of result.\")\n"%argstr
512 text+=" self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result.\")\n"%(argstr,str(sh))
513 return text
514
515 def mkCode(txt,args=[],intend=""):
516 s=txt.split("\n")
517 if len(s)>1:
518 out=""
519 for l in s:
520 out+=intend+l+"\n"
521 else:
522 out="%sreturn %s\n"%(intend,txt)
523 c=1
524 for r in args:
525 out=out.replace("%%a%s%%"%c,r)
526 return out
527
528 def innerTEST(arg0,arg1):
529 if isinstance(arg0,float):
530 out=numarray.array(arg0*arg1)
531 else:
532 out=(arg0*arg1).sum()
533 return out
534
535 def outerTEST(arg0,arg1):
536 if isinstance(arg0,float):
537 out=numarray.array(arg0*arg1)
538 elif isinstance(arg1,float):
539 out=numarray.array(arg0*arg1)
540 else:
541 out=numarray.outerproduct(arg0,arg1).resize(arg0.shape+arg1.shape)
542 return out
543
544 def tensorProductTest(arg0,arg1,sh_s):
545 if isinstance(arg0,float):
546 out=numarray.array(arg0*arg1)
547 elif isinstance(arg1,float):
548 out=numarray.array(arg0*arg1)
549 elif len(sh_s)==0:
550 out=numarray.outerproduct(arg0,arg1).resize(arg0.shape+arg1.shape)
551 else:
552 l=len(sh_s)
553 sh0=arg0.shape[:arg0.rank-l]
554 sh1=arg1.shape[l:]
555 ls,l0,l1=1,1,1
556 for i in sh_s: ls*=i
557 for i in sh0: l0*=i
558 for i in sh1: l1*=i
559 out1=numarray.outerproduct(arg0,arg1).resize((l0,ls,ls,l1))
560 out2=numarray.zeros((l0,l1),numarray.Float)
561 for i0 in range(l0):
562 for i1 in range(l1):
563 for i in range(ls): out2[i0,i1]+=out1[i0,i,i,i1]
564 out=out2.resize(sh0+sh1)
565 return out
566
567 def testMatrixMult(arg0,arg1,sh_s):
568 return numarray.matrixmultiply(arg0,arg1)
569
570
571 def testTensorMult(arg0,arg1,sh_s):
572 if len(arg0)==2:
573 return numarray.matrixmultiply(arg0,arg1)
574 else:
575 if arg1.rank==4:
576 out=numarray.zeros((arg0.shape[0],arg0.shape[1],arg1.shape[2],arg1.shape[3]),numarray.Float)
577 for i0 in range(arg0.shape[0]):
578 for i1 in range(arg0.shape[1]):
579 for i2 in range(arg0.shape[2]):
580 for i3 in range(arg0.shape[3]):
581 for j2 in range(arg1.shape[2]):
582 for j3 in range(arg1.shape[3]):
583 out[i0,i1,j2,j3]+=arg0[i0,i1,i2,i3]*arg1[i2,i3,j2,j3]
584 elif arg1.rank==3:
585 out=numarray.zeros((arg0.shape[0],arg0.shape[1],arg1.shape[2]),numarray.Float)
586 for i0 in range(arg0.shape[0]):
587 for i1 in range(arg0.shape[1]):
588 for i2 in range(arg0.shape[2]):
589 for i3 in range(arg0.shape[3]):
590 for j2 in range(arg1.shape[2]):
591 out[i0,i1,j2]+=arg0[i0,i1,i2,i3]*arg1[i2,i3,j2]
592 elif arg1.rank==2:
593 out=numarray.zeros((arg0.shape[0],arg0.shape[1]),numarray.Float)
594 for i0 in range(arg0.shape[0]):
595 for i1 in range(arg0.shape[1]):
596 for i2 in range(arg0.shape[2]):
597 for i3 in range(arg0.shape[3]):
598 out[i0,i1]+=arg0[i0,i1,i2,i3]*arg1[i2,i3]
599 return out
600
601 def testReduce(arg0,init_val,test_expr,post_expr):
602 out=init_val
603 if isinstance(arg0,float):
604 out=eval(test_expr.replace("%a1%","arg0"))
605 elif arg0.rank==0:
606 out=eval(test_expr.replace("%a1%","arg0"))
607 elif arg0.rank==1:
608 for i0 in range(arg0.shape[0]):
609 out=eval(test_expr.replace("%a1%","arg0[i0]"))
610 elif arg0.rank==2:
611 for i0 in range(arg0.shape[0]):
612 for i1 in range(arg0.shape[1]):
613 out=eval(test_expr.replace("%a1%","arg0[i0,i1]"))
614 elif arg0.rank==3:
615 for i0 in range(arg0.shape[0]):
616 for i1 in range(arg0.shape[1]):
617 for i2 in range(arg0.shape[2]):
618 out=eval(test_expr.replace("%a1%","arg0[i0,i1,i2]"))
619 elif arg0.rank==4:
620 for i0 in range(arg0.shape[0]):
621 for i1 in range(arg0.shape[1]):
622 for i2 in range(arg0.shape[2]):
623 for i3 in range(arg0.shape[3]):
624 out=eval(test_expr.replace("%a1%","arg0[i0,i1,i2,i3]"))
625 return eval(post_expr)
626
627 def clipTEST(arg0,mn,mx):
628 if isinstance(arg0,float):
629 return max(min(arg0,mx),mn)
630 out=numarray.zeros(arg0.shape,numarray.Float64)
631 if arg0.rank==1:
632 for i0 in range(arg0.shape[0]):
633 out[i0]=max(min(arg0[i0],mx),mn)
634 elif arg0.rank==2:
635 for i0 in range(arg0.shape[0]):
636 for i1 in range(arg0.shape[1]):
637 out[i0,i1]=max(min(arg0[i0,i1],mx),mn)
638 elif arg0.rank==3:
639 for i0 in range(arg0.shape[0]):
640 for i1 in range(arg0.shape[1]):
641 for i2 in range(arg0.shape[2]):
642 out[i0,i1,i2]=max(min(arg0[i0,i1,i2],mx),mn)
643 elif arg0.rank==4:
644 for i0 in range(arg0.shape[0]):
645 for i1 in range(arg0.shape[1]):
646 for i2 in range(arg0.shape[2]):
647 for i3 in range(arg0.shape[3]):
648 out[i0,i1,i2,i3]=max(min(arg0[i0,i1,i2,i3],mx),mn)
649 return out
650 def minimumTEST(arg0,arg1):
651 if isinstance(arg0,float):
652 if isinstance(arg1,float):
653 if arg0>arg1:
654 return arg1
655 else:
656 return arg0
657 else:
658 arg0=numarray.ones(arg1.shape)*arg0
659 else:
660 if isinstance(arg1,float):
661 arg1=numarray.ones(arg0.shape)*arg1
662 out=numarray.zeros(arg0.shape,numarray.Float64)
663 if arg0.rank==0:
664 if arg0>arg1:
665 out=arg1
666 else:
667 out=arg0
668 elif arg0.rank==1:
669 for i0 in range(arg0.shape[0]):
670 if arg0[i0]>arg1[i0]:
671 out[i0]=arg1[i0]
672 else:
673 out[i0]=arg0[i0]
674 elif arg0.rank==2:
675 for i0 in range(arg0.shape[0]):
676 for i1 in range(arg0.shape[1]):
677 if arg0[i0,i1]>arg1[i0,i1]:
678 out[i0,i1]=arg1[i0,i1]
679 else:
680 out[i0,i1]=arg0[i0,i1]
681 elif arg0.rank==3:
682 for i0 in range(arg0.shape[0]):
683 for i1 in range(arg0.shape[1]):
684 for i2 in range(arg0.shape[2]):
685 if arg0[i0,i1,i2]>arg1[i0,i1,i2]:
686 out[i0,i1,i2]=arg1[i0,i1,i2]
687 else:
688 out[i0,i1,i2]=arg0[i0,i1,i2]
689 elif arg0.rank==4:
690 for i0 in range(arg0.shape[0]):
691 for i1 in range(arg0.shape[1]):
692 for i2 in range(arg0.shape[2]):
693 for i3 in range(arg0.shape[3]):
694 if arg0[i0,i1,i2,i3]>arg1[i0,i1,i2,i3]:
695 out[i0,i1,i2,i3]=arg1[i0,i1,i2,i3]
696 else:
697 out[i0,i1,i2,i3]=arg0[i0,i1,i2,i3]
698 return out
699 def unrollLoops(a,b,o,arg,tap="",x="x"):
700 out=""
701 if a.rank==1:
702 z=""
703 for i99 in range(a.shape[0]):
704 if not z=="": z+="+"
705 if o=="1":
706 z+="(%s)*%s[%s]"%(a[i99]+b[i99],x,i99)
707 else:
708 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i99],x,i99,b[i99],x,i99)
709
710 out+=tap+"%s=%s\n"%(arg,z)
711
712 elif a.rank==2:
713 for i0 in range(a.shape[0]):
714 z=""
715 for i99 in range(a.shape[1]):
716 if not z=="": z+="+"
717 if o=="1":
718 z+="(%s)*x[%s]"%(a[i0,i99]+b[i0,i99],i99)
719 else:
720 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i99],x,i99,b[i0,i99],x,i99)
721
722 out+=tap+"%s[%s]=%s\n"%(arg,i0,z)
723 elif a.rank==3:
724 for i0 in range(a.shape[0]):
725 for i1 in range(a.shape[1]):
726 z=""
727 for i99 in range(a.shape[2]):
728 if not z=="": z+="+"
729 if o=="1":
730 z+="(%s)*%s[%s]"%(a[i0,i1,i99]+b[i0,i1,i99],x,i99)
731 else:
732 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i99],x,i99,b[i0,i1,i99],x,i99)
733
734 out+=tap+"%s[%s,%s]=%s\n"%(arg,i0,i1,z)
735 elif a.rank==4:
736 for i0 in range(a.shape[0]):
737 for i1 in range(a.shape[1]):
738 for i2 in range(a.shape[2]):
739 z=""
740 for i99 in range(a.shape[3]):
741 if not z=="": z+="+"
742 if o=="1":
743 z+="(%s)*%s[%s]"%(a[i0,i1,i2,i99]+b[i0,i1,i2,i99],x,i99)
744 else:
745 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i99],x,i99,b[i0,i1,i2,i99],x,i99)
746
747 out+=tap+"%s[%s,%s,%s]=%s\n"%(arg,i0,i1,i2,z)
748 elif a.rank==5:
749 for i0 in range(a.shape[0]):
750 for i1 in range(a.shape[1]):
751 for i2 in range(a.shape[2]):
752 for i3 in range(a.shape[3]):
753 z=""
754 for i99 in range(a.shape[4]):
755 if not z=="": z+="+"
756 if o=="1":
757 z+="(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99]+b[i0,i1,i2,i3,i99],x,i99)
758 else:
759 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99],x,i99,b[i0,i1,i2,i3,i99],x,i99)
760
761 out+=tap+"%s[%s,%s,%s,%s]=%s\n"%(arg,i0,i1,i2,i3,z)
762 return out
763
764 def unrollLoopsOfGrad(a,b,o,arg,tap=""):
765 out=""
766 if a.rank==1:
767 for i99 in range(a.shape[0]):
768 if o=="1":
769 out+=tap+"%s[%s]=(%s)\n"%(arg,i99,a[i99]+b[i99])
770 else:
771 out+=tap+"%s[%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i99,a[i99],i99,b[i99])
772
773 elif a.rank==2:
774 for i0 in range(a.shape[0]):
775 for i99 in range(a.shape[1]):
776 if o=="1":
777 out+=tap+"%s[%s,%s]=(%s)\n"%(arg,i0,i99,a[i0,i99]+b[i0,i99])
778 else:
779 out+=tap+"%s[%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i99,a[i0,i99],i99,b[i0,i99])
780
781 elif a.rank==3:
782 for i0 in range(a.shape[0]):
783 for i1 in range(a.shape[1]):
784 for i99 in range(a.shape[2]):
785 if o=="1":
786 out+=tap+"%s[%s,%s,%s]=(%s)\n"%(arg,i0,i1,i99,a[i0,i1,i99]+b[i0,i1,i99])
787 else:
788 out+=tap+"%s[%s,%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i1,i99,a[i0,i1,i99],i99,b[i0,i1,i99])
789
790 elif a.rank==4:
791 for i0 in range(a.shape[0]):
792 for i1 in range(a.shape[1]):
793 for i2 in range(a.shape[2]):
794 for i99 in range(a.shape[3]):
795 if o=="1":
796 out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99])
797 else:
798 out+=tap+"%s[%s,%s,%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99],i99,b[i0,i1,i2,i99])
799 return out
800 def unrollLoopsOfDiv(a,b,o,arg,tap=""):
801
802
803
804 out=tap+arg+"="
805 if o=="1":
806 z=0.
807 for i99 in range(a.shape[0]):
808 z+=b[i99,i99]+a[i99,i99]
809 out+="(%s)"%z
810 else:
811 z=0.
812 for i99 in range(a.shape[0]):
813 z+=b[i99,i99]
814 if i99>0: out+="+"
815 out+="o*(%s)*x_ref[%s]**(o-1)"%(a[i99,i99],i99)
816 out+="+(%s)"%z
817 return out
818
819 def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):
820 if where=="Function":
821 xfac_o=1.
822 xfac_op=0.
823 z_fac=1./2.
824 z_fac_s=""
825 zo_fac_s="/(o+1.)"
826 elif where=="FunctionOnBoundary":
827 xfac_o=1.
828 xfac_op=0.
829 z_fac=1.
830 z_fac_s="*dim"
831 zo_fac_s="*(1+2.*(dim-1.)/(o+1.))"
832 elif where in ["FunctionOnContactZero","FunctionOnContactOne"]:
833 xfac_o=0.
834 xfac_op=1.
835 z_fac=1./2.
836 z_fac_s=""
837 zo_fac_s="/(o+1.)"
838 out=""
839 if a.rank==1:
840 zo=0.
841 zop=0.
842 z=0.
843 for i99 in range(a.shape[0]):
844 if i99==0:
845 zo+= xfac_o*a[i99]
846 zop+= xfac_op*a[i99]
847 else:
848 zo+=a[i99]
849 z+=b[i99]
850
851 out+=tap+"%s=(%s)%s+(%s)%s"%(arg,zo,zo_fac_s,z*z_fac,z_fac_s)
852 if zop==0.:
853 out+="\n"
854 else:
855 out+="+(%s)*0.5**o\n"%zop
856 elif a.rank==2:
857 for i0 in range(a.shape[0]):
858 zo=0.
859 zop=0.
860 z=0.
861 for i99 in range(a.shape[1]):
862 if i99==0:
863 zo+= xfac_o*a[i0,i99]
864 zop+= xfac_op*a[i0,i99]
865 else:
866 zo+=a[i0,i99]
867 z+=b[i0,i99]
868
869 out+=tap+"%s[%s]=(%s)%s+(%s)%s"%(arg,i0,zo,zo_fac_s,z*z_fac,z_fac_s)
870 if zop==0.:
871 out+="\n"
872 else:
873 out+="+(%s)*0.5**o\n"%zop
874 elif a.rank==3:
875 for i0 in range(a.shape[0]):
876 for i1 in range(a.shape[1]):
877 zo=0.
878 zop=0.
879 z=0.
880 for i99 in range(a.shape[2]):
881 if i99==0:
882 zo+= xfac_o*a[i0,i1,i99]
883 zop+= xfac_op*a[i0,i1,i99]
884 else:
885 zo+=a[i0,i1,i99]
886 z+=b[i0,i1,i99]
887
888 out+=tap+"%s[%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,zo,zo_fac_s,z*z_fac,z_fac_s)
889 if zop==0.:
890 out+="\n"
891 else:
892 out+="+(%s)*0.5**o\n"%zop
893 elif a.rank==4:
894 for i0 in range(a.shape[0]):
895 for i1 in range(a.shape[1]):
896 for i2 in range(a.shape[2]):
897 zo=0.
898 zop=0.
899 z=0.
900 for i99 in range(a.shape[3]):
901 if i99==0:
902 zo+= xfac_o*a[i0,i1,i2,i99]
903 zop+= xfac_op*a[i0,i1,i2,i99]
904
905 else:
906 zo+=a[i0,i1,i2,i99]
907 z+=b[i0,i1,i2,i99]
908
909 out+=tap+"%s[%s,%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,i2,zo,zo_fac_s,z*z_fac,z_fac_s)
910 if zop==0.:
911 out+="\n"
912 else:
913 out+="+(%s)*0.5**o\n"%zop
914
915 elif a.rank==5:
916 for i0 in range(a.shape[0]):
917 for i1 in range(a.shape[1]):
918 for i2 in range(a.shape[2]):
919 for i3 in range(a.shape[3]):
920 zo=0.
921 zop=0.
922 z=0.
923 for i99 in range(a.shape[4]):
924 if i99==0:
925 zo+= xfac_o*a[i0,i1,i2,i3,i99]
926 zop+= xfac_op*a[i0,i1,i2,i3,i99]
927
928 else:
929 zo+=a[i0,i1,i2,i3,i99]
930 z+=b[i0,i1,i2,i3,i99]
931 out+=tap+"%s[%s,%s,%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,i2,i3,zo,zo_fac_s,z*z_fac,z_fac_s)
932 if zop==0.:
933 out+="\n"
934 else:
935 out+="+(%s)*0.5**o\n"%zop
936
937 return out
938 #=======================================================================================================
939 # div
940 #=======================================================================================================
941 for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
942 for data in ["Data","Symbol"]:
943 for case in ["ContinuousFunction","Solution","ReducedSolution"]:
944 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
945 tname="test_div_on%s_from%s_%s"%(where,data,case)
946 text+=" def %s(self):\n"%tname
947 text+=" \"\"\"\n"
948 text+=" tests divergence of %s on the %s\n\n"%(data,where)
949 text+=" assumptions: %s(self.domain) exists\n"%case
950 text+=" self.domain supports div on %s\n"%where
951 text+=" \"\"\"\n"
952 if case=="ReducedSolution":
953 text+=" o=1\n"
954 o="1"
955 else:
956 text+=" o=self.order\n"
957 o="o"
958 text+=" dim=self.domain.getDim()\n"
959 text+=" w_ref=%s(self.domain)\n"%where
960 text+=" x_ref=w_ref.getX()\n"
961 text+=" w=%s(self.domain)\n"%case
962 text+=" x=w.getX()\n"
963 a_2=makeArray((2,2),[-1.,1])
964 b_2=makeArray((2,2),[-1.,1])
965 a_3=makeArray((3,3),[-1.,1])
966 b_3=makeArray((3,3),[-1.,1])
967 if data=="Symbol":
968 text+=" arg=Symbol(shape=(dim,),dim=dim)\n"
969 val="s"
970 res="sub"
971 else:
972 val="arg"
973 res="res"
974 text+=" %s=Vector(0,w)\n"%val
975 text+=" if dim==2:\n"
976 text+=unrollLoops(a_2,b_2,o,val,tap=" ")
977 text+=unrollLoopsOfDiv(a_2,b_2,o,"ref",tap=" ")
978 text+="\n else:\n"
979
980 text+=unrollLoops(a_3,b_3,o,val,tap=" ")
981 text+=unrollLoopsOfDiv(a_3,b_3,o,"ref",tap=" ")
982 text+="\n res=div(arg,where=w_ref)\n"
983 if data=="Symbol":
984 text+=" sub=res.substitute({arg:s})\n"
985 text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
986 text+=" self.failUnlessEqual(res.getShape(),(),\"wrong shape of result.\")\n"
987 text+=" self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong function space of result.\")\n"%res
988 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
989
990
991 t_prog+=text
992 print t_prog
993 1/0
994
995 #=======================================================================================================
996 # interpolation
997 #=======================================================================================================
998 for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne","Solution","ReducedSolution"]:
999 for data in ["Data","Symbol"]:
1000 for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1001 for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
1002 if where==case or \
1003 ( case in ["ContinuousFunction","Solution","ReducedSolution"] and where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"] ) or \
1004 ( case in ["FunctionOnContactZero","FunctionOnContactOne"] and where in ["FunctionOnContactZero","FunctionOnContactOne"] ) or \
1005 (case=="ContinuousFunction" and where in ["Solution","ReducedSolution"]) or \
1006 (case=="Solution" and where=="ReducedSolution") :
1007
1008
1009 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1010 tname="test_interpolation_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1011 text+=" def %s(self):\n"%tname
1012 text+=" \"\"\"\n"
1013 text+=" tests interpolation for rank %s %s onto the %s\n\n"%(len(sh),data,where)
1014 text+=" assumptions: self.domain supports inpterpolation from %s onto %s\n"%(case,where)
1015 text+=" \"\"\"\n"
1016 if case=="ReducedSolution" or where=="ReducedSolution":
1017 text+=" o=1\n"
1018 o="1"
1019 else:
1020 text+=" o=self.order\n"
1021 o="o"
1022 text+=" dim=self.domain.getDim()\n"
1023 text+=" w_ref=%s(self.domain)\n"%where
1024 text+=" x_ref=w_ref.getX()\n"
1025 text+=" w=%s(self.domain)\n"%case
1026 text+=" x=w.getX()\n"
1027 a_2=makeArray(sh+(2,),[-1.,1])
1028 b_2=makeArray(sh+(2,),[-1.,1])
1029 a_3=makeArray(sh+(3,),[-1.,1])
1030 b_3=makeArray(sh+(3,),[-1.,1])
1031 if data=="Symbol":
1032 text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
1033 val="s"
1034 res="sub"
1035 else:
1036 val="arg"
1037 res="res"
1038 text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
1039 text+=" ref=Data(0,%s,w_ref)\n"%str(sh)
1040 text+=" if dim==2:\n"
1041 text+=unrollLoops(a_2,b_2,o,val,tap=" ")
1042 text+=unrollLoops(a_2,b_2,o,"ref",tap=" ",x="x_ref")
1043 text+=" else:\n"
1044
1045 text+=unrollLoops(a_3,b_3,o,val,tap=" ")
1046 text+=unrollLoops(a_3,b_3,o,"ref",tap=" ",x="x_ref")
1047 text+=" res=interpolate(arg,where=w_ref)\n"
1048 if data=="Symbol":
1049 text+=" sub=res.substitute({arg:s})\n"
1050 text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
1051 text+=" self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong functionspace of result.\")\n"%res
1052 text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
1053 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1054 t_prog+=text
1055 print test_header
1056 print t_prog
1057 print test_tail
1058 1/0
1059
1060 #=======================================================================================================
1061 # grad
1062 #=======================================================================================================
1063 for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1064 for data in ["Data","Symbol"]:
1065 for case in ["ContinuousFunction","Solution","ReducedSolution"]:
1066 for sh in [ (),(2,), (4,5), (6,2,2)]:
1067 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1068 tname="test_grad_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1069 text+=" def %s(self):\n"%tname
1070 text+=" \"\"\"\n"
1071 text+=" tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where)
1072 text+=" assumptions: %s(self.domain) exists\n"%case
1073 text+=" self.domain supports gardient on %s\n"%where
1074 text+=" \"\"\"\n"
1075 if case=="ReducedSolution":
1076 text+=" o=1\n"
1077 o="1"
1078 else:
1079 text+=" o=self.order\n"
1080 o="o"
1081 text+=" dim=self.domain.getDim()\n"
1082 text+=" w_ref=%s(self.domain)\n"%where
1083 text+=" x_ref=w_ref.getX()\n"
1084 text+=" w=%s(self.domain)\n"%case
1085 text+=" x=w.getX()\n"
1086 a_2=makeArray(sh+(2,),[-1.,1])
1087 b_2=makeArray(sh+(2,),[-1.,1])
1088 a_3=makeArray(sh+(3,),[-1.,1])
1089 b_3=makeArray(sh+(3,),[-1.,1])
1090 if data=="Symbol":
1091 text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
1092 val="s"
1093 res="sub"
1094 else:
1095 val="arg"
1096 res="res"
1097 text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
1098 text+=" ref=Data(0,%s+(dim,),w_ref)\n"%str(sh)
1099 text+=" if dim==2:\n"
1100 text+=unrollLoops(a_2,b_2,o,val,tap=" ")
1101 text+=unrollLoopsOfGrad(a_2,b_2,o,"ref",tap=" ")
1102 text+=" else:\n"
1103
1104 text+=unrollLoops(a_3,b_3,o,val,tap=" ")
1105 text+=unrollLoopsOfGrad(a_3,b_3,o,"ref",tap=" ")
1106 text+=" res=grad(arg,where=w_ref)\n"
1107 if data=="Symbol":
1108 text+=" sub=res.substitute({arg:s})\n"
1109 text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
1110 text+=" self.failUnlessEqual(res.getShape(),%s+(dim,),\"wrong shape of result.\")\n"%str(sh)
1111 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1112
1113
1114 t_prog+=text
1115 print test_header
1116 print t_prog
1117 print test_tail
1118 1/0
1119
1120
1121 #=======================================================================================================
1122 # integrate
1123 #=======================================================================================================
1124 for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1125 for data in ["Data","Symbol"]:
1126 for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1127 for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
1128 if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:
1129 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1130 tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1131 text+=" def %s(self):\n"%tname
1132 text+=" \"\"\"\n"
1133 text+=" tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where)
1134 text+=" assumptions: %s(self.domain) exists\n"%case
1135 text+=" self.domain supports integral on %s\n"%where
1136
1137 text+=" \"\"\"\n"
1138 if case=="ReducedSolution":
1139 text+=" o=1\n"
1140 o="1"
1141 else:
1142 text+=" o=self.order\n"
1143 o="o"
1144 text+=" dim=self.domain.getDim()\n"
1145 text+=" w_ref=%s(self.domain)\n"%where
1146 text+=" w=%s(self.domain)\n"%case
1147 text+=" x=w.getX()\n"
1148 a_2=makeArray(sh+(2,),[-1.,1])
1149 b_2=makeArray(sh+(2,),[-1.,1])
1150 a_3=makeArray(sh+(3,),[-1.,1])
1151 b_3=makeArray(sh+(3,),[-1.,1])
1152 if data=="Symbol":
1153 text+=" arg=Symbol(shape=%s)\n"%str(sh)
1154 val="s"
1155 res="sub"
1156 else:
1157 val="arg"
1158 res="res"
1159
1160 text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
1161 if not len(sh)==0:
1162 text+=" ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)
1163 text+=" if dim==2:\n"
1164 text+=unrollLoops(a_2,b_2,o,val,tap=" ")
1165 text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap=" ")
1166 text+=" else:\n"
1167
1168 text+=unrollLoops(a_3,b_3,o,val,tap=" ")
1169 text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap=" ")
1170 if case in ["ContinuousFunction","Solution","ReducedSolution"]:
1171 text+=" res=integrate(arg,where=w_ref)\n"
1172 else:
1173 text+=" res=integrate(arg)\n"
1174
1175 if data=="Symbol":
1176 text+=" sub=res.substitute({arg:s})\n"
1177 if len(sh)==0 and data=="Data":
1178 text+=" self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res
1179 else:
1180 if data=="Symbol":
1181 text+=" self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"
1182 text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
1183 else:
1184 text+=" self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n"
1185 text+=" self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)
1186 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1187
1188
1189 t_prog+=text
1190 print test_header
1191 print t_prog
1192 print test_tail
1193 1/0
1194 #=======================================================================================================
1195 # inverse
1196 #=======================================================================================================
1197 name="inverse"
1198 for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
1199 for sh0 in [ (1,1), (2,2), (3,3)]:
1200 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1201 tname="test_%s_%s_dim%s"%(name,case0,sh0[0])
1202 text+=" def %s(self):\n"%tname
1203 a_0=makeArray(sh0,[-1.,1])
1204 for i in range(sh0[0]): a_0[i,i]+=2
1205 if case0 in ["taggedData", "expandedData"]:
1206 a1_0=makeArray(sh0,[-1.,1])
1207 for i in range(sh0[0]): a1_0[i,i]+=3
1208 else:
1209 a1_0=a_0
1210
1211 text+=mkText(case0,"arg",a_0,a1_0)
1212 text+=" res=%s(arg)\n"%name
1213 if case0=="Symbol":
1214 text+=mkText("array","s",a_0,a1_0)
1215 text+=" sub=res.substitute({arg:s})\n"
1216 res="sub"
1217 ref="s"
1218 else:
1219 ref="arg"
1220 res="res"
1221 text+=mkTypeAndShapeTest(case0,sh0,"res")
1222 text+=" self.failUnless(Lsup(matrixmult(%s,%s)-kronecker(%s))<=self.RES_TOL,\"wrong result\")\n"%(res,ref,sh0[0])
1223
1224 if case0 == "taggedData" :
1225 t_prog_with_tags+=text
1226 else:
1227 t_prog+=text
1228
1229 print test_header
1230 # print t_prog
1231 print t_prog_with_tags
1232 print test_tail
1233 1/0
1234
1235 #=======================================================================================================
1236 # trace
1237 #=======================================================================================================
1238 def traceTest(r,offset):
1239 sh=r.shape
1240 r1=1
1241 for i in range(offset): r1*=sh[i]
1242 r2=1
1243 for i in range(offset+2,len(sh)): r2*=sh[i]
1244 r_s=numarray.reshape(r,(r1,sh[offset],sh[offset],r2))
1245 s=numarray.zeros([r1,r2],numarray.Float)
1246 for i1 in range(r1):
1247 for i2 in range(r2):
1248 for j in range(sh[offset]): s[i1,i2]+=r_s[i1,j,j,i2]
1249 return s.resize(sh[:offset]+sh[offset+2:])
1250 name,tt="trace",traceTest
1251 for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
1252 for sh0 in [ (4,5), (6,2,2),(3,2,3,4)]:
1253 for offset in range(len(sh0)-1):
1254 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1255 tname="test_%s_%s_rank%s_offset%s"%(name,case0,len(sh0),offset)
1256 text+=" def %s(self):\n"%tname
1257 sh_t=list(sh0)
1258 sh_t[offset+1]=sh_t[offset]
1259 sh_t=tuple(sh_t)
1260 sh_r=[]
1261 for i in range(offset): sh_r.append(sh0[i])
1262 for i in range(offset+2,len(sh0)): sh_r.append(sh0[i])
1263 sh_r=tuple(sh_r)
1264 a_0=makeArray(sh_t,[-1.,1])
1265 if case0 in ["taggedData", "expandedData"]:
1266 a1_0=makeArray(sh_t,[-1.,1])
1267 else:
1268 a1_0=a_0
1269 r=tt(a_0,offset)
1270 r1=tt(a1_0,offset)
1271 text+=mkText(case0,"arg",a_0,a1_0)
1272 text+=" res=%s(arg,%s)\n"%(name,offset)
1273 if case0=="Symbol":
1274 text+=mkText("array","s",a_0,a1_0)
1275 text+=" sub=res.substitute({arg:s})\n"
1276 res="sub"
1277 text+=mkText("array","ref",r,r1)
1278 else:
1279 res="res"
1280 text+=mkText(case0,"ref",r,r1)
1281 text+=mkTypeAndShapeTest(case0,sh_r,"res")
1282 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1283
1284 if case0 == "taggedData" :
1285 t_prog_with_tags+=text
1286 else:
1287 t_prog+=text
1288
1289 print test_header
1290 # print t_prog
1291 print t_prog_with_tags
1292 print test_tail
1293 1/0
1294
1295 #=======================================================================================================
1296 # clip
1297 #=======================================================================================================
1298 oper_L=[["clip",clipTEST]]
1299 for oper in oper_L:
1300 for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1301 for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1302 if len(sh0)==0 or not case0=="float":
1303 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1304 tname="test_%s_%s_rank%s"%(oper[0],case0,len(sh0))
1305 text+=" def %s(self):\n"%tname
1306 a_0=makeArray(sh0,[-1.,1])
1307 if case0 in ["taggedData", "expandedData"]:
1308 a1_0=makeArray(sh0,[-1.,1])
1309 else:
1310 a1_0=a_0
1311
1312 r=oper[1](a_0,-0.3,0.5)
1313 r1=oper[1](a1_0,-0.3,0.5)
1314 text+=mkText(case0,"arg",a_0,a1_0)
1315 text+=" res=%s(arg,-0.3,0.5)\n"%oper[0]
1316 if case0=="Symbol":
1317 text+=mkText("array","s",a_0,a1_0)
1318 text+=" sub=res.substitute({arg:s})\n"
1319 res="sub"
1320 text+=mkText("array","ref",r,r1)
1321 else:
1322 res="res"
1323 text+=mkText(case0,"ref",r,r1)
1324 text+=mkTypeAndShapeTest(case0,sh0,"res")
1325 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1326
1327 if case0 == "taggedData" :
1328 t_prog_with_tags+=text
1329 else:
1330 t_prog+=text
1331
1332 print test_header
1333 # print t_prog
1334 print t_prog_with_tags
1335 print test_tail
1336 1/0
1337
1338 #=======================================================================================================
1339 # maximum, minimum, clipping
1340 #=======================================================================================================
1341 oper_L=[ ["maximum",maximumTEST],
1342 ["minimum",minimumTEST]]
1343 for oper in oper_L:
1344 for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1345 for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1346 for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1347 for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1348 if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \
1349 and (sh0==sh1 or len(sh0)==0 or len(sh1)==0) :
1350 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1351
1352 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1353 tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1354 text+=" def %s(self):\n"%tname
1355 a_0=makeArray(sh0,[-1.,1])
1356 if case0 in ["taggedData", "expandedData"]:
1357 a1_0=makeArray(sh0,[-1.,1])
1358 else:
1359 a1_0=a_0
1360
1361 a_1=makeArray(sh1,[-1.,1])
1362 if case1 in ["taggedData", "expandedData"]:
1363 a1_1=makeArray(sh1,[-1.,1])
1364 else:
1365 a1_1=a_1
1366 r=oper[1](a_0,a_1)
1367 r1=oper[1](a1_0,a1_1)
1368 text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1369 text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1370 text+=" res=%s(arg0,arg1)\n"%oper[0]
1371 case=getResultCaseForBin(case0,case1)
1372 if case=="Symbol":
1373 c0_res,c1_res=case0,case1
1374 subs="{"
1375 if case0=="Symbol":
1376 text+=mkText("array","s0",a_0,a1_0)
1377 subs+="arg0:s0"
1378 c0_res="array"
1379 if case1=="Symbol":
1380 text+=mkText("array","s1",a_1,a1_1)
1381 if not subs.endswith("{"): subs+=","
1382 subs+="arg1:s1"
1383 c1_res="array"
1384 subs+="}"
1385 text+=" sub=res.substitute(%s)\n"%subs
1386 res="sub"
1387 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1388 else:
1389 res="res"
1390 text+=mkText(case,"ref",r,r1)
1391 if len(sh0)>len(sh1):
1392 text+=mkTypeAndShapeTest(case,sh0,"res")
1393 else:
1394 text+=mkTypeAndShapeTest(case,sh1,"res")
1395 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1396
1397 if case0 == "taggedData" or case1 == "taggedData":
1398 t_prog_with_tags+=text
1399 else:
1400 t_prog+=text
1401
1402 print test_header
1403 # print t_prog
1404 print t_prog_with_tags
1405 print test_tail
1406 1/0
1407
1408
1409 #=======================================================================================================
1410 # outer inner
1411 #=======================================================================================================
1412 oper=["outer",outerTEST]
1413 # oper=["inner",innerTEST]
1414 for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1415 for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1416 for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1417 for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1418 if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \
1419 and len(sh0+sh1)<5:
1420 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1421
1422 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1423 tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1424 text+=" def %s(self):\n"%tname
1425 a_0=makeArray(sh0,[-1.,1])
1426 if case0 in ["taggedData", "expandedData"]:
1427 a1_0=makeArray(sh0,[-1.,1])
1428 else:
1429 a1_0=a_0
1430
1431 a_1=makeArray(sh1,[-1.,1])
1432 if case1 in ["taggedData", "expandedData"]:
1433 a1_1=makeArray(sh1,[-1.,1])
1434 else:
1435 a1_1=a_1
1436 r=oper[1](a_0,a_1)
1437 r1=oper[1](a1_0,a1_1)
1438 text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1439 text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1440 text+=" res=%s(arg0,arg1)\n"%oper[0]
1441 case=getResultCaseForBin(case0,case1)
1442 if case=="Symbol":
1443 c0_res,c1_res=case0,case1
1444 subs="{"
1445 if case0=="Symbol":
1446 text+=mkText("array","s0",a_0,a1_0)
1447 subs+="arg0:s0"
1448 c0_res="array"
1449 if case1=="Symbol":
1450 text+=mkText("array","s1",a_1,a1_1)
1451 if not subs.endswith("{"): subs+=","
1452 subs+="arg1:s1"
1453 c1_res="array"
1454 subs+="}"
1455 text+=" sub=res.substitute(%s)\n"%subs
1456 res="sub"
1457 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1458 else:
1459 res="res"
1460 text+=mkText(case,"ref",r,r1)
1461 text+=mkTypeAndShapeTest(case,sh0+sh1,"res")
1462 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1463
1464 if case0 == "taggedData" or case1 == "taggedData":
1465 t_prog_with_tags+=text
1466 else:
1467 t_prog+=text
1468
1469 print test_header
1470 # print t_prog
1471 print t_prog_with_tags
1472 print test_tail
1473 1/0
1474
1475 #=======================================================================================================
1476 # local reduction
1477 #=======================================================================================================
1478 for oper in [["length",0.,"out+%a1%**2","math.sqrt(out)"],
1479 ["maxval",-1.e99,"max(out,%a1%)","out"],
1480 ["minval",1.e99,"min(out,%a1%)","out"] ]:
1481 for case in case_set:
1482 for sh in shape_set:
1483 if not case=="float" or len(sh)==0:
1484 text=""
1485 text+=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1486 tname="def test_%s_%s_rank%s"%(oper[0],case,len(sh))
1487 text+=" %s(self):\n"%tname
1488 a=makeArray(sh,[-1.,1.])
1489 a1=makeArray(sh,[-1.,1.])
1490 r1=testReduce(a1,oper[1],oper[2],oper[3])
1491 r=testReduce(a,oper[1],oper[2],oper[3])
1492
1493 text+=mkText(case,"arg",a,a1)
1494 text+=" res=%s(arg)\n"%oper[0]
1495 if case=="Symbol":
1496 text+=mkText("array","s",a,a1)
1497 text+=" sub=res.substitute({arg:s})\n"
1498 text+=mkText("array","ref",r,r1)
1499 res="sub"
1500 else:
1501 text+=mkText(case,"ref",r,r1)
1502 res="res"
1503 if oper[0]=="length":
1504 text+=mkTypeAndShapeTest(case,(),"res")
1505 else:
1506 if case=="float" or case=="array":
1507 text+=mkTypeAndShapeTest("float",(),"res")
1508 else:
1509 text+=mkTypeAndShapeTest(case,(),"res")
1510 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1511 if case == "taggedData":
1512 t_prog_with_tags+=text
1513 else:
1514 t_prog+=text
1515 print test_header
1516 # print t_prog
1517 print t_prog_with_tags
1518 print test_tail
1519 1/0
1520
1521 #=======================================================================================================
1522 # tensor multiply
1523 #=======================================================================================================
1524 # oper=["generalTensorProduct",tensorProductTest]
1525 # oper=["matrixmult",testMatrixMult]
1526 oper=["tensormult",testTensorMult]
1527
1528 for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1529 for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1530 for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1531 for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1532 for sh_s in [ (),(3,), (2,3), (2,4,3),(4,2,3,2)]:
1533 if (len(sh0+sh_s)==0 or not case0=="float") and (len(sh1+sh_s)==0 or not case1=="float") \
1534 and len(sh0+sh1)<5 and len(sh0+sh_s)<5 and len(sh1+sh_s)<5:
1535 # if len(sh_s)==1 and len(sh0+sh_s)==2 and (len(sh_s+sh1)==1 or len(sh_s+sh1)==2)): # test for matrixmult
1536 if ( len(sh_s)==1 and len(sh0+sh_s)==2 and ( len(sh1+sh_s)==2 or len(sh1+sh_s)==1 )) or (len(sh_s)==2 and len(sh0+sh_s)==4 and (len(sh1+sh_s)==2 or len(sh1+sh_s)==3 or len(sh1+sh_s)==4)): # test for tensormult
1537 case=getResultCaseForBin(case0,case1)
1538 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1539 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1540 # tname="test_generalTensorProduct_%s_rank%s_%s_rank%s_offset%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1),len(sh_s))
1541 #tname="test_matrixmult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
1542 tname="test_tensormult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
1543 # if tname=="test_generalTensorProduct_array_rank1_array_rank2_offset1":
1544 # print tnametest_generalTensorProduct_Symbol_rank1_Symbol_rank3_offset1
1545 text+=" def %s(self):\n"%tname
1546 a_0=makeArray(sh0+sh_s,[-1.,1])
1547 if case0 in ["taggedData", "expandedData"]:
1548 a1_0=makeArray(sh0+sh_s,[-1.,1])
1549 else:
1550 a1_0=a_0
1551
1552 a_1=makeArray(sh_s+sh1,[-1.,1])
1553 if case1 in ["taggedData", "expandedData"]:
1554 a1_1=makeArray(sh_s+sh1,[-1.,1])
1555 else:
1556 a1_1=a_1
1557 r=oper[1](a_0,a_1,sh_s)
1558 r1=oper[1](a1_0,a1_1,sh_s)
1559 text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1560 text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1561 #text+=" res=matrixmult(arg0,arg1)\n"
1562 text+=" res=tensormult(arg0,arg1)\n"
1563 #text+=" res=generalTensorProduct(arg0,arg1,offset=%s)\n"%(len(sh_s))
1564 if case=="Symbol":
1565 c0_res,c1_res=case0,case1
1566 subs="{"
1567 if case0=="Symbol":
1568 text+=mkText("array","s0",a_0,a1_0)
1569 subs+="arg0:s0"
1570 c0_res="array"
1571 if case1=="Symbol":
1572 text+=mkText("array","s1",a_1,a1_1)
1573 if not subs.endswith("{"): subs+=","
1574 subs+="arg1:s1"
1575 c1_res="array"
1576 subs+="}"
1577 text+=" sub=res.substitute(%s)\n"%subs
1578 res="sub"
1579 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1580 else:
1581 res="res"
1582 text+=mkText(case,"ref",r,r1)
1583 text+=mkTypeAndShapeTest(case,sh0+sh1,"res")
1584 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1585 if case0 == "taggedData" or case1 == "taggedData":
1586 t_prog_with_tags+=text
1587 else:
1588 t_prog+=text
1589 print test_header
1590 # print t_prog
1591 print t_prog_with_tags
1592 print test_tail
1593 1/0
1594 #=======================================================================================================
1595 # basic binary operation overloading (tests only!)
1596 #=======================================================================================================
1597 oper_range=[-5.,5.]
1598 for oper in [["add" ,"+",[-5.,5.]],
1599 ["sub" ,"-",[-5.,5.]],
1600 ["mult","*",[-5.,5.]],
1601 ["div" ,"/",[-5.,5.]],
1602 ["pow" ,"**",[0.01,5.]]]:
1603 for case0 in case_set:
1604 for sh0 in shape_set:
1605 for case1 in case_set:
1606 for sh1 in shape_set:
1607 if not case0=="array" and \
1608 (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1609 (sh0==() or sh1==() or sh1==sh0) and \
1610 not (case0 in ["float","array"] and case1 in ["float","array"]):
1611 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1612 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1613 tname="test_%s_overloaded_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1614 text+=" def %s(self):\n"%tname
1615 a_0=makeArray(sh0,oper[2])
1616 if case0 in ["taggedData", "expandedData"]:
1617 a1_0=makeArray(sh0,oper[2])
1618 else:
1619 a1_0=a_0
1620
1621 a_1=makeArray(sh1,oper[2])
1622 if case1 in ["taggedData", "expandedData"]:
1623 a1_1=makeArray(sh1,oper[2])
1624 else:
1625 a1_1=a_1
1626 r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1627 r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1628 text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1629 text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1630 text+=" res=arg0%sarg1\n"%oper[1]
1631
1632 case=getResultCaseForBin(case0,case1)
1633 if case=="Symbol":
1634 c0_res,c1_res=case0,case1
1635 subs="{"
1636 if case0=="Symbol":
1637 text+=mkText("array","s0",a_0,a1_0)
1638 subs+="arg0:s0"
1639 c0_res="array"
1640 if case1=="Symbol":
1641 text+=mkText("array","s1",a_1,a1_1)
1642 if not subs.endswith("{"): subs+=","
1643 subs+="arg1:s1"
1644 c1_res="array"
1645 subs+="}"
1646 text+=" sub=res.substitute(%s)\n"%subs
1647 res="sub"
1648 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1649 else:
1650 res="res"
1651 text+=mkText(case,"ref",r,r1)
1652 if isinstance(r,float):
1653 text+=mkTypeAndShapeTest(case,(),"res")
1654 else:
1655 text+=mkTypeAndShapeTest(case,r.shape,"res")
1656 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1657
1658 if case0 in [ "constData","taggedData","expandedData"] and case1 == "Symbol":
1659 t_prog_failing+=text
1660 else:
1661 if case0 == "taggedData" or case1 == "taggedData":
1662 t_prog_with_tags+=text
1663 else:
1664 t_prog+=text
1665
1666
1667 print test_header
1668 # print t_prog
1669 # print t_prog_with_tags
1670 print t_prog_failing
1671 print test_tail
1672 1/0
1673 #=======================================================================================================
1674 # basic binary operations (tests only!)
1675 #=======================================================================================================
1676 oper_range=[-5.,5.]
1677 for oper in [["add" ,"+",[-5.,5.]],
1678 ["mult","*",[-5.,5.]],
1679 ["quotient" ,"/",[-5.,5.]],
1680 ["power" ,"**",[0.01,5.]]]:
1681 for case0 in case_set:
1682 for case1 in case_set:
1683 for sh in shape_set:
1684 for sh_p in shape_set:
1685 if len(sh_p)>0:
1686 resource=[-1,1]
1687 else:
1688 resource=[1]
1689 for sh_d in resource:
1690 if sh_d>0:
1691 sh0=sh
1692 sh1=sh+sh_p
1693 else:
1694 sh1=sh
1695 sh0=sh+sh_p
1696
1697 if (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1698 len(sh0)<5 and len(sh1)<5:
1699 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1700 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1701 tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1702 text+=" def %s(self):\n"%tname
1703 a_0=makeArray(sh0,oper[2])
1704 if case0 in ["taggedData", "expandedData"]:
1705 a1_0=makeArray(sh0,oper[2])
1706 else:
1707 a1_0=a_0
1708
1709 a_1=makeArray(sh1,oper[2])
1710 if case1 in ["taggedData", "expandedData"]:
1711 a1_1=makeArray(sh1,oper[2])
1712 else:
1713 a1_1=a_1
1714 r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1715 r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1716 text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1717 text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1718 text+=" res=%s(arg0,arg1)\n"%oper[0]
1719
1720 case=getResultCaseForBin(case0,case1)
1721 if case=="Symbol":
1722 c0_res,c1_res=case0,case1
1723 subs="{"
1724 if case0=="Symbol":
1725 text+=mkText("array","s0",a_0,a1_0)
1726 subs+="arg0:s0"
1727 c0_res="array"
1728 if case1=="Symbol":
1729 text+=mkText("array","s1",a_1,a1_1)
1730 if not subs.endswith("{"): subs+=","
1731 subs+="arg1:s1"
1732 c1_res="array"
1733 subs+="}"
1734 text+=" sub=res.substitute(%s)\n"%subs
1735 res="sub"
1736 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1737 else:
1738 res="res"
1739 text+=mkText(case,"ref",r,r1)
1740 if isinstance(r,float):
1741 text+=mkTypeAndShapeTest(case,(),"res")
1742 else:
1743 text+=mkTypeAndShapeTest(case,r.shape,"res")
1744 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1745
1746 if case0 == "taggedData" or case1 == "taggedData":
1747 t_prog_with_tags+=text
1748 else:
1749 t_prog+=text
1750 print test_header
1751 # print t_prog
1752 print t_prog_with_tags
1753 print test_tail
1754 1/0
1755
1756 # print t_prog_with_tagsoper_range=[-5.,5.]
1757 for oper in [["add" ,"+",[-5.,5.]],
1758 ["sub" ,"-",[-5.,5.]],
1759 ["mult","*",[-5.,5.]],
1760 ["div" ,"/",[-5.,5.]],
1761 ["pow" ,"**",[0.01,5.]]]:
1762 for case0 in case_set:
1763 for sh0 in shape_set:
1764 for case1 in case_set:
1765 for sh1 in shape_set:
1766 if (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1767 (sh0==() or sh1==() or sh1==sh0) and \
1768 not (case0 in ["float","array"] and case1 in ["float","array"]):
1769 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1770 tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1771 text+=" def %s(self):\n"%tname
1772 a_0=makeArray(sh0,oper[2])
1773 if case0 in ["taggedData", "expandedData"]:
1774 a1_0=makeArray(sh0,oper[2])
1775 else:
1776 a1_0=a_0
1777
1778 a_1=makeArray(sh1,oper[2])
1779 if case1 in ["taggedData", "expandedData"]:
1780 a1_1=makeArray(sh1,oper[2])
1781 else:
1782 a1_1=a_1
1783 r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1784 r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1785 text+=mkText(case0,"arg0",a_0,a1_0)
1786 text+=mkText(case1,"arg1",a_1,a1_1)
1787 text+=" res=arg0%sarg1\n"%oper[1]
1788
1789 case=getResultCaseForBin(case0,case1)
1790 if case=="Symbol":
1791 c0_res,c1_res=case0,case1
1792 subs="{"
1793 if case0=="Symbol":
1794 text+=mkText("array","s0",a_0,a1_0)
1795 subs+="arg0:s0"
1796 c0_res="array"
1797 if case1=="Symbol":
1798 text+=mkText("array","s1",a_1,a1_1)
1799 if not subs.endswith("{"): subs+=","
1800 subs+="arg1:s1"
1801 c1_res="array"
1802 subs+="}"
1803 text+=" sub=res.substitute(%s)\n"%subs
1804 res="sub"
1805 text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1806 else:
1807 res="res"
1808 text+=mkText(case,"ref",r,r1)
1809 if isinstance(r,float):
1810 text+=mkTypeAndShapeTest(case,(),"res")
1811 else:
1812 text+=mkTypeAndShapeTest(case,r.shape,"res")
1813 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1814
1815 if case0 in [ "constData","taggedData","expandedData"] and case1 == "Symbol":
1816 t_prog_failing+=text
1817 else:
1818 if case0 == "taggedData" or case1 == "taggedData":
1819 t_prog_with_tags+=text
1820 else:
1821 t_prog+=text
1822
1823
1824 # print u_prog
1825 # 1/0
1826 print test_header
1827 print t_prog
1828 # print t_prog_with_tags
1829 # print t_prog_failing
1830 print test_tail
1831 # print t_prog_failing
1832 print test_tail
1833
1834 #=======================================================================================================
1835 # unary operations:
1836 #=======================================================================================================
1837 func= [
1838 OPERATOR(nickname="log10",\
1839 rng=[1.e-3,100.],\
1840 test_expr="math.log10(%a1%)",\
1841 math_expr="math.log10(%a1%)",\
1842 numarray_expr="numarray.log10(%a1%)",\
1843 symbol_expr="log(%a1%)/log(10.)",\
1844 name="base-10 logarithm"),
1845 OPERATOR(nickname="wherePositive",\
1846 rng=[-100.,100.],\
1847 test_expr="wherepos(%a1%)",\
1848 math_expr="if arg>0:\n return 1.\nelse:\n return 0.",
1849 numarray_expr="numarray.greater(arg,numarray.zeros(arg.shape,numarray.Float))",\
1850 name="mask of positive values"),
1851 OPERATOR(nickname="whereNegative",\
1852 rng=[-100.,100.],\
1853 test_expr="wherepos(-%a1%)",\
1854 math_expr="if arg<0:\n return 1.\nelse:\n return 0.",
1855 numarray_expr="numarray.less(arg,numarray.zeros(arg.shape,numarray.Float))",\
1856 name="mask of positive values"),
1857 OPERATOR(nickname="whereNonNegative",\
1858 rng=[-100.,100.],\
1859 test_expr="1-wherepos(-%a1%)", \
1860 math_expr="if arg<0:\n return 0.\nelse:\n return 1.",
1861 numarray_expr="numarray.greater_equal(arg,numarray.zeros(arg.shape,numarray.Float))",\
1862 symbol_expr="1-wherePositive(%a1%)",\
1863 name="mask of non-negative values"),
1864 OPERATOR(nickname="whereNonPositive",\
1865 rng=[-100.,100.],\
1866 test_expr="1-wherepos(%a1%)",\
1867 math_expr="if arg>0:\n return 0.\nelse:\n return 1.",
1868 numarray_expr="numarray.less_equal(arg,numarray.zeros(arg.shape,numarray.Float))",\
1869 symbol_expr="1-whereNegative(%a1%)",\
1870 name="mask of non-positive values"),
1871 OPERATOR(nickname="whereZero",\
1872 rng=[-100.,100.],\
1873 test_expr="1-wherepos(%a1%)-wherepos(-%a1%)",\
1874 math_expr="if abs(%a1%)<=tol:\n return 1.\nelse:\n return 0.",
1875 numarray_expr="numarray.less_equal(abs(%a1%)-tol,numarray.zeros(arg.shape,numarray.Float))",\
1876 name="mask of zero entries"),
1877 OPERATOR(nickname="whereNonZero",\
1878 rng=[-100.,100.],\
1879 test_expr="wherepos(%a1%)+wherepos(-%a1%)",\
1880 math_expr="if abs(%a1%)>tol:\n return 1.\nelse:\n return 0.",\
1881 numarray_expr="numarray.greater(abs(%a1%)-tol,numarray.zeros(arg.shape,numarray.Float))",\
1882 symbol_expr="1-whereZero(arg,tol)",\
1883 name="mask of values different from zero"),
1884 OPERATOR(nickname="sin",\
1885 rng=[-100.,100.],\
1886 test_expr="math.sin(%a1%)",
1887 numarray_expr="numarray.sin(%a1%)",\
1888 diff="cos(%a1%)",\
1889 name="sine"),
1890 OPERATOR(nickname="cos",\
1891 rng=[-100.,100.],\
1892 test_expr="math.cos(%a1%)",
1893 numarray_expr="numarray.cos(%a1%)",\
1894 diff="-sin(%a1%)",
1895 name="cosine"),
1896 OPERATOR(nickname="tan",\
1897 rng=[-100.,100.],\
1898 test_expr="math.tan(%a1%)",
1899 numarray_expr="numarray.tan(%a1%)",\
1900 diff="1./cos(%a1%)**2",
1901 name="tangent"),
1902 OPERATOR(nickname="asin",\
1903 rng=[-0.99,0.99],\
1904 test_expr="math.asin(%a1%)",
1905 numarray_expr="numarray.arcsin(%a1%)",
1906 diff="1./sqrt(1.-%a1%**2)",
1907 name="inverse sine"),
1908 OPERATOR(nickname="acos",\
1909 rng=[-0.99,0.99],\
1910 test_expr="math.acos(%a1%)",
1911 numarray_expr="numarray.arccos(%a1%)",
1912 diff="-1./sqrt(1.-%a1%**2)",
1913 name="inverse cosine"),
1914 OPERATOR(nickname="atan",\
1915 rng=[-100.,100.],\
1916 test_expr="math.atan(%a1%)",
1917 numarray_expr="numarray.arctan(%a1%)",
1918 diff="1./(1+%a1%**2)",
1919 name="inverse tangent"),
1920 OPERATOR(nickname="sinh",\
1921 rng=[-5,5],\
1922 test_expr="math.sinh(%a1%)",
1923 numarray_expr="numarray.sinh(%a1%)",
1924 diff="cosh(%a1%)",
1925 name="hyperbolic sine"),
1926 OPERATOR(nickname="cosh",\
1927 rng=[-5.,5.],
1928 test_expr="math.cosh(%a1%)",
1929 numarray_expr="numarray.cosh(%a1%)",
1930 diff="sinh(%a1%)",
1931 name="hyperbolic cosine"),
1932 OPERATOR(nickname="tanh",\
1933 rng=[-5.,5.],
1934 test_expr="math.tanh(%a1%)",
1935 numarray_expr="numarray.tanh(%a1%)",
1936 diff="1./cosh(%a1%)**2",
1937 name="hyperbolic tangent"),
1938 OPERATOR(nickname="asinh",\
1939 rng=[-100.,100.], \
1940 test_expr="numarray.arcsinh(%a1%)",
1941 math_expr="numarray.arcsinh(%a1%)",
1942 numarray_expr="numarray.arcsinh(%a1%)",
1943 diff="1./sqrt(%a1%**2+1)",
1944 name="inverse hyperbolic sine"),
1945 OPERATOR(nickname="acosh",\
1946 rng=[1.001,100.],\
1947 test_expr="numarray.arccosh(%a1%)",
1948 math_expr="numarray.arccosh(%a1%)",
1949 numarray_expr="numarray.arccosh(%a1%)",
1950 diff="1./sqrt(%a1%**2-1)",
1951 name="inverse hyperolic cosine"),
1952 OPERATOR(nickname="atanh",\
1953 rng=[-0.99,0.99], \
1954 test_expr="numarray.arctanh(%a1%)",
1955 math_expr="numarray.arctanh(%a1%)",
1956 numarray_expr="numarray.arctanh(%a1%)",
1957 diff="1./(1.-%a1%**2)",
1958 name="inverse hyperbolic tangent"),
1959 OPERATOR(nickname="exp",\
1960 rng=[-5.,5.],
1961 test_expr="math.exp(%a1%)",
1962 numarray_expr="numarray.exp(%a1%)",
1963 diff="self",
1964 name="exponential"),
1965 OPERATOR(nickname="sqrt",\
1966 rng=[1.e-3,100.],\
1967 test_expr="math.sqrt(%a1%)",
1968 numarray_expr="numarray.sqrt(%a1%)",
1969 diff="0.5/self",
1970 name="square root"),
1971 OPERATOR(nickname="log", \
1972 rng=[1.e-3,100.],\
1973 test_expr="math.log(%a1%)",
1974 numarray_expr="numarray.log(%a1%)",
1975 diff="1./arg",
1976 name="natural logarithm"),
1977 OPERATOR(nickname="sign",\
1978 rng=[-100.,100.], \
1979 math_expr="if %a1%>0:\n return 1.\nelif %a1%<0:\n return -1.\nelse:\n return 0.",
1980 test_expr="wherepos(%a1%)-wherepos(-%a1%)",
1981 numarray_expr="numarray.sign(%a1%)",
1982 symbol_expr="wherePositive(%a1%)-whereNegative(%a1%)",\
1983 name="sign"),
1984 OPERATOR(nickname="abs",\
1985 rng=[-100.,100.], \
1986 math_expr="if %a1%>0:\n return %a1% \nelif %a1%<0:\n return -(%a1%)\nelse:\n return 0.",
1987 test_expr="wherepos(%a1%)*(%a1%)-wherepos(-%a1%)*(%a1%)",
1988 numarray_expr="abs(%a1%)",
1989 diff="sign(%a1%)",
1990 name="absolute value")
1991
1992 ]
1993 for f in func:
1994 symbol_name=f.nickname[0].upper()+f.nickname[1:]
1995 if f.nickname!="abs":
1996 u_prog+="def %s(arg):\n"%f.nickname
1997 u_prog+=" \"\"\"\n"
1998 u_prog+=" returns %s of argument arg\n\n"%f.name
1999 u_prog+=" @param arg: argument\n"
2000 u_prog+=" @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.\n"
2001 u_prog+=" @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.\n"
2002 u_prog+=" @raises TypeError: if the type of the argument is not expected.\n"
2003 u_prog+=" \"\"\"\n"
2004 u_prog+=" if isinstance(arg,numarray.NumArray):\n"
2005 u_prog+=mkCode(f.numarray_expr,["arg"],2*" ")
2006 u_prog+=" elif isinstance(arg,escript.Data):\n"
2007 u_prog+=mkCode("arg._%s()"%f.nickname,[],2*" ")
2008 u_prog+=" elif isinstance(arg,float):\n"
2009 u_prog+=mkCode(f.math_expr,["arg"],2*" ")
2010 u_prog+=" elif isinstance(arg,int):\n"
2011 u_prog+=mkCode(f.math_expr,["float(arg)"],2*" ")
2012 u_prog+=" elif isinstance(arg,Symbol):\n"
2013 if f.symbol_expr==None:
2014 u_prog+=mkCode("%s_Symbol(arg)"%symbol_name,[],2*" ")
2015 else:
2016 u_prog+=mkCode(f.symbol_expr,["arg"],2*" ")
2017 u_prog+=" else:\n"
2018 u_prog+=" raise TypeError,\"%s: Unknown argument type.\"\n\n"%f.nickname
2019 if f.symbol_expr==None:
2020 u_prog+="class %s_Symbol(DependendSymbol):\n"%symbol_name
2021 u_prog+=" \"\"\"\n"
2022 u_prog+=" L{Symbol} representing the result of the %s function\n"%f.name
2023 u_prog+=" \"\"\"\n"
2024 u_prog+=" def __init__(self,arg):\n"
2025 u_prog+=" \"\"\"\n"
2026 u_prog+=" initialization of %s L{Symbol} with argument arg\n"%f.nickname
2027 u_prog+=" @param arg: argument of function\n"
2028 u_prog+=" @type arg: typically L{Symbol}.\n"
2029 u_prog+=" \"\"\"\n"
2030 u_prog+=" DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())\n"
2031 u_prog+="\n"
2032
2033 u_prog+=" def getMyCode(self,argstrs,format=\"escript\"):\n"
2034 u_prog+=" \"\"\"\n"
2035 u_prog+=" returns a program code that can be used to evaluate the symbol.\n\n"
2036
2037 u_prog+=" @param argstrs: gives for each argument a string representing the argument for the evaluation.\n"
2038 u_prog+=" @type argstrs: C{str} or a C{list} of length 1 of C{str}.\n"
2039 u_prog+=" @param format: specifies the format to be used. At the moment only \"escript\" ,\"text\" and \"str\" are supported.\n"
2040 u_prog+=" @type format: C{str}\n"
2041 u_prog+=" @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.\n"
2042 u_prog+=" @rtype: C{str}\n"
2043 u_prog+=" @raise: NotImplementedError: if the requested format is not available\n"
2044 u_prog+=" \"\"\"\n"
2045 u_prog+=" if isinstance(argstrs,list):\n"
2046 u_prog+=" argstrs=argstrs[0]\n"
2047 u_prog+=" if format==\"escript\" or format==\"str\" or format==\"text\":\n"
2048 u_prog+=" return \"%s(%%s)\"%%argstrs\n"%f.nickname
2049 u_prog+=" else:\n"
2050 u_prog+=" raise NotImplementedError,\"%s_Symbol does not provide program code for format %%s.\"%%format\n"%symbol_name
2051 u_prog+="\n"
2052
2053 u_prog+=" def substitute(self,argvals):\n"
2054 u_prog+=" \"\"\"\n"
2055 u_prog+=" assigns new values to symbols in the definition of the symbol.\n"
2056 u_prog+=" The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.\n"
2057 u_prog+="\n"
2058 u_prog+=" @param argvals: new values assigned to symbols\n"
2059 u_prog+=" @type argvals: C{dict} with keywords of type L{Symbol}.\n"
2060 u_prog+=" @return: result of the substitution process. Operations are executed as much as possible.\n"
2061 u_prog+=" @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution\n"
2062 u_prog+=" @raise TypeError: if a value for a L{Symbol} cannot be substituted.\n"
2063 u_prog+=" \"\"\"\n"
2064 u_prog+=" if argvals.has_key(self):\n"
2065 u_prog+=" arg=argvals[self]\n"
2066 u_prog+=" if self.isAppropriateValue(arg):\n"
2067 u_prog+=" return arg\n"
2068 u_prog+=" else:\n"
2069 u_prog+=" raise TypeError,\"%s: new value is not appropriate.\"%str(self)\n"
2070 u_prog+=" else:\n"
2071 u_prog+=" arg=self.getSubstitutedArguments(argvals)[0]\n"
2072 u_prog+=" return %s(arg)\n\n"%f.nickname
2073 if not f.diff==None:
2074 u_prog+=" def diff(self,arg):\n"
2075 u_prog+=" \"\"\"\n"
2076 u_prog+=" differential of this object\n"
2077 u_prog+="\n"
2078 u_prog+=" @param arg: the derivative is calculated with respect to arg\n"
2079 u_prog+=" @type arg: L{escript.Symbol}\n"
2080 u_prog+=" @return: derivative with respect to C{arg}\n"
2081 u_prog+=" @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray} are possible.\n"
2082 u_prog+=" \"\"\"\n"
2083 u_prog+=" if arg==self:\n"
2084 u_prog+=" return identity(self.getShape())\n"
2085 u_prog+=" else:\n"
2086 u_prog+=" myarg=self.getArgument()[0]\n"
2087 u_prog+=" val=matchShape(%s,self.getDifferentiatedArguments(arg)[0])\n"%f.diff.replace("%a1%","myarg")
2088 u_prog+=" return val[0]*val[1]\n\n"
2089
2090 for case in case_set:
2091 for sh in shape_set:
2092 if not case=="float" or len(sh)==0:
2093 text=""
2094 text+=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
2095 tname="def test_%s_%s_rank%s"%(f.nickname,case,len(sh))
2096 text+=" %s(self):\n"%tname
2097 a=makeArray(sh,f.rng)
2098 a1=makeArray(sh,f.rng)
2099 r1=makeResult(a1,f.test_expr)
2100 r=makeResult(a,f.test_expr)
2101
2102 text+=mkText(case,"arg",a,a1)
2103 text+=" res=%s(arg)\n"%f.nickname
2104 if case=="Symbol":
2105 text+=mkText("array","s",a,a1)
2106 text+=" sub=res.substitute({arg:s})\n"
2107 text+=mkText("array","ref",r,r1)
2108 res="sub"
2109 else:
2110 text+=mkText(case,"ref",r,r1)
2111 res="res"
2112 text+=mkTypeAndShapeTest(case,sh,"res")
2113 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
2114 if case == "taggedData":
2115 t_prog_with_tags+=text
2116 else:
2117 t_prog+=text
2118
2119 #=========== END OF GOOD CODE +++++++++++++++++++++++++++
2120
2121 1/0
2122
2123 def X():
2124 if args=="float":
2125 a=makeArray(sh,f[RANGE])
2126 r=makeResult(a,f)
2127 t_prog+=" arg=%s\n"%a[0]
2128 t_prog+=" ref=%s\n"%r[0]
2129 t_prog+=" res=%s(%a1%)\n"%f.nickname
2130 t_prog+=" self.failUnless(isinstance(res,float),\"wrong type of result.\")\n"
2131 t_prog+=" self.failUnless(Lsup(res-ref)<=self.tol*Lsup(ref),\"wrong result\")\n"
2132 elif args == "array":
2133 a=makeArray(sh,f[RANGE])
2134 r=makeResult(a,f)
2135 if len(sh)==0:
2136 t_prog+=" arg=numarray.array(%s)\n"%a[0]
2137 t_prog+=" ref=numarray.array(%s)\n"%r[0]
2138 else:
2139 t_prog+=" arg=numarray.array(%s)\n"%a.tolist()
2140 t_prog+=" ref=numarray.array(%s)\n"%r.tolist()
2141 t_prog+=" res=%s(%a1%)\n"%f.nickname
2142 t_prog+=" self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)
2143 t_prog+=" self.failUnless(Lsup(res-ref)<=self.tol*Lsup(ref),\"wrong result\")\n"
2144 elif args== "constData":
2145 a=makeArray(sh,f[RANGE])
2146 r=makeResult(a,f)
2147 if len(sh)==0:
2148 t_prog+=" arg=Data(%s,self.functionspace)\n"%(a)
2149 t_prog+=" ref=%s\n"%r
2150 else:
2151 t_prog+=" arg=Data(numarray.array(%s),self.functionspace)\n"%(a.tolist())
2152 t_prog+=" ref=numarray.array(%s)\n"%r.tolist()
2153 t_prog+=" res=%s(%a1%)\n"%f.nickname
2154 t_prog+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
2155 t_prog+=" self.failUnless(Lsup(res-ref)<=self.tol*Lsup(ref),\"wrong result\")\n"
2156 elif args in [ "taggedData","expandedData"]:
2157 a=makeArray(sh,f[RANGE])
2158 r=makeResult(a,f)
2159 a1=makeArray(sh,f[RANGE])
2160 r1=makeResult(a1,f)
2161 if len(sh)==0:
2162 if args=="expandedData":
2163 t_prog+=" arg=Data(%s,self.functionspace,True)\n"%(a)
2164 t_prog+=" ref=Data(%s,self.functionspace,True)\n"%(r)
2165 else:
2166 t_prog+=" arg=Data(%s,self.functionspace)\n"%(a)
2167 t_prog+=" ref=Data(%s,self.functionspace)\n"%(r)
2168 t_prog+=" arg.setTaggedValue(1,%s)\n"%a
2169 t_prog+=" ref.setTaggedValue(1,%s)\n"%r1
2170 else:
2171 if args=="expandedData":
2172 t_prog+=" arg=Data(numarray.array(%s),self.functionspace,True)\n"%(a.tolist())
2173 t_prog+=" ref=Data(numarray.array(%s),self.functionspace,True)\n"%(r.tolist())
2174 else:
2175 t_prog+=" arg=Data(numarray.array(%s),self.functionspace)\n"%(a.tolist())
2176 t_prog+=" ref=Data(numarray.array(%s),self.functionspace)\n"%(r.tolist())
2177 t_prog+=" arg.setTaggedValue(1,%s)\n"%a1.tolist()
2178 t_prog+=" ref.setTaggedValue(1,%s)\n"%r1.tolist()
2179 t_prog+=" res=%s(%a1%)\n"%f.nickname
2180 t_prog+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
2181 t_prog+=" self.failUnless(Lsup(res-ref)<=self.tol*Lsup(ref),\"wrong result\")\n"
2182 elif args=="Symbol":
2183 t_prog+=" arg=Symbol(shape=%s)\n"%str(sh)
2184 t_prog+=" v=%s(%a1%)\n"%f.nickname
2185 t_prog+=" self.failUnlessRaises(ValueError,v.substitute,Symbol(shape=(1,1)),\"illegal shape of substitute not identified.\")\n"
2186 a=makeArray(sh,f[RANGE])
2187 r=makeResult(a,f)
2188 if len(sh)==0:
2189 t_prog+=" res=v.substitute({arg : %s})\n"%a
2190 t_prog+=" ref=%s\n"%r
2191 t_prog+=" self.failUnless(isinstance(res,float),\"wrong type of result.\")\n"
2192 else:
2193 t_prog+=" res=v.substitute({arg : numarray.array(%s)})\n"%a.tolist()
2194 t_prog+=" ref=numarray.array(%s)\n"%r.tolist()
2195 t_prog+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of substitution result.\")\n"%str(sh)
2196 t_prog+=" self.failUnless(Lsup(res-ref)<=self.tol*Lsup(ref),\"wrong result\")\n"
2197
2198 if len(sh)==0:
2199 t_prog+=" # test derivative with respect to itself:\n"
2200 t_prog+=" dvdv=v.diff(v)\n"
2201 t_prog+=" self.failUnlessEqual(dvdv,1.,\"derivative with respect to self is not 1.\")\n"
2202 elif len(sh)==1:
2203 t_prog+=" # test derivative with respect to itself:\n"
2204 t_prog+=" dvdv=v.diff(v)\n"
2205 t_prog+=" self.failUnlessEqual(dvdv.shape,%s,\"shape of derivative with respect is wrong\")\n"%str(sh+sh)
2206 for i0_l in range(sh[0]):
2207 for i0_r in range(sh[0]):
2208 if i0_l == i0_r:
2209 v=1.
2210 else:
2211 v=0.
2212 t_prog+=" self.failUnlessEqual(dvdv[%s,%s],%s,\"derivative with respect to self: [%s,%s] is not %s\")\n"%(i0_l,i0_r,v,i0_l,i0_r,v)
2213 elif len(sh)==2:
2214 t_prog+=" # test derivative with respect to itself:\n"
2215 t_prog+=" dvdv=v.diff(v)\n"
2216 t_prog+=" self.failUnlessEqual(dvdv.shape,%s,\"shape of derivative with respect is wrong\")\n"%str(sh+sh)
2217 for i0_l in range(sh[0]):
2218 for i0_r in range(sh[0]):
2219 for i1_l in range(sh[1]):
2220 for i1_r in range(sh[1]):
2221 if i0_l == i0_r and i1_l == i1_r:
2222 v=1.
2223 else:
2224 v=0.
2225 t_prog+=" self.failUnlessEqual(dvdv[%s,%s,%s,%s],%s,\"derivative with respect to self: [%s,%s,%s,%s] is not %s\")\n"%(i0_l,i1_l,i0_r,i1_r,v,i0_l,i1_l,i0_r,i1_r,v)
2226
2227 for sh_in in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
2228 if len(sh_in)+len(sh)<=4:
2229
2230 t_prog+=" # test derivative with shape %s as argument\n"%str(sh_in)
2231 trafo=makeArray(sh+sh_in,[0,1.])
2232 a_in=makeArray(sh_in,f[RANGE])
2233 t_prog+=" arg_in=Symbol(shape=%s)\n"%str(sh_in)
2234 t_prog+=" arg2=Symbol(shape=%s)\n"%str(sh)
2235
2236 if len(sh)==0:
2237 t_prog+=" arg2="
2238 if len(sh_in)==0:
2239 t_prog+="%s*arg_in\n"%trafo
2240 elif len(sh_in)==1:
2241 for i0 in range(sh_in[0]):
2242 if i0>0: t_prog+="+"
2243 t_prog+="%s*arg_in[%s]"%(trafo[i0],i0)
2244 t_prog+="\n"
2245 elif len(sh_in)==2:
2246 for i0 in range(sh_in[0]):
2247 for i1 in range(sh_in[1]):
2248 if i0+i1>0: t_prog+="+"
2249 t_prog+="%s*arg_in[%s,%s]"%(trafo[i0,i1],i0,i1)
2250
2251 elif len(sh_in)==3:
2252 for i0 in range(sh_in[0]):
2253 for i1 in range(sh_in[1]):
2254 for i2 in range(sh_in[2]):
2255 if i0+i1+i2>0: t_prog+="+"
2256 t_prog+="%s*arg_in[%s,%s,%s]"%(trafo[i0,i1,i2],i0,i1,i2)
2257 elif len(sh_in)==4:
2258 for i0 in range(sh_in[0]):
2259 for i1 in range(sh_in[1]):
2260 for i2 in range(sh_in[2]):
2261 for i3 in range(sh_in[3]):
2262 if i0+i1+i2+i3>0: t_prog+="+"
2263 t_prog+="%s*arg_in[%s,%s,%s,%s]"%(trafo[i0,i1,i2,i3],i0,i1,i2,i3)
2264 t_prog+="\n"
2265 elif len(sh)==1:
2266 for j0 in range(sh[0]):
2267 t_prog+=" arg2[%s]="%j0
2268 if len(sh_in)==0:
2269 t_prog+="%s*arg_in"%trafo[j0]
2270 elif len(sh_in)==1:
2271 for i0 in range(sh_in[0]):
2272 if i0>0: t_prog+="+"
2273 t_prog+="%s*arg_in[%s]"%(trafo[j0,i0],i0)
2274 elif len(sh_in)==2: