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

Diff of /trunk/escript/py_src/generatediff

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

revision 291 by gross, Fri Dec 2 03:10:06 2005 UTC revision 588 by gross, Fri Mar 10 04:45:04 2006 UTC
# Line 33  def wherepos(arg): Line 33  def wherepos(arg):
33     else:     else:
34        return 0.        return 0.
35    
36    
37  class OPERATOR:  class OPERATOR:
38      def __init__(self,nickname,rng=[-1000.,1000],test_expr="",math_expr=None,      def __init__(self,nickname,rng=[-1000.,1000],test_expr="",math_expr=None,
39                  numarray_expr="",symbol_expr=None,diff=None,name=""):                  numarray_expr="",symbol_expr=None,diff=None,name=""):
# Line 172  def makeArray(shape,rng): Line 173  def makeArray(shape,rng):
173               for i2 in range(shape[2]):               for i2 in range(shape[2]):
174                  for i3 in range(shape[3]):                  for i3 in range(shape[3]):
175                     out[i0,i1,i2,i3]=l*random.random()+rng[0]                     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:     else:
184         raise SystemError,"rank is restricted to 4"         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 472  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 489  def mkText(case,name,a,a1=None,use_taggi Line 522  def mkText(case,name,a,a1=None,use_taggi
522    
523           return t_out           return t_out
524    
525  def mkTypeAndShapeTest(case,sh,argstr):  def mkTypeAndShapeTest(case,sh,argstr,name=""):
526      text=""      text=""
527      if case=="float":      if case=="float":
528           text+="      self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%argstr           text+="      self.failUnless(isinstance(%s,float),\"wrong type of result%s.\")\n"%(argstr,name)
529      elif case=="array":      elif case=="array":
530           text+="      self.failUnless(isinstance(%s,numarray.NumArray),\"wrong type of result.\")\n"%argstr           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.\")\n"%(argstr,str(sh))           text+="      self.failUnlessEqual(%s.shape,%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name)
532      elif case in ["constData","taggedData","expandedData"]:          elif case in ["constData","taggedData","expandedData"]:    
533           text+="      self.failUnless(isinstance(%s,Data),\"wrong type of result.\")\n"%argstr           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.\")\n"%(argstr,str(sh))           text+="      self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name)
535      elif case=="Symbol":      elif case=="Symbol":
536           text+="      self.failUnless(isinstance(%s,Symbol),\"wrong type of result.\")\n"%argstr           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.\")\n"%(argstr,str(sh))           text+="      self.failUnlessEqual(%s.getShape(),%s,\"wrong shape of result%s.\")\n"%(argstr,str(sh),name)
538      return text      return text
539    
540  def mkCode(txt,args=[],intend=""):  def mkCode(txt,args=[],intend=""):
# Line 517  def mkCode(txt,args=[],intend=""): Line 550  def mkCode(txt,args=[],intend=""):
550        out=out.replace("%%a%s%%"%c,r)        out=out.replace("%%a%s%%"%c,r)
551      return out        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):  def innerTEST(arg0,arg1):
1017      if isinstance(arg0,float):      if isinstance(arg0,float):
1018         out=numarray.array(arg0*arg1)         out=numarray.array(arg0*arg1)
# Line 589  def testTensorMult(arg0,arg1,sh_s): Line 1085  def testTensorMult(arg0,arg1,sh_s):
1085               for i3 in range(arg0.shape[3]):               for i3 in range(arg0.shape[3]):
1086                       out[i0,i1]+=arg0[i0,i1,i2,i3]*arg1[i2,i3]                       out[i0,i1]+=arg0[i0,i1,i2,i3]*arg1[i2,i3]
1087          return out          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  # tensor multiply  # transpose
1513  #=======================================================================================================  #=======================================================================================================
1514  # oper=["generalTensorProduct",tensorProductTest]  def transposeTest(r,offset):
1515  # oper=["matrixmult",testMatrixMult]      if isinstance(r,float): return r
1516  oper=["tensormult",testTensorMult]      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  for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:  #              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)]:    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"]:     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
2040       for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:       for sh0 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
2041         for sh_s in [ (),(3,), (2,3), (2,4,3),(4,2,3,2)]:          if (len(sh0)==0 or not case0=="float") and (len(sh1)==0 or not case1=="float") \
2042            if (len(sh0+sh_s)==0 or not case0=="float") and (len(sh1+sh_s)==0 or not case1=="float") \             and (sh0==sh1 or len(sh0)==0 or len(sh1)==0) :
                and len(sh0+sh1)<5 and len(sh0+sh_s)<5 and len(sh1+sh_s)<5:  
             # 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  
             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  
               case=getResultCaseForBin(case0,case1)    
2043                use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"                use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
2044    
2045                text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"                text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
2046                # tname="test_generalTensorProduct_%s_rank%s_%s_rank%s_offset%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1),len(sh_s))                tname="test_%s_%s_rank%s_%s_rank%s"%(oper[0],case0,len(sh0),case1,len(sh1))
               #tname="test_matrixmult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))  
               tname="test_tensormult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))  
               # if tname=="test_generalTensorProduct_array_rank1_array_rank2_offset1":  
               # print tnametest_generalTensorProduct_Symbol_rank1_Symbol_rank3_offset1  
2047                text+="   def %s(self):\n"%tname                text+="   def %s(self):\n"%tname
2048                a_0=makeArray(sh0+sh_s,[-1.,1])                a_0=makeArray(sh0,[-1.,1])
2049                if case0 in ["taggedData", "expandedData"]:                if case0 in ["taggedData", "expandedData"]:
2050                    a1_0=makeArray(sh0+sh_s,[-1.,1])                    a1_0=makeArray(sh0,[-1.,1])
2051                else:                else:
2052                    a1_0=a_0                    a1_0=a_0
2053    
2054                a_1=makeArray(sh_s+sh1,[-1.,1])                a_1=makeArray(sh1,[-1.,1])
2055                if case1 in ["taggedData", "expandedData"]:                if case1 in ["taggedData", "expandedData"]:
2056                    a1_1=makeArray(sh_s+sh1,[-1.,1])                    a1_1=makeArray(sh1,[-1.,1])
2057                else:                else:
2058                    a1_1=a_1                    a1_1=a_1
2059                r=oper[1](a_0,a_1,sh_s)                r=oper[1](a_0,a_1)
2060                r1=oper[1](a1_0,a1_1,sh_s)                r1=oper[1](a1_0,a1_1)
2061                text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)                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)                text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
2063                #text+="      res=matrixmult(arg0,arg1)\n"                text+="      res=%s(arg0,arg1)\n"%oper[0]
2064                text+="      res=tensormult(arg0,arg1)\n"                case=getResultCaseForBin(case0,case1)              
               #text+="      res=generalTensorProduct(arg0,arg1,offset=%s)\n"%(len(sh_s))  
2065                if case=="Symbol":                if case=="Symbol":
2066                   c0_res,c1_res=case0,case1                   c0_res,c1_res=case0,case1
2067                   subs="{"                   subs="{"
# Line 651  for case0 in ["float","array","Symbol"," Line 2081  for case0 in ["float","array","Symbol","
2081                else:                else:
2082                   res="res"                   res="res"
2083                   text+=mkText(case,"ref",r,r1)                   text+=mkText(case,"ref",r,r1)
2084                text+=mkTypeAndShapeTest(case,sh0+sh1,"res")                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                text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
2089                  
2090                if case0 == "taggedData" or case1 == "taggedData":                if case0 == "taggedData" or case1 == "taggedData":
2091                    t_prog_with_tags+=text                    t_prog_with_tags+=text
2092                else:                              else:              
2093                    t_prog+=text                    t_prog+=text
2094    
2095  print test_header  print test_header
2096  # print t_prog  # print t_prog
2097  print t_prog_with_tags  print t_prog_with_tags
2098  print test_tail            print test_tail          
2099  1/0  1/0
2100    
2101    
2102  #=======================================================================================================  #=======================================================================================================
2103  # outer/inner  # outer inner
2104  #=======================================================================================================  #=======================================================================================================
2105  oper=["inner",innerTEST]  oper=["outer",outerTEST]
2106  # oper=["outer",outerTEST]  # oper=["inner",innerTEST]
2107  for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:  for case0 in ["float","array","Symbol","constData","taggedData","expandedData"]:
2108    for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:    for sh1 in [ (),(2,), (4,5), (6,2,2),(3,2,3,4)]:
2109     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:     for case1 in ["float","array","Symbol","constData","taggedData","expandedData"]:
# Line 725  for case0 in ["float","array","Symbol"," Line 2162  for case0 in ["float","array","Symbol","
2162  print test_header  print test_header
2163  # print t_prog  # print t_prog
2164  print t_prog_with_tags  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                  use_tagging_for_expanded_data= case0=="taggedData" or case1=="taggedData"
2232                  text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
2233                  # tname="test_generalTensorProduct_%s_rank%s_%s_rank%s_offset%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1),len(sh_s))
2234                  #tname="test_matrixmult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
2235                  tname="test_tensormult_%s_rank%s_%s_rank%s"%(case0,len(sh0+sh_s),case1,len(sh_s+sh1))
2236                  # if tname=="test_generalTensorProduct_array_rank1_array_rank2_offset1":
2237                  # print tnametest_generalTensorProduct_Symbol_rank1_Symbol_rank3_offset1
2238                  text+="   def %s(self):\n"%tname
2239                  a_0=makeArray(sh0+sh_s,[-1.,1])
2240                  if case0 in ["taggedData", "expandedData"]:
2241                      a1_0=makeArray(sh0+sh_s,[-1.,1])
2242                  else:
2243                      a1_0=a_0
2244    
2245                  a_1=makeArray(sh_s+sh1,[-1.,1])
2246                  if case1 in ["taggedData", "expandedData"]:
2247                      a1_1=makeArray(sh_s+sh1,[-1.,1])
2248                  else:
2249                      a1_1=a_1
2250                  r=oper[1](a_0,a_1,sh_s)
2251                  r1=oper[1](a1_0,a1_1,sh_s)
2252                  text+=mkText(case0,"arg0",a_0,a1_0,use_tagging_for_expanded_data)
2253                  text+=mkText(case1,"arg1",a_1,a1_1,use_tagging_for_expanded_data)
2254                  #text+="      res=matrixmult(arg0,arg1)\n"
2255                  text+="      res=tensormult(arg0,arg1)\n"
2256                  #text+="      res=generalTensorProduct(arg0,arg1,offset=%s)\n"%(len(sh_s))
2257                  if case=="Symbol":
2258                     c0_res,c1_res=case0,case1
2259                     subs="{"
2260                     if case0=="Symbol":        
2261                        text+=mkText("array","s0",a_0,a1_0)
2262                        subs+="arg0:s0"
2263                        c0_res="array"
2264                     if case1=="Symbol":        
2265                        text+=mkText("array","s1",a_1,a1_1)
2266                        if not subs.endswith("{"): subs+=","
2267                        subs+="arg1:s1"
2268                        c1_res="array"
2269                     subs+="}"  
2270                     text+="      sub=res.substitute(%s)\n"%subs
2271                     res="sub"
2272                     text+=mkText(getResultCaseForBin(c0_res,c1_res),"ref",r,r1)
2273                  else:
2274                     res="res"
2275                     text+=mkText(case,"ref",r,r1)
2276                  text+=mkTypeAndShapeTest(case,sh0+sh1,"res")
2277                  text+="      self.failUnless(Lsup(%s-ref)<=self.RES_TOL*Lsup(ref),\"wrong result\")\n"%res
2278                  if case0 == "taggedData" or case1 == "taggedData":
2279                      t_prog_with_tags+=text
2280                  else:              
2281                      t_prog+=text
2282    print test_header
2283    # print t_prog
2284    print t_prog_with_tags
2285  print test_tail            print test_tail          
2286  1/0  1/0
2287  #=======================================================================================================  #=======================================================================================================

Legend:
Removed from v.291  
changed lines
  Added in v.588

  ViewVC Help
Powered by ViewVC 1.1.26