Parent Directory
|
Revision Log
1D and 3D tests for eigenvalues_and_eigenvector added.
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*ranm.random()+rng[0] |
183 | else: |
184 | raise SystemError,"rank is restricted to 5" |
185 | return out |
186 | |
187 | def makeNumberedArray(shape,s=1.): |
188 | out=numarray.zeros(shape,numarray.Float64) |
189 | if len(shape)==0: |
190 | out=s*1. |
191 | elif len(shape)==1: |
192 | for i0 in range(shape[0]): |
193 | out[i0]=s*int(8*random.random()+1) |
194 | elif len(shape)==2: |
195 | for i0 in range(shape[0]): |
196 | for i1 in range(shape[1]): |
197 | out[i0,i1]=s*int(8*random.random()+1) |
198 | elif len(shape)==3: |
199 | for i0 in range(shape[0]): |
200 | for i1 in range(shape[1]): |
201 | for i2 in range(shape[2]): |
202 | out[i0,i1,i2]=s*int(8*random.random()+1) |
203 | elif len(shape)==4: |
204 | for i0 in range(shape[0]): |
205 | for i1 in range(shape[1]): |
206 | for i2 in range(shape[2]): |
207 | for i3 in range(shape[3]): |
208 | out[i0,i1,i2,i3]=s*int(8*random.random()+1) |
209 | else: |
210 | raise SystemError,"rank is restricted to 4" |
211 | return out |
212 | |
213 | def makeResult(val,test_expr): |
214 | if isinstance(val,float): |
215 | out=eval(test_expr.replace("%a1%","val")) |
216 | elif len(val.shape)==0: |
217 | out=eval(test_expr.replace("%a1%","val")) |
218 | elif len(val.shape)==1: |
219 | out=numarray.zeros(val.shape,numarray.Float64) |
220 | for i0 in range(val.shape[0]): |
221 | out[i0]=eval(test_expr.replace("%a1%","val[i0]")) |
222 | elif len(val.shape)==2: |
223 | out=numarray.zeros(val.shape,numarray.Float64) |
224 | for i0 in range(val.shape[0]): |
225 | for i1 in range(val.shape[1]): |
226 | out[i0,i1]=eval(test_expr.replace("%a1%","val[i0,i1]")) |
227 | elif len(val.shape)==3: |
228 | out=numarray.zeros(val.shape,numarray.Float64) |
229 | for i0 in range(val.shape[0]): |
230 | for i1 in range(val.shape[1]): |
231 | for i2 in range(val.shape[2]): |
232 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val[i0,i1,i2]")) |
233 | elif len(val.shape)==4: |
234 | out=numarray.zeros(val.shape,numarray.Float64) |
235 | for i0 in range(val.shape[0]): |
236 | for i1 in range(val.shape[1]): |
237 | for i2 in range(val.shape[2]): |
238 | for i3 in range(val.shape[3]): |
239 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val[i0,i1,i2,i3]")) |
240 | else: |
241 | raise SystemError,"rank is restricted to 4" |
242 | return out |
243 | |
244 | def makeResult2(val0,val1,test_expr): |
245 | if isinstance(val0,float): |
246 | if isinstance(val1,float): |
247 | out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1")) |
248 | elif len(val1.shape)==0: |
249 | out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1")) |
250 | elif len(val1.shape)==1: |
251 | out=numarray.zeros(val1.shape,numarray.Float64) |
252 | for i0 in range(val1.shape[0]): |
253 | out[i0]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0]")) |
254 | elif len(val1.shape)==2: |
255 | out=numarray.zeros(val1.shape,numarray.Float64) |
256 | for i0 in range(val1.shape[0]): |
257 | for i1 in range(val1.shape[1]): |
258 | out[i0,i1]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1]")) |
259 | elif len(val1.shape)==3: |
260 | out=numarray.zeros(val1.shape,numarray.Float64) |
261 | for i0 in range(val1.shape[0]): |
262 | for i1 in range(val1.shape[1]): |
263 | for i2 in range(val1.shape[2]): |
264 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2]")) |
265 | elif len(val1.shape)==4: |
266 | out=numarray.zeros(val1.shape,numarray.Float64) |
267 | for i0 in range(val1.shape[0]): |
268 | for i1 in range(val1.shape[1]): |
269 | for i2 in range(val1.shape[2]): |
270 | for i3 in range(val1.shape[3]): |
271 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2,i3]")) |
272 | else: |
273 | raise SystemError,"rank of val1 is restricted to 4" |
274 | elif len(val0.shape)==0: |
275 | if isinstance(val1,float): |
276 | out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1")) |
277 | elif len(val1.shape)==0: |
278 | out=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1")) |
279 | elif len(val1.shape)==1: |
280 | out=numarray.zeros(val1.shape,numarray.Float64) |
281 | for i0 in range(val1.shape[0]): |
282 | out[i0]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0]")) |
283 | elif len(val1.shape)==2: |
284 | out=numarray.zeros(val1.shape,numarray.Float64) |
285 | for i0 in range(val1.shape[0]): |
286 | for i1 in range(val1.shape[1]): |
287 | out[i0,i1]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1]")) |
288 | elif len(val1.shape)==3: |
289 | out=numarray.zeros(val1.shape,numarray.Float64) |
290 | for i0 in range(val1.shape[0]): |
291 | for i1 in range(val1.shape[1]): |
292 | for i2 in range(val1.shape[2]): |
293 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2]")) |
294 | elif len(val1.shape)==4: |
295 | out=numarray.zeros(val1.shape,numarray.Float64) |
296 | for i0 in range(val1.shape[0]): |
297 | for i1 in range(val1.shape[1]): |
298 | for i2 in range(val1.shape[2]): |
299 | for i3 in range(val1.shape[3]): |
300 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0").replace("%a2%","val1[i0,i1,i2,i3]")) |
301 | else: |
302 | raise SystemError,"rank of val1 is restricted to 4" |
303 | elif len(val0.shape)==1: |
304 | if isinstance(val1,float): |
305 | out=numarray.zeros(val0.shape,numarray.Float64) |
306 | for i0 in range(val0.shape[0]): |
307 | out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1")) |
308 | elif len(val1.shape)==0: |
309 | out=numarray.zeros(val0.shape,numarray.Float64) |
310 | for i0 in range(val0.shape[0]): |
311 | out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1")) |
312 | elif len(val1.shape)==1: |
313 | out=numarray.zeros(val0.shape,numarray.Float64) |
314 | for i0 in range(val0.shape[0]): |
315 | out[i0]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0]")) |
316 | elif len(val1.shape)==2: |
317 | out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64) |
318 | for i0 in range(val0.shape[0]): |
319 | for j1 in range(val1.shape[1]): |
320 | out[i0,j1]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1]")) |
321 | elif len(val1.shape)==3: |
322 | out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64) |
323 | for i0 in range(val0.shape[0]): |
324 | for j1 in range(val1.shape[1]): |
325 | for j2 in range(val1.shape[2]): |
326 | out[i0,j1,j2]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1,j2]")) |
327 | elif len(val1.shape)==4: |
328 | out=numarray.zeros(val0.shape+val1.shape[1:],numarray.Float64) |
329 | for i0 in range(val0.shape[0]): |
330 | for j1 in range(val1.shape[1]): |
331 | for j2 in range(val1.shape[2]): |
332 | for j3 in range(val1.shape[3]): |
333 | out[i0,j1,j2,j3]=eval(test_expr.replace("%a1%","val0[i0]").replace("%a2%","val1[i0,j1,j2,j3]")) |
334 | else: |
335 | raise SystemError,"rank of val1 is restricted to 4" |
336 | elif len(val0.shape)==2: |
337 | if isinstance(val1,float): |
338 | out=numarray.zeros(val0.shape,numarray.Float64) |
339 | for i0 in range(val0.shape[0]): |
340 | for i1 in range(val0.shape[1]): |
341 | out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1")) |
342 | elif len(val1.shape)==0: |
343 | out=numarray.zeros(val0.shape,numarray.Float64) |
344 | for i0 in range(val0.shape[0]): |
345 | for i1 in range(val0.shape[1]): |
346 | out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1")) |
347 | elif len(val1.shape)==1: |
348 | out=numarray.zeros(val0.shape,numarray.Float64) |
349 | for i0 in range(val0.shape[0]): |
350 | for i1 in range(val0.shape[1]): |
351 | out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0]")) |
352 | elif len(val1.shape)==2: |
353 | out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64) |
354 | for i0 in range(val0.shape[0]): |
355 | for i1 in range(val0.shape[1]): |
356 | out[i0,i1]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1]")) |
357 | elif len(val1.shape)==3: |
358 | out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64) |
359 | for i0 in range(val0.shape[0]): |
360 | for i1 in range(val0.shape[1]): |
361 | for j2 in range(val1.shape[2]): |
362 | out[i0,i1,j2]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1,j2]")) |
363 | elif len(val1.shape)==4: |
364 | out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64) |
365 | for i0 in range(val0.shape[0]): |
366 | for i1 in range(val0.shape[1]): |
367 | for j2 in range(val1.shape[2]): |
368 | for j3 in range(val1.shape[3]): |
369 | out[i0,i1,j2,j3]=eval(test_expr.replace("%a1%","val0[i0,i1]").replace("%a2%","val1[i0,i1,j2,j3]")) |
370 | else: |
371 | raise SystemError,"rank of val1 is restricted to 4" |
372 | elif len(val0.shape)==3: |
373 | if isinstance(val1,float): |
374 | out=numarray.zeros(val0.shape,numarray.Float64) |
375 | for i0 in range(val0.shape[0]): |
376 | for i1 in range(val0.shape[1]): |
377 | for i2 in range(val0.shape[2]): |
378 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1")) |
379 | elif len(val1.shape)==0: |
380 | out=numarray.zeros(val0.shape,numarray.Float64) |
381 | for i0 in range(val0.shape[0]): |
382 | for i1 in range(val0.shape[1]): |
383 | for i2 in range(val0.shape[2]): |
384 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1")) |
385 | elif len(val1.shape)==1: |
386 | out=numarray.zeros(val0.shape,numarray.Float64) |
387 | for i0 in range(val0.shape[0]): |
388 | for i1 in range(val0.shape[1]): |
389 | for i2 in range(val0.shape[2]): |
390 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0]")) |
391 | elif len(val1.shape)==2: |
392 | out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64) |
393 | for i0 in range(val0.shape[0]): |
394 | for i1 in range(val0.shape[1]): |
395 | for i2 in range(val0.shape[2]): |
396 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1]")) |
397 | elif len(val1.shape)==3: |
398 | out=numarray.zeros(val0.shape,numarray.Float64) |
399 | for i0 in range(val0.shape[0]): |
400 | for i1 in range(val0.shape[1]): |
401 | for i2 in range(val0.shape[2]): |
402 | out[i0,i1,i2]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1,i2]")) |
403 | elif len(val1.shape)==4: |
404 | out=numarray.zeros(val0.shape+val1.shape[3:],numarray.Float64) |
405 | for i0 in range(val0.shape[0]): |
406 | for i1 in range(val0.shape[1]): |
407 | for i2 in range(val0.shape[2]): |
408 | for j3 in range(val1.shape[3]): |
409 | out[i0,i1,i2,j3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2]").replace("%a2%","val1[i0,i1,i2,j3]")) |
410 | else: |
411 | raise SystemError,"rank of val1 is rargs[1]estricted to 4" |
412 | elif len(val0.shape)==4: |
413 | if isinstance(val1,float): |
414 | out=numarray.zeros(val0.shape,numarray.Float64) |
415 | for i0 in range(val0.shape[0]): |
416 | for i1 in range(val0.shape[1]): |
417 | for i2 in range(val0.shape[2]): |
418 | for i3 in range(val0.shape[3]): |
419 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1")) |
420 | elif len(val1.shape)==0: |
421 | out=numarray.zeros(val0.shape,numarray.Float64) |
422 | for i0 in range(val0.shape[0]): |
423 | for i1 in range(val0.shape[1]): |
424 | for i2 in range(val0.shape[2]): |
425 | for i3 in range(val0.shape[3]): |
426 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1")) |
427 | elif len(val1.shape)==1: |
428 | out=numarray.zeros(val0.shape,numarray.Float64) |
429 | for i0 in range(val0.shape[0]): |
430 | for i1 in range(val0.shape[1]): |
431 | for i2 in range(val0.shape[2]): |
432 | for i3 in range(val0.shape[3]): |
433 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0]")) |
434 | elif len(val1.shape)==2: |
435 | out=numarray.zeros(val0.shape+val1.shape[2:],numarray.Float64) |
436 | for i0 in range(val0.shape[0]): |
437 | for i1 in range(val0.shape[1]): |
438 | for i2 in range(val0.shape[2]): |
439 | for i3 in range(val0.shape[3]): |
440 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1]")) |
441 | elif len(val1.shape)==3: |
442 | out=numarray.zeros(val0.shape,numarray.Float64) |
443 | for i0 in range(val0.shape[0]): |
444 | for i1 in range(val0.shape[1]): |
445 | for i2 in range(val0.shape[2]): |
446 | for i3 in range(val0.shape[3]): |
447 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1,i2]")) |
448 | elif len(val1.shape)==4: |
449 | out=numarray.zeros(val0.shape,numarray.Float64) |
450 | for i0 in range(val0.shape[0]): |
451 | for i1 in range(val0.shape[1]): |
452 | for i2 in range(val0.shape[2]): |
453 | for i3 in range(val0.shape[3]): |
454 | out[i0,i1,i2,i3]=eval(test_expr.replace("%a1%","val0[i0,i1,i2,i3]").replace("%a2%","val1[i0,i1,i2,i3]")) |
455 | else: |
456 | raise SystemError,"rank of val1 is restricted to 4" |
457 | else: |
458 | raise SystemError,"rank is restricted to 4" |
459 | return out |
460 | |
461 | |
462 | def mkText(case,name,a,a1=None,use_tagging_for_expanded_data=False): |
463 | t_out="" |
464 | if case=="float": |
465 | if isinstance(a,float): |
466 | t_out+=" %s=%s\n"%(name,a) |
467 | elif a.rank==0: |
468 | t_out+=" %s=%s\n"%(name,a) |
469 | else: |
470 | t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist()) |
471 | elif case=="array": |
472 | if isinstance(a,float): |
473 | t_out+=" %s=numarray.array(%s)\n"%(name,a) |
474 | elif a.rank==0: |
475 | t_out+=" %s=numarray.array(%s)\n"%(name,a) |
476 | else: |
477 | t_out+=" %s=numarray.array(%s)\n"%(name,a.tolist()) |
478 | elif case=="constData": |
479 | if isinstance(a,float): |
480 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
481 | elif a.rank==0: |
482 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
483 | else: |
484 | t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist()) |
485 | elif case=="taggedData": |
486 | if isinstance(a,float): |
487 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
488 | t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1) |
489 | elif a.rank==0: |
490 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
491 | t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1) |
492 | else: |
493 | t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist()) |
494 | t_out+=" %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist()) |
495 | elif case=="expandedData": |
496 | if use_tagging_for_expanded_data: |
497 | if isinstance(a,float): |
498 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
499 | t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1) |
500 | elif a.rank==0: |
501 | t_out+=" %s=Data(%s,self.functionspace)\n"%(name,a) |
502 | t_out+=" %s.setTaggedValue(1,%s)\n"%(name,a1) |
503 | else: |
504 | t_out+=" %s=Data(numarray.array(%s),self.functionspace)\n"%(name,a.tolist()) |
505 | t_out+=" %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist()) |
506 | t_out+=" %s.expand()\n"%name |
507 | else: |
508 | t_out+=" msk_%s=whereZero(self.functionspace.getX()[0],1.e-8)\n"%name |
509 | if isinstance(a,float): |
510 | t_out+=" %s=(1.-msk_%s)*(%s)+msk_%s*(%s)\n"%(name,name,a,name,a1) |
511 | elif a.rank==0: |
512 | t_out+=" %s=msk_%s*numarray.array(%s)+(1.-msk_%s)*numarray.array(%s)\n"%(name,name,a,name,a1) |
513 | else: |
514 | t_out+=" %s=msk_%s*numarray.array(%s)+(1.-msk_%s)*numarray.array(%s)\n"%(name,name,a.tolist(),name,a1.tolist()) |
515 | elif case=="Symbol": |
516 | if isinstance(a,float): |
517 | t_out+=" %s=Symbol(shape=())\n"%(name) |
518 | elif a.rank==0: |
519 | t_out+=" %s=Symbol(shape=())\n"%(name) |
520 | else: |
521 | t_out+=" %s=Symbol(shape=%s)\n"%(name,str(a.shape)) |
522 | |
523 | return t_out |
524 | |
525 | def mkTypeAndShapeTest(case,sh,argstr,name=""): |
526 | text="" |
527 | if case=="float": |
528 | text+=" self.failUnless(isinstance(%s,float),\"wrong type of result%s.\")\n"%(argstr,name) |
529 | elif case=="array": |
530 | text+=" self.failUnless(isinstance(%s,numarray.NumArray),\"wrong type of result%s.\")\n"%(argstr,name) |
531 | text+=" self.failUnlessEqual(%s.shape,%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name) |
532 | elif case in ["constData","taggedData","expandedData"]: |
533 | text+=" self.failUnless(isinstance(%s,Data),\"wrong type of result%s.\")\n"%(argstr,name) |
534 | text+=" self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name) |
535 | elif case=="Symbol": |
536 | text+=" self.failUnless(isinstance(%s,Symbol),\"wrong type of result%s.\")\n"%(argstr,name) |
537 | text+=" self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name) |
538 | return text |
539 | |
540 | def mkCode(txt,args=[],intend=""): |
541 | s=txt.split("\n") |
542 | if len(s)>1: |
543 | out="" |
544 | for l in s: |
545 | out+=intend+l+"\n" |
546 | else: |
547 | out="%sreturn %s\n"%(intend,txt) |
548 | c=1 |
549 | for r in args: |
550 | out=out.replace("%%a%s%%"%c,r) |
551 | return out |
552 | |
553 | |
554 | #======================================================================================================= |
555 | # eigenvalues and eigen vectors 2D: |
556 | #======================================================================================================= |
557 | alltests= \ |
558 | [ ("case0",[[0.0, 0.0], [0.0, 0.0]],[0.0, 0.0]) \ |
559 | , ("case3",[[-1.0, 0.0], [0.0, -1.0]],[-1.0, -1.0]) \ |
560 | , ("case5",[[-0.99999999999999967, -6.4606252205695602e-16], [-6.4606252205695602e-16, -0.99999999999999967]],[-1.0, -1.0]) \ |
561 | , ("case6",[[0.0, 0.0], [0.0, 0.0001]],[0.0, 0.0001]) \ |
562 | , ("case7",[[0.0001, 0.0], [0.0, 0.0]],[0.0, 0.0001]) \ |
563 | , ("case8",[[6.0598371831785722e-06, 2.3859213977648625e-05], [2.3859213977648629e-05, 9.3940162816821425e-05]],[0.0, 0.0001]) \ |
564 | , ("case9",[[1.0, 0.0], [0.0, 2.0]],[1.0, 2.0]) \ |
565 | , ("case10",[[2.0, 0.0], [0.0, 1.0]],[1.0, 2.0]) \ |
566 | , ("case11",[[1.0605983718317855, 0.23859213977648688], [0.23859213977648688, 1.9394016281682138]],[1.0, 2.0]) \ |
567 | , ("case12",[[1.0, 0.0], [0.0, 1000000.0]],[1.0, 1000000.0]) \ |
568 | , ("case13",[[1000000.0, 0.0], [0.0, 1.0]],[1.0, 1000000.0]) \ |
569 | , ("case14",[[60599.311233413886, 238591.90118434647], [238591.90118434647, 939401.68876658613]],[1.0, 1000000.0]) \ |
570 | ] |
571 | dim=2 |
572 | |
573 | |
574 | alltests= \ |
575 | [ ("case0",[[0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]],[0.0, 0.0, 0.0]) \ |
576 | , ("case5",[[10.0, 0.0, 0.0], [0.0, 10.0, 0.0], [0.0, 0.0, 10.0]],[10.0, 10.0, 10.0]) \ |
577 | , ("case10",[[0.9, 0.0, 0.0], [0.0, 0.9, 0.0], [0.0, 0.0, 1.0]],[0.9, 0.9, 1.0]) \ |
578 | , ("case11",[[0.9, 0.0, 0.0], [0.0, 0.97060899725040983, -0.045555123008643325], [0.0, -0.045555123008643339, 0.92939100274959041]],[0.9, 0.9, 1.0]) \ |
579 | , ("case12",[[0.92694799760252555, 0.0, 0.044368966468320177], [0.0, 0.9, 0.0], [0.044368966468320184, 0.0, 0.97305200239747425]],[0.9, 0.9, 1.0]) \ |
580 | , ("case13",[[1.0, 0.0, 0.0], [0.0, 0.9, 0.], [0.0, 0., 0.9]],[0.9, 0.9, 1.0]) \ |
581 | , ("case14",[[0.92379770619813639, 0.041031106298491521, -0.011396846732439278], [0.041031106298491535, 0.97074428392640366, -0.019650012730342326], [-0.011396846732439236, -0.019650012730342337, 0.90545800987545966]],[0.9, 0.9, 1.0]) \ |
582 | , ("case15",[[1.0, 0.0, 0.0], [0.0, 1.1, 0.0], [0.0, 0.0, 1.1]],[1.0, 1.1, 1.1]) \ |
583 | , ("case17",[[1.0269479976025255, 0.0, 0.044368966468320309], [0.0, 1.1, 0.0], [0.044368966468320295, 0.0, 1.0730520023974743]],[1.0, 1.1, 1.1]) \ |
584 | , ("case18",[[1.1, 0.0, 0.0], [0.0, 1.0153410887977139, -0.036038311201720394], [0.0, -0.036038311201720373, 1.084658911202286]],[1.0, 1.1, 1.1]) \ |
585 | , ("case19",[[1.035487967756175, 0.026317079185831614, -0.039960133424212368], [0.026317079185831618, 1.0892641940924184, 0.016301362071911414], [-0.039960133424212355, 0.016301362071911431, 1.0752478381514063]],[1.0, 1.1, 1.1]) \ |
586 | , ("case20",[[1.0, 0.0, 0.0], [0.0, 2.0, 0.0], [0.0, 0.0, 3.0]],[1.0, 2.0, 3.0]) \ |
587 | , ("case21",[[1.0, 0.0, 0.0], [0.0, 2.7060899725040968, -0.45555123008643206], [0.0, -0.45555123008643228, 2.2939100274959037]],[1.0, 2.0, 3.0]) \ |
588 | , ("case22",[[1.5389599520505153, 0.0, 0.88737932936638753], [0.0, 2.0, 0.0], [0.88737932936638753, 0.0, 2.4610400479494858]],[1.0, 2.0, 3.0]) \ |
589 | , ("case23",[[3.0, 0.0, 0.0], [0.0, 1.153410887977139, -0.36038311201720391], [0.0, -0.36038311201720391, 1.8465891120228608]],[1.0, 2.0, 3.0]) \ |
590 | , ("case24",[[1.5928567395431172, 0.67348185484323142, -0.51356980156651744], [0.67348185484323153, 2.6000847801882254, -0.033486506584313548], [-0.51356980156651744, -0.033486506584313541, 1.8070584802686565]],[1.0, 2.0, 3.0]) \ |
591 | , ("case25",[[1.0, 0.0, 0.0], [0.0, 2.0, 0.0], [0.0, 0.0, 30000.0]],[1.0, 2.0, 30000.0]) \ |
592 | , ("case26",[[1.0, 0.0, 0.0], [0.0, 21183.286995177881, -13665.625800132779], [0.0, -13665.625800132779, 8818.7130048221279]],[1.0, 2.0, 30000.0]) \ |
593 | , ("case27",[[8085.1298007817086, 0.0, 13310.246250831115], [0.0, 2.0, 0.0], [13310.246250831115, 0.0, 21915.870199218316]],[1.0, 2.0, 30000.0]) \ |
594 | , ("case28",[[30000.0, 0.0, 0.0], [0.0, 1.153410887977139, -0.36038311201720391], [0.0, -0.36038311201720391, 1.8465891120228608]],[1.0, 2.0, 30000.0]) \ |
595 | , ("case29",[[7140.1907849945546, 12308.774438213351, -3419.2256841313947], [12308.774438213351, 21223.762934183575, -5894.4478052274408], [-3419.2256841313947, -5894.4478052274408, 1639.0462808218595]],[1.0, 2.0, 30000.0]) \ |
596 | ] |
597 | |
598 | dim=3 |
599 | |
600 | alltests= \ |
601 | [ ("case0",[[0.0]],[0.0]) \ |
602 | , ("case1",[[1.0]],[1.0]) \ |
603 | ] |
604 | dim=1 |
605 | |
606 | for case in ["constData","taggedData","expandedData"]: |
607 | n=0 |
608 | while n<len(alltests): |
609 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
610 | tname="test_eigenvalues_and_eigenvectors_%s_dim%s_%s"%(case,dim,alltests[n][0]) |
611 | text+=" def %s(self):\n"%tname |
612 | a_0=numarray.array(alltests[n][1],numarray.Float64) |
613 | ev_0=numarray.array(alltests[n][2],numarray.Float64) |
614 | if case in ["taggedData", "expandedData"]: |
615 | if n+1<len(alltests): |
616 | a1_0=numarray.array(alltests[n+1][1],numarray.Float64) |
617 | ev1_0=numarray.array(alltests[n+1][2],numarray.Float64) |
618 | else: |
619 | a1_0=numarray.array(alltests[0][1],numarray.Float64) |
620 | ev1_0=numarray.array(alltests[0][2],numarray.Float64) |
621 | n+=2 |
622 | else: |
623 | a1_0=a_0 |
624 | ev1_0=ev_0 |
625 | n+=1 |
626 | text+=mkText(case,"arg",a_0,a1_0) |
627 | text+=" res=eigenvalues_and_eigenvectors(arg)\n" |
628 | text+=mkText(case,"ref_ev",ev_0,ev1_0) |
629 | text+=mkTypeAndShapeTest(case,(dim,),"res[0]"," for eigenvalues") |
630 | text+=mkTypeAndShapeTest(case,(dim,dim),"res[1]"," for eigenvectors") |
631 | text+=" self.failUnless(Lsup(res[0]-ref_ev)<=self.RES_TOL*Lsup(ref_ev),\"wrong eigenvalues\")\n" |
632 | for i in range(dim): |
633 | text+=" self.failUnless(Lsup(matrixmult(arg,res[1][:,%s])-res[0][%s]*res[1][:,%s])<=self.RES_TOL*Lsup(res[0]),\"wrong eigenvector %s\")\n"%(i,i,i,i) |
634 | text+=" self.failUnless(Lsup(matrixmult(transpose(res[1]),res[1])-kronecker(%s))<=self.RES_TOL,\"eigenvectors are not orthonormal\")\n"%dim |
635 | if case == "taggedData" : |
636 | t_prog_with_tags+=text |
637 | else: |
638 | t_prog+=text |
639 | print test_header |
640 | print t_prog |
641 | print t_prog_with_tags |
642 | print test_tail |
643 | 1/0 |
644 | #======================================================================================================= |
645 | # get slices |
646 | #======================================================================================================= |
647 | from esys.escript import * |
648 | for case0 in ["constData","taggedData","expandedData"]: |
649 | for sh0 in [ (3,), (3,4), (3,4,3) ,(4,3,3,3)]: |
650 | # get perm: |
651 | if len(sh0)==2: |
652 | check=[[1,0]] |
653 | elif len(sh0)==3: |
654 | check=[[1,0,2], |
655 | [1,2,0], |
656 | [2,1,0], |
657 | [2,0,2], |
658 | [0,2,1]] |
659 | elif len(sh0)==4: |
660 | check=[[0,1,3,2], |
661 | [0,2,1,3], |
662 | [0,2,3,1], |
663 | [0,3,2,1], |
664 | [0,3,1,2] , |
665 | [1,0,2,3], |
666 | [1,0,3,2], |
667 | [1,2,0,3], |
668 | [1,2,3,0], |
669 | [1,3,2,0], |
670 | [1,3,0,2], |
671 | [2,0,1,3], |
672 | [2,0,3,1], |
673 | [2,1,0,3], |
674 | [2,1,3,0], |
675 | [2,3,1,0], |
676 | [2,3,0,1], |
677 | [3,0,1,2], |
678 | [3,0,2,1], |
679 | [3,1,0,2], |
680 | [3,1,2,0], |
681 | [3,2,1,0], |
682 | [3,2,0,1]] |
683 | else: |
684 | check=[] |
685 | |
686 | # create the test cases: |
687 | processed=[] |
688 | l=["R","U","L","P","C","N"] |
689 | c=[""] |
690 | for i in range(len(sh0)): |
691 | tmp=[] |
692 | for ci in c: |
693 | tmp+=[ci+li for li in l] |
694 | c=tmp |
695 | # SHUFFLE |
696 | c2=[] |
697 | while len(c)>0: |
698 | i=int(random.random()*len(c)) |
699 | c2.append(c[i]) |
700 | del c[i] |
701 | c=c2 |
702 | for ci in c: |
703 | t="" |
704 | sh=() |
705 | sl=() |
706 | for i in range(len(ci)): |
707 | if ci[i]=="R": |
708 | s="%s:%s"%(1,sh0[i]-1) |
709 | sl=sl+(slice(1,sh0[i]-1),) |
710 | sh=sh+(sh0[i]-2,) |
711 | if ci[i]=="U": |
712 | s=":%s"%(sh0[i]-1) |
713 | sh=sh+(sh0[i]-1,) |
714 | sl=sl+(slice(0,sh0[i]-1),) |
715 | if ci[i]=="L": |
716 | s="2:" |
717 | sh=sh+(sh0[i]-2,) |
718 | sl=sl+(slice(2,sh0[i]),) |
719 | if ci[i]=="P": |
720 | s="%s"%(int(sh0[i]/2)) |
721 | sl=sl+(int(sh0[i]/2),) |
722 | if ci[i]=="C": |
723 | s=":" |
724 | sh=sh+(sh0[i],) |
725 | sl=sl+(slice(0,sh0[i]),) |
726 | if ci[i]=="N": |
727 | s="" |
728 | sh=sh+(sh0[i],) |
729 | if len(s)>0: |
730 | if not t=="": t+="," |
731 | t+=s |
732 | if len(sl)==1: sl=sl[0] |
733 | N_found=False |
734 | noN_found=False |
735 | process=len(t)>0 |
736 | for i in ci: |
737 | if i=="N": |
738 | if not noN_found and N_found: process=False |
739 | N_found=True |
740 | else: |
741 | if N_found: process=False |
742 | noNfound=True |
743 | # is there a similar one processed allready |
744 | if process and ci.find("N")==-1: |
745 | for ci2 in processed: |
746 | for chi in check: |
747 | is_perm=True |
748 | for i in range(len(chi)): |
749 | if not ci[i]==ci2[chi[i]]: is_perm=False |
750 | if is_perm: process=False |
751 | # if not process: print ci," rejected" |
752 | if process: |
753 | processed.append(ci) |
754 | for case1 in ["array","constData","taggedData","expandedData"]: |
755 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
756 | tname="test_setslice_%s_rank%s_%s_%s"%(case0,len(sh0),case1,ci) |
757 | text+=" def %s(self):\n"%tname |
758 | a_0=makeNumberedArray(sh0) |
759 | if case0 in ["taggedData", "expandedData"]: |
760 | a1_0=makeNumberedArray(sh0) |
761 | else: |
762 | a1_0=a_0*1. |
763 | |
764 | a_1=makeNumberedArray(sh) |
765 | if case1 in ["taggedData", "expandedData"]: |
766 | a1_1=makeNumberedArray(sh) |
767 | else: |
768 | a1_1=a_1*1. |
769 | |
770 | text+=mkText(case0,"arg",a_0,a1_0) |
771 | text+=mkText(case1,"val",a_1,a1_1) |
772 | text+=" arg[%s]=val\n"%t |
773 | a_0.__setitem__(sl,a_1) |
774 | a1_0.__setitem__(sl,a1_1) |
775 | if Lsup(a_0-a1_0)==0.: |
776 | text+=mkText("constData","ref",a_0,a1_0) |
777 | else: |
778 | text+=mkText("expandedData","ref",a_0,a1_0) |
779 | text+=" self.failUnless(Lsup(arg-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n" |
780 | |
781 | if case0 == "taggedData" or case1 == "taggedData": |
782 | t_prog_with_tags+=text |
783 | else: |
784 | t_prog+=text |
785 | |
786 | print test_header |
787 | # print t_prog |
788 | print t_prog_with_tags |
789 | print test_tail |
790 | 1/0 |
791 | |
792 | #======================================================================================================= |
793 | # (non)symmetric part |
794 | #======================================================================================================= |
795 | from esys.escript import * |
796 | for name in ["symmetric", "nonsymmetric"]: |
797 | f=1. |
798 | if name=="nonsymmetric": f=-1 |
799 | for case0 in ["array","Symbol","constData","taggedData","expandedData"]: |
800 | for sh0 in [ (3,3), (2,3,2,3)]: |
801 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
802 | tname="test_%s_%s_rank%s"%(name,case0,len(sh0)) |
803 | text+=" def %s(self):\n"%tname |
804 | a_0=makeNumberedArray(sh0,s=1.) |
805 | r_0=(a_0+f*transpose(a_0))/2. |
806 | if case0 in ["taggedData", "expandedData"]: |
807 | a1_0=makeNumberedArray(sh0,s=-1.) |
808 | r1_0=(a1_0+f*transpose(a1_0))/2. |
809 | else: |
810 | a1_0=a_0 |
811 | r1_0=r_0 |
812 | text+=mkText(case0,"arg",a_0,a1_0) |
813 | text+=" res=%s(arg)\n"%name |
814 | if case0=="Symbol": |
815 | text+=mkText("array","s",a_0,a1_0) |
816 | text+=" sub=res.substitute({arg:s})\n" |
817 | res="sub" |
818 | text+=mkText("array","ref",r_0,r1_0) |
819 | else: |
820 | res="res" |
821 | text+=mkText(case0,"ref",r_0,r1_0) |
822 | text+=mkTypeAndShapeTest(case0,sh0,"res") |
823 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
824 | |
825 | if case0 == "taggedData" : |
826 | t_prog_with_tags+=text |
827 | else: |
828 | t_prog+=text |
829 | print test_header |
830 | print t_prog |
831 | # print t_prog_with_tags |
832 | print test_tail |
833 | 1/0 |
834 | |
835 | #======================================================================================================= |
836 | # eigenvalues |
837 | #======================================================================================================= |
838 | import numarray.linear_algebra |
839 | name="eigenvalues" |
840 | for case0 in ["array","Symbol","constData","taggedData","expandedData"]: |
841 | for sh0 in [ (1,1), (2,2), (3,3)]: |
842 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
843 | tname="test_%s_%s_dim%s"%(name,case0,sh0[0]) |
844 | text+=" def %s(self):\n"%tname |
845 | a_0=makeArray(sh0,[-1.,1]) |
846 | a_0=(a_0+numarray.transpose(a_0))/2. |
847 | ev=numarray.linear_algebra.eigenvalues(a_0) |
848 | ev.sort() |
849 | if case0 in ["taggedData", "expandedData"]: |
850 | a1_0=makeArray(sh0,[-1.,1]) |
851 | a1_0=(a1_0+numarray.transpose(a1_0))/2. |
852 | ev1=numarray.linear_algebra.eigenvalues(a1_0) |
853 | ev1.sort() |
854 | else: |
855 | a1_0=a_0 |
856 | ev1=ev |
857 | text+=mkText(case0,"arg",a_0,a1_0) |
858 | text+=" res=%s(arg)\n"%name |
859 | if case0=="Symbol": |
860 | text+=mkText("array","s",a_0,a1_0) |
861 | text+=" sub=res.substitute({arg:s})\n" |
862 | res="sub" |
863 | text+=mkText("array","ref",ev,ev1) |
864 | else: |
865 | res="res" |
866 | text+=mkText(case0,"ref",ev,ev1) |
867 | text+=mkTypeAndShapeTest(case0,(sh0[0],),"res") |
868 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
869 | |
870 | if case0 == "taggedData" : |
871 | t_prog_with_tags+=text |
872 | else: |
873 | t_prog+=text |
874 | print test_header |
875 | # print t_prog |
876 | print t_prog_with_tags |
877 | print test_tail |
878 | 1/0 |
879 | |
880 | #======================================================================================================= |
881 | # get slices |
882 | #======================================================================================================= |
883 | for case0 in ["constData","taggedData","expandedData","Symbol"]: |
884 | for sh0 in [ (3,), (3,4), (3,4,3) ,(4,3,5,3)]: |
885 | # get perm: |
886 | if len(sh0)==2: |
887 | check=[[1,0]] |
888 | elif len(sh0)==3: |
889 | check=[[1,0,2], |
890 | [1,2,0], |
891 | [2,1,0], |
892 | [2,0,2], |
893 | [0,2,1]] |
894 | elif len(sh0)==4: |
895 | check=[[0,1,3,2], |
896 | [0,2,1,3], |
897 | [0,2,3,1], |
898 | [0,3,2,1], |
899 | [0,3,1,2] , |
900 | [1,0,2,3], |
901 | [1,0,3,2], |
902 | [1,2,0,3], |
903 | [1,2,3,0], |
904 | [1,3,2,0], |
905 | [1,3,0,2], |
906 | [2,0,1,3], |
907 | [2,0,3,1], |
908 | [2,1,0,3], |
909 | [2,1,3,0], |
910 | [2,3,1,0], |
911 | [2,3,0,1], |
912 | [3,0,1,2], |
913 | [3,0,2,1], |
914 | [3,1,0,2], |
915 | [3,1,2,0], |
916 | [3,2,1,0], |
917 | [3,2,0,1]] |
918 | else: |
919 | check=[] |
920 | |
921 | # create the test cases: |
922 | processed=[] |
923 | l=["R","U","L","P","C","N"] |
924 | c=[""] |
925 | for i in range(len(sh0)): |
926 | tmp=[] |
927 | for ci in c: |
928 | tmp+=[ci+li for li in l] |
929 | c=tmp |
930 | # SHUFFLE |
931 | c2=[] |
932 | while len(c)>0: |
933 | i=int(random.random()*len(c)) |
934 | c2.append(c[i]) |
935 | del c[i] |
936 | c=c2 |
937 | for ci in c: |
938 | t="" |
939 | sh=() |
940 | for i in range(len(ci)): |
941 | if ci[i]=="R": |
942 | s="%s:%s"%(1,sh0[i]-1) |
943 | sh=sh+(sh0[i]-2,) |
944 | if ci[i]=="U": |
945 | s=":%s"%(sh0[i]-1) |
946 | sh=sh+(sh0[i]-1,) |
947 | if ci[i]=="L": |
948 | s="2:" |
949 | sh=sh+(sh0[i]-2,) |
950 | if ci[i]=="P": |
951 | s="%s"%(int(sh0[i]/2)) |
952 | if ci[i]=="C": |
953 | s=":" |
954 | sh=sh+(sh0[i],) |
955 | if ci[i]=="N": |
956 | s="" |
957 | sh=sh+(sh0[i],) |
958 | if len(s)>0: |
959 | if not t=="": t+="," |
960 | t+=s |
961 | N_found=False |
962 | noN_found=False |
963 | process=len(t)>0 |
964 | for i in ci: |
965 | if i=="N": |
966 | if not noN_found and N_found: process=False |
967 | N_found=True |
968 | else: |
969 | if N_found: process=False |
970 | noNfound=True |
971 | # is there a similar one processed allready |
972 | if process and ci.find("N")==-1: |
973 | for ci2 in processed: |
974 | for chi in check: |
975 | is_perm=True |
976 | for i in range(len(chi)): |
977 | if not ci[i]==ci2[chi[i]]: is_perm=False |
978 | if is_perm: process=False |
979 | # if not process: print ci," rejected" |
980 | if process: |
981 | processed.append(ci) |
982 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
983 | tname="test_getslice_%s_rank%s_%s"%(case0,len(sh0),ci) |
984 | text+=" def %s(self):\n"%tname |
985 | a_0=makeNumberedArray(sh0,s=1) |
986 | if case0 in ["taggedData", "expandedData"]: |
987 | a1_0=makeNumberedArray(sh0,s=-1.) |
988 | else: |
989 | a1_0=a_0 |
990 | r=eval("a_0[%s]"%t) |
991 | r1=eval("a1_0[%s]"%t) |
992 | text+=mkText(case0,"arg",a_0,a1_0) |
993 | text+=" res=arg[%s]\n"%t |
994 | if case0=="Symbol": |
995 | text+=mkText("array","s",a_0,a1_0) |
996 | text+=" sub=res.substitute({arg:s})\n" |
997 | res="sub" |
998 | text+=mkText("array","ref",r,r1) |
999 | else: |
1000 | res="res" |
1001 | text+=mkText(case0,"ref",r,r1) |
1002 | text+=mkTypeAndShapeTest(case0,sh,"res") |
1003 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1004 | |
1005 | if case0 == "taggedData" : |
1006 | t_prog_with_tags+=text |
1007 | else: |
1008 | t_prog+=text |
1009 | |
1010 | print test_header |
1011 | # print t_prog |
1012 | print t_prog_with_tags |
1013 | print test_tail |
1014 | 1/0 |
1015 | #============================================================================================ |
1016 | def innerTEST(arg0,arg1): |
1017 | if isinstance(arg0,float): |
1018 | out=numarray.array(arg0*arg1) |
1019 | else: |
1020 | out=(arg0*arg1).sum() |
1021 | return out |
1022 | |
1023 | def outerTEST(arg0,arg1): |
1024 | if isinstance(arg0,float): |
1025 | out=numarray.array(arg0*arg1) |
1026 | elif isinstance(arg1,float): |
1027 | out=numarray.array(arg0*arg1) |
1028 | else: |
1029 | out=numarray.outerproduct(arg0,arg1).resize(arg0.shape+arg1.shape) |
1030 | return out |
1031 | |
1032 | def tensorProductTest(arg0,arg1,sh_s): |
1033 | if isinstance(arg0,float): |
1034 | out=numarray.array(arg0*arg1) |
1035 | elif isinstance(arg1,float): |
1036 | out=numarray.array(arg0*arg1) |
1037 | elif len(sh_s)==0: |
1038 | out=numarray.outerproduct(arg0,arg1).resize(arg0.shape+arg1.shape) |
1039 | else: |
1040 | l=len(sh_s) |
1041 | sh0=arg0.shape[:arg0.rank-l] |
1042 | sh1=arg1.shape[l:] |
1043 | ls,l0,l1=1,1,1 |
1044 | for i in sh_s: ls*=i |
1045 | for i in sh0: l0*=i |
1046 | for i in sh1: l1*=i |
1047 | out1=numarray.outerproduct(arg0,arg1).resize((l0,ls,ls,l1)) |
1048 | out2=numarray.zeros((l0,l1),numarray.Float) |
1049 | for i0 in range(l0): |
1050 | for i1 in range(l1): |
1051 | for i in range(ls): out2[i0,i1]+=out1[i0,i,i,i1] |
1052 | out=out2.resize(sh0+sh1) |
1053 | return out |
1054 | |
1055 | def testMatrixMult(arg0,arg1,sh_s): |
1056 | return numarray.matrixmultiply(arg0,arg1) |
1057 | |
1058 | |
1059 | def testTensorMult(arg0,arg1,sh_s): |
1060 | if len(arg0)==2: |
1061 | return numarray.matrixmultiply(arg0,arg1) |
1062 | else: |
1063 | if arg1.rank==4: |
1064 | out=numarray.zeros((arg0.shape[0],arg0.shape[1],arg1.shape[2],arg1.shape[3]),numarray.Float) |
1065 | for i0 in range(arg0.shape[0]): |
1066 | for i1 in range(arg0.shape[1]): |
1067 | for i2 in range(arg0.shape[2]): |
1068 | for i3 in range(arg0.shape[3]): |
1069 | for j2 in range(arg1.shape[2]): |
1070 | for j3 in range(arg1.shape[3]): |
1071 | out[i0,i1,j2,j3]+=arg0[i0,i1,i2,i3]*arg1[i2,i3,j2,j3] |
1072 | elif arg1.rank==3: |
1073 | out=numarray.zeros((arg0.shape[0],arg0.shape[1],arg1.shape[2]),numarray.Float) |
1074 | for i0 in range(arg0.shape[0]): |
1075 | for i1 in range(arg0.shape[1]): |
1076 | for i2 in range(arg0.shape[2]): |
1077 | for i3 in range(arg0.shape[3]): |
1078 | for j2 in range(arg1.shape[2]): |
1079 | out[i0,i1,j2]+=arg0[i0,i1,i2,i3]*arg1[i2,i3,j2] |
1080 | elif arg1.rank==2: |
1081 | out=numarray.zeros((arg0.shape[0],arg0.shape[1]),numarray.Float) |
1082 | for i0 in range(arg0.shape[0]): |
1083 | for i1 in range(arg0.shape[1]): |
1084 | for i2 in range(arg0.shape[2]): |
1085 | for i3 in range(arg0.shape[3]): |
1086 | out[i0,i1]+=arg0[i0,i1,i2,i3]*arg1[i2,i3] |
1087 | return out |
1088 | |
1089 | def testReduce(arg0,init_val,test_expr,post_expr): |
1090 | out=init_val |
1091 | if isinstance(arg0,float): |
1092 | out=eval(test_expr.replace("%a1%","arg0")) |
1093 | elif arg0.rank==0: |
1094 | out=eval(test_expr.replace("%a1%","arg0")) |
1095 | elif arg0.rank==1: |
1096 | for i0 in range(arg0.shape[0]): |
1097 | out=eval(test_expr.replace("%a1%","arg0[i0]")) |
1098 | elif arg0.rank==2: |
1099 | for i0 in range(arg0.shape[0]): |
1100 | for i1 in range(arg0.shape[1]): |
1101 | out=eval(test_expr.replace("%a1%","arg0[i0,i1]")) |
1102 | elif arg0.rank==3: |
1103 | for i0 in range(arg0.shape[0]): |
1104 | for i1 in range(arg0.shape[1]): |
1105 | for i2 in range(arg0.shape[2]): |
1106 | out=eval(test_expr.replace("%a1%","arg0[i0,i1,i2]")) |
1107 | elif arg0.rank==4: |
1108 | for i0 in range(arg0.shape[0]): |
1109 | for i1 in range(arg0.shape[1]): |
1110 | for i2 in range(arg0.shape[2]): |
1111 | for i3 in range(arg0.shape[3]): |
1112 | out=eval(test_expr.replace("%a1%","arg0[i0,i1,i2,i3]")) |
1113 | return eval(post_expr) |
1114 | |
1115 | def clipTEST(arg0,mn,mx): |
1116 | if isinstance(arg0,float): |
1117 | return max(min(arg0,mx),mn) |
1118 | out=numarray.zeros(arg0.shape,numarray.Float64) |
1119 | if arg0.rank==1: |
1120 | for i0 in range(arg0.shape[0]): |
1121 | out[i0]=max(min(arg0[i0],mx),mn) |
1122 | elif arg0.rank==2: |
1123 | for i0 in range(arg0.shape[0]): |
1124 | for i1 in range(arg0.shape[1]): |
1125 | out[i0,i1]=max(min(arg0[i0,i1],mx),mn) |
1126 | elif arg0.rank==3: |
1127 | for i0 in range(arg0.shape[0]): |
1128 | for i1 in range(arg0.shape[1]): |
1129 | for i2 in range(arg0.shape[2]): |
1130 | out[i0,i1,i2]=max(min(arg0[i0,i1,i2],mx),mn) |
1131 | elif arg0.rank==4: |
1132 | for i0 in range(arg0.shape[0]): |
1133 | for i1 in range(arg0.shape[1]): |
1134 | for i2 in range(arg0.shape[2]): |
1135 | for i3 in range(arg0.shape[3]): |
1136 | out[i0,i1,i2,i3]=max(min(arg0[i0,i1,i2,i3],mx),mn) |
1137 | return out |
1138 | def minimumTEST(arg0,arg1): |
1139 | if isinstance(arg0,float): |
1140 | if isinstance(arg1,float): |
1141 | if arg0>arg1: |
1142 | return arg1 |
1143 | else: |
1144 | return arg0 |
1145 | else: |
1146 | arg0=numarray.ones(arg1.shape)*arg0 |
1147 | else: |
1148 | if isinstance(arg1,float): |
1149 | arg1=numarray.ones(arg0.shape)*arg1 |
1150 | out=numarray.zeros(arg0.shape,numarray.Float64) |
1151 | if arg0.rank==0: |
1152 | if arg0>arg1: |
1153 | out=arg1 |
1154 | else: |
1155 | out=arg0 |
1156 | elif arg0.rank==1: |
1157 | for i0 in range(arg0.shape[0]): |
1158 | if arg0[i0]>arg1[i0]: |
1159 | out[i0]=arg1[i0] |
1160 | else: |
1161 | out[i0]=arg0[i0] |
1162 | elif arg0.rank==2: |
1163 | for i0 in range(arg0.shape[0]): |
1164 | for i1 in range(arg0.shape[1]): |
1165 | if arg0[i0,i1]>arg1[i0,i1]: |
1166 | out[i0,i1]=arg1[i0,i1] |
1167 | else: |
1168 | out[i0,i1]=arg0[i0,i1] |
1169 | elif arg0.rank==3: |
1170 | for i0 in range(arg0.shape[0]): |
1171 | for i1 in range(arg0.shape[1]): |
1172 | for i2 in range(arg0.shape[2]): |
1173 | if arg0[i0,i1,i2]>arg1[i0,i1,i2]: |
1174 | out[i0,i1,i2]=arg1[i0,i1,i2] |
1175 | else: |
1176 | out[i0,i1,i2]=arg0[i0,i1,i2] |
1177 | elif arg0.rank==4: |
1178 | for i0 in range(arg0.shape[0]): |
1179 | for i1 in range(arg0.shape[1]): |
1180 | for i2 in range(arg0.shape[2]): |
1181 | for i3 in range(arg0.shape[3]): |
1182 | if arg0[i0,i1,i2,i3]>arg1[i0,i1,i2,i3]: |
1183 | out[i0,i1,i2,i3]=arg1[i0,i1,i2,i3] |
1184 | else: |
1185 | out[i0,i1,i2,i3]=arg0[i0,i1,i2,i3] |
1186 | return out |
1187 | |
1188 | def unrollLoops(a,b,o,arg,tap="",x="x"): |
1189 | out="" |
1190 | if a.rank==1: |
1191 | z="" |
1192 | for i99 in range(a.shape[0]): |
1193 | if not z=="": z+="+" |
1194 | if o=="1": |
1195 | z+="(%s)*%s[%s]"%(a[i99]+b[i99],x,i99) |
1196 | else: |
1197 | z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i99],x,i99,b[i99],x,i99) |
1198 | |
1199 | out+=tap+"%s=%s\n"%(arg,z) |
1200 | |
1201 | elif a.rank==2: |
1202 | for i0 in range(a.shape[0]): |
1203 | z="" |
1204 | for i99 in range(a.shape[1]): |
1205 | if not z=="": z+="+" |
1206 | if o=="1": |
1207 | z+="(%s)*x[%s]"%(a[i0,i99]+b[i0,i99],i99) |
1208 | else: |
1209 | z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i99],x,i99,b[i0,i99],x,i99) |
1210 | |
1211 | out+=tap+"%s[%s]=%s\n"%(arg,i0,z) |
1212 | elif a.rank==3: |
1213 | for i0 in range(a.shape[0]): |
1214 | for i1 in range(a.shape[1]): |
1215 | z="" |
1216 | for i99 in range(a.shape[2]): |
1217 | if not z=="": z+="+" |
1218 | if o=="1": |
1219 | z+="(%s)*%s[%s]"%(a[i0,i1,i99]+b[i0,i1,i99],x,i99) |
1220 | else: |
1221 | z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i99],x,i99,b[i0,i1,i99],x,i99) |
1222 | |
1223 | out+=tap+"%s[%s,%s]=%s\n"%(arg,i0,i1,z) |
1224 | elif a.rank==4: |
1225 | for i0 in range(a.shape[0]): |
1226 | for i1 in range(a.shape[1]): |
1227 | for i2 in range(a.shape[2]): |
1228 | z="" |
1229 | for i99 in range(a.shape[3]): |
1230 | if not z=="": z+="+" |
1231 | if o=="1": |
1232 | z+="(%s)*%s[%s]"%(a[i0,i1,i2,i99]+b[i0,i1,i2,i99],x,i99) |
1233 | else: |
1234 | z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i99],x,i99,b[i0,i1,i2,i99],x,i99) |
1235 | |
1236 | out+=tap+"%s[%s,%s,%s]=%s\n"%(arg,i0,i1,i2,z) |
1237 | elif a.rank==5: |
1238 | for i0 in range(a.shape[0]): |
1239 | for i1 in range(a.shape[1]): |
1240 | for i2 in range(a.shape[2]): |
1241 | for i3 in range(a.shape[3]): |
1242 | z="" |
1243 | for i99 in range(a.shape[4]): |
1244 | if not z=="": z+="+" |
1245 | if o=="1": |
1246 | z+="(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99]+b[i0,i1,i2,i3,i99],x,i99) |
1247 | else: |
1248 | z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99],x,i99,b[i0,i1,i2,i3,i99],x,i99) |
1249 | |
1250 | out+=tap+"%s[%s,%s,%s,%s]=%s\n"%(arg,i0,i1,i2,i3,z) |
1251 | return out |
1252 | |
1253 | def unrollLoopsOfGrad(a,b,o,arg,tap=""): |
1254 | out="" |
1255 | if a.rank==1: |
1256 | for i99 in range(a.shape[0]): |
1257 | if o=="1": |
1258 | out+=tap+"%s[%s]=(%s)\n"%(arg,i99,a[i99]+b[i99]) |
1259 | else: |
1260 | out+=tap+"%s[%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i99,a[i99],i99,b[i99]) |
1261 | |
1262 | elif a.rank==2: |
1263 | for i0 in range(a.shape[0]): |
1264 | for i99 in range(a.shape[1]): |
1265 | if o=="1": |
1266 | out+=tap+"%s[%s,%s]=(%s)\n"%(arg,i0,i99,a[i0,i99]+b[i0,i99]) |
1267 | else: |
1268 | out+=tap+"%s[%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i99,a[i0,i99],i99,b[i0,i99]) |
1269 | |
1270 | elif a.rank==3: |
1271 | for i0 in range(a.shape[0]): |
1272 | for i1 in range(a.shape[1]): |
1273 | for i99 in range(a.shape[2]): |
1274 | if o=="1": |
1275 | out+=tap+"%s[%s,%s,%s]=(%s)\n"%(arg,i0,i1,i99,a[i0,i1,i99]+b[i0,i1,i99]) |
1276 | else: |
1277 | 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]) |
1278 | |
1279 | elif a.rank==4: |
1280 | for i0 in range(a.shape[0]): |
1281 | for i1 in range(a.shape[1]): |
1282 | for i2 in range(a.shape[2]): |
1283 | for i99 in range(a.shape[3]): |
1284 | if o=="1": |
1285 | out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99]) |
1286 | else: |
1287 | 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]) |
1288 | return out |
1289 | def unrollLoopsOfDiv(a,b,o,arg,tap=""): |
1290 | out=tap+arg+"=" |
1291 | if o=="1": |
1292 | z=0. |
1293 | for i99 in range(a.shape[0]): |
1294 | z+=b[i99,i99]+a[i99,i99] |
1295 | out+="(%s)"%z |
1296 | else: |
1297 | z=0. |
1298 | for i99 in range(a.shape[0]): |
1299 | z+=b[i99,i99] |
1300 | if i99>0: out+="+" |
1301 | out+="o*(%s)*x_ref[%s]**(o-1)"%(a[i99,i99],i99) |
1302 | out+="+(%s)"%z |
1303 | return out |
1304 | |
1305 | def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""): |
1306 | if where=="Function": |
1307 | xfac_o=1. |
1308 | xfac_op=0. |
1309 | z_fac=1./2. |
1310 | z_fac_s="" |
1311 | zo_fac_s="/(o+1.)" |
1312 | elif where=="FunctionOnBoundary": |
1313 | xfac_o=1. |
1314 | xfac_op=0. |
1315 | z_fac=1. |
1316 | z_fac_s="*dim" |
1317 | zo_fac_s="*(1+2.*(dim-1.)/(o+1.))" |
1318 | elif where in ["FunctionOnContactZero","FunctionOnContactOne"]: |
1319 | xfac_o=0. |
1320 | xfac_op=1. |
1321 | z_fac=1./2. |
1322 | z_fac_s="" |
1323 | zo_fac_s="/(o+1.)" |
1324 | out="" |
1325 | if a.rank==1: |
1326 | zo=0. |
1327 | zop=0. |
1328 | z=0. |
1329 | for i99 in range(a.shape[0]): |
1330 | if i99==0: |
1331 | zo+= xfac_o*a[i99] |
1332 | zop+= xfac_op*a[i99] |
1333 | else: |
1334 | zo+=a[i99] |
1335 | z+=b[i99] |
1336 | |
1337 | out+=tap+"%s=(%s)%s+(%s)%s"%(arg,zo,zo_fac_s,z*z_fac,z_fac_s) |
1338 | if zop==0.: |
1339 | out+="\n" |
1340 | else: |
1341 | out+="+(%s)*0.5**o\n"%zop |
1342 | elif a.rank==2: |
1343 | for i0 in range(a.shape[0]): |
1344 | zo=0. |
1345 | zop=0. |
1346 | z=0. |
1347 | for i99 in range(a.shape[1]): |
1348 | if i99==0: |
1349 | zo+= xfac_o*a[i0,i99] |
1350 | zop+= xfac_op*a[i0,i99] |
1351 | else: |
1352 | zo+=a[i0,i99] |
1353 | z+=b[i0,i99] |
1354 | |
1355 | out+=tap+"%s[%s]=(%s)%s+(%s)%s"%(arg,i0,zo,zo_fac_s,z*z_fac,z_fac_s) |
1356 | if zop==0.: |
1357 | out+="\n" |
1358 | else: |
1359 | out+="+(%s)*0.5**o\n"%zop |
1360 | elif a.rank==3: |
1361 | for i0 in range(a.shape[0]): |
1362 | for i1 in range(a.shape[1]): |
1363 | zo=0. |
1364 | zop=0. |
1365 | z=0. |
1366 | for i99 in range(a.shape[2]): |
1367 | if i99==0: |
1368 | zo+= xfac_o*a[i0,i1,i99] |
1369 | zop+= xfac_op*a[i0,i1,i99] |
1370 | else: |
1371 | zo+=a[i0,i1,i99] |
1372 | z+=b[i0,i1,i99] |
1373 | |
1374 | out+=tap+"%s[%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,zo,zo_fac_s,z*z_fac,z_fac_s) |
1375 | if zop==0.: |
1376 | out+="\n" |
1377 | else: |
1378 | out+="+(%s)*0.5**o\n"%zop |
1379 | elif a.rank==4: |
1380 | for i0 in range(a.shape[0]): |
1381 | for i1 in range(a.shape[1]): |
1382 | for i2 in range(a.shape[2]): |
1383 | zo=0. |
1384 | zop=0. |
1385 | z=0. |
1386 | for i99 in range(a.shape[3]): |
1387 | if i99==0: |
1388 | zo+= xfac_o*a[i0,i1,i2,i99] |
1389 | zop+= xfac_op*a[i0,i1,i2,i99] |
1390 | |
1391 | else: |
1392 | zo+=a[i0,i1,i2,i99] |
1393 | z+=b[i0,i1,i2,i99] |
1394 | |
1395 | out+=tap+"%s[%s,%s,%s]=(%s)%s+(%s)%s"%(arg,i0,i1,i2,zo,zo_fac_s,z*z_fac,z_fac_s) |
1396 | if zop==0.: |
1397 | out+="\n" |
1398 | else: |
1399 | out+="+(%s)*0.5**o\n"%zop |
1400 | |
1401 | elif a.rank==5: |
1402 | for i0 in range(a.shape[0]): |
1403 | for i1 in range(a.shape[1]): |
1404 | for i2 in range(a.shape[2]): |
1405 | for i3 in range(a.shape[3]): |
1406 | zo=0. |
1407 | zop=0. |
1408 | z=0. |
1409 | for i99 in range(a.shape[4]): |
1410 | if i99==0: |
1411 | zo+= xfac_o*a[i0,i1,i2,i3,i99] |
1412 | zop+= xfac_op*a[i0,i1,i2,i3,i99] |
1413 | |
1414 | else: |
1415 | zo+=a[i0,i1,i2,i3,i99] |
1416 | z+=b[i0,i1,i2,i3,i99] |
1417 | 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) |
1418 | if zop==0.: |
1419 | out+="\n" |
1420 | else: |
1421 | out+="+(%s)*0.5**o\n"%zop |
1422 | |
1423 | return out |
1424 | def unrollLoopsSimplified(b,arg,tap=""): |
1425 | out="" |
1426 | if isinstance(b,float) or b.rank==0: |
1427 | out+=tap+"%s=(%s)*x[0]\n"%(arg,str(b)) |
1428 | |
1429 | elif b.rank==1: |
1430 | for i0 in range(b.shape[0]): |
1431 | out+=tap+"%s[%s]=(%s)*x[%s]\n"%(arg,i0,b[i0],i0) |
1432 | elif b.rank==2: |
1433 | for i0 in range(b.shape[0]): |
1434 | for i1 in range(b.shape[1]): |
1435 | out+=tap+"%s[%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,b[i0,i1],i1) |
1436 | elif b.rank==3: |
1437 | for i0 in range(b.shape[0]): |
1438 | for i1 in range(b.shape[1]): |
1439 | for i2 in range(b.shape[2]): |
1440 | out+=tap+"%s[%s,%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,i2,b[i0,i1,i2],i2) |
1441 | elif b.rank==4: |
1442 | for i0 in range(b.shape[0]): |
1443 | for i1 in range(b.shape[1]): |
1444 | for i2 in range(b.shape[2]): |
1445 | for i3 in range(b.shape[3]): |
1446 | out+=tap+"%s[%s,%s,%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,i2,i3,b[i0,i1,i2,i3],i3) |
1447 | return out |
1448 | |
1449 | def unrollLoopsOfL2(b,where,arg,tap=""): |
1450 | out="" |
1451 | z=[] |
1452 | if isinstance(b,float) or b.rank==0: |
1453 | z.append(b**2) |
1454 | elif b.rank==1: |
1455 | for i0 in range(b.shape[0]): |
1456 | z.append(b[i0]**2) |
1457 | elif b.rank==2: |
1458 | for i1 in range(b.shape[1]): |
1459 | s=0 |
1460 | for i0 in range(b.shape[0]): |
1461 | s+=b[i0,i1]**2 |
1462 | z.append(s) |
1463 | elif b.rank==3: |
1464 | for i2 in range(b.shape[2]): |
1465 | s=0 |
1466 | for i0 in range(b.shape[0]): |
1467 | for i1 in range(b.shape[1]): |
1468 | s+=b[i0,i1,i2]**2 |
1469 | z.append(s) |
1470 | |
1471 | elif b.rank==4: |
1472 | for i3 in range(b.shape[3]): |
1473 | s=0 |
1474 | for i0 in range(b.shape[0]): |
1475 | for i1 in range(b.shape[1]): |
1476 | for i2 in range(b.shape[2]): |
1477 | s+=b[i0,i1,i2,i3]**2 |
1478 | z.append(s) |
1479 | if where=="Function": |
1480 | xfac_o=1. |
1481 | xfac_op=0. |
1482 | z_fac_s="" |
1483 | zo_fac_s="" |
1484 | zo_fac=1./3. |
1485 | elif where=="FunctionOnBoundary": |
1486 | xfac_o=1. |
1487 | xfac_op=0. |
1488 | z_fac_s="*dim" |
1489 | zo_fac_s="*(2.*dim+1.)/3." |
1490 | zo_fac=1. |
1491 | elif where in ["FunctionOnContactZero","FunctionOnContactOne"]: |
1492 | xfac_o=0. |
1493 | xfac_op=1. |
1494 | z_fac_s="" |
1495 | zo_fac_s="" |
1496 | zo_fac=1./3. |
1497 | zo=0. |
1498 | zop=0. |
1499 | for i99 in range(len(z)): |
1500 | if i99==0: |
1501 | zo+=xfac_o*z[i99] |
1502 | zop+=xfac_op*z[i99] |
1503 | else: |
1504 | zo+=z[i99] |
1505 | out+=tap+"%s=sqrt((%s)%s"%(arg,zo*zo_fac,zo_fac_s) |
1506 | if zop==0.: |
1507 | out+=")\n" |
1508 | else: |
1509 | out+="+(%s))\n"%(zop*0.5**2) |
1510 | return out |
1511 | #======================================================================================================= |
1512 | # transpose |
1513 | #======================================================================================================= |
1514 | def transposeTest(r,offset): |
1515 | if isinstance(r,float): return r |
1516 | s=r.shape |
1517 | s1=1 |
1518 | for i in s[:offset]: s1*=i |
1519 | s2=1 |
1520 | for i in s[offset:]: s2*=i |
1521 | out=numarray.reshape(r,(s1,s2)) |
1522 | out.transpose() |
1523 | return numarray.resize(out,s[offset:]+s[:offset]) |
1524 | |
1525 | name,tt="transpose",transposeTest |
1526 | for case0 in ["array","Symbol","constData","taggedData","expandedData"]: |
1527 | for sh0 in [ (), (3,), (4,5), (6,2,2),(3,2,3,4)]: |
1528 | for offset in range(len(sh0)+1): |
1529 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1530 | tname="test_%s_%s_rank%s_offset%s"%(name,case0,len(sh0),offset) |
1531 | text+=" def %s(self):\n"%tname |
1532 | sh_t=sh0[offset:]+sh0[:offset] |
1533 | |
1534 | # sh_t=list(sh0) |
1535 | # sh_t[offset+1]=sh_t[offset] |
1536 | # sh_t=tuple(sh_t) |
1537 | # sh_r=[] |
1538 | # for i in range(offset): sh_r.append(sh0[i]) |
1539 | # for i in range(offset+2,len(sh0)): sh_r.append(sh0[i]) |
1540 | # sh_r=tuple(sh_r) |
1541 | |
1542 | a_0=makeArray(sh0,[-1.,1]) |
1543 | if case0 in ["taggedData", "expandedData"]: |
1544 | a1_0=makeArray(sh0,[-1.,1]) |
1545 | else: |
1546 | a1_0=a_0 |
1547 | r=tt(a_0,offset) |
1548 | r1=tt(a1_0,offset) |
1549 | text+=mkText(case0,"arg",a_0,a1_0) |
1550 | text+=" res=%s(arg,%s)\n"%(name,offset) |
1551 | if case0=="Symbol": |
1552 | text+=mkText("array","s",a_0,a1_0) |
1553 | text+=" sub=res.substitute({arg:s})\n" |
1554 | res="sub" |
1555 | text+=mkText("array","ref",r,r1) |
1556 | else: |
1557 | res="res" |
1558 | text+=mkText(case0,"ref",r,r1) |
1559 | text+=mkTypeAndShapeTest(case0,sh_t,"res") |
1560 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1561 | |
1562 | if case0 == "taggedData" : |
1563 | t_prog_with_tags+=text |
1564 | else: |
1565 | t_prog+=text |
1566 | |
1567 | print test_header |
1568 | # print t_prog |
1569 | print t_prog_with_tags |
1570 | print test_tail |
1571 | 1/0 |
1572 | #======================================================================================================= |
1573 | # L2 |
1574 | #======================================================================================================= |
1575 | for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1576 | for data in ["Data","Symbol"]: |
1577 | for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]: |
1578 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1579 | tname="test_L2_on%s_from%s_rank%s"%(where,data,len(sh)) |
1580 | text+=" def %s(self):\n"%tname |
1581 | text+=" \"\"\"\n" |
1582 | text+=" tests L2-norm of %s on the %s\n\n"%(data,where) |
1583 | text+=" assumptions: self.domain supports integration on %s\n"%where |
1584 | text+=" \"\"\"\n" |
1585 | text+=" dim=self.domain.getDim()\n" |
1586 | text+=" w=%s(self.domain)\n"%where |
1587 | text+=" x=w.getX()\n" |
1588 | o="1" |
1589 | if len(sh)>0: |
1590 | sh_2=sh[:len(sh)-1]+(2,) |
1591 | sh_3=sh[:len(sh)-1]+(3,) |
1592 | b_2=makeArray(sh[:len(sh)-1]+(2,),[-1.,1]) |
1593 | b_3=makeArray(sh[:len(sh)-1]+(3,),[-1.,1]) |
1594 | else: |
1595 | sh_2=() |
1596 | sh_3=() |
1597 | b_2=makeArray(sh,[-1.,1]) |
1598 | b_3=makeArray(sh,[-1.,1]) |
1599 | |
1600 | if data=="Symbol": |
1601 | val="s" |
1602 | res="sub" |
1603 | else: |
1604 | val="arg" |
1605 | res="res" |
1606 | text+=" if dim==2:\n" |
1607 | if data=="Symbol": |
1608 | text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh_2) |
1609 | |
1610 | text+=" %s=Data(0,%s,w)\n"%(val,sh_2) |
1611 | text+=unrollLoopsSimplified(b_2,val,tap=" ") |
1612 | text+=unrollLoopsOfL2(b_2,where,"ref",tap=" ") |
1613 | text+="\n else:\n" |
1614 | if data=="Symbol": |
1615 | text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh_3) |
1616 | text+=" %s=Data(0,%s,w)\n"%(val,sh_3) |
1617 | text+=unrollLoopsSimplified(b_3,val,tap=" ") |
1618 | text+=unrollLoopsOfL2(b_3,where,"ref",tap=" ") |
1619 | text+="\n res=L2(arg)\n" |
1620 | if data=="Symbol": |
1621 | text+=" sub=res.substitute({arg:s})\n" |
1622 | text+=" self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n" |
1623 | text+=" self.failUnlessEqual(res.getShape(),(),\"wrong shape of result.\")\n" |
1624 | else: |
1625 | text+=" self.failUnless(isinstance(res,float),\"wrong type of result.\")\n" |
1626 | text+=" self.failUnlessAlmostEqual(%s,ref,int(-log10(self.RES_TOL)),\"wrong result\")\n"%res |
1627 | t_prog+=text |
1628 | print t_prog |
1629 | 1/0 |
1630 | |
1631 | #======================================================================================================= |
1632 | # div |
1633 | #======================================================================================================= |
1634 | for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1635 | for data in ["Data","Symbol"]: |
1636 | for case in ["ContinuousFunction","Solution","ReducedSolution"]: |
1637 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1638 | tname="test_div_on%s_from%s_%s"%(where,data,case) |
1639 | text+=" def %s(self):\n"%tname |
1640 | text+=" \"\"\"\n" |
1641 | text+=" tests divergence of %s on the %s\n\n"%(data,where) |
1642 | text+=" assumptions: %s(self.domain) exists\n"%case |
1643 | text+=" self.domain supports div on %s\n"%where |
1644 | text+=" \"\"\"\n" |
1645 | if case=="ReducedSolution": |
1646 | text+=" o=1\n" |
1647 | o="1" |
1648 | else: |
1649 | text+=" o=self.order\n" |
1650 | o="o" |
1651 | text+=" dim=self.domain.getDim()\n" |
1652 | text+=" w_ref=%s(self.domain)\n"%where |
1653 | text+=" x_ref=w_ref.getX()\n" |
1654 | text+=" w=%s(self.domain)\n"%case |
1655 | text+=" x=w.getX()\n" |
1656 | a_2=makeArray((2,2),[-1.,1]) |
1657 | b_2=makeArray((2,2),[-1.,1]) |
1658 | a_3=makeArray((3,3),[-1.,1]) |
1659 | b_3=makeArray((3,3),[-1.,1]) |
1660 | if data=="Symbol": |
1661 | text+=" arg=Symbol(shape=(dim,),dim=dim)\n" |
1662 | val="s" |
1663 | res="sub" |
1664 | else: |
1665 | val="arg" |
1666 | res="res" |
1667 | text+=" %s=Vector(0,w)\n"%val |
1668 | text+=" if dim==2:\n" |
1669 | text+=unrollLoops(a_2,b_2,o,val,tap=" ") |
1670 | text+=unrollLoopsOfDiv(a_2,b_2,o,"ref",tap=" ") |
1671 | text+="\n else:\n" |
1672 | |
1673 | text+=unrollLoops(a_3,b_3,o,val,tap=" ") |
1674 | text+=unrollLoopsOfDiv(a_3,b_3,o,"ref",tap=" ") |
1675 | text+="\n res=div(arg,where=w_ref)\n" |
1676 | if data=="Symbol": |
1677 | text+=" sub=res.substitute({arg:s})\n" |
1678 | text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data |
1679 | text+=" self.failUnlessEqual(res.getShape(),(),\"wrong shape of result.\")\n" |
1680 | text+=" self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong function space of result.\")\n"%res |
1681 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1682 | |
1683 | |
1684 | t_prog+=text |
1685 | print t_prog |
1686 | 1/0 |
1687 | |
1688 | #======================================================================================================= |
1689 | # interpolation |
1690 | #======================================================================================================= |
1691 | for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne","Solution","ReducedSolution"]: |
1692 | for data in ["Data","Symbol"]: |
1693 | for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1694 | for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]: |
1695 | if where==case or \ |
1696 | ( case in ["ContinuousFunction","Solution","ReducedSolution"] and where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"] ) or \ |
1697 | ( case in ["FunctionOnContactZero","FunctionOnContactOne"] and where in ["FunctionOnContactZero","FunctionOnContactOne"] ) or \ |
1698 | (case=="ContinuousFunction" and where in ["Solution","ReducedSolution"]) or \ |
1699 | (case=="Solution" and where=="ReducedSolution") : |
1700 | |
1701 | |
1702 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1703 | tname="test_interpolation_on%s_from%s_%s_rank%s"%(where,data,case,len(sh)) |
1704 | text+=" def %s(self):\n"%tname |
1705 | text+=" \"\"\"\n" |
1706 | text+=" tests interpolation for rank %s %s onto the %s\n\n"%(len(sh),data,where) |
1707 | text+=" assumptions: self.domain supports inpterpolation from %s onto %s\n"%(case,where) |
1708 | text+=" \"\"\"\n" |
1709 | if case=="ReducedSolution" or where=="ReducedSolution": |
1710 | text+=" o=1\n" |
1711 | o="1" |
1712 | else: |
1713 | text+=" o=self.order\n" |
1714 | o="o" |
1715 | text+=" dim=self.domain.getDim()\n" |
1716 | text+=" w_ref=%s(self.domain)\n"%where |
1717 | text+=" x_ref=w_ref.getX()\n" |
1718 | text+=" w=%s(self.domain)\n"%case |
1719 | text+=" x=w.getX()\n" |
1720 | a_2=makeArray(sh+(2,),[-1.,1]) |
1721 | b_2=makeArray(sh+(2,),[-1.,1]) |
1722 | a_3=makeArray(sh+(3,),[-1.,1]) |
1723 | b_3=makeArray(sh+(3,),[-1.,1]) |
1724 | if data=="Symbol": |
1725 | text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh) |
1726 | val="s" |
1727 | res="sub" |
1728 | else: |
1729 | val="arg" |
1730 | res="res" |
1731 | text+=" %s=Data(0,%s,w)\n"%(val,str(sh)) |
1732 | text+=" ref=Data(0,%s,w_ref)\n"%str(sh) |
1733 | text+=" if dim==2:\n" |
1734 | text+=unrollLoops(a_2,b_2,o,val,tap=" ") |
1735 | text+=unrollLoops(a_2,b_2,o,"ref",tap=" ",x="x_ref") |
1736 | text+=" else:\n" |
1737 | |
1738 | text+=unrollLoops(a_3,b_3,o,val,tap=" ") |
1739 | text+=unrollLoops(a_3,b_3,o,"ref",tap=" ",x="x_ref") |
1740 | text+=" res=interpolate(arg,where=w_ref)\n" |
1741 | if data=="Symbol": |
1742 | text+=" sub=res.substitute({arg:s})\n" |
1743 | text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data |
1744 | text+=" self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong functionspace of result.\")\n"%res |
1745 | text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh) |
1746 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1747 | t_prog+=text |
1748 | print test_header |
1749 | print t_prog |
1750 | print test_tail |
1751 | 1/0 |
1752 | |
1753 | #======================================================================================================= |
1754 | # grad |
1755 | #======================================================================================================= |
1756 | for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1757 | for data in ["Data","Symbol"]: |
1758 | for case in ["ContinuousFunction","Solution","ReducedSolution"]: |
1759 | for sh in [ (),(2,), (4,5), (6,2,2)]: |
1760 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1761 | tname="test_grad_on%s_from%s_%s_rank%s"%(where,data,case,len(sh)) |
1762 | text+=" def %s(self):\n"%tname |
1763 | text+=" \"\"\"\n" |
1764 | text+=" tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where) |
1765 | text+=" assumptions: %s(self.domain) exists\n"%case |
1766 | text+=" self.domain supports gardient on %s\n"%where |
1767 | text+=" \"\"\"\n" |
1768 | if case=="ReducedSolution": |
1769 | text+=" o=1\n" |
1770 | o="1" |
1771 | else: |
1772 | text+=" o=self.order\n" |
1773 | o="o" |
1774 | text+=" dim=self.domain.getDim()\n" |
1775 | text+=" w_ref=%s(self.domain)\n"%where |
1776 | text+=" x_ref=w_ref.getX()\n" |
1777 | text+=" w=%s(self.domain)\n"%case |
1778 | text+=" x=w.getX()\n" |
1779 | a_2=makeArray(sh+(2,),[-1.,1]) |
1780 | b_2=makeArray(sh+(2,),[-1.,1]) |
1781 | a_3=makeArray(sh+(3,),[-1.,1]) |
1782 | b_3=makeArray(sh+(3,),[-1.,1]) |
1783 | if data=="Symbol": |
1784 | text+=" arg=Symbol(shape=%s,dim=dim)\n"%str(sh) |
1785 | val="s" |
1786 | res="sub" |
1787 | else: |
1788 | val="arg" |
1789 | res="res" |
1790 | text+=" %s=Data(0,%s,w)\n"%(val,str(sh)) |
1791 | text+=" ref=Data(0,%s+(dim,),w_ref)\n"%str(sh) |
1792 | text+=" if dim==2:\n" |
1793 | text+=unrollLoops(a_2,b_2,o,val,tap=" ") |
1794 | text+=unrollLoopsOfGrad(a_2,b_2,o,"ref",tap=" ") |
1795 | text+=" else:\n" |
1796 | |
1797 | text+=unrollLoops(a_3,b_3,o,val,tap=" ") |
1798 | text+=unrollLoopsOfGrad(a_3,b_3,o,"ref",tap=" ") |
1799 | text+=" res=grad(arg,where=w_ref)\n" |
1800 | if data=="Symbol": |
1801 | text+=" sub=res.substitute({arg:s})\n" |
1802 | text+=" self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data |
1803 | text+=" self.failUnlessEqual(res.getShape(),%s+(dim,),\"wrong shape of result.\")\n"%str(sh) |
1804 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1805 | |
1806 | |
1807 | t_prog+=text |
1808 | print test_header |
1809 | print t_prog |
1810 | print test_tail |
1811 | 1/0 |
1812 | |
1813 | |
1814 | #======================================================================================================= |
1815 | # integrate |
1816 | #======================================================================================================= |
1817 | for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1818 | for data in ["Data","Symbol"]: |
1819 | for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]: |
1820 | for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]: |
1821 | if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case: |
1822 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1823 | tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh)) |
1824 | text+=" def %s(self):\n"%tname |
1825 | text+=" \"\"\"\n" |
1826 | text+=" tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where) |
1827 | text+=" assumptions: %s(self.domain) exists\n"%case |
1828 | text+=" self.domain supports integral on %s\n"%where |
1829 | |
1830 | text+=" \"\"\"\n" |
1831 | if case=="ReducedSolution": |
1832 | text+=" o=1\n" |
1833 | o="1" |
1834 | else: |
1835 | text+=" o=self.order\n" |
1836 | o="o" |
1837 | text+=" dim=self.domain.getDim()\n" |
1838 | text+=" w_ref=%s(self.domain)\n"%where |
1839 | text+=" w=%s(self.domain)\n"%case |
1840 | text+=" x=w.getX()\n" |
1841 | a_2=makeArray(sh+(2,),[-1.,1]) |
1842 | b_2=makeArray(sh+(2,),[-1.,1]) |
1843 | a_3=makeArray(sh+(3,),[-1.,1]) |
1844 | b_3=makeArray(sh+(3,),[-1.,1]) |
1845 | if data=="Symbol": |
1846 | text+=" arg=Symbol(shape=%s)\n"%str(sh) |
1847 | val="s" |
1848 | res="sub" |
1849 | else: |
1850 | val="arg" |
1851 | res="res" |
1852 | |
1853 | text+=" %s=Data(0,%s,w)\n"%(val,str(sh)) |
1854 | if not len(sh)==0: |
1855 | text+=" ref=numarray.zeros(%s,numarray.Float)\n"%str(sh) |
1856 | text+=" if dim==2:\n" |
1857 | text+=unrollLoops(a_2,b_2,o,val,tap=" ") |
1858 | text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap=" ") |
1859 | text+=" else:\n" |
1860 | |
1861 | text+=unrollLoops(a_3,b_3,o,val,tap=" ") |
1862 | text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap=" ") |
1863 | if case in ["ContinuousFunction","Solution","ReducedSolution"]: |
1864 | text+=" res=integrate(arg,where=w_ref)\n" |
1865 | else: |
1866 | text+=" res=integrate(arg)\n" |
1867 | |
1868 | if data=="Symbol": |
1869 | text+=" sub=res.substitute({arg:s})\n" |
1870 | if len(sh)==0 and data=="Data": |
1871 | text+=" self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res |
1872 | else: |
1873 | if data=="Symbol": |
1874 | text+=" self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n" |
1875 | text+=" self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh) |
1876 | else: |
1877 | text+=" self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n" |
1878 | text+=" self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh) |
1879 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1880 | |
1881 | |
1882 | t_prog+=text |
1883 | print test_header |
1884 | print t_prog |
1885 | print test_tail |
1886 | 1/0 |
1887 | #======================================================================================================= |
1888 | # inverse |
1889 | #======================================================================================================= |
1890 | name="inverse" |
1891 | for case0 in ["array","Symbol","constData","taggedData","expandedData"]: |
1892 | for sh0 in [ (1,1), (2,2), (3,3)]: |
1893 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1894 | tname="test_%s_%s_dim%s"%(name,case0,sh0[0]) |
1895 | text+=" def %s(self):\n"%tname |
1896 | a_0=makeArray(sh0,[-1.,1]) |
1897 | for i in range(sh0[0]): a_0[i,i]+=2 |
1898 | if case0 in ["taggedData", "expandedData"]: |
1899 | a1_0=makeArray(sh0,[-1.,1]) |
1900 | for i in range(sh0[0]): a1_0[i,i]+=3 |
1901 | else: |
1902 | a1_0=a_0 |
1903 | |
1904 | text+=mkText(case0,"arg",a_0,a1_0) |
1905 | text+=" res=%s(arg)\n"%name |
1906 | if case0=="Symbol": |
1907 | text+=mkText("array","s",a_0,a1_0) |
1908 | text+=" sub=res.substitute({arg:s})\n" |
1909 | res="sub" |
1910 | ref="s" |
1911 | else: |
1912 | ref="arg" |
1913 | res="res" |
1914 | text+=mkTypeAndShapeTest(case0,sh0,"res") |
1915 | text+=" self.failUnless(Lsup(matrixmult(%s,%s)-kronecker(%s))<=self.RES_TOL,\"wrong result\")\n"%(res,ref,sh0[0]) |
1916 | |
1917 | if case0 == "taggedData" : |
1918 | t_prog_with_tags+=text |
1919 | else: |
1920 | t_prog+=text |
1921 | |
1922 | print test_header |
1923 | # print t_prog |
1924 | print t_prog_with_tags |
1925 | print test_tail |
1926 | 1/0 |
1927 | |
1928 | #======================================================================================================= |
1929 | # trace |
1930 | #======================================================================================================= |
1931 | def traceTest(r,offset): |
1932 | sh=r.shape |
1933 | r1=1 |
1934 | for i in range(offset): r1*=sh[i] |
1935 | r2=1 |
1936 | for i in range(offset+2,len(sh)): r2*=sh[i] |
1937 | r_s=numarray.reshape(r,(r1,sh[offset],sh[offset],r2)) |
1938 | s=numarray.zeros([r1,r2],numarray.Float) |
1939 | for i1 in range(r1): |
1940 | for i2 in range(r2): |
1941 | for j in range(sh[offset]): s[i1,i2]+=r_s[i1,j,j,i2] |
1942 | return s.resize(sh[:offset]+sh[offset+2:]) |
1943 | name,tt="trace",traceTest |
1944 | for case0 in ["array","Symbol","constData","taggedData","expandedData"]: |
1945 | for sh0 in [ (4,5), (6,2,2),(3,2,3,4)]: |
1946 | for offset in range(len(sh0)-1): |
1947 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1948 | tname="test_%s_%s_rank%s_offset%s"%(name,case0,len(sh0),offset) |
1949 | text+=" def %s(self):\n"%tname |
1950 | sh_t=list(sh0) |
1951 | sh_t[offset+1]=sh_t[offset] |
1952 | sh_t=tuple(sh_t) |
1953 | sh_r=[] |
1954 | for i in range(offset): sh_r.append(sh0[i]) |
1955 | for i in range(offset+2,len(sh0)): sh_r.append(sh0[i]) |
1956 | sh_r=tuple(sh_r) |
1957 | a_0=makeArray(sh_t,[-1.,1]) |
1958 | if case0 in ["taggedData", "expandedData"]: |
1959 | a1_0=makeArray(sh_t,[-1.,1]) |
1960 | else: |
1961 | a1_0=a_0 |
1962 | r=tt(a_0,offset) |
1963 | r1=tt(a1_0,offset) |
1964 | text+=mkText(case0,"arg",a_0,a1_0) |
1965 | text+=" res=%s(arg,%s)\n"%(name,offset) |
1966 | if case0=="Symbol": |
1967 | text+=mkText("array","s",a_0,a1_0) |
1968 | text+=" sub=res.substitute({arg:s})\n" |
1969 | res="sub" |
1970 | text+=mkText("array","ref",r,r1) |
1971 | else: |
1972 | res="res" |
1973 | text+=mkText(case0,"ref",r,r1) |
1974 | text+=mkTypeAndShapeTest(case0,sh_r,"res") |
1975 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
1976 | |
1977 | if case0 == "taggedData" : |
1978 | t_prog_with_tags+=text |
1979 | else: |
1980 | t_prog+=text |
1981 | |
1982 | print test_header |
1983 | # print t_prog |
1984 | print t_prog_with_tags |
1985 | print test_tail |
1986 | 1/0 |
1987 | |
1988 | #======================================================================================================= |
1989 | # clip |
1990 | #======================================================================================================= |
1991 | oper_L=[["clip",clipTEST]] |
1992 | for oper in oper_L: |
1993 | for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
1994 | for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
1995 | if len(sh0)==0 or not case0=="float": |
1996 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
1997 | tname="test_%s_%s_rank%s"%(oper[0],case0,len(sh0)) |
1998 | text+=" def %s(self):\n"%tname |
1999 | a_0=makeArray(sh0,[-1.,1]) |
2000 | if case0 in ["taggedData", "expandedData"]: |
2001 | a1_0=makeArray(sh0,[-1.,1]) |
2002 | else: |
2003 | a1_0=a_0 |
2004 | |
2005 | r=oper[1](a_0,-0.3,0.5) |
2006 | r1=oper[1](a1_0,-0.3,0.5) |
2007 | text+=mkText(case0,"arg",a_0,a1_0) |
2008 | text+=" res=%s(arg,-0.3,0.5)\n"%oper[0] |
2009 | if case0=="Symbol": |
2010 | text+=mkText("array","s",a_0,a1_0) |
2011 | text+=" sub=res.substitute({arg:s})\n" |
2012 | res="sub" |
2013 | text+=mkText("array","ref",r,r1) |
2014 | else: |
2015 | res="res" |
2016 | text+=mkText(case0,"ref",r,r1) |
2017 | text+=mkTypeAndShapeTest(case0,sh0,"res") |
2018 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
2019 | |
2020 | if case0 == "taggedData" : |
2021 | t_prog_with_tags+=text |
2022 | else: |
2023 | t_prog+=text |
2024 | |
2025 | print test_header |
2026 | # print t_prog |
2027 | print t_prog_with_tags |
2028 | print test_tail |
2029 | 1/0 |
2030 | |
2031 | #======================================================================================================= |
2032 | # maximum, minimum, clipping |
2033 | #======================================================================================================= |
2034 | oper_L=[ ["maximum",maximumTEST], |
2035 | ["minimum",minimumTEST]] |
2036 | for oper in oper_L: |
2037 | for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2038 | for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2039 | for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2040 | for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2041 | if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \ |
2042 | and (sh0==sh1 or len(sh0)==0 or len(sh1)==0) : |
2043 | use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData" |
2044 | |
2045 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
2046 | tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1)) |
2047 | text+=" def %s(self):\n"%tname |
2048 | a_0=makeArray(sh0,[-1.,1]) |
2049 | if case0 in ["taggedData", "expandedData"]: |
2050 | a1_0=makeArray(sh0,[-1.,1]) |
2051 | else: |
2052 | a1_0=a_0 |
2053 | |
2054 | a_1=makeArray(sh1,[-1.,1]) |
2055 | if case1 in ["taggedData", "expandedData"]: |
2056 | a1_1=makeArray(sh1,[-1.,1]) |
2057 | else: |
2058 | a1_1=a_1 |
2059 | r=oper[1](a_0,a_1) |
2060 | r1=oper[1](a1_0,a1_1) |
2061 | text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data) |
2062 | text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data) |
2063 | text+=" res=%s(arg0,arg1)\n"%oper[0] |
2064 | case=getResultCaseForBin(case0,case1) |
2065 | if case=="Symbol": |
2066 | c0_res,c1_res=case0,case1 |
2067 | subs="{" |
2068 | if case0=="Symbol": |
2069 | text+=mkText("array","s0",a_0,a1_0) |
2070 | subs+="arg0:s0" |
2071 | c0_res="array" |
2072 | if case1=="Symbol": |
2073 | text+=mkText("array","s1",a_1,a1_1) |
2074 | if not subs.endswith("{"): subs+="," |
2075 | subs+="arg1:s1" |
2076 | c1_res="array" |
2077 | subs+="}" |
2078 | text+=" sub=res.substitute(%s)\n"%subs |
2079 | res="sub" |
2080 | text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1) |
2081 | else: |
2082 | res="res" |
2083 | text+=mkText(case,"ref",r,r1) |
2084 | if len(sh0)>len(sh1): |
2085 | text+=mkTypeAndShapeTest(case,sh0,"res") |
2086 | else: |
2087 | text+=mkTypeAndShapeTest(case,sh1,"res") |
2088 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
2089 | |
2090 | if case0 == "taggedData" or case1 == "taggedData": |
2091 | t_prog_with_tags+=text |
2092 | else: |
2093 | t_prog+=text |
2094 | |
2095 | print test_header |
2096 | # print t_prog |
2097 | print t_prog_with_tags |
2098 | print test_tail |
2099 | 1/0 |
2100 | |
2101 | |
2102 | #======================================================================================================= |
2103 | # outer inner |
2104 | #======================================================================================================= |
2105 | oper=["outer",outerTEST] |
2106 | # oper=["inner",innerTEST] |
2107 | for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2108 | for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2109 | for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2110 | for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2111 | if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \ |
2112 | and len(sh0+sh1)<5: |
2113 | use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData" |
2114 | |
2115 | text=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
2116 | tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1)) |
2117 | text+=" def %s(self):\n"%tname |
2118 | a_0=makeArray(sh0,[-1.,1]) |
2119 | if case0 in ["taggedData", "expandedData"]: |
2120 | a1_0=makeArray(sh0,[-1.,1]) |
2121 | else: |
2122 | a1_0=a_0 |
2123 | |
2124 | a_1=makeArray(sh1,[-1.,1]) |
2125 | if case1 in ["taggedData", "expandedData"]: |
2126 | a1_1=makeArray(sh1,[-1.,1]) |
2127 | else: |
2128 | a1_1=a_1 |
2129 | r=oper[1](a_0,a_1) |
2130 | r1=oper[1](a1_0,a1_1) |
2131 | text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data) |
2132 | text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data) |
2133 | text+=" res=%s(arg0,arg1)\n"%oper[0] |
2134 | case=getResultCaseForBin(case0,case1) |
2135 | if case=="Symbol": |
2136 | c0_res,c1_res=case0,case1 |
2137 | subs="{" |
2138 | if case0=="Symbol": |
2139 | text+=mkText("array","s0",a_0,a1_0) |
2140 | subs+="arg0:s0" |
2141 | c0_res="array" |
2142 | if case1=="Symbol": |
2143 | text+=mkText("array","s1",a_1,a1_1) |
2144 | if not subs.endswith("{"): subs+="," |
2145 | subs+="arg1:s1" |
2146 | c1_res="array" |
2147 | subs+="}" |
2148 | text+=" sub=res.substitute(%s)\n"%subs |
2149 | res="sub" |
2150 | text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1) |
2151 | else: |
2152 | res="res" |
2153 | text+=mkText(case,"ref",r,r1) |
2154 | text+=mkTypeAndShapeTest(case,sh0+sh1,"res") |
2155 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
2156 | |
2157 | if case0 == "taggedData" or case1 == "taggedData": |
2158 | t_prog_with_tags+=text |
2159 | else: |
2160 | t_prog+=text |
2161 | |
2162 | print test_header |
2163 | # print t_prog |
2164 | print t_prog_with_tags |
2165 | print test_tail |
2166 | 1/0 |
2167 | |
2168 | #======================================================================================================= |
2169 | # local reduction |
2170 | #======================================================================================================= |
2171 | for oper in [["length",0.,"out+%a1%**2","math.sqrt(out)"], |
2172 | ["maxval",-1.e99,"max(out,%a1%)","out"], |
2173 | ["minval",1.e99,"min(out,%a1%)","out"] ]: |
2174 | for case in case_set: |
2175 | for sh in shape_set: |
2176 | if not case=="float" or len(sh)==0: |
2177 | text="" |
2178 | text+=" #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n" |
2179 | tname="def test_%s_%s_rank%s"%(oper[0],case,len(sh)) |
2180 | text+=" %s(self):\n"%tname |
2181 | a=makeArray(sh,[-1.,1.]) |
2182 | a1=makeArray(sh,[-1.,1.]) |
2183 | r1=testReduce(a1,oper[1],oper[2],oper[3]) |
2184 | r=testReduce(a,oper[1],oper[2],oper[3]) |
2185 | |
2186 | text+=mkText(case,"arg",a,a1) |
2187 | text+=" res=%s(arg)\n"%oper[0] |
2188 | if case=="Symbol": |
2189 | text+=mkText("array","s",a,a1) |
2190 | text+=" sub=res.substitute({arg:s})\n" |
2191 | text+=mkText("array","ref",r,r1) |
2192 | res="sub" |
2193 | else: |
2194 | text+=mkText(case,"ref",r,r1) |
2195 | res="res" |
2196 | if oper[0]=="length": |
2197 | text+=mkTypeAndShapeTest(case,(),"res") |
2198 | else: |
2199 | if case=="float" or case=="array": |
2200 | text+=mkTypeAndShapeTest("float",(),"res") |
2201 | else: |
2202 | text+=mkTypeAndShapeTest(case,(),"res") |
2203 | text+=" self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res |
2204 | if case == "taggedData": |
2205 | t_prog_with_tags+=text |
2206 | else: |
2207 | t_prog+=text |
2208 | print test_header |
2209 | # print t_prog |
2210 | print t_prog_with_tags |
2211 | print test_tail |
2212 | 1/0 |
2213 | |
2214 | #======================================================================================================= |
2215 | # tensor multiply |
2216 | #======================================================================================================= |
2217 | # oper=["generalTensorProduct",tensorProductTest] |
2218 | # oper=["matrixmult",testMatrixMult] |
2219 | oper=["tensormult",testTensorMult] |
2220 | |
2221 | for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2222 | for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2223 | for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]: |
2224 | for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]: |
2225 | for sh_s in [ (),(3,), (2,3), (2,4,3),(4,2,3,2)]: |
2226 | if (len(sh0+sh_s)==0 or not case0=="float") and (len(sh1+sh_s)==0 or not case1=="float") \ |
2227 | and len(sh0+sh1)<5 and len(sh0+sh_s)<5 and len(sh1+sh_s)<5: |
2228 | # 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 |
2229 | 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 |
2230 | case=getResultCaseForBin(case0,case1) |
2231 |