# Diff of /trunk/escript/py_src/generatediff

revision 437 by gross, Fri Jan 20 00:16:58 2006 UTC revision 438 by gross, Fri Jan 20 00:39:00 2006 UTC
# Line 796  def unrollLoopsOfGrad(a,b,o,arg,tap=""): Line 796  def unrollLoopsOfGrad(a,b,o,arg,tap=""):
796                   out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99])                   out+=tap+"%s[%s,%s,%s,%s]=(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99]+b[i0,i1,i2,i99])
797                 else:                 else:
798                   out+=tap+"%s[%s,%s,%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99],i99,b[i0,i1,i2,i99])                   out+=tap+"%s[%s,%s,%s,%s]=o*(%s)*x_ref[%s]**(o-1)+(%s)\n"%(arg,i0,i1,i2,i99,a[i0,i1,i2,i99],i99,b[i0,i1,i2,i99])
799        return out
800  def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):  def unrollLoopsOfInteriorIntegral(a,b,where,arg,tap=""):
801      if where=="Function":      if where=="Function":
802         xfac_o=1.         xfac_o=1.
# Line 918  def unrollLoopsOfInteriorIntegral(a,b,wh Line 918  def unrollLoopsOfInteriorIntegral(a,b,wh
918      return out      return out
919
920

921  #=======================================================================================================  #=======================================================================================================
922  # integrate  # grad
923  #=======================================================================================================  #=======================================================================================================
924  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
925    for data in ["Data","Symbol"]:    for data in ["Data","Symbol"]:
926      for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:       for case in ["ContinuousFunction","Solution","ReducedSolution"]:
927        for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:         for sh in [ (),(2,), (4,5), (6,2,2)]:
if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:
928           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
929           tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))           tname="test_grad_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
930           text+="   def %s(self):\n"%tname           text+="   def %s(self):\n"%tname
931           text+="      \"\"\"\n"           text+="      \"\"\"\n"
932           text+="      tests integral of 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)
933           text+="      assumptions: %s(self.domain) exists\n"%case           text+="      assumptions: %s(self.domain) exists\n"%case
934           text+="                   self.domain supports integral on %s\n"%where           text+="                   self.domain supports gardient on %s\n"%where

935           text+="      \"\"\"\n"           text+="      \"\"\"\n"
936           if case=="ReducedSolution":           if case=="ReducedSolution":
937              text+="      o=1\n"              text+="      o=1\n"
# Line 944  for where in ["Function","FunctionOnBoun Line 941  for where in ["Function","FunctionOnBoun
941              o="o"              o="o"
942           text+="      dim=self.domain.getDim()\n"           text+="      dim=self.domain.getDim()\n"
943           text+="      w_ref=%s(self.domain)\n"%where           text+="      w_ref=%s(self.domain)\n"%where
944             text+="      x_ref=w_ref.getX()\n"
945           text+="      w=%s(self.domain)\n"%case           text+="      w=%s(self.domain)\n"%case
946           text+="      x=w.getX()\n"           text+="      x=w.getX()\n"
947           a_2=makeArray(sh+(2,),[-1.,1])           a_2=makeArray(sh+(2,),[-1.,1])
# Line 951  for where in ["Function","FunctionOnBoun Line 949  for where in ["Function","FunctionOnBoun
949           a_3=makeArray(sh+(3,),[-1.,1])           a_3=makeArray(sh+(3,),[-1.,1])
950           b_3=makeArray(sh+(3,),[-1.,1])           b_3=makeArray(sh+(3,),[-1.,1])
951           if data=="Symbol":           if data=="Symbol":
952              text+="      arg=Symbol(shape=%s)\n"%str(sh)              text+="      arg=Symbol(shape=%s,dim=dim)\n"%str(sh)
953              val="s"              val="s"
954              res="sub"              res="sub"
955           else:           else:
956              val="arg"              val="arg"
957              res="res"              res="res"

958           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))
959           if not len(sh)==0:           text+="      ref=Data(0,%s+(dim,),w_ref)\n"%str(sh)
text+="      ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)
960           text+="      if dim==2:\n"           text+="      if dim==2:\n"
961           text+=unrollLoops(a_2,b_2,o,val,tap="        ")           text+=unrollLoops(a_2,b_2,o,val,tap="        ")
962           text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap="        ")           text+=unrollLoopsOfGrad(a_2,b_2,o,"ref",tap="        ")
963           text+="      else:\n"           text+="      else:\n"
964
965           text+=unrollLoops(a_3,b_3,o,val,tap="        ")           text+=unrollLoops(a_3,b_3,o,val,tap="        ")
966           text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap="        ")           text+=unrollLoopsOfGrad(a_3,b_3,o,"ref",tap="        ")
967           if case in ["ContinuousFunction","Solution","ReducedSolution"]:           text+="      res=grad(arg,where=w_ref)\n"
text+="      res=integrate(arg,where=w_ref)\n"
else:
text+="      res=integrate(arg)\n"

968           if data=="Symbol":           if data=="Symbol":
969              text+="      sub=res.substitute({arg:s})\n"              text+="      sub=res.substitute({arg:s})\n"
970           if len(sh)==0 and data=="Data":           text+="      self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data
971              text+="      self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res           text+="      self.failUnlessEqual(res.getShape(),%s+(dim,),\"wrong shape of result.\")\n"%str(sh)
else:
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)
972           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
973
974
# Line 992  print test_header Line 977  print test_header
977  print t_prog  print t_prog
978  print test_tail            print test_tail
979  1/0  1/0
980
981
982  #=======================================================================================================  #=======================================================================================================
983  # grad  # integrate
984  #=======================================================================================================  #=======================================================================================================
985  for where in ["Function","FunctionOnBoundary"]:  for where in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
986    for data in ["Data","Symbol"]:    for data in ["Data","Symbol"]:
987       for case in ["ContinuousFunction","Solution","ReducedSolution"]:      for case in ["ContinuousFunction","Solution","ReducedSolution","Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]:
988         for sh in [ (),(2,), (4,5), (6,2,2)]:        for sh in [ (),(2,), (4,5), (6,2,2),(4,5,3,2)]:
989            if (not case in ["Function","FunctionOnBoundary","FunctionOnContactZero","FunctionOnContactOne"]) or where==case:
990           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"           text="   #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
991           tname="test_grad_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))           tname="test_integrate_on%s_from%s_%s_rank%s"%(where,data,case,len(sh))
992           text+="   def %s(self):\n"%tname           text+="   def %s(self):\n"%tname
993           text+="      \"\"\"\n"           text+="      \"\"\"\n"
994           text+="      tests gradient for rank %s %s on the %s\n\n"%(len(sh),data,where)           text+="      tests integral of rank %s %s on the %s\n\n"%(len(sh),data,where)
995           text+="      assumptions: %s(self.domain) exists\n"%where           text+="      assumptions: %s(self.domain) exists\n"%case
996           if where=="FunctionOnBoundary": text+="                   self.domain supports gradient on the boundary\n"           text+="                   self.domain supports integral on %s\n"%where
997
998           text+="      \"\"\"\n"           text+="      \"\"\"\n"
999           if case=="ReducedSolution":           if case=="ReducedSolution":
1000              text+="      o=1\n"              text+="      o=1\n"
# Line 1015  for where in ["Function","FunctionOnBoun Line 1004  for where in ["Function","FunctionOnBoun
1004              o="o"              o="o"
1005           text+="      dim=self.domain.getDim()\n"           text+="      dim=self.domain.getDim()\n"
1006           text+="      w_ref=%s(self.domain)\n"%where           text+="      w_ref=%s(self.domain)\n"%where
text+="      x_ref=w_ref.getX()\n"
1007           text+="      w=%s(self.domain)\n"%case           text+="      w=%s(self.domain)\n"%case
1008           text+="      x=w.getX()\n"           text+="      x=w.getX()\n"
1009           a_2=makeArray(sh+(2,),[-1.,1])           a_2=makeArray(sh+(2,),[-1.,1])
# Line 1023  for where in ["Function","FunctionOnBoun Line 1011  for where in ["Function","FunctionOnBoun
1011           a_3=makeArray(sh+(3,),[-1.,1])           a_3=makeArray(sh+(3,),[-1.,1])
1012           b_3=makeArray(sh+(3,),[-1.,1])           b_3=makeArray(sh+(3,),[-1.,1])
1013           if data=="Symbol":           if data=="Symbol":
1014              text+="      arg=Symbol(shape=%s,dim=dim)\n"%str(sh)              text+="      arg=Symbol(shape=%s)\n"%str(sh)
1015              val="s"              val="s"
1016              res="sub"              res="sub"
1017           else:           else:
1018              val="arg"              val="arg"
1019              res="res"              res="res"
1020
1021           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))           text+="      %s=Data(0,%s,w)\n"%(val,str(sh))
1022           text+="      ref=Data(0,%s+(dim,),w_ref)\n"%str(sh)           if not len(sh)==0:
1023                text+="      ref=numarray.zeros(%s,numarray.Float)\n"%str(sh)
1024           text+="      if dim==2:\n"           text+="      if dim==2:\n"
1025           text+=unrollLoops(a_2,b_2,o,val,tap="        ")           text+=unrollLoops(a_2,b_2,o,val,tap="        ")
1026           text+=unrollLoopsOfGrad(a_2,b_2,o,"ref",tap="        ")           text+=unrollLoopsOfInteriorIntegral(a_2,b_2,where,"ref",tap="        ")
1027           text+="      else:\n"           text+="      else:\n"
1028
1029           text+=unrollLoops(a_3,b_3,o,val,tap="        ")           text+=unrollLoops(a_3,b_3,o,val,tap="        ")
1030           text+=unrollLoopsOfGrad(a_3,b_3,o,"ref",tap="        ")           text+=unrollLoopsOfInteriorIntegral(a_3,b_3,where,"ref",tap="        ")
1031           text+="      res=grad(arg,where=w_ref)\n"           if case in ["ContinuousFunction","Solution","ReducedSolution"]:
1032                 text+="      res=integrate(arg,where=w_ref)\n"
1033             else:
1034                 text+="      res=integrate(arg)\n"
1035
1036           if data=="Symbol":           if data=="Symbol":
1037              text+="      sub=res.substitute({arg:s})\n"              text+="      sub=res.substitute({arg:s})\n"
1038           text+="      self.failUnless(isinstance(res,%s),\"wrong type of result.\")\n"%data           if len(sh)==0 and data=="Data":
1039           text+="      self.failUnlessEqual(res.getShape(),%s+(dim,),\"wrong shape of result.\")\n"%str(sh)              text+="      self.failUnless(isinstance(%s,float),\"wrong type of result.\")\n"%res
1040             else:
1041                if data=="Symbol":
1042                   text+="      self.failUnless(isinstance(res,Symbol),\"wrong type of result.\")\n"
1043                   text+="      self.failUnlessEqual(res.getShape(),%s,\"wrong shape of result.\")\n"%str(sh)
1044                else:
1045                   text+="      self.failUnless(isinstance(res,numarray.NumArray),\"wrong type of result.\")\n"
1046                   text+="      self.failUnlessEqual(res.shape,%s,\"wrong shape of result.\")\n"%str(sh)
1047           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
1048
1049
# Line 1051  print test_header Line 1052  print test_header
1052  print t_prog  print t_prog
1053  print test_tail            print test_tail
1054  1/0  1/0

1055  #=======================================================================================================  #=======================================================================================================
1056  # inverse  # inverse
1057  #=======================================================================================================  #=======================================================================================================

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

 ViewVC Help Powered by ViewVC 1.1.26