/[escript]/trunk-mpi-branch/escript/py_src/util.py
ViewVC logotype

Diff of /trunk-mpi-branch/escript/py_src/util.py

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

revision 142 by jgs, Mon Jul 25 05:28:20 2005 UTC revision 147 by jgs, Fri Aug 12 01:45:47 2005 UTC
# Line 747  class WhereNegative_Symbol(Symbol): Line 747  class WhereNegative_Symbol(Symbol):
747     def eval(self,argval):     def eval(self,argval):
748         return whereNegative(self.getEvaluatedArguments(argval)[0])         return whereNegative(self.getEvaluatedArguments(argval)[0])
749    
750    def maximum(arg0,arg1):
751        """return arg1 where arg1 is bigger then arg0 otherwise arg0 is returned"""
752        m=whereNegative(arg0-arg1)
753        return m*arg1+(1.-m)*arg0
754      
755    def minimum(arg0,arg1):
756        """return arg0 where arg1 is bigger then arg0 otherwise arg1 is returned"""
757        m=whereNegative(arg0-arg1)
758        return m*arg0+(1.-m)*arg1
759      
760  def outer(arg0,arg1):  def outer(arg0,arg1):
761     if _testForZero(arg0) or _testForZero(arg1):     if _testForZero(arg0) or _testForZero(arg1):
762        return 0        return 0
# Line 809  def Interpolated_Symbol(Symbol): Line 819  def Interpolated_Symbol(Symbol):
819         a=self.getDifferentiatedArguments(arg)         a=self.getDifferentiatedArguments(arg)
820         return integrate(a[0],where=self.getArgument(1))         return integrate(a[0],where=self.getArgument(1))
821    
822    def div(arg,where=None):
823        """
824        @brief returns the divergence of arg at where.
825    
826        @param arg:   Data object representing the function which gradient to be calculated.
827        @param where: FunctionSpace in which the gradient will be calculated. If not present or
828                      None an appropriate default is used.
829        """
830        return trace(grad(arg,where))
831    
832  def grad(arg,where=None):  def grad(arg,where=None):
833      """      """
834      @brief returns the spatial gradient of arg at where.      @brief returns the spatial gradient of arg at where.
# Line 931  def trace(arg,axis0=0,axis1=1): Line 951  def trace(arg,axis0=0,axis1=1):
951         if s[axis0]!=s[axis1]:         if s[axis0]!=s[axis1]:
952             raise ValueError,"illegal axis in trace"             raise ValueError,"illegal axis in trace"
953         out=escript.Scalar(0.,arg.getFunctionSpace())         out=escript.Scalar(0.,arg.getFunctionSpace())
954         for i in range(s[0]):         for i in range(s[axis0]):
955            for j in range(s[1]):            out+=arg[i,i]
              out+=arg[i,j]  
956         return out         return out
957         # end hack for trace         # end hack for trace
        return arg.transpose(axis0=axis0,axis1=axis1)  
958      else:      else:
959         return numarray.trace(arg,axis0=axis0,axis1=axis1)         return numarray.trace(arg,axis0=axis0,axis1=axis1)
960    
# Line 954  def length(arg): Line 972  def length(arg):
972         if arg.getRank()==0:         if arg.getRank()==0:
973            return abs(arg)            return abs(arg)
974         elif arg.getRank()==1:         elif arg.getRank()==1:
975            sum=escript.Scalar(0,arg.getFunctionSpace())            out=escript.Scalar(0,arg.getFunctionSpace())
976            for i in range(arg.getShape()[0]):            for i in range(arg.getShape()[0]):
977               sum+=arg[i]**2               out+=arg[i]**2
978            return sqrt(sum)            return sqrt(out)
979         elif arg.getRank()==2:         elif arg.getRank()==2:
980            sum=escript.Scalar(0,arg.getFunctionSpace())            out=escript.Scalar(0,arg.getFunctionSpace())
981            for i in range(arg.getShape()[0]):            for i in range(arg.getShape()[0]):
982               for j in range(arg.getShape()[1]):               for j in range(arg.getShape()[1]):
983                  sum+=arg[i,j]**2                  out+=arg[i,j]**2
984            return sqrt(sum)            return sqrt(out)
985         elif arg.getRank()==3:         elif arg.getRank()==3:
986            sum=escript.Scalar(0,arg.getFunctionSpace())            out=escript.Scalar(0,arg.getFunctionSpace())
987            for i in range(arg.getShape()[0]):            for i in range(arg.getShape()[0]):
988               for j in range(arg.getShape()[1]):               for j in range(arg.getShape()[1]):
989                  for k in range(arg.getShape()[2]):                  for k in range(arg.getShape()[2]):
990                     sum+=arg[i,j,k]**2                     out+=arg[i,j,k]**2
991            return sqrt(sum)            return sqrt(out)
992         elif arg.getRank()==4:         elif arg.getRank()==4:
993            sum=escript.Scalar(0,arg.getFunctionSpace())            out=escript.Scalar(0,arg.getFunctionSpace())
994            for i in range(arg.getShape()[0]):            for i in range(arg.getShape()[0]):
995               for j in range(arg.getShape()[1]):               for j in range(arg.getShape()[1]):
996                  for k in range(arg.getShape()[2]):                  for k in range(arg.getShape()[2]):
997                     for l in range(arg.getShape()[3]):                     for l in range(arg.getShape()[3]):
998                        sum+=arg[i,j,k,l]**2                        out+=arg[i,j,k,l]**2
999            return sqrt(sum)            return sqrt(out)
1000         else:         else:
1001            raise SystemError,"length is not been fully implemented yet"            raise SystemError,"length is not been fully implemented yet"
1002            # return arg.length()            # return arg.length()
1003        elif isinstance(arg,float):
1004           return abs(arg)
1005      else:      else:
1006         return sqrt((arg**2).sum())         return sqrt((arg**2).sum())
1007    
# Line 1007  def inner(arg0,arg1): Line 1027  def inner(arg0,arg1):
1027    
1028      @param arg0, arg1      @param arg0, arg1
1029      """      """
1030      sum=escript.Scalar(0,arg0.getFunctionSpace())      if isinstance(arg0,escript.Data):
1031      if arg0.getRank()==0:         arg=arg0
1032        else:
1033           arg=arg1
1034    
1035        out=escript.Scalar(0,arg.getFunctionSpace())
1036        if arg.getRank()==0:
1037            return arg0*arg1            return arg0*arg1
1038      elif arg0.getRank()==1:      elif arg.getRank()==1:
1039           sum=escript.Scalar(0,arg0.getFunctionSpace())           out=escript.Scalar(0,arg.getFunctionSpace())
1040           for i in range(arg0.getShape()[0]):           for i in range(arg.getShape()[0]):
1041              sum+=arg0[i]*arg1[i]              out+=arg0[i]*arg1[i]
1042      elif arg0.getRank()==2:      elif arg.getRank()==2:
1043          sum=escript.Scalar(0,arg0.getFunctionSpace())          out=escript.Scalar(0,arg.getFunctionSpace())
1044          for i in range(arg0.getShape()[0]):          for i in range(arg.getShape()[0]):
1045             for j in range(arg0.getShape()[1]):             for j in range(arg.getShape()[1]):
1046                sum+=arg0[i,j]*arg1[i,j]                out+=arg0[i,j]*arg1[i,j]
1047      elif arg0.getRank()==3:      elif arg.getRank()==3:
1048          sum=escript.Scalar(0,arg0.getFunctionSpace())          out=escript.Scalar(0,arg.getFunctionSpace())
1049          for i in range(arg0.getShape()[0]):          for i in range(arg.getShape()[0]):
1050              for j in range(arg0.getShape()[1]):              for j in range(arg.getShape()[1]):
1051                 for k in range(arg0.getShape()[2]):                 for k in range(arg.getShape()[2]):
1052                    sum+=arg0[i,j,k]*arg1[i,j,k]                    out+=arg0[i,j,k]*arg1[i,j,k]
1053      elif arg0.getRank()==4:      elif arg.getRank()==4:
1054          sum=escript.Scalar(0,arg0.getFunctionSpace())          out=escript.Scalar(0,arg.getFunctionSpace())
1055          for i in range(arg0.getShape()[0]):          for i in range(arg.getShape()[0]):
1056             for j in range(arg0.getShape()[1]):             for j in range(arg.getShape()[1]):
1057                for k in range(arg0.getShape()[2]):                for k in range(arg.getShape()[2]):
1058                   for l in range(arg0.getShape()[3]):                   for l in range(arg.getShape()[3]):
1059                      sum+=arg0[i,j,k,l]*arg1[i,j,k,l]                      out+=arg0[i,j,k,l]*arg1[i,j,k,l]
1060      else:      else:
1061            raise SystemError,"inner is not been implemented yet"            raise SystemError,"inner is not been implemented yet"
1062      return sum      return out
1063    
1064  def matrixmult(arg0,arg1):  def matrixmult(arg0,arg1):
1065    
# Line 1053  def matrixmult(arg0,arg1): Line 1078  def matrixmult(arg0,arg1):
1078                 # uses Data object slicing, plus Data * and += operators                 # uses Data object slicing, plus Data * and += operators
1079                 out[i]+=arg0[i,j]*arg1[j]                 out[i]+=arg0[i,j]*arg1[j]
1080            return out            return out
1081          elif arg0.getRank()==1 and arg1.getRank()==1:
1082              return inner(arg0,arg1)
1083        else:        else:
1084            raise SystemError,"matrixmult is not fully implemented yet!"            raise SystemError,"matrixmult is not fully implemented yet!"
1085    
# Line 1134  def dot(arg0,arg1): Line 1161  def dot(arg0,arg1):
1161    
1162  def kronecker(d):  def kronecker(d):
1163     if hasattr(d,"getDim"):     if hasattr(d,"getDim"):
1164        return numarray.identity(d.getDim())        return numarray.identity(d.getDim())*1.
1165     else:     else:
1166        return numarray.identity(d)        return numarray.identity(d)*1.
1167    
1168  def unit(i,d):  def unit(i,d):
1169     """     """
# Line 1202  if __name__=="__main__": Line 1229  if __name__=="__main__":
1229    
1230  #  #
1231  # $Log$  # $Log$
1232  # Revision 1.14  2005/07/25 05:28:13  jgs  # Revision 1.15  2005/08/12 01:45:36  jgs
1233  # Merge of development branch back to main trunk on 2005-07-25  # erge of development branch dev-02 back to main trunk on 2005-08-12
1234  #  #
1235  # Revision 1.13  2005/07/22 03:53:01  jgs  # Revision 1.14.2.3  2005/08/03 09:55:33  gross
1236  # Merge of development branch back to main trunk on 2005-07-22  # ContactTest is passing now./mk install!
1237  #  #
1238  # Revision 1.12  2005/07/20 06:14:58  jgs  # Revision 1.14.2.2  2005/08/02 03:15:14  gross
1239  # added ln(data) style wrapper for data.ln() - also added corresponding  # bug inb trace fixed!
 # implementation of Ln_Symbol class (not sure if this is right though)  
1240  #  #
1241  # Revision 1.11  2005/07/08 04:07:35  jgs  # Revision 1.14.2.1  2005/07/29 07:10:28  gross
1242  # Merge of development branch back to main trunk on 2005-07-08  # new functions in util and a new pde type in linearPDEs
1243  #  #
1244  # Revision 1.10  2005/06/09 05:37:59  jgs  # Revision 1.2.2.21  2005/07/28 04:19:23  gross
1245  # Merge of development branch back to main trunk on 2005-06-09  # new functions maximum and minimum introduced.
1246  #  #
1247  # Revision 1.2.2.20  2005/07/25 01:26:27  gross  # Revision 1.2.2.20  2005/07/25 01:26:27  gross
1248  # bug in inner fixed  # bug in inner fixed

Legend:
Removed from v.142  
changed lines
  Added in v.147

  ViewVC Help
Powered by ViewVC 1.1.26