/[escript]/trunk/finley/test/python/gentest
ViewVC logotype

Annotation of /trunk/finley/test/python/gentest

Parent Directory Parent Directory | Revision Log Revision Log


Revision 766 - (hide annotations)
Fri Jun 30 06:44:14 2006 UTC (13 years, 7 months ago) by gross
File size: 31719 byte(s)
modified version of a test generator (private use)
1 gross 589 #!/usr/bin/python
2 gross 766 import numarray
3     import random
4 gross 589 t_prog={}
5 gross 766 t_prog["2Do1_strong"]=""
6     t_prog["2Do2_strong"]=""
7     t_prog["3Do1_strong"]=""
8     t_prog["3Do2_strong"]=""
9     t_prog["2Do1_weak"]=""
10     t_prog["2Do2_weak"]=""
11     t_prog["3Do1_weak"]=""
12     t_prog["3Do2_weak"]=""
13 gross 589 t_prog["2Do1_contact"]=""
14     t_prog["2Do2_contact"]=""
15     t_prog["3Do1_contact"]=""
16     t_prog["3Do2_contact"]=""
17     intend=" "
18 gross 766 max_order=3
19 gross 589
20 gross 766 def makeTitle(d,coeffo,o,s,coef,case,typ,body="",mark="",pdeargs="",add_jump=False):
21 gross 589 if mark=="":
22 gross 766 name="test_assemblage_%sD_solO%s_coeffO%s_neq%s_%s_%s_type%s"%(d,o,coeffo,s,coef,case,typ)
23 gross 589 else:
24 gross 766 name="test_assemblage_%sD_solO%s_coeffO%s_neq%s_%s_%s_type%s_comp%s"%(d,o,coeffo,s,coef,case,typ,mark)
25     key="%sDo%s_%s"%(d,o,typ.lower())
26     t_prog[key]+=intend+"#"+50*"="+"\n"+intend+"def %s(self):\n"%name
27     t_prog[key]+=2*intend+"x=self.domain.getX()\n"
28     if add_jump:
29     t_prog[key]+=2*intend+"jump=Data(0.,(),ContinuousFunction(self.domain))\n"
30     t_prog[key]+=2*intend+"jump.setTaggedValue(2,1.)\n"
31     t_prog[key]+=body
32     t_prog[key]+=2*intend+"pde=LinearPDE(self.domain)\n"
33     t_prog[key]+=2*intend+"pde.setValue(%s)\n"%pdeargs
34     t_prog[key]+=2*intend+"r=pde.getResidual(u)\n"
35     t_prog[key]+=2*intend+"rhs=pde.getRightHandSide()\n"
36     t_prog[key]+=2*intend+"self.failUnless(Lsup(rhs)>0,\"right hand side is zero\")\n"
37     t_prog[key]+=2*intend+"self.failUnless(Lsup(r)<=self.RES_TOL*Lsup(rhs),\"residual is too big\")\n"
38    
39     def getRandom():
40     out=0
41     while out==0:
42     out = int(18*random.random())-9
43     return out
44    
45     def makeTestSolution(order,ncomps,dim):
46     if dim == 2:
47     out=numarray.zeros((ncomps,max_order+1,max_order+1))
48     for d1 in range(order+1):
49     for d2 in range(order+1):
50     for k in range(ncomps):
51     out[k,d1,d2]=getRandom()
52     if dim == 3:
53     out=numarray.zeros((ncomps,max_order+1,max_order+1,max_order+1))
54     for d1 in range(order+1):
55     for d2 in range(order+1):
56     for d3 in range(order+1):
57     for k in range(ncomps):
58     out[k,d1,d2,d3]=getRandom()
59     return out
60    
61     def makeGradient(func,dim):
62     sh=func.shape
63     if len(sh)==dim+1:
64     sh_out=sh[:1]+(dim,)+sh[1:]
65     out=numarray.zeros(sh_out)
66     if dim==2:
67     for d1 in range(sh[1]-1):
68     for d2 in range(sh[2]-1):
69     for k in range(sh[0]):
70     out[k,0,d1,d2]=(d1+1)*func[k,d1+1,d2]
71     out[k,1,d1,d2]=(d2+1)*func[k,d1,d2+1]
72     else:
73     for d1 in range(sh[1]-1):
74     for d2 in range(sh[2]-1):
75     for d3 in range(sh[3]-1):
76     for k in range(sh[0]):
77     out[k,0,d1,d2,d3]=(d1+1)*func[k,d1+1,d2,d3]
78     out[k,1,d1,d2,d3]=(d2+1)*func[k,d1,d2+1,d3]
79     out[k,2,d1,d2,d3]=(d3+1)*func[k,d1,d2,d3+1]
80     elif len(sh)==dim+2:
81     raise RunTimeError,"sdfadsa"
82     sh_out=sh[:2]+(dim,)+sh[2:]
83     out=numarray.zeros(sh_out)
84 gross 589 else:
85 gross 766 raise RunTimeError,"sdfadsa"
86     return out
87 gross 589
88 gross 766 def makeDiv(func,dim):
89     sh=func.shape
90     if len(sh)==dim+1:
91     sh_out=(1,)+sh[1:]
92     out=numarray.zeros(sh_out)
93     if dim==2:
94     for d1 in range(sh[1]-1):
95     for d2 in range(sh[2]-1):
96     for k in range(sh[0]):
97     out[0,d1,d2]=(d1+1)*func[0,d1+1,d2]+(d2+1)*func[1,d1,d2+1]
98     else:
99     for d1 in range(sh[1]-1):
100     for d2 in range(sh[2]-1):
101     for d3 in range(sh[3]-1):
102     for k in range(sh[0]):
103     out[0,d1,d2,d3]=(d1+1)*func[0,d1+1,d2,d3]+(d2+1)*func[1,d1,d2+1,d3]+(d3+1)*func[2,d1,d2,d3+1]
104     elif len(sh)==dim+2:
105     sh_out=(sh[0],)+sh[2:]
106     out=numarray.zeros(sh_out)
107     if dim==2:
108     for d1 in range(sh[2]-1):
109     for d2 in range(sh[3]-1):
110     for k in range(sh[0]):
111     out[k,d1,d2]=(d1+1)*func[k,0,d1+1,d2]+(d2+1)*func[k,1,d1,d2+1]
112     else:
113     for d1 in range(sh[2]-1):
114     for d2 in range(sh[3]-1):
115     for d3 in range(sh[4]-1):
116     for k in range(sh[0]):
117     out[k,d1,d2,d3]=(d1+1)*func[k,0,d1+1,d2,d3]+(d2+1)*func[k,1,d1,d2+1,d3]+(d3+1)*func[k,2,d1,d2,d3+1]
118     else:
119     raise RunTimeError,"sdfadsa"
120     return out
121    
122     def multiplyByFandX(f,n,func,in_index,out_shape,out_index):
123     if n==0:
124     n1=1
125     else:
126     n1=0
127     if n==1:
128     n2=1
129     else:
130     n2=0
131     if n==2:
132     n3=1
133     else:
134     n3=0
135    
136     sh=func.shape
137     dim=len(sh)-len(in_index)
138     sh_out=out_shape+sh[len(in_index):]
139     out=numarray.zeros(sh_out)
140     if len(in_index)==1:
141     if len(out_index)==1:
142     if dim==2:
143     for d1 in range(n1,sh[1]):
144     for d2 in range(n2,sh[2]):
145     out[out_index[0],d1,d2]=f*func[in_index[0],d1-n1,d2-n2]
146    
147     else:
148     for d1 in range(n1,sh[1]):
149     for d2 in range(n2,sh[2]):
150     for d3 in range(n3,sh[3]):
151     out[out_index[0],d1,d2,d3]=f*func[in_index[0],d1-n1,d2-n2,d3-n3]
152    
153     elif len(out_index)==2:
154     if dim==2:
155     for d1 in range(n1,sh[1]):
156     for d2 in range(n2,sh[2]):
157     out[out_index[0],out_index[1],d1,d2]=f*func[in_index[0],d1-n1,d2-n2]
158     else:
159     for d1 in range(n1,sh[1]):
160     for d2 in range(n2,sh[2]):
161     for d3 in range(n3,sh[3]):
162     out[out_index[0],out_index[1],d1,d2,d3]=f*func[in_index[0],d1-n1,d2-n2,d3-n3]
163    
164     else:
165     raise RunTimeError,"sdfadsa"
166    
167     elif len(in_index)==2:
168     if len(out_index)==1:
169     if dim==2:
170     for d1 in range(n1,sh[2]):
171     for d2 in range(n2,sh[3]):
172     out[out_index[0],d1,d2]=f*func[in_index[0],in_index[1],d1-n1,d2-n2]
173     else:
174     for d1 in range(n1,sh[2]):
175     for d2 in range(n2,sh[3]):
176     for d3 in range(n3,sh[4]):
177     out[out_index[0],d1,d2,d3]=f*func[in_index[0],in_index[1],d1-n1,d2-n2,d3-n3]
178    
179     elif len(out_index)==2:
180     if dim==2:
181     for d1 in range(n1,sh[2]):
182     for d2 in range(n2,sh[3]):
183     out[out_index[0],out_index[1],d1,d2]=f*func[in_index[0],in_index[1],d1-n1,d2-n2]
184     else:
185     for d1 in range(n1,sh[2]):
186     for d2 in range(n2,sh[3]):
187     for d3 in range(n3,sh[4]):
188     out[out_index[0],out_index[1],d1,d2,d3]=f*func[in_index[0],in_index[1],d1-n1,d2-n2,d3-n3]
189    
190     else:
191     raise RunTimeError,"sdfadsa"
192     else:
193     raise RunTimeError,"sdfadsa"
194     return out
195    
196     def powTextforX(c,e):
197     if e==0:
198     out=""
199     elif e==1:
200     out="x[%s]"%c
201     else:
202     out="x[%s]**%s"%(c,e)
203     return out
204     def multTextFactors(*txt):
205     out=""
206     for t in txt:
207     if len(t)>0:
208     if len(out)==0:
209     out=t
210     else:
211     out+="*"+t
212     return out
213    
214     def addSumTermText(txt,f,fac):
215     if len(fac)>0:
216     fac2="*"+fac
217     if f>0:
218     if f==1:
219     if len(txt)==0:
220     out="%s"%fac
221     else:
222     out="+%s"%fac
223     else:
224     if len(txt)==0:
225     out="%s%s"%(f,fac2)
226     else:
227     out="+%s%s"%(f,fac2)
228     elif f<0:
229     if f==1:
230     out="-%s"%fac
231     else:
232     if len(txt)==0:
233     out="(-%s)%s"%(-f,fac2)
234     else:
235     out="-%s%s"%(-f,fac2)
236     else:
237     out=""
238     else:
239     if f>0:
240     if f==1:
241     if len(txt)==0:
242     out="1"
243     else:
244     out="+1"
245     else:
246     if len(txt)==0:
247     out="%s"%f
248     else:
249     out="+%s"%f
250     elif f<0:
251     if f==1:
252     out="-1"
253     else:
254     if len(txt)==0:
255     out="(-%s)"%(-f,)
256     else:
257     out="-%s"%(-f,)
258     else:
259     out=""
260    
261     return txt+out
262    
263     def makeTermText(func):
264     sh=func.shape
265     out=""
266     if len(sh)==2:
267     for d1 in range(sh[0]):
268     for d2 in range(sh[1]):
269     fac=multTextFactors(powTextforX(0,d1),powTextforX(1,d2))
270     out=addSumTermText(out,func[d1,d2],fac)
271     else:
272     for d1 in range(sh[0]):
273     for d2 in range(sh[1]):
274     for d3 in range(sh[2]):
275     fac=multTextFactors(powTextforX(0,d1),powTextforX(1,d2),powTextforX(2,d3))
276     out=addSumTermText(out,func[d1,d2,d3],fac)
277     if len(out)==0:
278     return "0"
279     else:
280     return out
281     def makeContactDerivativeText(func,dim,func_name):
282     sh=func.shape
283     out=""
284     if len(sh)==dim+2:
285     if sh[0]==1:
286     out2=""
287     for j in range(sh[1]):
288     t=makeTermText(func[0,j])
289     if not t == "0":
290     if len(out2)>0: out2+="+"
291     if j==0: out2+=t
292     if len(out2)==0:
293     out+=2*intend+"%s=%s\n"%(func_name,0)
294     else:
295     out+=2*intend+"%s=%s\n"%(func_name,out2)
296     else:
297     out+=2*intend+"%s=Data(0.,%s,FunctionOnContactZero(self.domain))\n"%(func_name,str(tuple(sh[:1])))
298     for i in range(sh[0]):
299     out2=""
300     for j in range(sh[1]):
301     t=makeTermText(func[i,j])
302     if not t == "0":
303     if len(out2)>0: out2+="+"
304     out2+="%s"%t
305     if len(out2)>0:
306     out+=2*intend+"%s[%s]=%s\n"%(func_name,i,out2)
307     else:
308     raise RunTimeError,"KKK"
309     return out
310    
311     def makeNormalDerivativeText(func,dim,func_name,add_jump=False):
312     sh=func.shape
313     out=""
314     if add_jump:
315     extra=")*jump"
316     extra_pre="("
317     else:
318     extra=""
319     extra_pre=""
320    
321     out=2*intend+"x_boundary=FunctionOnBoundary(self.domain).getX()\n"
322     if dim==3:
323     out+=2*intend+"n=%swhereZero(x_boundary[0] ,self.ABS_TOL)*numarray.array([-1., 0., 0.])"%extra_pre
324     out+="+whereZero(x_boundary[0]-1.,self.ABS_TOL)*numarray.array([ 1., 0., 0.])"
325     out+="+whereZero(x_boundary[1] ,self.ABS_TOL)*numarray.array([ 0.,-1., 0.])"
326     out+="+whereZero(x_boundary[1]-1.,self.ABS_TOL)*numarray.array([ 0., 1., 0.])"
327     out+="+whereZero(x_boundary[2] ,self.ABS_TOL)*numarray.array([ 0., 0.,-1.])"
328     out+="+whereZero(x_boundary[2]-1.,self.ABS_TOL)*numarray.array([ 0., 0., 1.])%s\n"%extra
329     else:
330     out+=2*intend+"n=%swhereZero(x_boundary[0] ,self.ABS_TOL)*numarray.array([-1., 0.])"%extra_pre
331     out+="+whereZero(x_boundary[0]-1.,self.ABS_TOL)*numarray.array([ 1., 0.])"
332     out+="+whereZero(x_boundary[1] ,self.ABS_TOL)*numarray.array([ 0.,-1.])"
333     out+="+whereZero(x_boundary[1]-1.,self.ABS_TOL)*numarray.array([ 0., 1.])%s\n"%extra
334     if len(sh)==dim+2:
335     if sh[0]==1:
336     out2=""
337     for j in range(sh[1]):
338     t=makeTermText(func[0,j])
339     if not t == "0":
340     if len(out2)>0: out2+="+"
341     out2+="n[%s]*(%s)"%(j,t)
342     if len(out2)==0:
343     out+=2*intend+"%s=%s\n"%(func_name,0)
344     else:
345     out+=2*intend+"%s=%s\n"%(func_name,out2)
346     else:
347     out+=2*intend+"%s=Data(0.,%s,FunctionOnBoundary(self.domain))\n"%(func_name,str(tuple(sh[:1])))
348     for i in range(sh[0]):
349     out2=""
350     for j in range(sh[1]):
351     t=makeTermText(func[i,j])
352     if not t == "0":
353     if len(out2)>0: out2+="+"
354     out2+="n[%s]*(%s)"%(j,t)
355     if len(out2)>0:
356     out+=2*intend+"%s[%s]=%s\n"%(func_name,i,out2)
357     else:
358     raise RunTimeError,"KKK"
359     return out
360    
361     def makeFunctionText(func,dim,func_name,add_jump=False):
362     sh=func.shape
363     out=""
364     if add_jump:
365     extra=")*jump"
366     extra_pre="("
367     else:
368     extra=""
369     extra_pre=""
370    
371     if len(sh)==dim:
372     out+=2*intend+"%s=%s%s%s\n"%(func_name,extra_pre,makeTermText(func),extra)
373     elif len(sh)==dim+1:
374     if sh[0]==1:
375     out+=2*intend+"%s=%s%s%s\n"%(func_name,extra_pre,makeTermText(func[0]),extra)
376     else:
377     out+=2*intend+"%s=Data(0.,%s,ContinuousFunction(self.domain))\n"%(func_name,str(tuple(sh[:1])))
378     for i in range(sh[0]):
379     t=makeTermText(func[i])
380     if not t=="0": out+=2*intend+"%s[%s]=%s%s%s\n"%(func_name,i,extra_pre,t,extra)
381     elif len(sh)==dim+2:
382     if sh[0]==1:
383     if sh[1]==1:
384     out+=2*intend+"%s=%s%s%s\n"%(func_name,extra_pre,makeTermText(func[0,0]),extra)
385     else:
386     out+=2*intend+"%s=Data(0.,%s,ContinuousFunction(self.domain))\n"%(func_name,str(tuple(sh[1:2])))
387     for j in range(sh[1]):
388     t=makeTermText(func[0,j])
389     if not t=="0": out+=2*intend+"%s[%s]=%s%s%s\n"%(func_name,j,extra_pre,t,extra)
390     else:
391     if sh[1]==1:
392     out+=2*intend+"%s=Data(0.,%s,ContinuousFunction(self.domain))\n"%(func_name,str(tuple(sh[:1])))
393     for i in range(sh[0]):
394     t=makeTermText(func[i,0])
395     if not t=="0": out+=2*intend+"%s[%s]=%s%s%s\n"%(func_name,i,extra_pre,t,extra)
396     else:
397     out+=2*intend+"%s=Data(0.,%s,ContinuousFunction(self.domain))\n"%(func_name,str(tuple(sh[:2])))
398     for i in range(sh[0]):
399     for j in range(sh[1]):
400     t=makeTermText(func[i,j])
401     if not t=="0": out+=2*intend+"%s[%s,%s]=%s%s%s\n"%(func_name,i,j,extra_pre,t,extra)
402     return out
403     for s in [1,2]:
404 gross 589 for d in [2,3]:
405 gross 766 for typ in ["Strong","Weak", "Contact"]:
406     for coffo in ["Full"]: # ["Full", "Reduced" ]:
407     if coffo=="Full":
408     func_i="Function"
409     else:
410     func_i="ReducedFunction"
411     for case in ["Const", "Vario" ]:
412     for solo in [1, 2 ]:
413     order=solo
414     if typ in ["Strong","Weak"]:
415     # coefficient A:
416     if s==1:
417 gross 589 for i in range(d):
418     for j in range(d):
419 gross 766 test_func=makeTestSolution(order,s,d)
420     grad_test_func=makeGradient(test_func,d)
421     body2=makeFunctionText(test_func,d,"u")
422     body2+=2*intend+"A_test=Data(0.,(%d,%d),%s(self.domain))\n"%(d,d,func_i)
423     if case == "Const" :
424     f=int(8*random.random())+1
425     body2+=2*intend+"A_test[%d,%d]=%s\n"%(i,j,f)
426     A_x_grad=multiplyByFandX(f,-1,grad_test_func,(0,j),(1,d),(0,i))
427     else:
428     body2+=2*intend+"A_test[%d,%d]=x[%s]\n"%(i,j,i)
429     A_x_grad=multiplyByFandX(1.,i,grad_test_func,(0,j),(1,d),(0,i))
430     if typ=="Weak":
431     body2+=makeFunctionText(A_x_grad,d,"X_test")
432     args="A=A_test, X=X_test"
433     else:
434     div_A_x_grad=makeDiv(A_x_grad,d)
435     body2+=makeFunctionText(-div_A_x_grad,d,"Y_test")
436     body2+=makeNormalDerivativeText(A_x_grad,d,"y_test")
437     args="A=A_test, Y=Y_test, y=y_test"
438    
439     makeTitle(d,coffo,solo,s,"A",case,typ,body2,mark="%s%s"%(i,j),pdeargs=args)
440     else:
441 gross 589 for p in range(s):
442     for i in range(d):
443     for q in range(s):
444     for j in range(d):
445 gross 766 test_func=makeTestSolution(order,s,d)
446     grad_test_func=makeGradient(test_func,d)
447     body2=makeFunctionText(test_func,d,"u")
448     body2+=2*intend+"A_test=Data(0.,(%d,%d,%d,%d),%s(self.domain))\n"%(s,d,s,d,func_i)
449     if case == "Const" :
450     f=int(8*random.random())+1
451     body2+=2*intend+"A_test[%d,%d,%d,%d]=%s\n"%(p,i,q,j,f)
452     A_x_grad=multiplyByFandX(f,-1,grad_test_func,(q,j),(s,d),(p,i))
453     else:
454     body2+=2*intend+"A_test[%d,%d,%d,%d]=x[%s]\n"%(p,i,q,j,i)
455     A_x_grad=multiplyByFandX(1.,i,grad_test_func,(q,j),(s,d),(p,i))
456     if typ=="Weak":
457     body2+=makeFunctionText(A_x_grad,d,"X_test")
458     args="A=A_test, X=X_test"
459     else:
460     div_A_x_grad=makeDiv(A_x_grad,d)
461     body2+=makeFunctionText(-div_A_x_grad,d,"Y_test")
462     body2+=makeNormalDerivativeText(A_x_grad,d,"y_test")
463     args="A=A_test, Y=Y_test, y=y_test"
464     makeTitle(d,coffo,solo,s,"A",case,typ,body2,mark="%s%s%s%s"%(p,i,q,j),pdeargs=args)
465 gross 589 # coefficient B:
466     if s==1:
467     for i in range(d):
468 gross 766 test_func=makeTestSolution(order,s,d)
469     body2=makeFunctionText(test_func,d,"u",typ=="Contact")
470     body2+=2*intend+"B_test=Data(0.,(%d,),%s(self.domain))\n"%(d,func_i)
471     if case == "Const" :
472     f=int(8*random.random())+1
473     body2+=2*intend+"B_test[%d]=%s\n"%(i,f)
474     B=multiplyByFandX(f,-1,test_func,(0,),(1,d),(0,i))
475     else:
476     body2+=2*intend+"B_test[%d]=x[%i]\n"%(i,i)
477     B=multiplyByFandX(1.,i,test_func,(0,),(1,d),(0,i))
478     if typ=="Weak":
479     body2+=makeFunctionText(B,d,"X_test")
480     args="B=B_test, X=X_test"
481     elif typ=="Strong":
482     div_B=makeDiv(B,d)
483     body2+=makeFunctionText(-div_B,d,"Y_test")
484     body2+=makeNormalDerivativeText(B,d,"y_test")
485     args="B=B_test, Y=Y_test, y=y_test"
486     else:
487     div_B=makeDiv(B,d)
488     body2+=makeFunctionText(-div_B,d,"Y_test",typ=="Contact")
489     body2+=makeNormalDerivativeText(B,d,"y_test",typ=="Contact")
490     body2+=makeContactDerivativeText(-B,d,"y_contact_test")
491     args="B=B_test, Y=Y_test, y=y_test, y_contact=y_contact_test"
492     makeTitle(d,coffo,solo,s,"B",case,typ,body2,mark="%s"%(i),pdeargs=args,add_jump=(typ=="Contact"))
493 gross 589 else:
494     for p in range(s):
495     for i in range(d):
496     for q in range(s):
497 gross 766 test_func=makeTestSolution(order,s,d)
498     grad_test_func=makeGradient(test_func,d)
499     body2=makeFunctionText(test_func,d,"u",typ=="Contact")
500     body2+=2*intend+"B_test=Data(0.,(%d,%d,%d),%s(self.domain))\n"%(s,d,s,func_i)
501     if case == "Const" :
502     f=int(8*random.random())+1
503     body2+=2*intend+"B_test[%d,%d,%d]=%s\n"%(p,i,q,f)
504     B=multiplyByFandX(f,-1,test_func,(q,),(s,d),(p,i))
505     else:
506     body2+=2*intend+"B_test[%d,%d,%d]=x[%i]\n"%(p,i,q,i)
507     B=multiplyByFandX(1.,i,test_func,(q,),(s,d),(p,i))
508     if typ=="Weak":
509     body2+=makeFunctionText(B,d,"X_test")
510     args="B=B_test, X=X_test"
511     elif typ=="Strong":
512     div_B=makeDiv(B,d)
513     body2+=makeFunctionText(-div_B,d,"Y_test")
514     body2+=makeNormalDerivativeText(B,d,"y_test")
515     args="B=B_test, Y=Y_test, y=y_test"
516     else:
517     div_B=makeDiv(B,d)
518     body2+=makeFunctionText(-div_B,d,"Y_test",typ=="Contact")
519     body2+=makeNormalDerivativeText(B,d,"y_test",typ=="Contact")
520     body2+=makeContactDerivativeText(B,d,"y_contact_test")
521     args="B=B_test, Y=Y_test, y=y_test, y_contact=y_contact_test"
522     makeTitle(d,coffo,solo,s,"B",case,typ,body2,mark="%s%s%s"%(p,i,q),pdeargs=args,add_jump=(typ=="Contact"))
523     if typ=="Strong":
524 gross 589 # coefficient C:
525     if s==1:
526     for j in range(d):
527 gross 766 test_func=makeTestSolution(order,s,d)
528     grad_test_func=makeGradient(test_func,d)
529     body2=makeFunctionText(test_func,d,"u")
530     body2+=2*intend+"C_test=Data(0.,(%d,),%s(self.domain))\n"%(d,func_i)
531     if case == "Const" :
532     f=int(8*random.random())+1
533     body2+=2*intend+"C_test[%d]=%s\n"%(j,f)
534     C_x_grad=multiplyByFandX(f,-1,grad_test_func,(0,j),(1,),(0,))
535     else:
536     body2+=2*intend+"C_test[%d]=x[%i]\n"%(j,j)
537     C_x_grad=multiplyByFandX(1,j,grad_test_func,(0,j),(1,),(0,))
538     body2+=makeFunctionText(C_x_grad,d,"Y_test")
539     args="C=C_test, Y=Y_test"
540     makeTitle(d,coffo,solo,s,"C",case,typ,body2,mark="%s"%(j),pdeargs=args)
541 gross 589 else:
542     for p in range(s):
543 gross 766 for q in range(s):
544     for j in range(d):
545     test_func=makeTestSolution(order,s,d)
546     grad_test_func=makeGradient(test_func,d)
547     body2=makeFunctionText(test_func,d,"u")
548     body2+=2*intend+"C_test=Data(0.,(%d,%d,%d),%s(self.domain))\n"%(s,s,d,func_i)
549     if case == "Const" :
550     f=int(8*random.random())+1
551     body2+=2*intend+"C_test[%d,%d,%d]=%s\n"%(p,q,j,f)
552     C_x_grad=multiplyByFandX(f,-1,grad_test_func,(q,j),(s,),(p,))
553     else:
554     body2+=2*intend+"C_test[%d,%d,%d]=x[%i]\n"%(p,q,j,j)
555     C_x_grad=multiplyByFandX(1.,j,grad_test_func,(q,j),(s,),(p,))
556     args="C=C_test, Y=Y_test"
557     body2+=makeFunctionText(C_x_grad,d,"Y_test")
558     makeTitle(d,coffo,solo,s,"C",case,typ,body2,mark="%s%s%s"%(p,q,j),pdeargs=args)
559 gross 589 # coefficient D:
560     if s==1:
561 gross 766 if case == "Const" :
562     test_func=makeTestSolution(order,s,d)
563     else:
564     test_func=makeTestSolution(order-1,s,d)
565     body2=makeFunctionText(test_func,d,"u")
566     if case == "Const" :
567     f=int(8*random.random())+1
568     body2+=2*intend+"D_test=Data(%s,(),%s(self.domain))\n"%(f,func_i)
569     D=multiplyByFandX(f,-1,test_func,(0,),(1,),(0,))
570     else:
571     body2+=2*intend+"D_test=%s(self.domain).getX()[0]\n"%func_i
572     D=multiplyByFandX(1,0,test_func,(0,),(1,),(0,))
573     args="D=D_test, Y=Y_test"
574     body2+=makeFunctionText(D,d,"Y_test")
575     makeTitle(d,coffo,solo,s,"D",case,typ,body2,mark="",pdeargs=args)
576 gross 589 else:
577     for p in range(s):
578     for q in range(s):
579 gross 766 if case == "Const" :
580     test_func=makeTestSolution(order,s,d)
581     else:
582     test_func=makeTestSolution(order-1,s,d)
583     body2=makeFunctionText(test_func,d,"u")
584     body2+=2*intend+"D_test=Data(0.,(%d,%d),%s(self.domain))\n"%(s,s,func_i)
585     if case == "Const" :
586     f=int(8*random.random())+1
587     body2+=2*intend+"D_test[%d,%d]=%s\n"%(p,q,f)
588     D=multiplyByFandX(f,-1,test_func,(q,),(s,),(p,))
589     else:
590     body2+=2*intend+"D_test[%d,%d]=x[0]\n"%(p,q)
591     D=multiplyByFandX(1.,0,test_func,(q,),(s,),(p,))
592     args="D=D_test, Y=Y_test"
593     body2+=makeFunctionText(D,d,"Y_test")
594     makeTitle(d,coffo,solo,s,"D",case,typ,body2,mark="%s%s"%(p,q),pdeargs=args)
595     if coffo=="Full":
596     # coefficient d:
597     if typ == "Strong":
598     if s==1:
599     if case == "Const" :
600     test_func=makeTestSolution(order,s,d)
601     else:
602     test_func=makeTestSolution(order-1,s,d)
603     body2=makeFunctionText(test_func,d,"u")
604     if case == "Const" :
605     f=int(8*random.random())+1
606     body2+=2*intend+"d_test=Data(%s,(),FunctionOnBoundary(self.domain))\n"%(f,)
607     D=multiplyByFandX(f,-1,test_func,(0,),(1,),(0,))
608     else:
609     body2+=2*intend+"d_test=interpolate(x[0],FunctionOnBoundary(self.domain))\n"
610     D=multiplyByFandX(1.,0,test_func,(0,),(1,),(0,))
611     body2+=makeFunctionText(D,d,"y_test")
612     args="d=d_test, y=y_test"
613     makeTitle(d,coffo,solo,s,"d",case,typ,body2,mark="",pdeargs=args)
614     else:
615     for p in range(s):
616     for q in range(s):
617     if case == "Const" :
618     test_func=makeTestSolution(order,s,d)
619     else:
620     test_func=makeTestSolution(order-1,s,d)
621     body2=makeFunctionText(test_func,d,"u")
622     body2+=2*intend+"d_test=Data(0.,(%d,%d),FunctionOnBoundary(self.domain))\n"%(s,s)
623     if case == "Const" :
624     f=int(8*random.random())+1
625     body2+=2*intend+"d_test[%d,%d]=%s\n"%(p,q,f)
626     D=multiplyByFandX(f,-1,test_func,(q,),(s,),(p,))
627     else:
628     body2+=2*intend+"d_test[%d,%d]=x[0]\n"%(p,q)
629    
630     D=multiplyByFandX(1.,0,test_func,(q,),(s,),(p,))
631     args="d=d_test, y=y_test"
632     body2+=makeFunctionText(D,d,"y_test")
633     makeTitle(d,coffo,solo,s,"d",case,typ,body2,mark="%s%s"%(p,q),pdeargs=args)
634     # coefficient d_contact:
635     if typ == "Contact":
636     if s==1:
637     if case == "Const" :
638     test_func=makeTestSolution(order,s,d)
639     else:
640     test_func=makeTestSolution(order-1,s,d)
641    
642     body2=makeFunctionText(test_func,d,"u",add_jump=True)
643     if case == "Const" :
644     f=int(8*random.random())+1
645     body2+=2*intend+"d_contact_test=Data(%s,(),FunctionOnContactZero(self.domain))\n"%(f,)
646     D=multiplyByFandX(f,-1,test_func,(0,),(1,),(0,))
647     else:
648     body2+=2*intend+"d_contact_test=interpolate(x[0],FunctionOnContactZero(self.domain))\n"
649     D=multiplyByFandX(1.,0,test_func,(0,),(1,),(0,))
650     body2+=makeFunctionText(D,d,"y_contact_test")
651     args="d_contact=d_contact_test, y_contact=y_contact_test"
652     makeTitle(d,coffo,solo,s,"d_contact",case,typ,body2,mark="",pdeargs=args,add_jump=True)
653     else:
654     for p in range(s):
655     for q in range(s):
656     if case == "Const" :
657     test_func=makeTestSolution(order,s,d)
658     else:
659     test_func=makeTestSolution(order-1,s,d)
660    
661     body2=makeFunctionText(test_func,d,"u",add_jump=True)
662     body2+=2*intend+"d_contact_test=Data(0.,(%d,%d),FunctionOnContactZero(self.domain))\n"%(s,s)
663     if case == "Const" :
664     f=int(8*random.random())+1
665     body2+=2*intend+"d_contact_test[%d,%d]=%s\n"%(p,q,f)
666     D=multiplyByFandX(f,-1,test_func,(q,),(s,),(p,))
667     else:
668     body2+=2*intend+"d_contact_test[%d,%d]=x[0]\n"%(p,q)
669     D=multiplyByFandX(1.,0,test_func,(q,),(s,),(p,))
670     body2+=makeFunctionText(D,d,"y_contact_test")
671     args="d_contact=d_contact_test, y_contact=y_contact_test"
672     makeTitle(d,coffo,solo,s,"d_contact",case,typ,body2,mark="%s%s"%(p,q),pdeargs=args,add_jump=True
673     )
674     #================================
675 gross 589 print "import unittest"
676     print "import numarray"
677     print "from esys.escript import *"
678     print "from esys.finley import Rectangle,Brick"
679 gross 766 #print "class Test_assemblage_2Do1_Strong(unittest.TestCase):"
680     #print t_prog["2Do1_strong"]
681     #print "class Test_assemblage_2Do2_Strong(unittest.TestCase):"
682     #print t_prog["2Do2_strong"]
683     #print "class Test_assemblage_3Do1_Strong(unittest.TestCase):"
684     #print t_prog["3Do1_strong"]
685     #print "class Test_assemblage_3Do2_Strong(unittest.TestCase):"
686     #print t_prog["3Do2_strong"]
687     #print "class Test_assemblage_2Do1_Weak(unittest.TestCase):"
688     #print t_prog["2Do1_weak"]
689     #print "class Test_assemblage_2Do2_Weak(unittest.TestCase):"
690     #print t_prog["2Do2_weak"]
691     #print "class Test_assemblage_3Do1_Weak(unittest.TestCase):"
692     #print t_prog["3Do1_weak"]
693     #print "class Test_assemblage_3Do2_Weak(unittest.TestCase):"
694     #print t_prog["3Do2_weak"]
695     print "class Test_assemblage_2Do1_Contact(unittest.TestCase):"
696 gross 589 print t_prog["2Do1_contact"]
697 gross 766 #print "class Test_assemblage_2Do2_Contact(unittest.TestCase):"
698     #print t_prog["2Do2_contact"]
699     #print "class Test_assemblage_3Do1_Contact(unittest.TestCase):"
700     #print t_prog["3Do1_contact"]
701     #print "class Test_assemblage_3Do2_contact(unittest.TestCase):"
702     #print t_prog["3Do2_contact"]
703 gross 589
704 gross 766
705     print "from esys.escript import *\n"
706     print "from esys.finley import *\n"
707     print "from esys.escript.linearPDEs import LinearPDE\n"
708     print "import numarray"
709     print "import unittest"
710    
711     # print "class Test_Finley_assemblage_2Do1(Test_assemblage_2Do1_Strong,Test_assemblage_2Do1_Weak):"
712     # print " RES_TOL=1.e-7"
713     # print " ABS_TOL=1.e-8"
714     # print " def setUp(self):"
715     # print " self.domain =Rectangle(2,2,1)"
716     print "class Test_Finley_assemblage_2Do1_Contact(Test_assemblage_2Do1_Contact):"
717 gross 589 print " RES_TOL=1.e-7"
718 gross 766 print " ABS_TOL=1.e-8"
719 gross 589 print " def setUp(self):"
720 gross 766 print " d1 = Rectangle(n0=2,n1=1,l0=0.5,order=1,useElementsOnFace=True)"
721     print " x1 = ContinuousFunction(d1).getX()"
722     print " ContinuousFunction(d1).setTags(1,Scalar(1,ContinuousFunction(d1)))"
723     print " d2 = Rectangle(n0=2,n1=1,l0=0.5,order=1,useElementsOnFace=True)"
724     print " ContinuousFunction(d2).setTags(2,Scalar(1,ContinuousFunction(d2)))"
725     print " d2.setX(d2.getX()+[0.5,0.])"
726     print " self.domain = JoinFaces([d1,d2])"
727 gross 589 print "suite = unittest.TestSuite()"
728 gross 766 print "suite.addTest(unittest.makeSuite(Test_Finley_assemblage_2Do1_Contact))"
729 gross 589 print "unittest.TextTestRunner(verbosity=2).run(suite)"

  ViewVC Help
Powered by ViewVC 1.1.26