/[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

trunk/escript/py_src/util.py revision 804 by gross, Thu Aug 10 01:12:16 2006 UTC trunk-mpi-branch/escript/py_src/util.py revision 1295 by ksteube, Mon Sep 10 06:07:09 2007 UTC
# Line 26  import math Line 26  import math
26  import numarray  import numarray
27  import escript  import escript
28  import os  import os
29    from esys.escript import C_GeneralTensorProduct
30    from esys.escript import getVersion
31    
32  #=========================================================  #=========================================================
33  #   some helpers:  #   some helpers:
34  #=========================================================  #=========================================================
35    def getTagNames(domain):
36        """
37        returns a list of the tag names used by the domain
38    
39        
40        @param domain: a domain object
41        @type domain: L{escript.Domain}
42        @return: a list of the tag name used by the domain.
43        @rtype: C{list} of C{str}
44        """
45        return [n.strip() for n in domain.showTagNames().split(",") ]
46    
47    def insertTagNames(domain,**kwargs):
48        """
49        inserts tag names into the domain
50    
51        @param domain: a domain object
52        @type domain: C{escript.Domain}
53        @keyword <tag name>: tag key assigned to <tag name>
54        @type <tag name>: C{int}
55        """
56        for  k in kwargs:
57             domain.setTagMap(k,kwargs[k])
58    
59    def insertTaggedValues(target,**kwargs):
60        """
61        inserts tagged values into the tagged using tag names
62    
63        @param target: data to be filled by tagged values
64        @type target: L{escript.Data}
65        @keyword <tag name>: value to be used for <tag name>
66        @type <tag name>: C{float} or {numarray.NumArray}
67        @return: C{target}
68        @rtype: L{escript.Data}
69        """
70        for k in kwargs:
71            target.setTaggedValue(k,kwargs[k])
72        return target
73    
74  def saveVTK(filename,domain=None,**data):  def saveVTK(filename,domain=None,**data):
75      """      """
76      writes a L{Data} objects into a files using the the VTK XML file format.      writes a L{Data} objects into a files using the the VTK XML file format.
# Line 50  def saveVTK(filename,domain=None,**data) Line 91  def saveVTK(filename,domain=None,**data)
91      @type <name>: L{Data} object.      @type <name>: L{Data} object.
92      @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.      @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.
93      """      """
94      if domain==None:      new_data={}
95         for i in data.keys():      for n,d in data.items():
96            if not data[i].isEmpty(): domain=data[i].getFunctionSpace().getDomain()            if not d.isEmpty():
97                fs=d.getFunctionSpace()
98                domain2=fs.getDomain()
99                if fs == escript.Solution(domain2):
100                   new_data[n]=interpolate(d,escript.ContinuousFunction(domain2))
101                elif fs == escript.ReducedSolution(domain2):
102                   new_data[n]=interpolate(d,escript.ReducedContinuousFunction(domain2))
103                else:
104                   new_data[n]=d
105                if domain==None: domain=domain2
106      if domain==None:      if domain==None:
107          raise ValueError,"no domain detected."          raise ValueError,"no domain detected."
108      else:      domain.saveVTK(filename,new_data)
         domain.saveVTK(filename,data)  
109    
110  def saveDX(filename,domain=None,**data):  def saveDX(filename,domain=None,**data):
111      """      """
# Line 78  def saveDX(filename,domain=None,**data): Line 127  def saveDX(filename,domain=None,**data):
127      @type <name>: L{Data} object.      @type <name>: L{Data} object.
128      @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.      @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.
129      """      """
130      if domain==None:      new_data={}
131         for i in data.keys():      for n,d in data.items():
132            if not data[i].isEmpty(): domain=data[i].getFunctionSpace().getDomain()            if not d.isEmpty():
133                fs=d.getFunctionSpace()
134                domain2=fs.getDomain()
135                if fs == escript.Solution(domain2):
136                   new_data[n]=interpolate(d,escript.ReducedContinuousFunction(domain2))
137                elif fs == escript.ReducedSolution(domain2):
138                   new_data[n]=interpolate(d,escript.ReducedContinuousFunction(domain2))
139                elif fs == escript.ContinuousFunction(domain2):
140                   new_data[n]=interpolate(d,escript.ReducedContinuousFunction(domain2))
141                else:
142                   new_data[n]=d
143                if domain==None: domain=domain2
144      if domain==None:      if domain==None:
145          raise ValueError,"no domain detected."          raise ValueError,"no domain detected."
146      else:      domain.saveDX(filename,new_data)
         domain.saveDX(filename,data)  
147    
148  def kronecker(d=3):  def kronecker(d=3):
149     """     """
# Line 238  def inf(arg): Line 297  def inf(arg):
297  #=========================================================================  #=========================================================================
298  #   some little helpers  #   some little helpers
299  #=========================================================================  #=========================================================================
300  def pokeShape(arg):  def getRank(arg):
301        """
302        identifies the rank of its argument
303    
304        @param arg: a given object
305        @type arg: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, C{Symbol}
306        @return: the rank of the argument
307        @rtype: C{int}
308        @raise TypeError: if type of arg cannot be processed
309        """
310    
311        if isinstance(arg,numarray.NumArray):
312            return arg.rank
313        elif isinstance(arg,escript.Data):
314            return arg.getRank()
315        elif isinstance(arg,float):
316            return 0
317        elif isinstance(arg,int):
318            return 0
319        elif isinstance(arg,Symbol):
320            return arg.getRank()
321        else:
322          raise TypeError,"getShape: cannot identify shape"
323    def getShape(arg):
324      """      """
325      identifies the shape of its argument      identifies the shape of its argument
326    
# Line 260  def pokeShape(arg): Line 342  def pokeShape(arg):
342      elif isinstance(arg,Symbol):      elif isinstance(arg,Symbol):
343          return arg.getShape()          return arg.getShape()
344      else:      else:
345        raise TypeError,"pokeShape: cannot identify shape"        raise TypeError,"getShape: cannot identify shape"
346    
347  def pokeDim(arg):  def pokeDim(arg):
348      """      """
# Line 283  def commonShape(arg0,arg1): Line 365  def commonShape(arg0,arg1):
365      """      """
366      returns a shape to which arg0 can be extendent from the right and arg1 can be extended from the left.      returns a shape to which arg0 can be extendent from the right and arg1 can be extended from the left.
367    
368      @param arg0: an object with a shape (see L{pokeShape})      @param arg0: an object with a shape (see L{getShape})
369      @param arg1: an object with a shape (see L{pokeShape})      @param arg1: an object with a shape (see L{getShape})
370      @return: the shape of arg0 or arg1 such that the left port equals the shape of arg0 and the right end equals the shape of arg1.      @return: the shape of arg0 or arg1 such that the left port equals the shape of arg0 and the right end equals the shape of arg1.
371      @rtype: C{tuple} of C{int}      @rtype: C{tuple} of C{int}
372      @raise ValueError: if no shape can be found.      @raise ValueError: if no shape can be found.
373      """      """
374      sh0=pokeShape(arg0)      sh0=getShape(arg0)
375      sh1=pokeShape(arg1)      sh1=getShape(arg1)
376      if len(sh0)<len(sh1):      if len(sh0)<len(sh1):
377         if not sh0==sh1[:len(sh0)]:         if not sh0==sh1[:len(sh0)]:
378               raise ValueError,"argument 0 cannot be extended to the shape of argument 1"               raise ValueError,"argument 0 cannot be extended to the shape of argument 1"
# Line 444  def matchShape(arg0,arg1): Line 526  def matchShape(arg0,arg1):
526      @rtype: C{tuple}      @rtype: C{tuple}
527      """      """
528      sh=commonShape(arg0,arg1)      sh=commonShape(arg0,arg1)
529      sh0=pokeShape(arg0)      sh0=getShape(arg0)
530      sh1=pokeShape(arg1)      sh1=getShape(arg1)
531      if len(sh0)<len(sh):      if len(sh0)<len(sh):
532         return outer(arg0,numarray.ones(sh[len(sh0):],numarray.Float64)),arg1         return outer(arg0,numarray.ones(sh[len(sh0):],numarray.Float64)),arg1
533      elif len(sh1)<len(sh):      elif len(sh1)<len(sh):
# Line 573  class Symbol(object): Line 655  class Symbol(object):
655            if isinstance(a,Symbol):            if isinstance(a,Symbol):
656               out.append(a.substitute(argvals))               out.append(a.substitute(argvals))
657            else:            else:
658                s=pokeShape(s)+arg.getShape()                s=getShape(s)+arg.getShape()
659                if len(s)>0:                if len(s)>0:
660                   out.append(numarray.zeros(s),numarray.Float64)                   out.append(numarray.zeros(s),numarray.Float64)
661                else:                else:
# Line 1310  def whereNonZero(arg,tol=0.): Line 1392  def whereNonZero(arg,tol=0.):
1392     else:     else:
1393        raise TypeError,"whereNonZero: Unknown argument type."        raise TypeError,"whereNonZero: Unknown argument type."
1394    
1395    def erf(arg):
1396       """
1397       returns erf of argument arg
1398    
1399       @param arg: argument
1400       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1401       @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1402       @raises TypeError: if the type of the argument is not expected.
1403       """
1404       if isinstance(arg,escript.Data):
1405          return arg._erf()
1406       else:
1407          raise TypeError,"erf: Unknown argument type."
1408    
1409  def sin(arg):  def sin(arg):
1410     """     """
1411     returns sine of argument arg     returns sine of argument arg
# Line 3645  class Add_Symbol(DependendSymbol): Line 3741  class Add_Symbol(DependendSymbol):
3741         @raise ValueError: if both arguments do not have the same shape.         @raise ValueError: if both arguments do not have the same shape.
3742         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
3743         """         """
3744         sh0=pokeShape(arg0)         sh0=getShape(arg0)
3745         sh1=pokeShape(arg1)         sh1=getShape(arg1)
3746         if not sh0==sh1:         if not sh0==sh1:
3747            raise ValueError,"Add_Symbol: shape of arguments must match"            raise ValueError,"Add_Symbol: shape of arguments must match"
3748         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])
# Line 3720  def mult(arg0,arg1): Line 3816  def mult(arg0,arg1):
3816         """         """
3817         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
3818         if testForZero(args[0]) or testForZero(args[1]):         if testForZero(args[0]) or testForZero(args[1]):
3819            return numarray.zeros(pokeShape(args[0]),numarray.Float64)            return numarray.zeros(getShape(args[0]),numarray.Float64)
3820         else:         else:
3821            if isinstance(args[0],Symbol) or isinstance(args[1],Symbol) :            if isinstance(args[0],Symbol) or isinstance(args[1],Symbol) :
3822                return Mult_Symbol(args[0],args[1])                return Mult_Symbol(args[0],args[1])
# Line 3744  class Mult_Symbol(DependendSymbol): Line 3840  class Mult_Symbol(DependendSymbol):
3840         @raise ValueError: if both arguments do not have the same shape.         @raise ValueError: if both arguments do not have the same shape.
3841         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
3842         """         """
3843         sh0=pokeShape(arg0)         sh0=getShape(arg0)
3844         sh1=pokeShape(arg1)         sh1=getShape(arg1)
3845         if not sh0==sh1:         if not sh0==sh1:
3846            raise ValueError,"Mult_Symbol: shape of arguments must match"            raise ValueError,"Mult_Symbol: shape of arguments must match"
3847         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])
# Line 3820  def quotient(arg0,arg1): Line 3916  def quotient(arg0,arg1):
3916         """         """
3917         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
3918         if testForZero(args[0]):         if testForZero(args[0]):
3919            return numarray.zeros(pokeShape(args[0]),numarray.Float64)            return numarray.zeros(getShape(args[0]),numarray.Float64)
3920         elif isinstance(args[0],Symbol):         elif isinstance(args[0],Symbol):
3921            if isinstance(args[1],Symbol):            if isinstance(args[1],Symbol):
3922               return Quotient_Symbol(args[0],args[1])               return Quotient_Symbol(args[0],args[1])
# Line 3849  class Quotient_Symbol(DependendSymbol): Line 3945  class Quotient_Symbol(DependendSymbol):
3945         @raise ValueError: if both arguments do not have the same shape.         @raise ValueError: if both arguments do not have the same shape.
3946         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
3947         """         """
3948         sh0=pokeShape(arg0)         sh0=getShape(arg0)
3949         sh1=pokeShape(arg1)         sh1=getShape(arg1)
3950         if not sh0==sh1:         if not sh0==sh1:
3951            raise ValueError,"Quotient_Symbol: shape of arguments must match"            raise ValueError,"Quotient_Symbol: shape of arguments must match"
3952         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0,args=[arg0,arg1])
# Line 3926  def power(arg0,arg1): Line 4022  def power(arg0,arg1):
4022         """         """
4023         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
4024         if testForZero(args[0]):         if testForZero(args[0]):
4025            return numarray.zeros(pokeShape(args[0]),numarray.Float64)            return numarray.zeros(getShape(args[0]),numarray.Float64)
4026         elif testForZero(args[1]):         elif testForZero(args[1]):
4027            return numarray.ones(pokeShape(args[1]),numarray.Float64)            return numarray.ones(getShape(args[1]),numarray.Float64)
4028         elif isinstance(args[0],Symbol) or isinstance(args[1],Symbol):         elif isinstance(args[0],Symbol) or isinstance(args[1],Symbol):
4029            return Power_Symbol(args[0],args[1])            return Power_Symbol(args[0],args[1])
4030         elif isinstance(args[0],numarray.NumArray) and not isinstance(args[1],numarray.NumArray):         elif isinstance(args[0],numarray.NumArray) and not isinstance(args[1],numarray.NumArray):
# Line 3951  class Power_Symbol(DependendSymbol): Line 4047  class Power_Symbol(DependendSymbol):
4047         @raise ValueError: if both arguments do not have the same shape.         @raise ValueError: if both arguments do not have the same shape.
4048         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
4049         """         """
4050         sh0=pokeShape(arg0)         sh0=getShape(arg0)
4051         sh1=pokeShape(arg1)         sh1=getShape(arg1)
4052         if not sh0==sh1:         if not sh0==sh1:
4053            raise ValueError,"Power_Symbol: shape of arguments must match"            raise ValueError,"Power_Symbol: shape of arguments must match"
4054         d0=pokeDim(arg0)         d0=pokeDim(arg0)
# Line 4051  def minimum(*args): Line 4147  def minimum(*args):
4147            out=add(out,mult(whereNegative(diff),diff))            out=add(out,mult(whereNegative(diff),diff))
4148      return out      return out
4149    
4150  def clip(arg,minval=0.,maxval=1.):  def clip(arg,minval=None,maxval=None):
4151      """      """
4152      cuts the values of arg between minval and maxval      cuts the values of arg between minval and maxval
4153    
4154      @param arg: argument      @param arg: argument
4155      @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float}      @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float}
4156      @param minval: lower range      @param minval: lower range. If None no lower range is applied
4157      @type minval: C{float}      @type minval: C{float} or C{None}
4158      @param maxval: upper range      @param maxval: upper range. If None no upper range is applied
4159      @type maxval: C{float}      @type maxval: C{float} or C{None}
4160      @return: is on object with all its value between minval and maxval. value of the argument that greater then minval and      @return: is on object with all its value between minval and maxval. value of the argument that greater then minval and
4161               less then maxval are unchanged.               less then maxval are unchanged.
4162      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float} depending on the input
4163      @raise ValueError: if minval>maxval      @raise ValueError: if minval>maxval
4164      """      """
4165      if minval>maxval:      if not minval==None and not maxval==None:
4166         raise ValueError,"minval = %s must be less then maxval %s"%(minval,maxval)         if minval>maxval:
4167      return minimum(maximum(minval,arg),maxval)            raise ValueError,"minval = %s must be less then maxval %s"%(minval,maxval)
4168        if minval == None:
4169            tmp=arg
4170        else:
4171            tmp=maximum(minval,arg)
4172        if maxval == None:
4173            return tmp
4174        else:
4175            return minimum(tmp,maxval)
4176    
4177        
4178  def inner(arg0,arg1):  def inner(arg0,arg1):
# Line 4089  def inner(arg0,arg1): Line 4193  def inner(arg0,arg1):
4193      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float} depending on the input
4194      @raise ValueError: if the shapes of the arguments are not identical      @raise ValueError: if the shapes of the arguments are not identical
4195      """      """
4196      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4197      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4198      if not sh0==sh1:      if not sh0==sh1:
4199          raise ValueError,"inner: shape of arguments does not match"          raise ValueError,"inner: shape of arguments does not match"
4200      return generalTensorProduct(arg0,arg1,axis_offset=len(sh0))      return generalTensorProduct(arg0,arg1,axis_offset=len(sh0))
# Line 4141  def matrix_mult(arg0,arg1): Line 4245  def matrix_mult(arg0,arg1):
4245      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4246      @raise ValueError: if the shapes of the arguments are not appropriate      @raise ValueError: if the shapes of the arguments are not appropriate
4247      """      """
4248      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4249      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4250      if not len(sh0)==2 :      if not len(sh0)==2 :
4251          raise ValueError,"first argument must have rank 2"          raise ValueError,"first argument must have rank 2"
4252      if not len(sh1)==2 and not len(sh1)==1:      if not len(sh1)==2 and not len(sh1)==1:
# Line 4190  def tensor_mult(arg0,arg1): Line 4294  def tensor_mult(arg0,arg1):
4294      @return: the tensor product of arg0 and arg1 at each data point      @return: the tensor product of arg0 and arg1 at each data point
4295      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4296      """      """
4297      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4298      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4299      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4300         return generalTensorProduct(arg0,arg1,axis_offset=1)         return generalTensorProduct(arg0,arg1,axis_offset=1)
4301      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):
# Line 4265  class GeneralTensorProduct_Symbol(Depend Line 4369  class GeneralTensorProduct_Symbol(Depend
4369         @raise ValueError: illegal dimension         @raise ValueError: illegal dimension
4370         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
4371         """         """
4372         sh_arg0=pokeShape(arg0)         sh_arg0=getShape(arg0)
4373         sh_arg1=pokeShape(arg1)         sh_arg1=getShape(arg1)
4374         sh0=sh_arg0[:len(sh_arg0)-axis_offset]         sh0=sh_arg0[:len(sh_arg0)-axis_offset]
4375         sh01=sh_arg0[len(sh_arg0)-axis_offset:]         sh01=sh_arg0[len(sh_arg0)-axis_offset:]
4376         sh10=sh_arg1[:axis_offset]         sh10=sh_arg1[:axis_offset]
# Line 4313  class GeneralTensorProduct_Symbol(Depend Line 4417  class GeneralTensorProduct_Symbol(Depend
4417           args=self.getSubstitutedArguments(argvals)           args=self.getSubstitutedArguments(argvals)
4418           return generalTensorProduct(args[0],args[1],args[2])           return generalTensorProduct(args[0],args[1],args[2])
4419    
4420  def escript_generalTensorProductNew(arg0,arg1,axis_offset):  def escript_generalTensorProduct(arg0,arg1,axis_offset,transpose=0):
4421      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"
4422      # calculate the return shape:      return C_GeneralTensorProduct(arg0, arg1, axis_offset, transpose)
     shape0=arg0.getShape()[:arg0.getRank()-axis_offset]  
     shape01=arg0.getShape()[arg0.getRank()-axis_offset:]  
     shape10=arg1.getShape()[:axis_offset]  
     shape1=arg1.getShape()[axis_offset:]  
     if not shape01==shape10:  
         raise ValueError,"dimensions of last %s components in left argument don't match the first %s components in the right argument."%(axis_offset,axis_offset)  
     # Figure out which functionspace to use (look at where operator+ is defined maybe in BinaryOp.h to get the logic for this)  
     # fs=(escript.Scalar(0.,arg0.getFunctionSpace())+escript.Scalar(0.,arg1.getFunctionSpace())).getFunctionSpace()  
     out=GeneralTensorProduct(arg0, arg1, axis_offset)  
     return out  
   
 def escript_generalTensorProduct(arg0,arg1,axis_offset): # this should be escript._generalTensorProduct  
     "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"  
     # calculate the return shape:  
     shape0=arg0.getShape()[:arg0.getRank()-axis_offset]  
     shape01=arg0.getShape()[arg0.getRank()-axis_offset:]  
     shape10=arg1.getShape()[:axis_offset]  
     shape1=arg1.getShape()[axis_offset:]  
     if not shape01==shape10:  
         raise ValueError,"dimensions of last %s components in left argument don't match the first %s components in the right argument."%(axis_offset,axis_offset)  
   
     # whatr function space should be used? (this here is not good!)  
     fs=(escript.Scalar(0.,arg0.getFunctionSpace())+escript.Scalar(0.,arg1.getFunctionSpace())).getFunctionSpace()  
     # create return value:  
     out=escript.Data(0.,tuple(shape0+shape1),fs)  
     #  
     s0=[[]]  
     for k in shape0:  
           s=[]  
           for j in s0:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s0=s  
     s1=[[]]  
     for k in shape1:  
           s=[]  
           for j in s1:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s1=s  
     s01=[[]]  
     for k in shape01:  
           s=[]  
           for j in s01:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s01=s  
   
     for i0 in s0:  
        for i1 in s1:  
          s=escript.Scalar(0.,fs)  
          for i01 in s01:  
             s+=arg0.__getitem__(tuple(i0+i01))*arg1.__getitem__(tuple(i01+i1))  
          out.__setitem__(tuple(i0+i1),s)  
     return out  
   
4423    
4424  def transposed_matrix_mult(arg0,arg1):  def transposed_matrix_mult(arg0,arg1):
4425      """      """
# Line 4392  def transposed_matrix_mult(arg0,arg1): Line 4443  def transposed_matrix_mult(arg0,arg1):
4443      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4444      @raise ValueError: if the shapes of the arguments are not appropriate      @raise ValueError: if the shapes of the arguments are not appropriate
4445      """      """
4446      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4447      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4448      if not len(sh0)==2 :      if not len(sh0)==2 :
4449          raise ValueError,"first argument must have rank 2"          raise ValueError,"first argument must have rank 2"
4450      if not len(sh1)==2 and not len(sh1)==1:      if not len(sh1)==2 and not len(sh1)==1:
# Line 4437  def transposed_tensor_mult(arg0,arg1): Line 4488  def transposed_tensor_mult(arg0,arg1):
4488      @return: the tensor product of tarnsposed of arg0 and arg1 at each data point      @return: the tensor product of tarnsposed of arg0 and arg1 at each data point
4489      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4490      """      """
4491      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4492      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4493      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4494         return generalTransposedTensorProduct(arg0,arg1,axis_offset=1)         return generalTransposedTensorProduct(arg0,arg1,axis_offset=1)
4495      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):
# Line 4515  class GeneralTransposedTensorProduct_Sym Line 4566  class GeneralTransposedTensorProduct_Sym
4566         @raise ValueError: inconsistent dimensions of arguments.         @raise ValueError: inconsistent dimensions of arguments.
4567         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
4568         """         """
4569         sh_arg0=pokeShape(arg0)         sh_arg0=getShape(arg0)
4570         sh_arg1=pokeShape(arg1)         sh_arg1=getShape(arg1)
4571         sh01=sh_arg0[:axis_offset]         sh01=sh_arg0[:axis_offset]
4572         sh10=sh_arg1[:axis_offset]         sh10=sh_arg1[:axis_offset]
4573         sh0=sh_arg0[axis_offset:]         sh0=sh_arg0[axis_offset:]
# Line 4565  class GeneralTransposedTensorProduct_Sym Line 4616  class GeneralTransposedTensorProduct_Sym
4616    
4617  def escript_generalTransposedTensorProduct(arg0,arg1,axis_offset): # this should be escript._generalTransposedTensorProduct  def escript_generalTransposedTensorProduct(arg0,arg1,axis_offset): # this should be escript._generalTransposedTensorProduct
4618      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"
4619      # calculate the return shape:      return C_GeneralTensorProduct(arg0, arg1, axis_offset, 1)
     shape01=arg0.getShape()[:axis_offset]  
     shape10=arg1.getShape()[:axis_offset]  
     shape0=arg0.getShape()[axis_offset:]  
     shape1=arg1.getShape()[axis_offset:]  
     if not shape01==shape10:  
         raise ValueError,"dimensions of first %s components in left argument don't match the first %s components in the right argument."%(axis_offset,axis_offset)  
   
     # whatr function space should be used? (this here is not good!)  
     fs=(escript.Scalar(0.,arg0.getFunctionSpace())+escript.Scalar(0.,arg1.getFunctionSpace())).getFunctionSpace()  
     # create return value:  
     out=escript.Data(0.,tuple(shape0+shape1),fs)  
     #  
     s0=[[]]  
     for k in shape0:  
           s=[]  
           for j in s0:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s0=s  
     s1=[[]]  
     for k in shape1:  
           s=[]  
           for j in s1:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s1=s  
     s01=[[]]  
     for k in shape01:  
           s=[]  
           for j in s01:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s01=s  
   
     for i0 in s0:  
        for i1 in s1:  
          s=escript.Scalar(0.,fs)  
          for i01 in s01:  
             s+=arg0.__getitem__(tuple(i01+i0))*arg1.__getitem__(tuple(i01+i1))  
          out.__setitem__(tuple(i0+i1),s)  
     return out  
   
4620    
4621  def matrix_transposed_mult(arg0,arg1):  def matrix_transposed_mult(arg0,arg1):
4622      """      """
# Line 4624  def matrix_transposed_mult(arg0,arg1): Line 4636  def matrix_transposed_mult(arg0,arg1):
4636      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4637      @raise ValueError: if the shapes of the arguments are not appropriate      @raise ValueError: if the shapes of the arguments are not appropriate
4638      """      """
4639      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4640      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4641      if not len(sh0)==2 :      if not len(sh0)==2 :
4642          raise ValueError,"first argument must have rank 2"          raise ValueError,"first argument must have rank 2"
4643      if not len(sh1)==2 and not len(sh1)==1:      if not len(sh1)==2 and not len(sh1)==1:
# Line 4660  def tensor_transposed_mult(arg0,arg1): Line 4672  def tensor_transposed_mult(arg0,arg1):
4672      @return: the tensor product of tarnsposed of arg0 and arg1 at each data point      @return: the tensor product of tarnsposed of arg0 and arg1 at each data point
4673      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4674      """      """
4675      sh0=pokeShape(arg0)      sh0=getShape(arg0)
4676      sh1=pokeShape(arg1)      sh1=getShape(arg1)
4677      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4678         return generalTensorTransposedProduct(arg0,arg1,axis_offset=1)         return generalTensorTransposedProduct(arg0,arg1,axis_offset=1)
4679      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):      elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):
# Line 4738  class GeneralTensorTransposedProduct_Sym Line 4750  class GeneralTensorTransposedProduct_Sym
4750         @raise ValueError: inconsistent dimensions of arguments.         @raise ValueError: inconsistent dimensions of arguments.
4751         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
4752         """         """
4753         sh_arg0=pokeShape(arg0)         sh_arg0=getShape(arg0)
4754         sh_arg1=pokeShape(arg1)         sh_arg1=getShape(arg1)
4755         sh0=sh_arg0[:len(sh_arg0)-axis_offset]         sh0=sh_arg0[:len(sh_arg0)-axis_offset]
4756         sh01=sh_arg0[len(sh_arg0)-axis_offset:]         sh01=sh_arg0[len(sh_arg0)-axis_offset:]
4757         sh10=sh_arg1[len(sh_arg1)-axis_offset:]         sh10=sh_arg1[len(sh_arg1)-axis_offset:]
# Line 4788  class GeneralTensorTransposedProduct_Sym Line 4800  class GeneralTensorTransposedProduct_Sym
4800    
4801  def escript_generalTensorTransposedProduct(arg0,arg1,axis_offset): # this should be escript._generalTensorTransposedProduct  def escript_generalTensorTransposedProduct(arg0,arg1,axis_offset): # this should be escript._generalTensorTransposedProduct
4802      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"      "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"
4803      # calculate the return shape:      return C_GeneralTensorProduct(arg0, arg1, axis_offset, 2)
     shape01=arg0.getShape()[arg0.getRank()-axis_offset:]  
     shape0=arg0.getShape()[:arg0.getRank()-axis_offset]  
     shape10=arg1.getShape()[arg1.getRank()-axis_offset:]  
     shape1=arg1.getShape()[:arg1.getRank()-axis_offset]  
     if not shape01==shape10:  
         raise ValueError,"dimensions of first %s components in left argument don't match the first %s components in the right argument."%(axis_offset,axis_offset)  
   
     # whatr function space should be used? (this here is not good!)  
     fs=(escript.Scalar(0.,arg0.getFunctionSpace())+escript.Scalar(0.,arg1.getFunctionSpace())).getFunctionSpace()  
     # create return value:  
     out=escript.Data(0.,tuple(shape0+shape1),fs)  
     #  
     s0=[[]]  
     for k in shape0:  
           s=[]  
           for j in s0:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s0=s  
     s1=[[]]  
     for k in shape1:  
           s=[]  
           for j in s1:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s1=s  
     s01=[[]]  
     for k in shape01:  
           s=[]  
           for j in s01:  
                 for i in range(k): s.append(j+[slice(i,i)])  
           s01=s  
   
     for i0 in s0:  
        for i1 in s1:  
          s=escript.Scalar(0.,fs)  
          for i01 in s01:  
             s+=arg0.__getitem__(tuple(i0+i01))*arg1.__getitem__(tuple(i1+i01))  
          out.__setitem__(tuple(i0+i1),s)  
     return out  
   
4804    
4805  #=========================================================  #=========================================================
4806  #  functions dealing with spatial dependency  #  functions dealing with spatial dependency
# Line 4953  def integrate(arg,where=None): Line 4926  def integrate(arg,where=None):
4926         else:         else:
4927            return arg._integrate()            return arg._integrate()
4928      else:      else:
4929        raise TypeError,"integrate: Unknown argument type."         arg2=escript.Data(arg,where)
4930           if arg2.getRank()==0:
4931              return arg2._integrate()[0]
4932           else:
4933              return arg2._integrate()
4934    
4935  class Integrate_Symbol(DependendSymbol):  class Integrate_Symbol(DependendSymbol):
4936     """     """

Legend:
Removed from v.804  
changed lines
  Added in v.1295

  ViewVC Help
Powered by ViewVC 1.1.26