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

Annotation of /trunk/escript/py_src/generatediff

Parent Directory Parent Directory | Revision Log Revision Log


Revision 439 - (hide annotations)
Fri Jan 20 01:46:22 2006 UTC (13 years, 9 months ago) by gross
Original Path: trunk/escript/py_src/generateutil
File size: 120373 byte(s)
tests for interpolation added
1 jgs 154 #!/usr/bin/python
2 gross 283 # $Id$
3 jgs 154
4     """
5     program generates parts of the util.py and the test_util.py script
6     """
7 gross 157 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 jgs 154 t_prog=""
26 gross 157 t_prog_with_tags=""
27     t_prog_failing=""
28     u_prog=""
29 jgs 154
30 gross 157 def wherepos(arg):
31     if arg>0.:
32     return 1.
33     else:
34     return 0.
35    
36 gross 313
37 gross 157 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 jgs 154 import random
53     import numarray
54     import math
55     finc=1.e-6
56    
57 gross 157 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 jgs 154
152 gross 157
153 jgs 154 def makeArray(shape,rng):
154     l=rng[1]-rng[0]
155     out=numarray.zeros(shape,numarray.Float64)
156     if len(shape)==0:
157 gross 157 out=l*random.random()+rng[0]
158 jgs 154 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 gross 437 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 jgs 154 else:
184 gross 437 raise SystemError,"rank is restricted to 5"
185 jgs 154 return out
186    
187    
188 gross 157 def makeResult(val,test_expr):
189 jgs 154 if isinstance(val,float):
190 gross 157 out=eval(test_expr.replace("%a1%","val"))
191 jgs 154 elif len(val.shape)==0:
192 gross 157 out=eval(test_expr.replace("%a1%","val"))
193 jgs 154 elif len(val.shape)==1:
194     out=numarray.zeros(val.shape,numarray.Float64)
195     for i0 in range(val.shape[0]):
196 gross 157 out[i0]=eval(test_expr.replace("%a1%","val[i0]"))
197 jgs 154 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 gross 157 out[i0,i1]=eval(test_expr.replace("%a1%","val[i0,i1]"))
202 jgs 154 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 gross 157 out[i0,i1,i2]=eval(test_expr.replace("%a1%","val[i0,i1,i2]"))
208 jgs 154 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 gross 157 out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val[i0,i1,i2,i3]"))
215 jgs 154 else:
216     raise SystemError,"rank is restricted to 4"
217     return out
218    
219 gross 157 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 jgs 154
436 gross 157
437     def mkText(case,name,a,a1=None,use_tagging_for_expanded_data=False):
438 jgs 154 t_out=""
439 gross 157 if case=="float":
440 jgs 154 if isinstance(a,float):
441     t_out+=" %s=%s\n"%(name,a)
442 gross 291 elif a.rank==0:
443 jgs 154 t_out+=" %s=%s\n"%(name,a)
444     else:
445     t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist())
446 gross 157 elif case=="array":
447     if isinstance(a,float):
448     t_out+=" %s=numarray.array(%s)\n"%(name,a)
449 gross 291 elif a.rank==0:
450 gross 157 t_out+=" %s=numarray.array(%s)\n"%(name,a)
451     else:
452     t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist())
453 jgs 154 elif case=="constData":
454     if isinstance(a,float):
455     t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a)
456 gross 291 elif a.rank==0:
457 jgs 154 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 gross 291 elif a.rank==0:
465 jgs 154 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 gross 157 t_out+=" %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist())
470 jgs 154 elif case=="expandedData":
471 gross 157 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 gross 291 elif a.rank==0:
476 gross 157 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 gross 291 elif a.rank==0:
487 gross 157 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 jgs 154 elif case=="Symbol":
491     if isinstance(a,float):
492     t_out+=" %s=Symbol(shape=())\n"%(name)
493 gross 291 elif a.rank==0:
494 jgs 154 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 gross 157 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 jgs 154
515 gross 157 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 jgs 154
528 gross 291 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 gross 313
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 gross 396
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 gross 439 def unrollLoops(a,b,o,arg,tap="",x="x"):
700 gross 437 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 gross 439 z+="(%s)*%s[%s]"%(a[i99]+b[i99],x,i99)
707 gross 437 else:
708 gross 439 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i99],x,i99,b[i99],x,i99)
709 gross 437
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 gross 439 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i99],x,i99,b[i0,i99],x,i99)
721 gross 437
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 gross 439 z+="(%s)*%s[%s]"%(a[i0,i1,i99]+b[i0,i1,i99],x,i99)
731 gross 437 else:
732 gross 439 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i99],x,i99,b[i0,i1,i99],x,i99)
733 gross 437
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 gross 439 z+="(%s)*%s[%s]"%(a[i0,i1,i2,i99]+b[i0,i1,i2,i99],x,i99)
744 gross 437 else:
745 gross 439 z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i99],x,i99,b[i0,i1,i2,i99],x,i99)
746 gross 437
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 gross 439 z+="(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99]+b[i0,i1,i2,i3,i99],x,i99)
758 gross 437 else:
759 gross 439 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 gross 437
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 gross 438 return out
800 gross 437 def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):
801     if where=="Function":
802     xfac_o=1.
803     xfac_op=0.
804     z_fac=1./2.
805     z_fac_s=""
806     zo_fac_s="/(o+1.)"
807     elif where=="FunctionOnBoundary":
808     xfac_o=1.
809     xfac_op=0.
810     z_fac=1.
811     z_fac_s="*dim"
812     zo_fac_s="*(1+2.*(dim-1.)/(o+1.))"
813     elif where in ["FunctionOnContactZero","FunctionOnContactOne"]:
814     xfac_o=0.
815     xfac_op=1.
816     z_fac=1./2.
817     z_fac_s=""
818     zo_fac_s="/(o+1.)"
819     out=""
820     if a.rank==1:
821     zo=0.
822     zop=0.
823     z=0.
824     for i99 in range(a.shape[0]):
825     if i99==0:
826     zo+= xfac_o*a[i99]
827     zop+= xfac_op*a[i99]
828     else:
829     zo+=a[i99]
830     z+=b[i99]
831    
832     out+=tap+"%s=(%s)%s+(%s)%s"%(arg,zo,zo_fac_s,z*z_fac,z_fac_s)
833     if zop==0.:
834     out+="\n"
835     else:
836     out+="+(%s)*0.5**o\n"%zop
837     elif a.rank==2:
838     for i0 in range(a.shape[0]):
839     zo=0.
840     zop=0.
841     z=0.
842     for i99 in range(a.shape[1]):
843     if i99==0:
844     zo+= xfac_o*a[i0,i99]
845     zop+= xfac_op*a[i0,i99]
846     else:
847     zo+=a[i0,i99]
848     z+=b[i0,i99]
849    
850     out+=tap+"%s[%s]=(%s)%s+(%s)%s"%(arg,i0,zo,zo_fac_s,z*z_fac,z_fac_s)
851     if zop==0.:
852     out+="\n"
853     else:
854     out+="+(%s)*0.5**o\n"%zop
855     elif a.rank==3:
856     for i0 in range(a.shape[0]):
857     for i1 in range(a.shape[1]):
858     zo=0.
859     zop=0.
860     z=0.
861     for i99 in range(a.shape[2]):
862     if i99==0:
863     zo+= xfac_o*a[i0,i1,i99]
864     zop+= xfac_op*a[i0,i1,i99]
865     else:
866     zo+=a[i0,i1,i99]
867     z+=b[i0,i1,i99]
868    
869     out+=tap+"%s[%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,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==4:
875     for i0 in range(a.shape[0]):
876     for i1 in range(a.shape[1]):
877     for i2 in range(a.shape[2]):
878     zo=0.
879     zop=0.
880     z=0.
881     for i99 in range(a.shape[3]):
882     if i99==0:
883     zo+= xfac_o*a[i0,i1,i2,i99]
884     zop+= xfac_op*a[i0,i1,i2,i99]
885    
886     else:
887     zo+=a[i0,i1,i2,i99]
888     z+=b[i0,i1,i2,i99]
889    
890     out+=tap+"%s[%s,%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,i2,zo,zo_fac_s,z*z_fac,z_fac_s)
891     if zop==0.:
892     out+="\n"
893     else:
894     out+="+(%s)*0.5**o\n"%zop
895    
896     elif a.rank==5:
897     for i0 in range(a.shape[0]):
898     for i1 in range(a.shape[1]):
899     for i2 in range(a.shape[2]):
900     for i3 in range(a.shape[3]):
901     zo=0.
902     zop=0.
903     z=0.
904     for i99 in range(a.shape[4]):
905     if i99==0:
906     zo+= xfac_o*a[i0,i1,i2,i3,i99]
907     zop+= xfac_op*a[i0,i1,i2,i3,i99]
908    
909     else:
910     zo+=a[i0,i1,i2,i3,i99]
911     z+=b[i0,i1,i2,i3,i99]
912     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)
913     if zop==0.:
914     out+="\n"
915     else:
916     out+="+(%s)*0.5**o\n"%zop
917    
918     return out
919    
920 gross 439 #=======================================================================================================
921     # interpolation
922     #=======================================================================================================
923     for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne","Solution","ReducedSolution"]:
924     for data in ["Data","Symbol"]:
925     for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
926     for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
927     if where==case or \
928     ( case in ["ContinuousFunction","Solution","ReducedSolution"] and where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"] ) or \
929     ( case in ["FunctionOnContactZero","FunctionOnContactOne"] and where in ["FunctionOnContactZero","FunctionOnContactOne"] ) or \
930     (case=="ContinuousFunction" and where in ["Solution","ReducedSolution"]) or \
931     (case=="Solution" and where=="ReducedSolution") :
932    
933    
934     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
935     tname="test_interpolation_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
936     text+=" def %s(self):\n"%tname
937     text+=" \"\"\"\n"
938     text+=" tests interpolation for rank %s %s onto the %s\n\n"%(len(sh),data,where)
939     text+=" assumptions: self.domain supports inpterpolation from %s onto %s\n"%(case,where)
940     text+=" \"\"\"\n"
941     if case=="ReducedSolution" or where=="ReducedSolution":
942     text+=" o=1\n"
943     o="1"
944     else:
945     text+=" o=self.order\n"
946     o="o"
947     text+=" dim=self.domain.getDim()\n"
948     text+=" w_ref=%s(self.domain)\n"%where
949     text+=" x_ref=w_ref.getX()\n"
950     text+=" w=%s(self.domain)\n"%case
951     text+=" x=w.getX()\n"
952     a_2=makeArray(sh+(2,),[-1.,1])
953     b_2=makeArray(sh+(2,),[-1.,1])
954     a_3=makeArray(sh+(3,),[-1.,1])
955     b_3=makeArray(sh+(3,),[-1.,1])
956     if data=="Symbol":
957     text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
958     val="s"
959     res="sub"
960     else:
961     val="arg"
962     res="res"
963     text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
964     text+=" ref=Data(0,%s,w_ref)\n"%str(sh)
965     text+=" if dim==2:\n"
966     text+=unrollLoops(a_2,b_2,o,val,tap=" ")
967     text+=unrollLoops(a_2,b_2,o,"ref",tap=" ",x="x_ref")
968     text+=" else:\n"
969    
970     text+=unrollLoops(a_3,b_3,o,val,tap=" ")
971     text+=unrollLoops(a_3,b_3,o,"ref",tap=" ",x="x_ref")
972     text+=" res=interpolate(arg,where=w_ref)\n"
973     if data=="Symbol":
974     text+=" sub=res.substitute({arg:s})\n"
975     text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
976     text+=" self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong functionspace of result.\")\n"%res
977     text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
978     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
979     t_prog+=text
980     print test_header
981     print t_prog
982     print test_tail
983     1/0
984 gross 437
985 gross 429 #=======================================================================================================
986 gross 438 # grad
987 gross 437 #=======================================================================================================
988     for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
989     for data in ["Data","Symbol"]:
990 gross 438 for case in ["ContinuousFunction","Solution","ReducedSolution"]:
991     for sh in [ (),(2,), (4,5), (6,2,2)]:
992 gross 437 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
993 gross 438 tname="test_grad_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
994 gross 437 text+=" def %s(self):\n"%tname
995     text+=" \"\"\"\n"
996 gross 438 text+=" tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where)
997 gross 437 text+=" assumptions: %s(self.domain) exists\n"%case
998 gross 438 text+=" self.domain supports gardient on %s\n"%where
999 gross 437 text+=" \"\"\"\n"
1000     if case=="ReducedSolution":
1001     text+=" o=1\n"
1002     o="1"
1003     else:
1004     text+=" o=self.order\n"
1005     o="o"
1006     text+=" dim=self.domain.getDim()\n"
1007     text+=" w_ref=%s(self.domain)\n"%where
1008 gross 438 text+=" x_ref=w_ref.getX()\n"
1009 gross 437 text+=" w=%s(self.domain)\n"%case
1010     text+=" x=w.getX()\n"
1011     a_2=makeArray(sh+(2,),[-1.,1])
1012     b_2=makeArray(sh+(2,),[-1.,1])
1013     a_3=makeArray(sh+(3,),[-1.,1])
1014     b_3=makeArray(sh+(3,),[-1.,1])
1015     if data=="Symbol":
1016 gross 438 text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
1017 gross 437 val="s"
1018     res="sub"
1019     else:
1020     val="arg"
1021     res="res"
1022     text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
1023 gross 438 text+=" ref=Data(0,%s+(dim,),w_ref)\n"%str(sh)
1024 gross 437 text+=" if dim==2:\n"
1025     text+=unrollLoops(a_2,b_2,o,val,tap=" ")
1026 gross 438 text+=unrollLoopsOfGrad(a_2,b_2,o,"ref",tap=" ")
1027 gross 437 text+=" else:\n"
1028    
1029     text+=unrollLoops(a_3,b_3,o,val,tap=" ")
1030 gross 438 text+=unrollLoopsOfGrad(a_3,b_3,o,"ref",tap=" ")
1031     text+=" res=grad(arg,where=w_ref)\n"
1032 gross 437 if data=="Symbol":
1033     text+=" sub=res.substitute({arg:s})\n"
1034 gross 438 text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
1035     text+=" self.failUnlessEqual(res.getShape(),%s+(dim,),\"wrong shape of result.\")\n"%str(sh)
1036 gross 437 text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1037    
1038    
1039     t_prog+=text
1040     print test_header
1041     print t_prog
1042     print test_tail
1043     1/0
1044 gross 438
1045    
1046 gross 437 #=======================================================================================================
1047 gross 438 # integrate
1048 gross 437 #=======================================================================================================
1049 gross 438 for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1050 gross 437 for data in ["Data","Symbol"]:
1051 gross 438 for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1052     for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
1053     if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:
1054 gross 437 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1055 gross 438 tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1056 gross 437 text+=" def %s(self):\n"%tname
1057     text+=" \"\"\"\n"
1058 gross 438 text+=" tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where)
1059     text+=" assumptions: %s(self.domain) exists\n"%case
1060     text+=" self.domain supports integral on %s\n"%where
1061    
1062 gross 437 text+=" \"\"\"\n"
1063     if case=="ReducedSolution":
1064     text+=" o=1\n"
1065     o="1"
1066     else:
1067     text+=" o=self.order\n"
1068     o="o"
1069     text+=" dim=self.domain.getDim()\n"
1070     text+=" w_ref=%s(self.domain)\n"%where
1071     text+=" w=%s(self.domain)\n"%case
1072     text+=" x=w.getX()\n"
1073     a_2=makeArray(sh+(2,),[-1.,1])
1074     b_2=makeArray(sh+(2,),[-1.,1])
1075     a_3=makeArray(sh+(3,),[-1.,1])
1076     b_3=makeArray(sh+(3,),[-1.,1])
1077     if data=="Symbol":
1078 gross 438 text+=" arg=Symbol(shape=%s)\n"%str(sh)
1079 gross 437 val="s"
1080     res="sub"
1081     else:
1082     val="arg"
1083     res="res"
1084 gross 438
1085 gross 437 text+=" %s=Data(0,%s,w)\n"%(val,str(sh))
1086 gross 438 if not len(sh)==0:
1087     text+=" ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)
1088 gross 437 text+=" if dim==2:\n"
1089     text+=unrollLoops(a_2,b_2,o,val,tap=" ")
1090 gross 438 text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap=" ")
1091 gross 437 text+=" else:\n"
1092    
1093     text+=unrollLoops(a_3,b_3,o,val,tap=" ")
1094 gross 438 text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap=" ")
1095     if case in ["ContinuousFunction","Solution","ReducedSolution"]:
1096     text+=" res=integrate(arg,where=w_ref)\n"
1097     else:
1098     text+=" res=integrate(arg)\n"
1099    
1100 gross 437 if data=="Symbol":
1101     text+=" sub=res.substitute({arg:s})\n"
1102 gross 438 if len(sh)==0 and data=="Data":
1103     text+=" self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res
1104     else:
1105     if data=="Symbol":
1106     text+=" self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"
1107     text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
1108     else:
1109     text+=" self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n"
1110     text+=" self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)
1111 gross 437 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 gross 433 # inverse
1121     #=======================================================================================================
1122     name="inverse"
1123     for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
1124     for sh0 in [ (1,1), (2,2), (3,3)]:
1125     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1126     tname="test_%s_%s_dim%s"%(name,case0,sh0[0])
1127     text+=" def %s(self):\n"%tname
1128     a_0=makeArray(sh0,[-1.,1])
1129     for i in range(sh0[0]): a_0[i,i]+=2
1130     if case0 in ["taggedData", "expandedData"]:
1131     a1_0=makeArray(sh0,[-1.,1])
1132     for i in range(sh0[0]): a1_0[i,i]+=3
1133     else:
1134     a1_0=a_0
1135    
1136     text+=mkText(case0,"arg",a_0,a1_0)
1137     text+=" res=%s(arg)\n"%name
1138     if case0=="Symbol":
1139     text+=mkText("array","s",a_0,a1_0)
1140     text+=" sub=res.substitute({arg:s})\n"
1141     res="sub"
1142     ref="s"
1143     else:
1144     ref="arg"
1145     res="res"
1146     text+=mkTypeAndShapeTest(case0,sh0,"res")
1147     text+=" self.failUnless(Lsup(matrixmult(%s,%s)-kronecker(%s))<=self.RES_TOL,\"wrong result\")\n"%(res,ref,sh0[0])
1148    
1149     if case0 == "taggedData" :
1150     t_prog_with_tags+=text
1151     else:
1152     t_prog+=text
1153    
1154     print test_header
1155     # print t_prog
1156     print t_prog_with_tags
1157     print test_tail
1158     1/0
1159    
1160     #=======================================================================================================
1161 gross 429 # trace
1162     #=======================================================================================================
1163     def traceTest(r,offset):
1164     sh=r.shape
1165     r1=1
1166     for i in range(offset): r1*=sh[i]
1167     r2=1
1168     for i in range(offset+2,len(sh)): r2*=sh[i]
1169     r_s=numarray.reshape(r,(r1,sh[offset],sh[offset],r2))
1170     s=numarray.zeros([r1,r2],numarray.Float)
1171     for i1 in range(r1):
1172     for i2 in range(r2):
1173     for j in range(sh[offset]): s[i1,i2]+=r_s[i1,j,j,i2]
1174     return s.resize(sh[:offset]+sh[offset+2:])
1175     name,tt="trace",traceTest
1176     for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
1177     for sh0 in [ (4,5), (6,2,2),(3,2,3,4)]:
1178     for offset in range(len(sh0)-1):
1179     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1180     tname="test_%s_%s_rank%s_offset%s"%(name,case0,len(sh0),offset)
1181     text+=" def %s(self):\n"%tname
1182     sh_t=list(sh0)
1183     sh_t[offset+1]=sh_t[offset]
1184     sh_t=tuple(sh_t)
1185     sh_r=[]
1186     for i in range(offset): sh_r.append(sh0[i])
1187     for i in range(offset+2,len(sh0)): sh_r.append(sh0[i])
1188     sh_r=tuple(sh_r)
1189     a_0=makeArray(sh_t,[-1.,1])
1190     if case0 in ["taggedData", "expandedData"]:
1191     a1_0=makeArray(sh_t,[-1.,1])
1192     else:
1193     a1_0=a_0
1194     r=tt(a_0,offset)
1195     r1=tt(a1_0,offset)
1196     text+=mkText(case0,"arg",a_0,a1_0)
1197     text+=" res=%s(arg,%s)\n"%(name,offset)
1198     if case0=="Symbol":
1199     text+=mkText("array","s",a_0,a1_0)
1200     text+=" sub=res.substitute({arg:s})\n"
1201     res="sub"
1202     text+=mkText("array","ref",r,r1)
1203     else:
1204     res="res"
1205     text+=mkText(case0,"ref",r,r1)
1206     text+=mkTypeAndShapeTest(case0,sh_r,"res")
1207     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1208    
1209     if case0 == "taggedData" :
1210     t_prog_with_tags+=text
1211     else:
1212     t_prog+=text
1213 gross 396
1214 gross 429 print test_header
1215     # print t_prog
1216     print t_prog_with_tags
1217     print test_tail
1218     1/0
1219 gross 396
1220 gross 157 #=======================================================================================================
1221 gross 396 # clip
1222     #=======================================================================================================
1223     oper_L=[["clip",clipTEST]]
1224     for oper in oper_L:
1225     for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1226     for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1227     if len(sh0)==0 or not case0=="float":
1228     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1229     tname="test_%s_%s_rank%s"%(oper[0],case0,len(sh0))
1230     text+=" def %s(self):\n"%tname
1231     a_0=makeArray(sh0,[-1.,1])
1232     if case0 in ["taggedData", "expandedData"]:
1233     a1_0=makeArray(sh0,[-1.,1])
1234     else:
1235     a1_0=a_0
1236    
1237     r=oper[1](a_0,-0.3,0.5)
1238     r1=oper[1](a1_0,-0.3,0.5)
1239     text+=mkText(case0,"arg",a_0,a1_0)
1240     text+=" res=%s(arg,-0.3,0.5)\n"%oper[0]
1241     if case0=="Symbol":
1242     text+=mkText("array","s",a_0,a1_0)
1243     text+=" sub=res.substitute({arg:s})\n"
1244     res="sub"
1245     text+=mkText("array","ref",r,r1)
1246     else:
1247     res="res"
1248     text+=mkText(case0,"ref",r,r1)
1249     text+=mkTypeAndShapeTest(case0,sh0,"res")
1250     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1251    
1252     if case0 == "taggedData" :
1253     t_prog_with_tags+=text
1254     else:
1255     t_prog+=text
1256    
1257     print test_header
1258     # print t_prog
1259     print t_prog_with_tags
1260     print test_tail
1261     1/0
1262 gross 429
1263 gross 396 #=======================================================================================================
1264     # maximum, minimum, clipping
1265     #=======================================================================================================
1266     oper_L=[ ["maximum",maximumTEST],
1267     ["minimum",minimumTEST]]
1268     for oper in oper_L:
1269     for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1270     for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1271     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1272     for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1273     if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \
1274     and (sh0==sh1 or len(sh0)==0 or len(sh1)==0) :
1275     use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1276    
1277     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1278     tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1279     text+=" def %s(self):\n"%tname
1280     a_0=makeArray(sh0,[-1.,1])
1281     if case0 in ["taggedData", "expandedData"]:
1282     a1_0=makeArray(sh0,[-1.,1])
1283     else:
1284     a1_0=a_0
1285    
1286     a_1=makeArray(sh1,[-1.,1])
1287     if case1 in ["taggedData", "expandedData"]:
1288     a1_1=makeArray(sh1,[-1.,1])
1289     else:
1290     a1_1=a_1
1291     r=oper[1](a_0,a_1)
1292     r1=oper[1](a1_0,a1_1)
1293     text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1294     text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1295     text+=" res=%s(arg0,arg1)\n"%oper[0]
1296     case=getResultCaseForBin(case0,case1)
1297     if case=="Symbol":
1298     c0_res,c1_res=case0,case1
1299     subs="{"
1300     if case0=="Symbol":
1301     text+=mkText("array","s0",a_0,a1_0)
1302     subs+="arg0:s0"
1303     c0_res="array"
1304     if case1=="Symbol":
1305     text+=mkText("array","s1",a_1,a1_1)
1306     if not subs.endswith("{"): subs+=","
1307     subs+="arg1:s1"
1308     c1_res="array"
1309     subs+="}"
1310     text+=" sub=res.substitute(%s)\n"%subs
1311     res="sub"
1312     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1313     else:
1314     res="res"
1315     text+=mkText(case,"ref",r,r1)
1316     if len(sh0)>len(sh1):
1317     text+=mkTypeAndShapeTest(case,sh0,"res")
1318     else:
1319     text+=mkTypeAndShapeTest(case,sh1,"res")
1320     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1321    
1322     if case0 == "taggedData" or case1 == "taggedData":
1323     t_prog_with_tags+=text
1324     else:
1325     t_prog+=text
1326    
1327     print test_header
1328     # print t_prog
1329     print t_prog_with_tags
1330     print test_tail
1331     1/0
1332    
1333    
1334     #=======================================================================================================
1335     # outer inner
1336     #=======================================================================================================
1337     oper=["outer",outerTEST]
1338     # oper=["inner",innerTEST]
1339     for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1340     for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1341     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1342     for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1343     if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \
1344     and len(sh0+sh1)<5:
1345     use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1346    
1347     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1348     tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1349     text+=" def %s(self):\n"%tname
1350     a_0=makeArray(sh0,[-1.,1])
1351     if case0 in ["taggedData", "expandedData"]:
1352     a1_0=makeArray(sh0,[-1.,1])
1353     else:
1354     a1_0=a_0
1355    
1356     a_1=makeArray(sh1,[-1.,1])
1357     if case1 in ["taggedData", "expandedData"]:
1358     a1_1=makeArray(sh1,[-1.,1])
1359     else:
1360     a1_1=a_1
1361     r=oper[1](a_0,a_1)
1362     r1=oper[1](a1_0,a1_1)
1363     text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1364     text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1365     text+=" res=%s(arg0,arg1)\n"%oper[0]
1366     case=getResultCaseForBin(case0,case1)
1367     if case=="Symbol":
1368     c0_res,c1_res=case0,case1
1369     subs="{"
1370     if case0=="Symbol":
1371     text+=mkText("array","s0",a_0,a1_0)
1372     subs+="arg0:s0"
1373     c0_res="array"
1374     if case1=="Symbol":
1375     text+=mkText("array","s1",a_1,a1_1)
1376     if not subs.endswith("{"): subs+=","
1377     subs+="arg1:s1"
1378     c1_res="array"
1379     subs+="}"
1380     text+=" sub=res.substitute(%s)\n"%subs
1381     res="sub"
1382     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1383     else:
1384     res="res"
1385     text+=mkText(case,"ref",r,r1)
1386     text+=mkTypeAndShapeTest(case,sh0+sh1,"res")
1387     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1388    
1389     if case0 == "taggedData" or case1 == "taggedData":
1390     t_prog_with_tags+=text
1391     else:
1392     t_prog+=text
1393    
1394     print test_header
1395     # print t_prog
1396     print t_prog_with_tags
1397     print test_tail
1398     1/0
1399    
1400     #=======================================================================================================
1401 gross 313 # local reduction
1402     #=======================================================================================================
1403     for oper in [["length",0.,"out+%a1%**2","math.sqrt(out)"],
1404     ["maxval",-1.e99,"max(out,%a1%)","out"],
1405     ["minval",1.e99,"min(out,%a1%)","out"] ]:
1406     for case in case_set:
1407     for sh in shape_set:
1408     if not case=="float" or len(sh)==0:
1409     text=""
1410     text+=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1411     tname="def test_%s_%s_rank%s"%(oper[0],case,len(sh))
1412     text+=" %s(self):\n"%tname
1413     a=makeArray(sh,[-1.,1.])
1414     a1=makeArray(sh,[-1.,1.])
1415     r1=testReduce(a1,oper[1],oper[2],oper[3])
1416     r=testReduce(a,oper[1],oper[2],oper[3])
1417    
1418     text+=mkText(case,"arg",a,a1)
1419     text+=" res=%s(arg)\n"%oper[0]
1420     if case=="Symbol":
1421     text+=mkText("array","s",a,a1)
1422     text+=" sub=res.substitute({arg:s})\n"
1423     text+=mkText("array","ref",r,r1)
1424     res="sub"
1425     else:
1426     text+=mkText(case,"ref",r,r1)
1427     res="res"
1428     if oper[0]=="length":
1429     text+=mkTypeAndShapeTest(case,(),"res")
1430     else:
1431     if case=="float" or case=="array":
1432     text+=mkTypeAndShapeTest("float",(),"res")
1433     else:
1434     text+=mkTypeAndShapeTest(case,(),"res")
1435     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1436     if case == "taggedData":
1437     t_prog_with_tags+=text
1438     else:
1439     t_prog+=text
1440     print test_header
1441     # print t_prog
1442     print t_prog_with_tags
1443     print test_tail
1444     1/0
1445    
1446     #=======================================================================================================
1447 gross 291 # tensor multiply
1448 gross 157 #=======================================================================================================
1449 gross 291 # oper=["generalTensorProduct",tensorProductTest]
1450     # oper=["matrixmult",testMatrixMult]
1451     oper=["tensormult",testTensorMult]
1452    
1453     for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1454     for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1455     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
1456     for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
1457     for sh_s in [ (),(3,), (2,3), (2,4,3),(4,2,3,2)]:
1458     if (len(sh0+sh_s)==0 or not case0=="float") and (len(sh1+sh_s)==0 or not case1=="float") \
1459     and len(sh0+sh1)<5 and len(sh0+sh_s)<5 and len(sh1+sh_s)<5:
1460     # 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
1461     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
1462     case=getResultCaseForBin(case0,case1)
1463 gross 157 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1464     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1465 gross 291 # tname="test_generalTensorProduct_%s_rank%s_%s_rank%s_offset%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1),len(sh_s))
1466     #tname="test_matrixmult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
1467     tname="test_tensormult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
1468     # if tname=="test_generalTensorProduct_array_rank1_array_rank2_offset1":
1469     # print tnametest_generalTensorProduct_Symbol_rank1_Symbol_rank3_offset1
1470     text+=" def %s(self):\n"%tname
1471     a_0=makeArray(sh0+sh_s,[-1.,1])
1472     if case0 in ["taggedData", "expandedData"]:
1473     a1_0=makeArray(sh0+sh_s,[-1.,1])
1474     else:
1475     a1_0=a_0
1476    
1477     a_1=makeArray(sh_s+sh1,[-1.,1])
1478     if case1 in ["taggedData", "expandedData"]:
1479     a1_1=makeArray(sh_s+sh1,[-1.,1])
1480     else:
1481     a1_1=a_1
1482     r=oper[1](a_0,a_1,sh_s)
1483     r1=oper[1](a1_0,a1_1,sh_s)
1484     text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1485     text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1486     #text+=" res=matrixmult(arg0,arg1)\n"
1487     text+=" res=tensormult(arg0,arg1)\n"
1488     #text+=" res=generalTensorProduct(arg0,arg1,offset=%s)\n"%(len(sh_s))
1489     if case=="Symbol":
1490     c0_res,c1_res=case0,case1
1491     subs="{"
1492     if case0=="Symbol":
1493     text+=mkText("array","s0",a_0,a1_0)
1494     subs+="arg0:s0"
1495     c0_res="array"
1496     if case1=="Symbol":
1497     text+=mkText("array","s1",a_1,a1_1)
1498     if not subs.endswith("{"): subs+=","
1499     subs+="arg1:s1"
1500     c1_res="array"
1501     subs+="}"
1502     text+=" sub=res.substitute(%s)\n"%subs
1503     res="sub"
1504     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1505     else:
1506     res="res"
1507     text+=mkText(case,"ref",r,r1)
1508     text+=mkTypeAndShapeTest(case,sh0+sh1,"res")
1509     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1510     if case0 == "taggedData" or case1 == "taggedData":
1511     t_prog_with_tags+=text
1512     else:
1513     t_prog+=text
1514     print test_header
1515     # print t_prog
1516     print t_prog_with_tags
1517     print test_tail
1518     1/0
1519     #=======================================================================================================
1520 gross 157 # basic binary operation overloading (tests only!)
1521     #=======================================================================================================
1522     oper_range=[-5.,5.]
1523     for oper in [["add" ,"+",[-5.,5.]],
1524     ["sub" ,"-",[-5.,5.]],
1525     ["mult","*",[-5.,5.]],
1526     ["div" ,"/",[-5.,5.]],
1527     ["pow" ,"**",[0.01,5.]]]:
1528     for case0 in case_set:
1529     for sh0 in shape_set:
1530     for case1 in case_set:
1531     for sh1 in shape_set:
1532 gross 291 if not case0=="array" and \
1533     (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1534 gross 157 (sh0==() or sh1==() or sh1==sh0) and \
1535     not (case0 in ["float","array"] and case1 in ["float","array"]):
1536 gross 291 use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1537 gross 157 text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1538     tname="test_%s_overloaded_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1539     text+=" def %s(self):\n"%tname
1540     a_0=makeArray(sh0,oper[2])
1541     if case0 in ["taggedData", "expandedData"]:
1542     a1_0=makeArray(sh0,oper[2])
1543     else:
1544     a1_0=a_0
1545 jgs 154
1546 gross 157 a_1=makeArray(sh1,oper[2])
1547     if case1 in ["taggedData", "expandedData"]:
1548     a1_1=makeArray(sh1,oper[2])
1549 jgs 154 else:
1550 gross 157 a1_1=a_1
1551     r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1552     r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1553     text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1554     text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1555     text+=" res=arg0%sarg1\n"%oper[1]
1556    
1557     case=getResultCaseForBin(case0,case1)
1558     if case=="Symbol":
1559     c0_res,c1_res=case0,case1
1560     subs="{"
1561     if case0=="Symbol":
1562     text+=mkText("array","s0",a_0,a1_0)
1563     subs+="arg0:s0"
1564     c0_res="array"
1565     if case1=="Symbol":
1566     text+=mkText("array","s1",a_1,a1_1)
1567     if not subs.endswith("{"): subs+=","
1568     subs+="arg1:s1"
1569     c1_res="array"
1570     subs+="}"
1571     text+=" sub=res.substitute(%s)\n"%subs
1572     res="sub"
1573     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1574 jgs 154 else:
1575 gross 157 res="res"
1576     text+=mkText(case,"ref",r,r1)
1577     if isinstance(r,float):
1578     text+=mkTypeAndShapeTest(case,(),"res")
1579     else:
1580     text+=mkTypeAndShapeTest(case,r.shape,"res")
1581     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1582    
1583     if case0 in [ "constData","taggedData","expandedData"] and case1 == "Symbol":
1584     t_prog_failing+=text
1585     else:
1586     if case0 == "taggedData" or case1 == "taggedData":
1587     t_prog_with_tags+=text
1588     else:
1589     t_prog+=text
1590 jgs 154
1591 gross 157
1592     print test_header
1593 gross 291 # print t_prog
1594     # print t_prog_with_tags
1595     print t_prog_failing
1596     print test_tail
1597 jgs 154 1/0
1598 gross 291 #=======================================================================================================
1599     # basic binary operations (tests only!)
1600     #=======================================================================================================
1601     oper_range=[-5.,5.]
1602     for oper in [["add" ,"+",[-5.,5.]],
1603     ["mult","*",[-5.,5.]],
1604     ["quotient" ,"/",[-5.,5.]],
1605     ["power" ,"**",[0.01,5.]]]:
1606     for case0 in case_set:
1607     for case1 in case_set:
1608     for sh in shape_set:
1609     for sh_p in shape_set:
1610     if len(sh_p)>0:
1611     resource=[-1,1]
1612     else:
1613     resource=[1]
1614     for sh_d in resource:
1615     if sh_d>0:
1616     sh0=sh
1617     sh1=sh+sh_p
1618     else:
1619     sh1=sh
1620     sh0=sh+sh_p
1621    
1622     if (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1623     len(sh0)<5 and len(sh1)<5:
1624     use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
1625     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1626     tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1627     text+=" def %s(self):\n"%tname
1628     a_0=makeArray(sh0,oper[2])
1629     if case0 in ["taggedData", "expandedData"]:
1630     a1_0=makeArray(sh0,oper[2])
1631     else:
1632     a1_0=a_0
1633    
1634     a_1=makeArray(sh1,oper[2])
1635     if case1 in ["taggedData", "expandedData"]:
1636     a1_1=makeArray(sh1,oper[2])
1637     else:
1638     a1_1=a_1
1639     r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1640     r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1641     text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
1642     text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
1643     text+=" res=%s(arg0,arg1)\n"%oper[0]
1644    
1645     case=getResultCaseForBin(case0,case1)
1646     if case=="Symbol":
1647     c0_res,c1_res=case0,case1
1648     subs="{"
1649     if case0=="Symbol":
1650     text+=mkText("array","s0",a_0,a1_0)
1651     subs+="arg0:s0"
1652     c0_res="array"
1653     if case1=="Symbol":
1654     text+=mkText("array","s1",a_1,a1_1)
1655     if not subs.endswith("{"): subs+=","
1656     subs+="arg1:s1"
1657     c1_res="array"
1658     subs+="}"
1659     text+=" sub=res.substitute(%s)\n"%subs
1660     res="sub"
1661     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1662     else:
1663     res="res"
1664     text+=mkText(case,"ref",r,r1)
1665     if isinstance(r,float):
1666     text+=mkTypeAndShapeTest(case,(),"res")
1667     else:
1668     text+=mkTypeAndShapeTest(case,r.shape,"res")
1669     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1670    
1671     if case0 == "taggedData" or case1 == "taggedData":
1672     t_prog_with_tags+=text
1673     else:
1674     t_prog+=text
1675     print test_header
1676     # print t_prog
1677     print t_prog_with_tags
1678     print test_tail
1679     1/0
1680    
1681 gross 157 # print t_prog_with_tagsoper_range=[-5.,5.]
1682     for oper in [["add" ,"+",[-5.,5.]],
1683     ["sub" ,"-",[-5.,5.]],
1684     ["mult","*",[-5.,5.]],
1685     ["div" ,"/",[-5.,5.]],
1686     ["pow" ,"**",[0.01,5.]]]:
1687     for case0 in case_set:
1688     for sh0 in shape_set:
1689     for case1 in case_set:
1690     for sh1 in shape_set:
1691     if (not case0=="float" or len(sh0)==0) and (not case1=="float" or len(sh1)==0) and \
1692     (sh0==() or sh1==() or sh1==sh0) and \
1693     not (case0 in ["float","array"] and case1 in ["float","array"]):
1694     text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1695     tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
1696     text+=" def %s(self):\n"%tname
1697     a_0=makeArray(sh0,oper[2])
1698     if case0 in ["taggedData", "expandedData"]:
1699     a1_0=makeArray(sh0,oper[2])
1700     else:
1701     a1_0=a_0
1702 jgs 154
1703 gross 157 a_1=makeArray(sh1,oper[2])
1704     if case1 in ["taggedData", "expandedData"]:
1705     a1_1=makeArray(sh1,oper[2])
1706 jgs 154 else:
1707 gross 157 a1_1=a_1
1708     r1=makeResult2(a1_0,a1_1,"%a1%"+oper[1]+"%a2%")
1709     r=makeResult2(a_0,a_1,"%a1%"+oper[1]+"%a2%")
1710     text+=mkText(case0,"arg0",a_0,a1_0)
1711     text+=mkText(case1,"arg1",a_1,a1_1)
1712     text+=" res=arg0%sarg1\n"%oper[1]
1713    
1714     case=getResultCaseForBin(case0,case1)
1715     if case=="Symbol":
1716     c0_res,c1_res=case0,case1
1717     subs="{"
1718     if case0=="Symbol":
1719     text+=mkText("array","s0",a_0,a1_0)
1720     subs+="arg0:s0"
1721     c0_res="array"
1722     if case1=="Symbol":
1723     text+=mkText("array","s1",a_1,a1_1)
1724     if not subs.endswith("{"): subs+=","
1725     subs+="arg1:s1"
1726     c1_res="array"
1727     subs+="}"
1728     text+=" sub=res.substitute(%s)\n"%subs
1729     res="sub"
1730     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
1731 jgs 154 else:
1732 gross 157 res="res"
1733     text+=mkText(case,"ref",r,r1)
1734     if isinstance(r,float):
1735     text+=mkTypeAndShapeTest(case,(),"res")
1736     else:
1737     text+=mkTypeAndShapeTest(case,r.shape,"res")
1738     text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1739    
1740     if case0 in [ "constData","taggedData","expandedData"] and case1 == "Symbol":
1741     t_prog_failing+=text
1742     else:
1743     if case0 == "taggedData" or case1 == "taggedData":
1744     t_prog_with_tags+=text
1745     else:
1746     t_prog+=text
1747 jgs 154
1748 gross 157
1749     # print u_prog
1750     # 1/0
1751     print test_header
1752 jgs 154 print t_prog
1753 gross 157 # print t_prog_with_tags
1754     # print t_prog_failing
1755     print test_tail
1756     # print t_prog_failing
1757     print test_tail
1758 jgs 154
1759 gross 157 #=======================================================================================================
1760     # unary operations:
1761     #=======================================================================================================
1762     func= [
1763     OPERATOR(nickname="log10",\
1764     rng=[1.e-3,100.],\
1765     test_expr="math.log10(%a1%)",\
1766     math_expr="math.log10(%a1%)",\
1767     numarray_expr="numarray.log10(%a1%)",\
1768     symbol_expr="log(%a1%)/log(10.)",\
1769     name="base-10 logarithm"),
1770     OPERATOR(nickname="wherePositive",\
1771     rng=[-100.,100.],\
1772     test_expr="wherepos(%a1%)",\
1773     math_expr="if arg>0:\n return 1.\nelse:\n return 0.",
1774     numarray_expr="numarray.greater(arg,numarray.zeros(arg.shape,numarray.Float))",\
1775     name="mask of positive values"),
1776     OPERATOR(nickname="whereNegative",\
1777     rng=[-100.,100.],\
1778     test_expr="wherepos(-%a1%)",\
1779     math_expr="if arg<0:\n return 1.\nelse:\n return 0.",
1780     numarray_expr="numarray.less(arg,numarray.zeros(arg.shape,numarray.Float))",\
1781     name="mask of positive values"),
1782     OPERATOR(nickname="whereNonNegative",\
1783     rng=[-100.,100.],\
1784     test_expr="1-wherepos(-%a1%)", \
1785     math_expr="if arg<0:\n return 0.\nelse:\n return 1.",
1786     numarray_expr="numarray.greater_equal(arg,numarray.zeros(arg.shape,numarray.Float))",\
1787     symbol_expr="1-wherePositive(%a1%)",\
1788     name="mask of non-negative values"),
1789     OPERATOR(nickname="whereNonPositive",\
1790     rng=[-100.,100.],\
1791     test_expr="1-wherepos(%a1%)",\
1792     math_expr="if arg>0:\n return 0.\nelse:\n return 1.",
1793     numarray_expr="numarray.less_equal(arg,numarray.zeros(arg.shape,numarray.Float))",\
1794     symbol_expr="1-whereNegative(%a1%)",\
1795     name="mask of non-positive values"),
1796     OPERATOR(nickname="whereZero",\
1797     rng=[-100.,100.],\
1798     test_expr="1-wherepos(%a1%)-wherepos(-%a1%)",\
1799     math_expr="if abs(%a1%)<=tol:\n return 1.\nelse:\n return 0.",
1800     numarray_expr="numarray.less_equal(abs(%a1%)-tol,numarray.zeros(arg.shape,numarray.Float))",\
1801     name="mask of zero entries"),
1802     OPERATOR(nickname="whereNonZero",\
1803     rng=[-100.,100.],\
1804     test_expr="wherepos(%a1%)+wherepos(-%a1%)",\
1805     math_expr="if abs(%a1%)>tol:\n return 1.\nelse:\n return 0.",\
1806     numarray_expr="numarray.greater(abs(%a1%)-tol,numarray.zeros(arg.shape,numarray.Float))",\
1807     symbol_expr="1-whereZero(arg,tol)",\
1808     name="mask of values different from zero"),
1809     OPERATOR(nickname="sin",\
1810     rng=[-100.,100.],\
1811     test_expr="math.sin(%a1%)",
1812     numarray_expr="numarray.sin(%a1%)",\
1813     diff="cos(%a1%)",\
1814     name="sine"),
1815     OPERATOR(nickname="cos",\
1816     rng=[-100.,100.],\
1817     test_expr="math.cos(%a1%)",
1818     numarray_expr="numarray.cos(%a1%)",\
1819     diff="-sin(%a1%)",
1820     name="cosine"),
1821     OPERATOR(nickname="tan",\
1822     rng=[-100.,100.],\
1823     test_expr="math.tan(%a1%)",
1824     numarray_expr="numarray.tan(%a1%)",\
1825     diff="1./cos(%a1%)**2",
1826     name="tangent"),
1827     OPERATOR(nickname="asin",\
1828     rng=[-0.99,0.99],\
1829     test_expr="math.asin(%a1%)",
1830     numarray_expr="numarray.arcsin(%a1%)",
1831     diff="1./sqrt(1.-%a1%**2)",
1832     name="inverse sine"),
1833     OPERATOR(nickname="acos",\
1834     rng=[-0.99,0.99],\
1835     test_expr="math.acos(%a1%)",
1836     numarray_expr="numarray.arccos(%a1%)",
1837     diff="-1./sqrt(1.-%a1%**2)",
1838     name="inverse cosine"),
1839     OPERATOR(nickname="atan",\
1840     rng=[-100.,100.],\
1841     test_expr="math.atan(%a1%)",
1842     numarray_expr="numarray.arctan(%a1%)",
1843     diff="1./(1+%a1%**2)",
1844     name="inverse tangent"),
1845     OPERATOR(nickname="sinh",\
1846     rng=[-5,5],\
1847     test_expr="math.sinh(%a1%)",
1848     numarray_expr="numarray.sinh(%a1%)",
1849     diff="cosh(%a1%)",
1850     name="hyperbolic sine"),
1851     OPERATOR(nickname="cosh",\
1852     rng=[-5.,5.],
1853     test_expr="math.cosh(%a1%)",
1854     numarray_expr="numarray.cosh(%a1%)",
1855     diff="sinh(%a1%)",
1856     name="hyperbolic cosine"),
1857     OPERATOR(nickname="tanh",\
1858     rng=[-5.,5.],
1859     test_expr="math.tanh(%a1%)",
1860     numarray_expr="numarray.tanh(%a1%)",
1861     diff="1./cosh(%a1%)**2",
1862     name="hyperbolic tangent"),
1863     OPERATOR(nickname="asinh",\
1864     rng=[-100.,100.], \
1865     test_expr="numarray.arcsinh(%a1%)",
1866     math_expr="numarray.arcsinh(%a1%)",
1867     numarray_expr="numarray.arcsinh(%a1%)",
1868     diff="1./sqrt(%a1%**2+1)",
1869     name="inverse hyperbolic sine"),
1870     OPERATOR(nickname="acosh",\
1871     rng=[1.001,100.],\
1872     test_expr="numarray.arccosh(%a1%)",
1873     math_expr="numarray.arccosh(%a1%)",
1874     numarray_expr="numarray.arccosh(%a1%)",
1875     diff="1./sqrt(%a1%**2-1)",
1876     name<