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

Diff of /trunk/escript/py_src/generateutil

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 437 by gross, Fri Jan 20 00:16:58 2006 UTC revision 550 by gross, Wed Feb 22 02:14:38 2006 UTC
# Line 179  def makeArray(shape,rng): Line 179  def makeArray(shape,rng):
179               for i2 in range(shape[2]):               for i2 in range(shape[2]):
180                  for i3 in range(shape[3]):                  for i3 in range(shape[3]):
181                    for i4 in range(shape[4]):                    for i4 in range(shape[4]):
182                     out[i0,i1,i2,i3,i4]=l*random.random()+rng[0]                     out[i0,i1,i2,i3,i4]=l*ranm.random()+rng[0]
183     else:     else:
184         raise SystemError,"rank is restricted to 5"         raise SystemError,"rank is restricted to 5"
185     return out             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):  def makeResult(val,test_expr):
214     if isinstance(val,float):     if isinstance(val,float):
# Line 480  def mkText(case,name,a,a1=None,use_taggi Line 505  def mkText(case,name,a,a1=None,use_taggi
505                   t_out+="      %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist())                   t_out+="      %s.setTaggedValue(1,numarray.array(%s))\n"%(name,a1.tolist())
506                t_out+="      %s.expand()\n"%name                          t_out+="      %s.expand()\n"%name          
507             else:             else:
508                t_out+="      msk_%s=whereNegative(self.functionspace.getX()[0]-0.5)\n"%name                t_out+="      msk_%s=whereZero(self.functionspace.getX()[0],1.e-8)\n"%name
509                if isinstance(a,float):                if isinstance(a,float):
510                     t_out+="      %s=msk_%s*(%s)+(1.-msk_%s)*(%s)\n"%(name,name,a,name,a1)                     t_out+="      %s=(1.-msk_%s)*(%s)+msk_%s*(%s)\n"%(name,name,a,name,a1)
511                elif a.rank==0:                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)                     t_out+="      %s=msk_%s*numarray.array(%s)+(1.-msk_%s)*numarray.array(%s)\n"%(name,name,a,name,a1)
513                else:                else:
# Line 524  def mkCode(txt,args=[],intend=""): Line 549  def mkCode(txt,args=[],intend=""):
549      for r in args:      for r in args:
550        out=out.replace("%%a%s%%"%c,r)        out=out.replace("%%a%s%%"%c,r)
551      return out        return out  
552    #=======================================================================================================
553    # get slices
554    #=======================================================================================================
555    from esys.escript import *
556    for case0 in ["constData","taggedData","expandedData"]:
557       for sh0 in [ (3,), (3,4), (3,4,3) ,(4,3,3,3)]:
558        # get perm:
559        if len(sh0)==2:
560            check=[[1,0]]
561        elif len(sh0)==3:
562            check=[[1,0,2],
563                   [1,2,0],
564                   [2,1,0],
565                   [2,0,2],
566                   [0,2,1]]
567        elif len(sh0)==4:
568            check=[[0,1,3,2],
569                   [0,2,1,3],
570                   [0,2,3,1],
571                   [0,3,2,1],
572                   [0,3,1,2] ,          
573                   [1,0,2,3],
574                   [1,0,3,2],
575                   [1,2,0,3],
576                   [1,2,3,0],
577                   [1,3,2,0],
578                   [1,3,0,2],
579                   [2,0,1,3],
580                   [2,0,3,1],
581                   [2,1,0,3],
582                   [2,1,3,0],
583                   [2,3,1,0],
584                   [2,3,0,1],
585                   [3,0,1,2],
586                   [3,0,2,1],
587                   [3,1,0,2],
588                   [3,1,2,0],
589                   [3,2,1,0],
590                   [3,2,0,1]]
591        else:
592             check=[]
593        
594        # create the test cases:
595        processed=[]
596        l=["R","U","L","P","C","N"]
597        c=[""]
598        for i in range(len(sh0)):
599           tmp=[]
600           for ci in c:
601              tmp+=[ci+li for li in l]
602           c=tmp
603        # SHUFFLE
604        c2=[]
605        while len(c)>0:
606            i=int(random.random()*len(c))
607            c2.append(c[i])
608            del c[i]
609        c=c2
610        for ci in c:
611          t=""
612          sh=()
613          sl=()
614          for i in range(len(ci)):
615              if ci[i]=="R":
616                 s="%s:%s"%(1,sh0[i]-1)
617                 sl=sl+(slice(1,sh0[i]-1),)
618                 sh=sh+(sh0[i]-2,)            
619              if ci[i]=="U":
620                  s=":%s"%(sh0[i]-1)
621                  sh=sh+(sh0[i]-1,)
622                  sl=sl+(slice(0,sh0[i]-1),)            
623              if ci[i]=="L":
624                  s="2:"
625                  sh=sh+(sh0[i]-2,)            
626                  sl=sl+(slice(2,sh0[i]),)            
627              if ci[i]=="P":
628                  s="%s"%(int(sh0[i]/2))
629                  sl=sl+(int(sh0[i]/2),)            
630              if ci[i]=="C":
631                  s=":"
632                  sh=sh+(sh0[i],)            
633                  sl=sl+(slice(0,sh0[i]),)            
634              if ci[i]=="N":
635                  s=""
636                  sh=sh+(sh0[i],)
637              if len(s)>0:
638                 if not t=="": t+=","
639                 t+=s
640          if len(sl)==1: sl=sl[0]
641          N_found=False
642          noN_found=False
643          process=len(t)>0
644          for i in ci:
645             if i=="N":
646                if not noN_found and N_found: process=False
647                N_found=True
648             else:
649               if N_found: process=False
650               noNfound=True
651          # is there a similar one processed allready
652          if process and ci.find("N")==-1:
653             for ci2 in processed:
654               for chi in check:
655                   is_perm=True
656                   for i in range(len(chi)):
657                       if not ci[i]==ci2[chi[i]]: is_perm=False
658                   if is_perm: process=False
659          # if not process: print ci," rejected"
660          if process:
661           processed.append(ci)
662           for case1 in ["array","constData","taggedData","expandedData"]:
663              text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
664              tname="test_setslice_%s_rank%s_%s_%s"%(case0,len(sh0),case1,ci)
665              text+="   def %s(self):\n"%tname
666              a_0=makeNumberedArray(sh0)
667              if case0 in ["taggedData", "expandedData"]:
668                  a1_0=makeNumberedArray(sh0)
669              else:
670                  a1_0=a_0*1.
671    
672              a_1=makeNumberedArray(sh)
673              if case1 in ["taggedData", "expandedData"]:
674                  a1_1=makeNumberedArray(sh)
675              else:
676                  a1_1=a_1*1.
677    
678              text+=mkText(case0,"arg",a_0,a1_0)                                  
679              text+=mkText(case1,"val",a_1,a1_1)                                  
680              text+="      arg[%s]=val\n"%t
681              a_0.__setitem__(sl,a_1)
682              a1_0.__setitem__(sl,a1_1)
683              if Lsup(a_0-a1_0)==0.:
684                 text+=mkText("constData","ref",a_0,a1_0)
685              else:
686                 text+=mkText("expandedData","ref",a_0,a1_0)
687              text+="      self.failUnless(Lsup(arg-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"
688                  
689              if case0 == "taggedData" or case1 == "taggedData":
690                t_prog_with_tags+=text
691              else:              
692                t_prog+=text
693    
694    print test_header
695    # print t_prog
696    print t_prog_with_tags
697    print test_tail          
698    1/0
699    
700    #=======================================================================================================
701    # (non)symmetric part
702    #=======================================================================================================
703    from esys.escript import *
704    for name in ["symmetric", "nonsymmetric"]:
705     f=1.
706     if name=="nonsymmetric": f=-1
707     for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
708      for sh0 in [ (3,3), (2,3,2,3)]:
709                  text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
710                  tname="test_%s_%s_rank%s"%(name,case0,len(sh0))
711                  text+="   def %s(self):\n"%tname
712                  a_0=makeNumberedArray(sh0,s=1.)
713                  r_0=(a_0+f*transpose(a_0))/2.
714                  if case0 in ["taggedData", "expandedData"]:
715                     a1_0=makeNumberedArray(sh0,s=-1.)
716                     r1_0=(a1_0+f*transpose(a1_0))/2.
717                  else:
718                      a1_0=a_0                  
719                      r1_0=r_0
720                  text+=mkText(case0,"arg",a_0,a1_0)
721                  text+="      res=%s(arg)\n"%name
722                  if case0=="Symbol":
723                     text+=mkText("array","s",a_0,a1_0)
724                     text+="      sub=res.substitute({arg:s})\n"
725                     res="sub"
726                     text+=mkText("array","ref",r_0,r1_0)
727                  else:
728                     res="res"
729                     text+=mkText(case0,"ref",r_0,r1_0)  
730                  text+=mkTypeAndShapeTest(case0,sh0,"res")
731                  text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
732                  
733                  if case0 == "taggedData" :
734                      t_prog_with_tags+=text
735                  else:              
736                      t_prog+=text
737    print test_header
738    print t_prog
739    # print t_prog_with_tags
740    print test_tail          
741    1/0
742    
743    #=======================================================================================================
744    # eigenvalues
745    #=======================================================================================================
746    import numarray.linear_algebra
747    name="eigenvalues"
748    for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
749      for sh0 in [ (1,1), (2,2), (3,3)]:
750                  text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
751                  tname="test_%s_%s_dim%s"%(name,case0,sh0[0])
752                  text+="   def %s(self):\n"%tname
753                  a_0=makeArray(sh0,[-1.,1])
754                  a_0=(a_0+numarray.transpose(a_0))/2.
755                  ev=numarray.linear_algebra.eigenvalues(a_0)
756                  ev.sort()
757                  if case0 in ["taggedData", "expandedData"]:
758                      a1_0=makeArray(sh0,[-1.,1])
759                      a1_0=(a1_0+numarray.transpose(a1_0))/2.
760                      ev1=numarray.linear_algebra.eigenvalues(a1_0)
761                      ev1.sort()
762                  else:
763                      a1_0=a_0                  
764                      ev1=ev
765                  text+=mkText(case0,"arg",a_0,a1_0)
766                  text+="      res=%s(arg)\n"%name
767                  if case0=="Symbol":
768                     text+=mkText("array","s",a_0,a1_0)
769                     text+="      sub=res.substitute({arg:s})\n"
770                     res="sub"
771                     text+=mkText("array","ref",ev,ev1)
772                  else:
773                     res="res"
774                     text+=mkText(case0,"ref",ev,ev1)  
775                  text+=mkTypeAndShapeTest(case0,(sh0[0],),"res")
776                  text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
777                  
778                  if case0 == "taggedData" :
779                      t_prog_with_tags+=text
780                  else:              
781                      t_prog+=text
782    print test_header
783    # print t_prog
784    print t_prog_with_tags
785    print test_tail          
786    1/0
787    
788    #=======================================================================================================
789    # get slices
790    #=======================================================================================================
791    for case0 in ["constData","taggedData","expandedData","Symbol"]:
792      for sh0 in [ (3,), (3,4), (3,4,3) ,(4,3,5,3)]:
793        # get perm:
794        if len(sh0)==2:
795            check=[[1,0]]
796        elif len(sh0)==3:
797            check=[[1,0,2],
798                   [1,2,0],
799                   [2,1,0],
800                   [2,0,2],
801                   [0,2,1]]
802        elif len(sh0)==4:
803            check=[[0,1,3,2],
804                   [0,2,1,3],
805                   [0,2,3,1],
806                   [0,3,2,1],
807                   [0,3,1,2] ,          
808                   [1,0,2,3],
809                   [1,0,3,2],
810                   [1,2,0,3],
811                   [1,2,3,0],
812                   [1,3,2,0],
813                   [1,3,0,2],
814                   [2,0,1,3],
815                   [2,0,3,1],
816                   [2,1,0,3],
817                   [2,1,3,0],
818                   [2,3,1,0],
819                   [2,3,0,1],
820                   [3,0,1,2],
821                   [3,0,2,1],
822                   [3,1,0,2],
823                   [3,1,2,0],
824                   [3,2,1,0],
825                   [3,2,0,1]]
826        else:
827             check=[]
828        
829        # create the test cases:
830        processed=[]
831        l=["R","U","L","P","C","N"]
832        c=[""]
833        for i in range(len(sh0)):
834           tmp=[]
835           for ci in c:
836              tmp+=[ci+li for li in l]
837           c=tmp
838        # SHUFFLE
839        c2=[]
840        while len(c)>0:
841            i=int(random.random()*len(c))
842            c2.append(c[i])
843            del c[i]
844        c=c2
845        for ci in c:
846          t=""
847          sh=()
848          for i in range(len(ci)):
849              if ci[i]=="R":
850                 s="%s:%s"%(1,sh0[i]-1)
851                 sh=sh+(sh0[i]-2,)            
852              if ci[i]=="U":
853                  s=":%s"%(sh0[i]-1)
854                  sh=sh+(sh0[i]-1,)            
855              if ci[i]=="L":
856                  s="2:"
857                  sh=sh+(sh0[i]-2,)            
858              if ci[i]=="P":
859                  s="%s"%(int(sh0[i]/2))
860              if ci[i]=="C":
861                  s=":"
862                  sh=sh+(sh0[i],)            
863              if ci[i]=="N":
864                  s=""
865                  sh=sh+(sh0[i],)
866              if len(s)>0:
867                 if not t=="": t+=","
868                 t+=s
869          N_found=False
870          noN_found=False
871          process=len(t)>0
872          for i in ci:
873             if i=="N":
874                if not noN_found and N_found: process=False
875                N_found=True
876             else:
877               if N_found: process=False
878               noNfound=True
879          # is there a similar one processed allready
880          if process and ci.find("N")==-1:
881             for ci2 in processed:
882               for chi in check:
883                   is_perm=True
884                   for i in range(len(chi)):
885                       if not ci[i]==ci2[chi[i]]: is_perm=False
886                   if is_perm: process=False
887          # if not process: print ci," rejected"
888          if process:
889           processed.append(ci)
890           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
891           tname="test_getslice_%s_rank%s_%s"%(case0,len(sh0),ci)
892           text+="   def %s(self):\n"%tname
893           a_0=makeNumberedArray(sh0,s=1)
894           if case0 in ["taggedData", "expandedData"]:
895                a1_0=makeNumberedArray(sh0,s=-1.)
896           else:
897                a1_0=a_0
898           r=eval("a_0[%s]"%t)
899           r1=eval("a1_0[%s]"%t)
900           text+=mkText(case0,"arg",a_0,a1_0)
901           text+="      res=arg[%s]\n"%t
902           if case0=="Symbol":
903               text+=mkText("array","s",a_0,a1_0)
904               text+="      sub=res.substitute({arg:s})\n"
905               res="sub"
906               text+=mkText("array","ref",r,r1)
907           else:
908               res="res"
909               text+=mkText(case0,"ref",r,r1)
910           text+=mkTypeAndShapeTest(case0,sh,"res")
911           text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
912                  
913           if case0 == "taggedData" :
914                t_prog_with_tags+=text
915           else:              
916                t_prog+=text
917    
918    print test_header
919    # print t_prog
920    print t_prog_with_tags
921    print test_tail          
922    1/0
923    #============================================================================================
924  def innerTEST(arg0,arg1):  def innerTEST(arg0,arg1):
925      if isinstance(arg0,float):      if isinstance(arg0,float):
926         out=numarray.array(arg0*arg1)         out=numarray.array(arg0*arg1)
# Line 696  def minimumTEST(arg0,arg1): Line 1092  def minimumTEST(arg0,arg1):
1092                else:                else:
1093                 out[i0,i1,i2,i3]=arg0[i0,i1,i2,i3]                 out[i0,i1,i2,i3]=arg0[i0,i1,i2,i3]
1094       return out       return out
1095  def unrollLoops(a,b,o,arg,tap=""):      
1096    def unrollLoops(a,b,o,arg,tap="",x="x"):
1097      out=""      out=""
1098      if a.rank==1:      if a.rank==1:
1099               z=""               z=""
1100               for i99 in range(a.shape[0]):               for i99 in range(a.shape[0]):
1101                 if not z=="": z+="+"                 if not z=="": z+="+"
1102                 if o=="1":                 if o=="1":
1103                    z+="(%s)*x[%s]"%(a[i99]+b[i99],i99)                    z+="(%s)*%s[%s]"%(a[i99]+b[i99],x,i99)
1104                 else:                 else:
1105                    z+="(%s)*x[%s]**o+(%s)*x[%s]"%(a[i99],i99,b[i99],i99)                    z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i99],x,i99,b[i99],x,i99)
1106    
1107               out+=tap+"%s=%s\n"%(arg,z)               out+=tap+"%s=%s\n"%(arg,z)
1108    
# Line 717  def unrollLoops(a,b,o,arg,tap=""): Line 1114  def unrollLoops(a,b,o,arg,tap=""):
1114                 if o=="1":                 if o=="1":
1115                    z+="(%s)*x[%s]"%(a[i0,i99]+b[i0,i99],i99)                    z+="(%s)*x[%s]"%(a[i0,i99]+b[i0,i99],i99)
1116                 else:                 else:
1117                    z+="(%s)*x[%s]**o+(%s)*x[%s]"%(a[i0,i99],i99,b[i0,i99],i99)                    z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i99],x,i99,b[i0,i99],x,i99)
1118    
1119               out+=tap+"%s[%s]=%s\n"%(arg,i0,z)               out+=tap+"%s[%s]=%s\n"%(arg,i0,z)
1120      elif a.rank==3:      elif a.rank==3:
# Line 727  def unrollLoops(a,b,o,arg,tap=""): Line 1124  def unrollLoops(a,b,o,arg,tap=""):
1124               for i99 in range(a.shape[2]):               for i99 in range(a.shape[2]):
1125                 if not z=="": z+="+"                 if not z=="": z+="+"
1126                 if o=="1":                 if o=="1":
1127                    z+="(%s)*x[%s]"%(a[i0,i1,i99]+b[i0,i1,i99],i99)                    z+="(%s)*%s[%s]"%(a[i0,i1,i99]+b[i0,i1,i99],x,i99)
1128                 else:                 else:
1129                    z+="(%s)*x[%s]**o+(%s)*x[%s]"%(a[i0,i1,i99],i99,b[i0,i1,i99],i99)                    z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i99],x,i99,b[i0,i1,i99],x,i99)
1130    
1131               out+=tap+"%s[%s,%s]=%s\n"%(arg,i0,i1,z)               out+=tap+"%s[%s,%s]=%s\n"%(arg,i0,i1,z)
1132      elif a.rank==4:      elif a.rank==4:
# Line 740  def unrollLoops(a,b,o,arg,tap=""): Line 1137  def unrollLoops(a,b,o,arg,tap=""):
1137               for i99 in range(a.shape[3]):               for i99 in range(a.shape[3]):
1138                 if not z=="": z+="+"                 if not z=="": z+="+"
1139                 if o=="1":                 if o=="1":
1140                    z+="(%s)*x[%s]"%(a[i0,i1,i2,i99]+b[i0,i1,i2,i99],i99)                    z+="(%s)*%s[%s]"%(a[i0,i1,i2,i99]+b[i0,i1,i2,i99],x,i99)
1141                 else:                 else:
1142                    z+="(%s)*x[%s]**o+(%s)*x[%s]"%(a[i0,i1,i2,i99],i99,b[i0,i1,i2,i99],i99)                    z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i99],x,i99,b[i0,i1,i2,i99],x,i99)
1143    
1144               out+=tap+"%s[%s,%s,%s]=%s\n"%(arg,i0,i1,i2,z)               out+=tap+"%s[%s,%s,%s]=%s\n"%(arg,i0,i1,i2,z)
1145      elif a.rank==5:      elif a.rank==5:
# Line 754  def unrollLoops(a,b,o,arg,tap=""): Line 1151  def unrollLoops(a,b,o,arg,tap=""):
1151               for i99 in range(a.shape[4]):               for i99 in range(a.shape[4]):
1152                 if not z=="": z+="+"                 if not z=="": z+="+"
1153                 if o=="1":                 if o=="1":
1154                    z+="(%s)*x[%s]"%(a[i0,i1,i2,i3,i99]+b[i0,i1,i2,i3,i99],i99)                    z+="(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99]+b[i0,i1,i2,i3,i99],x,i99)
1155                 else:                 else:
1156                    z+="(%s)*x[%s]**o+(%s)*x[%s]"%(a[i0,i1,i2,i3,i99],i99,b[i0,i1,i2,i3,i99],i99)                    z+="(%s)*%s[%s]**o+(%s)*%s[%s]"%(a[i0,i1,i2,i3,i99],x,i99,b[i0,i1,i2,i3,i99],x,i99)
1157    
1158               out+=tap+"%s[%s,%s,%s,%s]=%s\n"%(arg,i0,i1,i2,i3,z)               out+=tap+"%s[%s,%s,%s,%s]=%s\n"%(arg,i0,i1,i2,i3,z)
1159      return out      return out
# Line 796  def unrollLoopsOfGrad(a,b,o,arg,tap=""): Line 1193  def unrollLoopsOfGrad(a,b,o,arg,tap=""):
1193                   out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99])                   out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99])
1194                 else:                 else:
1195                   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])                   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])
1196        return out
1197    def unrollLoopsOfDiv(a,b,o,arg,tap=""):
1198        out=tap+arg+"="
1199        if o=="1":
1200           z=0.
1201           for i99 in range(a.shape[0]):
1202                z+=b[i99,i99]+a[i99,i99]
1203           out+="(%s)"%z    
1204        else:
1205           z=0.
1206           for i99 in range(a.shape[0]):
1207                z+=b[i99,i99]
1208                if i99>0: out+="+"
1209                out+="o*(%s)*x_ref[%s]**(o-1)"%(a[i99,i99],i99)
1210           out+="+(%s)"%z    
1211        return out
1212    
1213  def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):  def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):
1214      if where=="Function":      if where=="Function":
# Line 916  def unrollLoopsOfInteriorIntegral(a,b,wh Line 1329  def unrollLoopsOfInteriorIntegral(a,b,wh
1329                 out+="+(%s)*0.5**o\n"%zop                 out+="+(%s)*0.5**o\n"%zop
1330    
1331      return out      return out
1332    def unrollLoopsSimplified(b,arg,tap=""):
1333        out=""
1334        if isinstance(b,float) or b.rank==0:
1335                 out+=tap+"%s=(%s)*x[0]\n"%(arg,str(b))
1336    
1337        elif b.rank==1:
1338            for i0 in range(b.shape[0]):
1339                 out+=tap+"%s[%s]=(%s)*x[%s]\n"%(arg,i0,b[i0],i0)
1340        elif b.rank==2:
1341            for i0 in range(b.shape[0]):
1342             for i1 in range(b.shape[1]):
1343                 out+=tap+"%s[%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,b[i0,i1],i1)
1344        elif b.rank==3:
1345            for i0 in range(b.shape[0]):
1346             for i1 in range(b.shape[1]):
1347               for i2 in range(b.shape[2]):
1348                 out+=tap+"%s[%s,%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,i2,b[i0,i1,i2],i2)
1349        elif b.rank==4:
1350            for i0 in range(b.shape[0]):
1351             for i1 in range(b.shape[1]):
1352               for i2 in range(b.shape[2]):
1353                for i3 in range(b.shape[3]):
1354                 out+=tap+"%s[%s,%s,%s,%s]=(%s)*x[%s]\n"%(arg,i0,i1,i2,i3,b[i0,i1,i2,i3],i3)
1355        return out
1356    
1357    def unrollLoopsOfL2(b,where,arg,tap=""):
1358        out=""
1359        z=[]
1360        if isinstance(b,float) or b.rank==0:
1361           z.append(b**2)
1362        elif b.rank==1:
1363            for i0 in range(b.shape[0]):
1364                 z.append(b[i0]**2)
1365        elif b.rank==2:
1366            for i1 in range(b.shape[1]):
1367               s=0
1368               for i0 in range(b.shape[0]):
1369                  s+=b[i0,i1]**2
1370               z.append(s)
1371        elif b.rank==3:
1372            for i2 in range(b.shape[2]):
1373              s=0
1374              for i0 in range(b.shape[0]):
1375                 for i1 in range(b.shape[1]):
1376                    s+=b[i0,i1,i2]**2
1377              z.append(s)
1378    
1379        elif b.rank==4:
1380          for i3 in range(b.shape[3]):
1381             s=0
1382             for i0 in range(b.shape[0]):
1383               for i1 in range(b.shape[1]):
1384                  for i2 in range(b.shape[2]):
1385                     s+=b[i0,i1,i2,i3]**2
1386             z.append(s)        
1387        if where=="Function":
1388           xfac_o=1.
1389           xfac_op=0.
1390           z_fac_s=""
1391           zo_fac_s=""
1392           zo_fac=1./3.
1393        elif where=="FunctionOnBoundary":
1394           xfac_o=1.
1395           xfac_op=0.
1396           z_fac_s="*dim"
1397           zo_fac_s="*(2.*dim+1.)/3."
1398           zo_fac=1.
1399        elif where in ["FunctionOnContactZero","FunctionOnContactOne"]:
1400           xfac_o=0.
1401           xfac_op=1.
1402           z_fac_s=""
1403           zo_fac_s=""    
1404           zo_fac=1./3.    
1405        zo=0.
1406        zop=0.
1407        for i99 in range(len(z)):
1408               if i99==0:
1409                   zo+=xfac_o*z[i99]
1410                   zop+=xfac_op*z[i99]
1411               else:
1412                   zo+=z[i99]
1413        out+=tap+"%s=sqrt((%s)%s"%(arg,zo*zo_fac,zo_fac_s)
1414        if zop==0.:
1415           out+=")\n"
1416        else:
1417           out+="+(%s))\n"%(zop*0.5**2)
1418        return out
1419    #=======================================================================================================
1420    # transpose
1421    #=======================================================================================================
1422    def transposeTest(r,offset):
1423        if isinstance(r,float): return r
1424        s=r.shape
1425        s1=1
1426        for i in s[:offset]: s1*=i
1427        s2=1
1428        for i in s[offset:]: s2*=i
1429        out=numarray.reshape(r,(s1,s2))
1430        out.transpose()
1431        return numarray.resize(out,s[offset:]+s[:offset])
1432    
1433    name,tt="transpose",transposeTest
1434    for case0 in ["array","Symbol","constData","taggedData","expandedData"]:
1435      for sh0 in [ (), (3,), (4,5), (6,2,2),(3,2,3,4)]:
1436        for offset in range(len(sh0)+1):
1437                  text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1438                  tname="test_%s_%s_rank%s_offset%s"%(name,case0,len(sh0),offset)
1439                  text+="   def %s(self):\n"%tname
1440                  sh_t=sh0[offset:]+sh0[:offset]
1441    
1442    #              sh_t=list(sh0)
1443    #              sh_t[offset+1]=sh_t[offset]
1444    #              sh_t=tuple(sh_t)
1445    #              sh_r=[]
1446    #              for i in range(offset): sh_r.append(sh0[i])
1447    #              for i in range(offset+2,len(sh0)): sh_r.append(sh0[i])              
1448    #              sh_r=tuple(sh_r)
1449    
1450                  a_0=makeArray(sh0,[-1.,1])
1451                  if case0 in ["taggedData", "expandedData"]:
1452                      a1_0=makeArray(sh0,[-1.,1])
1453                  else:
1454                      a1_0=a_0
1455                  r=tt(a_0,offset)
1456                  r1=tt(a1_0,offset)
1457                  text+=mkText(case0,"arg",a_0,a1_0)
1458                  text+="      res=%s(arg,%s)\n"%(name,offset)
1459                  if case0=="Symbol":
1460                     text+=mkText("array","s",a_0,a1_0)
1461                     text+="      sub=res.substitute({arg:s})\n"
1462                     res="sub"
1463                     text+=mkText("array","ref",r,r1)
1464                  else:
1465                     res="res"
1466                     text+=mkText(case0,"ref",r,r1)
1467                  text+=mkTypeAndShapeTest(case0,sh_t,"res")
1468                  text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1469                  
1470                  if case0 == "taggedData" :
1471                      t_prog_with_tags+=text
1472                  else:              
1473                      t_prog+=text
1474    
1475    print test_header
1476    # print t_prog
1477    print t_prog_with_tags
1478    print test_tail          
1479    1/0
1480  #=======================================================================================================  #=======================================================================================================
1481  # integrate  # L2
1482  #=======================================================================================================  #=======================================================================================================
1483  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1484    for data in ["Data","Symbol"]:    for data in ["Data","Symbol"]:
1485      for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:      for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
       for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:  
         if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:    
1486           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1487           tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))           tname="test_L2_on%s_from%s_rank%s"%(where,data,len(sh))
1488           text+="   def %s(self):\n"%tname           text+="   def %s(self):\n"%tname
1489           text+="      \"\"\"\n"           text+="      \"\"\"\n"
1490           text+="      tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where)           text+="      tests L2-norm of %s on the %s\n\n"%(data,where)
1491           text+="      assumptions: %s(self.domain) exists\n"%case           text+="      assumptions: self.domain supports integration on %s\n"%where
1492           text+="                   self.domain supports integral on %s\n"%where           text+="      \"\"\"\n"
1493             text+="      dim=self.domain.getDim()\n"
1494             text+="      w=%s(self.domain)\n"%where
1495             text+="      x=w.getX()\n"
1496             o="1"
1497             if len(sh)>0:
1498                sh_2=sh[:len(sh)-1]+(2,)
1499                sh_3=sh[:len(sh)-1]+(3,)            
1500                b_2=makeArray(sh[:len(sh)-1]+(2,),[-1.,1])
1501                b_3=makeArray(sh[:len(sh)-1]+(3,),[-1.,1])
1502             else:
1503                sh_2=()
1504                sh_3=()
1505                b_2=makeArray(sh,[-1.,1])
1506                b_3=makeArray(sh,[-1.,1])
1507    
1508             if data=="Symbol":
1509                val="s"
1510                res="sub"
1511             else:
1512                val="arg"
1513                res="res"
1514             text+="      if dim==2:\n"
1515             if data=="Symbol":
1516                   text+="        arg=Symbol(shape=%s,dim=dim)\n"%str(sh_2)
1517    
1518             text+="        %s=Data(0,%s,w)\n"%(val,sh_2)
1519             text+=unrollLoopsSimplified(b_2,val,tap="        ")
1520             text+=unrollLoopsOfL2(b_2,where,"ref",tap="        ")
1521             text+="\n      else:\n"
1522             if data=="Symbol":
1523                   text+="        arg=Symbol(shape=%s,dim=dim)\n"%str(sh_3)
1524             text+="        %s=Data(0,%s,w)\n"%(val,sh_3)        
1525             text+=unrollLoopsSimplified(b_3,val,tap="        ")
1526             text+=unrollLoopsOfL2(b_3,where,"ref",tap="        ")
1527             text+="\n      res=L2(arg)\n"
1528             if data=="Symbol":
1529                text+="      sub=res.substitute({arg:s})\n"
1530                text+="      self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"
1531                text+="      self.failUnlessEqual(res.getShape(),(),\"wrong shape of result.\")\n"
1532             else:
1533                text+="      self.failUnless(isinstance(res,float),\"wrong type of result.\")\n"
1534             text+="      self.failUnlessAlmostEqual(%s,ref,int(-log10(self.RES_TOL)),\"wrong result\")\n"%res
1535             t_prog+=text
1536    print t_prog
1537    1/0
1538    
1539    #=======================================================================================================
1540    # div
1541    #=======================================================================================================
1542    for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1543      for data in ["Data","Symbol"]:
1544         for case in ["ContinuousFunction","Solution","ReducedSolution"]:
1545             text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1546             tname="test_div_on%s_from%s_%s"%(where,data,case)
1547             text+="   def %s(self):\n"%tname
1548             text+="      \"\"\"\n"
1549             text+="      tests divergence of %s on the %s\n\n"%(data,where)
1550             text+="      assumptions: %s(self.domain) exists\n"%case
1551             text+="                   self.domain supports div on %s\n"%where
1552           text+="      \"\"\"\n"           text+="      \"\"\"\n"
1553           if case=="ReducedSolution":           if case=="ReducedSolution":
1554              text+="      o=1\n"              text+="      o=1\n"
# Line 944  for where in ["Function","FunctionOnBoun Line 1558  for where in ["Function","FunctionOnBoun
1558              o="o"              o="o"
1559           text+="      dim=self.domain.getDim()\n"           text+="      dim=self.domain.getDim()\n"
1560           text+="      w_ref=%s(self.domain)\n"%where           text+="      w_ref=%s(self.domain)\n"%where
1561             text+="      x_ref=w_ref.getX()\n"
1562             text+="      w=%s(self.domain)\n"%case
1563             text+="      x=w.getX()\n"
1564             a_2=makeArray((2,2),[-1.,1])
1565             b_2=makeArray((2,2),[-1.,1])
1566             a_3=makeArray((3,3),[-1.,1])
1567             b_3=makeArray((3,3),[-1.,1])
1568             if data=="Symbol":
1569                text+="      arg=Symbol(shape=(dim,),dim=dim)\n"
1570                val="s"
1571                res="sub"
1572             else:
1573                val="arg"
1574                res="res"
1575             text+="      %s=Vector(0,w)\n"%val
1576             text+="      if dim==2:\n"
1577             text+=unrollLoops(a_2,b_2,o,val,tap="        ")
1578             text+=unrollLoopsOfDiv(a_2,b_2,o,"ref",tap="        ")
1579             text+="\n      else:\n"
1580            
1581             text+=unrollLoops(a_3,b_3,o,val,tap="        ")
1582             text+=unrollLoopsOfDiv(a_3,b_3,o,"ref",tap="        ")
1583             text+="\n      res=div(arg,where=w_ref)\n"
1584             if data=="Symbol":
1585                text+="      sub=res.substitute({arg:s})\n"
1586             text+="      self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
1587             text+="      self.failUnlessEqual(res.getShape(),(),\"wrong shape of result.\")\n"
1588             text+="      self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong function space of result.\")\n"%res
1589             text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1590    
1591    
1592             t_prog+=text
1593    print t_prog
1594    1/0
1595    
1596    #=======================================================================================================
1597    # interpolation
1598    #=======================================================================================================
1599    for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne","Solution","ReducedSolution"]:
1600      for data in ["Data","Symbol"]:
1601         for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1602          for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
1603            if  where==case or \
1604                ( case in ["ContinuousFunction","Solution","ReducedSolution"] and where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"] ) or \
1605                ( case in ["FunctionOnContactZero","FunctionOnContactOne"] and where in ["FunctionOnContactZero","FunctionOnContactOne"] ) or \
1606                (case=="ContinuousFunction" and  where in ["Solution","ReducedSolution"]) or \
1607                (case=="Solution" and  where=="ReducedSolution") :
1608                
1609    
1610             text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1611             tname="test_interpolation_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1612             text+="   def %s(self):\n"%tname
1613             text+="      \"\"\"\n"
1614             text+="      tests interpolation for rank %s %s onto the %s\n\n"%(len(sh),data,where)
1615             text+="      assumptions: self.domain supports inpterpolation from %s onto %s\n"%(case,where)
1616             text+="      \"\"\"\n"
1617             if case=="ReducedSolution" or where=="ReducedSolution":
1618                text+="      o=1\n"
1619                o="1"
1620             else:
1621                text+="      o=self.order\n"
1622                o="o"
1623             text+="      dim=self.domain.getDim()\n"
1624             text+="      w_ref=%s(self.domain)\n"%where
1625             text+="      x_ref=w_ref.getX()\n"
1626           text+="      w=%s(self.domain)\n"%case           text+="      w=%s(self.domain)\n"%case
1627           text+="      x=w.getX()\n"           text+="      x=w.getX()\n"
1628           a_2=makeArray(sh+(2,),[-1.,1])           a_2=makeArray(sh+(2,),[-1.,1])
# Line 951  for where in ["Function","FunctionOnBoun Line 1630  for where in ["Function","FunctionOnBoun
1630           a_3=makeArray(sh+(3,),[-1.,1])           a_3=makeArray(sh+(3,),[-1.,1])
1631           b_3=makeArray(sh+(3,),[-1.,1])           b_3=makeArray(sh+(3,),[-1.,1])
1632           if data=="Symbol":           if data=="Symbol":
1633              text+="      arg=Symbol(shape=%s)\n"%str(sh)              text+="      arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
1634              val="s"              val="s"
1635              res="sub"              res="sub"
1636           else:           else:
1637              val="arg"              val="arg"
1638              res="res"              res="res"
               
1639           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))
1640           if not len(sh)==0:           text+="      ref=Data(0,%s,w_ref)\n"%str(sh)
             text+="      ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)  
1641           text+="      if dim==2:\n"           text+="      if dim==2:\n"
1642           text+=unrollLoops(a_2,b_2,o,val,tap="        ")           text+=unrollLoops(a_2,b_2,o,val,tap="        ")
1643           text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap="        ")           text+=unrollLoops(a_2,b_2,o,"ref",tap="        ",x="x_ref")
1644           text+="      else:\n"           text+="      else:\n"
1645                    
1646           text+=unrollLoops(a_3,b_3,o,val,tap="        ")           text+=unrollLoops(a_3,b_3,o,val,tap="        ")
1647           text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap="        ")           text+=unrollLoops(a_3,b_3,o,"ref",tap="        ",x="x_ref")
1648           if case in ["ContinuousFunction","Solution","ReducedSolution"]:           text+="      res=interpolate(arg,where=w_ref)\n"
              text+="      res=integrate(arg,where=w_ref)\n"  
          else:  
              text+="      res=integrate(arg)\n"  
   
1649           if data=="Symbol":           if data=="Symbol":
1650              text+="      sub=res.substitute({arg:s})\n"              text+="      sub=res.substitute({arg:s})\n"
1651           if len(sh)==0 and data=="Data":           text+="      self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
1652              text+="      self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res           text+="      self.failUnlessEqual(%s.getFunctionSpace(),w_ref,\"wrong functionspace of result.\")\n"%res
1653           else:           text+="      self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
             if data=="Symbol":  
                text+="      self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"  
                text+="      self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)  
             else:  
                text+="      self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n"  
                text+="      self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)  
1654           text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res           text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
   
   
1655           t_prog+=text           t_prog+=text
1656  print test_header  print test_header
1657  print t_prog  print t_prog
1658  print test_tail            print test_tail          
1659  1/0  1/0
1660    
1661  #=======================================================================================================  #=======================================================================================================
1662  # grad  # grad
1663  #=======================================================================================================  #=======================================================================================================
1664  for where in ["Function","FunctionOnBoundary"]:  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1665    for data in ["Data","Symbol"]:    for data in ["Data","Symbol"]:
1666       for case in ["ContinuousFunction","Solution","ReducedSolution"]:       for case in ["ContinuousFunction","Solution","ReducedSolution"]:
1667         for sh in [ (),(2,), (4,5), (6,2,2)]:         for sh in [ (),(2,), (4,5), (6,2,2)]:
# Line 1004  for where in ["Function","FunctionOnBoun Line 1670  for where in ["Function","FunctionOnBoun
1670           text+="   def %s(self):\n"%tname           text+="   def %s(self):\n"%tname
1671           text+="      \"\"\"\n"           text+="      \"\"\"\n"
1672           text+="      tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where)           text+="      tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where)
1673           text+="      assumptions: %s(self.domain) exists\n"%where           text+="      assumptions: %s(self.domain) exists\n"%case
1674           if where=="FunctionOnBoundary": text+="                   self.domain supports gradient on the boundary\n"           text+="                   self.domain supports gardient on %s\n"%where
1675           text+="      \"\"\"\n"           text+="      \"\"\"\n"
1676           if case=="ReducedSolution":           if case=="ReducedSolution":
1677              text+="      o=1\n"              text+="      o=1\n"
# Line 1052  print t_prog Line 1718  print t_prog
1718  print test_tail            print test_tail          
1719  1/0  1/0
1720    
1721    
1722    #=======================================================================================================
1723    # integrate
1724    #=======================================================================================================
1725    for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1726      for data in ["Data","Symbol"]:
1727        for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
1728          for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
1729            if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:  
1730             text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
1731             tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
1732             text+="   def %s(self):\n"%tname
1733             text+="      \"\"\"\n"
1734             text+="      tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where)
1735             text+="      assumptions: %s(self.domain) exists\n"%case
1736             text+="                   self.domain supports integral on %s\n"%where
1737    
1738             text+="      \"\"\"\n"
1739             if case=="ReducedSolution":
1740                text+="      o=1\n"
1741                o="1"
1742             else:
1743                text+="      o=self.order\n"
1744                o="o"
1745             text+="      dim=self.domain.getDim()\n"
1746             text+="      w_ref=%s(self.domain)\n"%where
1747             text+="      w=%s(self.domain)\n"%case
1748             text+="      x=w.getX()\n"
1749             a_2=makeArray(sh+(2,),[-1.,1])
1750             b_2=makeArray(sh+(2,),[-1.,1])
1751             a_3=makeArray(sh+(3,),[-1.,1])
1752             b_3=makeArray(sh+(3,),[-1.,1])
1753             if data=="Symbol":
1754                text+="      arg=Symbol(shape=%s)\n"%str(sh)
1755                val="s"
1756                res="sub"
1757             else:
1758                val="arg"
1759                res="res"
1760                
1761             text+="      %s=Data(0,%s,w)\n"%(val,str(sh))
1762             if not len(sh)==0:
1763                text+="      ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)
1764             text+="      if dim==2:\n"
1765             text+=unrollLoops(a_2,b_2,o,val,tap="        ")
1766             text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap="        ")
1767             text+="      else:\n"
1768            
1769             text+=unrollLoops(a_3,b_3,o,val,tap="        ")
1770             text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap="        ")
1771             if case in ["ContinuousFunction","Solution","ReducedSolution"]:
1772                 text+="      res=integrate(arg,where=w_ref)\n"
1773             else:
1774                 text+="      res=integrate(arg)\n"
1775    
1776             if data=="Symbol":
1777                text+="      sub=res.substitute({arg:s})\n"
1778             if len(sh)==0 and data=="Data":
1779                text+="      self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res
1780             else:
1781                if data=="Symbol":
1782                   text+="      self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"
1783                   text+="      self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
1784                else:
1785                   text+="      self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n"
1786                   text+="      self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)
1787             text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
1788    
1789    
1790             t_prog+=text
1791    print test_header
1792    print t_prog
1793    print test_tail          
1794    1/0
1795  #=======================================================================================================  #=======================================================================================================
1796  # inverse  # inverse
1797  #=======================================================================================================  #=======================================================================================================

Legend:
Removed from v.437  
changed lines
  Added in v.550

  ViewVC Help
Powered by ViewVC 1.1.26