/[escript]/trunk/escript/py_src/util.py
ViewVC logotype

Diff of /trunk/escript/py_src/util.py

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

revision 341 by gross, Mon Dec 12 05:26:10 2005 UTC revision 813 by ksteube, Mon Aug 21 02:08:47 2006 UTC
# Line 1  Line 1 
1  # $Id$  # $Id$
 #  
 #      COPYRIGHT ACcESS 2004 -  All Rights Reserved  
 #  
 #   This software is the property of ACcESS.  No part of this code  
 #   may be copied in any form or by any means without the expressed written  
 #   consent of ACcESS.  Copying, use or modification of this software  
 #   by any unauthorised person is illegal unless that  
 #   person has a software license agreement with ACcESS.  
 #  
2    
3  """  """
4  Utility functions for escript  Utility functions for escript
5    
 @remark:  This module is under construction and is still tested!!!  
   
6  @var __author__: name of author  @var __author__: name of author
7  @var __licence__: licence agreement  @var __copyright__: copyrights
8    @var __license__: licence agreement
9  @var __url__: url entry point on documentation  @var __url__: url entry point on documentation
10  @var __version__: version  @var __version__: version
11  @var __date__: date of the version  @var __date__: date of the version
12  """  """
13                                                                                                                                                                                                                                                                                                                                                                                                            
14  __author__="Lutz Gross, l.gross@uq.edu.au"  __author__="Lutz Gross, l.gross@uq.edu.au"
15  __licence__="contact: esys@access.uq.edu.au"  __copyright__="""  Copyright (c) 2006 by ACcESS MNRF
16                        http://www.access.edu.au
17                    Primary Business: Queensland, Australia"""
18    __license__="""Licensed under the Open Software License version 3.0
19                 http://www.opensource.org/licenses/osl-3.0.php"""
20  __url__="http://www.iservo.edu.au/esys/escript"  __url__="http://www.iservo.edu.au/esys/escript"
21  __version__="$Revision: 329 $"  __version__="$Revision$"
22  __date__="$Date$"  __date__="$Date$"
23    
24    
# Line 32  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
 # missing tests:  
   
 # def pokeShape(arg):  
 # def pokeDim(arg):  
 # def commonShape(arg0,arg1):  
 # def commonDim(*args):  
 # def testForZero(arg):  
 # def matchType(arg0=0.,arg1=0.):  
 # def matchShape(arg0,arg1):  
   
 # def maximum(arg0,arg1):  
 # def minimum(arg0,arg1):  
   
 # def transpose(arg,axis=None):  
 # def trace(arg,axis0=0,axis1=1):  
 # def reorderComponents(arg,index):  
   
 # def integrate(arg,where=None):  
 # def interpolate(arg,where):  
 # def div(arg,where=None):  
 # def grad(arg,where=None):  
   
 #  
 # slicing: get  
 #          set  
 #  
 # and derivatives  
30    
31  #=========================================================  #=========================================================
32  #   some helpers:  #   some helpers:
# Line 68  def saveVTK(filename,domain=None,**data) Line 35  def saveVTK(filename,domain=None,**data)
35      """      """
36      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.
37    
38      Example:      Example::
39    
40         tmp=Scalar(..)         tmp=Scalar(..)
41         v=Vector(..)         v=Vector(..)
# Line 96  def saveDX(filename,domain=None,**data): Line 63  def saveDX(filename,domain=None,**data):
63      """      """
64      writes a L{Data} objects into a files using the the DX file format.      writes a L{Data} objects into a files using the the DX file format.
65    
66      Example:      Example::
67    
68         tmp=Scalar(..)         tmp=Scalar(..)
69         v=Vector(..)         v=Vector(..)
# Line 125  def kronecker(d=3): Line 92  def kronecker(d=3):
92     return the kronecker S{delta}-symbol     return the kronecker S{delta}-symbol
93    
94     @param d: dimension or an object that has the C{getDim} method defining the dimension     @param d: dimension or an object that has the C{getDim} method defining the dimension
95     @type d: C{int} or any object with a C{getDim} method     @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
96     @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise     @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise
97     @rtype d: L{numarray.NumArray} of rank 2.     @rtype: L{numarray.NumArray} or L{escript.Data} of rank 2.
    @remark: the function is identical L{identity}  
98     """     """
99     return identityTensor(d)     return identityTensor(d)
100    
# Line 143  def identity(shape=()): Line 109  def identity(shape=()):
109     @raise ValueError: if len(shape)>2.     @raise ValueError: if len(shape)>2.
110     """     """
111     if len(shape)>0:     if len(shape)>0:
112        out=numarray.zeros(shape+shape,numarray.Float)        out=numarray.zeros(shape+shape,numarray.Float64)
113        if len(shape)==1:        if len(shape)==1:
114            for i0 in range(shape[0]):            for i0 in range(shape[0]):
115               out[i0,i0]=1.               out[i0,i0]=1.
   
116        elif len(shape)==2:        elif len(shape)==2:
117            for i0 in range(shape[0]):            for i0 in range(shape[0]):
118               for i1 in range(shape[1]):               for i1 in range(shape[1]):
# Line 163  def identityTensor(d=3): Line 128  def identityTensor(d=3):
128     return the dxd identity matrix     return the dxd identity matrix
129    
130     @param d: dimension or an object that has the C{getDim} method defining the dimension     @param d: dimension or an object that has the C{getDim} method defining the dimension
131     @type d: C{int} or any object with a C{getDim} method     @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
132     @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise     @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise
133     @rtype: L{numarray.NumArray} of rank 2.     @rtype: L{numarray.NumArray} or L{escript.Data} of rank 2
134     """     """
135     if hasattr(d,"getDim"):     if isinstance(d,escript.FunctionSpace):
136        d=d.getDim()         return escript.Data(identity((d.getDim(),)),d)
137     return identity(shape=(d,))     elif isinstance(d,escript.Domain):
138           return identity((d.getDim(),))
139       else:
140           return identity((d,))
141    
142  def identityTensor4(d=3):  def identityTensor4(d=3):
143     """     """
# Line 178  def identityTensor4(d=3): Line 146  def identityTensor4(d=3):
146     @param d: dimension or an object that has the C{getDim} method defining the dimension     @param d: dimension or an object that has the C{getDim} method defining the dimension
147     @type d: C{int} or any object with a C{getDim} method     @type d: C{int} or any object with a C{getDim} method
148     @return: the object u of rank 4 with M{u[i,j,k,l]=1} for M{i=k and j=l} and M{u[i,j,k,l]=0} otherwise     @return: the object u of rank 4 with M{u[i,j,k,l]=1} for M{i=k and j=l} and M{u[i,j,k,l]=0} otherwise
149     @rtype: L{numarray.NumArray} of rank 4.     @rtype: L{numarray.NumArray} or L{escript.Data} of rank 4.
150     """     """
151     if hasattr(d,"getDim"):     if isinstance(d,escript.FunctionSpace):
152        d=d.getDim()         return escript.Data(identity((d.getDim(),d.getDim())),d)
153     return identity((d,d))     elif isinstance(d,escript.Domain):
154           return identity((d.getDim(),d.getDim()))
155       else:
156           return identity((d,d))
157    
158  def unitVector(i=0,d=3):  def unitVector(i=0,d=3):
159     """     """
# Line 191  def unitVector(i=0,d=3): Line 162  def unitVector(i=0,d=3):
162     @param i: index     @param i: index
163     @type i: C{int}     @type i: C{int}
164     @param d: dimension or an object that has the C{getDim} method defining the dimension     @param d: dimension or an object that has the C{getDim} method defining the dimension
165     @type d: C{int} or any object with a C{getDim} method     @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
166     @return: the object u of rank 1 with M{u[j]=1} for M{j=i} and M{u[i]=0} otherwise     @return: the object u of rank 1 with M{u[j]=1} for M{j=i} and M{u[i]=0} otherwise
167     @rtype: L{numarray.NumArray} of rank 1.     @rtype: L{numarray.NumArray} or L{escript.Data} of rank 1
168     """     """
169     return kronecker(d)[i]     return kronecker(d)[i]
170    
# Line 249  def inf(arg): Line 220  def inf(arg):
220    
221      @param arg: argument      @param arg: argument
222      @type arg: C{float}, C{int}, L{escript.Data}, L{numarray.NumArray}.      @type arg: C{float}, C{int}, L{escript.Data}, L{numarray.NumArray}.
223      @return : minimum value of arg over all components and all data points      @return: minimum value of arg over all components and all data points
224      @rtype: C{float}      @rtype: C{float}
225      @raise TypeError: if type of arg cannot be processed      @raise TypeError: if type of arg cannot be processed
226      """      """
# Line 338  def commonDim(*args): Line 309  def commonDim(*args):
309      """      """
310      identifies, if possible, the spatial dimension across a set of objects which may or my not have a spatial dimension.      identifies, if possible, the spatial dimension across a set of objects which may or my not have a spatial dimension.
311    
312      @param *args: given objects      @param args: given objects
313      @return: the spatial dimension of the objects with identifiable dimension (see L{pokeDim}). If none the objects has      @return: the spatial dimension of the objects with identifiable dimension (see L{pokeDim}). If none the objects has
314               a spatial dimension C{None} is returned.               a spatial dimension C{None} is returned.
315      @rtype: C{int} or C{None}      @rtype: C{int} or C{None}
# Line 360  def testForZero(arg): Line 331  def testForZero(arg):
331    
332      @param arg: a given object      @param arg: a given object
333      @type arg: typically L{numarray.NumArray},L{escript.Data},C{float}, C{int}      @type arg: typically L{numarray.NumArray},L{escript.Data},C{float}, C{int}
334      @return : True if the argument is identical to zero.      @return: True if the argument is identical to zero.
335      @rtype : C{bool}      @rtype: C{bool}
336      """      """
337      try:      if isinstance(arg,numarray.NumArray):
338           return not Lsup(arg)>0.
339        elif isinstance(arg,escript.Data):
340           return False
341        elif isinstance(arg,float):
342         return not Lsup(arg)>0.         return not Lsup(arg)>0.
343      except TypeError:      elif isinstance(arg,int):
344           return not Lsup(arg)>0.
345        elif isinstance(arg,Symbol):
346           return False
347        else:
348         return False         return False
349    
350  def matchType(arg0=0.,arg1=0.):  def matchType(arg0=0.,arg1=0.):
# Line 386  def matchType(arg0=0.,arg1=0.): Line 365  def matchType(arg0=0.,arg1=0.):
365         elif isinstance(arg1,escript.Data):         elif isinstance(arg1,escript.Data):
366            arg0=escript.Data(arg0,arg1.getFunctionSpace())            arg0=escript.Data(arg0,arg1.getFunctionSpace())
367         elif isinstance(arg1,float):         elif isinstance(arg1,float):
368            arg1=numarray.array(arg1)            arg1=numarray.array(arg1,type=numarray.Float64)
369         elif isinstance(arg1,int):         elif isinstance(arg1,int):
370            arg1=numarray.array(float(arg1))            arg1=numarray.array(float(arg1),type=numarray.Float64)
371         elif isinstance(arg1,Symbol):         elif isinstance(arg1,Symbol):
372            pass            pass
373         else:         else:
# Line 412  def matchType(arg0=0.,arg1=0.): Line 391  def matchType(arg0=0.,arg1=0.):
391         elif isinstance(arg1,escript.Data):         elif isinstance(arg1,escript.Data):
392            pass            pass
393         elif isinstance(arg1,float):         elif isinstance(arg1,float):
394            arg1=numarray.array(arg1)            arg1=numarray.array(arg1,type=numarray.Float64)
395         elif isinstance(arg1,int):         elif isinstance(arg1,int):
396            arg1=numarray.array(float(arg1))            arg1=numarray.array(float(arg1),type=numarray.Float64)
397         elif isinstance(arg1,Symbol):         elif isinstance(arg1,Symbol):
398            pass            pass
399         else:         else:
400            raise TypeError,"function: Unknown type of second argument."                raise TypeError,"function: Unknown type of second argument."    
401      elif isinstance(arg0,float):      elif isinstance(arg0,float):
402         if isinstance(arg1,numarray.NumArray):         if isinstance(arg1,numarray.NumArray):
403            arg0=numarray.array(arg0)            arg0=numarray.array(arg0,type=numarray.Float64)
404         elif isinstance(arg1,escript.Data):         elif isinstance(arg1,escript.Data):
405            arg0=escript.Data(arg0,arg1.getFunctionSpace())            arg0=escript.Data(arg0,arg1.getFunctionSpace())
406         elif isinstance(arg1,float):         elif isinstance(arg1,float):
407            arg0=numarray.array(arg0)            arg0=numarray.array(arg0,type=numarray.Float64)
408            arg1=numarray.array(arg1)            arg1=numarray.array(arg1,type=numarray.Float64)
409         elif isinstance(arg1,int):         elif isinstance(arg1,int):
410            arg0=numarray.array(arg0)            arg0=numarray.array(arg0,type=numarray.Float64)
411            arg1=numarray.array(float(arg1))            arg1=numarray.array(float(arg1),type=numarray.Float64)
412         elif isinstance(arg1,Symbol):         elif isinstance(arg1,Symbol):
413            arg0=numarray.array(arg0)            arg0=numarray.array(arg0,type=numarray.Float64)
414         else:         else:
415            raise TypeError,"function: Unknown type of second argument."                raise TypeError,"function: Unknown type of second argument."    
416      elif isinstance(arg0,int):      elif isinstance(arg0,int):
417         if isinstance(arg1,numarray.NumArray):         if isinstance(arg1,numarray.NumArray):
418            arg0=numarray.array(float(arg0))            arg0=numarray.array(float(arg0),type=numarray.Float64)
419         elif isinstance(arg1,escript.Data):         elif isinstance(arg1,escript.Data):
420            arg0=escript.Data(float(arg0),arg1.getFunctionSpace())            arg0=escript.Data(float(arg0),arg1.getFunctionSpace())
421         elif isinstance(arg1,float):         elif isinstance(arg1,float):
422            arg0=numarray.array(float(arg0))            arg0=numarray.array(float(arg0),type=numarray.Float64)
423            arg1=numarray.array(arg1)            arg1=numarray.array(arg1,type=numarray.Float64)
424         elif isinstance(arg1,int):         elif isinstance(arg1,int):
425            arg0=numarray.array(float(arg0))            arg0=numarray.array(float(arg0),type=numarray.Float64)
426            arg1=numarray.array(float(arg1))            arg1=numarray.array(float(arg1),type=numarray.Float64)
427         elif isinstance(arg1,Symbol):         elif isinstance(arg1,Symbol):
428            arg0=numarray.array(float(arg0))            arg0=numarray.array(float(arg0),type=numarray.Float64)
429         else:         else:
430            raise TypeError,"function: Unknown type of second argument."                raise TypeError,"function: Unknown type of second argument."    
431      else:      else:
# Line 456  def matchType(arg0=0.,arg1=0.): Line 435  def matchType(arg0=0.,arg1=0.):
435    
436  def matchShape(arg0,arg1):  def matchShape(arg0,arg1):
437      """      """
438            return representations of arg0 amd arg1 which ahve the same shape
   
     If shape is not given the shape "largest" shape of args is used.  
439    
440      @param args: a given ob      @param arg0: a given object
441      @type arg: typically L{numarray.NumArray},L{escript.Data},C{float}, C{int}      @type arg0: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, L{Symbol}
442      @return: True if the argument is identical to zero.      @param arg1: a given object
443      @rtype: C{list} of C{int}      @type arg1: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, L{Symbol}
444        @return: C{arg0} and C{arg1} where copies are returned when the shape has to be changed.
445        @rtype: C{tuple}
446      """      """
447      sh=commonShape(arg0,arg1)      sh=commonShape(arg0,arg1)
448      sh0=pokeShape(arg0)      sh0=pokeShape(arg0)
449      sh1=pokeShape(arg1)      sh1=pokeShape(arg1)
450      if len(sh0)<len(sh):      if len(sh0)<len(sh):
451         return outer(arg0,numarray.ones(sh[len(sh0):],numarray.Float)),arg1         return outer(arg0,numarray.ones(sh[len(sh0):],numarray.Float64)),arg1
452      elif len(sh1)<len(sh):      elif len(sh1)<len(sh):
453         return arg0,outer(arg1,numarray.ones(sh[len(sh1):],numarray.Float))         return arg0,outer(arg1,numarray.ones(sh[len(sh1):],numarray.Float64))
454      else:      else:
455         return arg0,arg1         return arg0,arg1
456  #=========================================================  #=========================================================
# Line 491  class Symbol(object): Line 470  class Symbol(object):
470         Creates an instance of a symbol of a given shape. The symbol may depending on a list of arguments args which may be         Creates an instance of a symbol of a given shape. The symbol may depending on a list of arguments args which may be
471         symbols or any other object.         symbols or any other object.
472    
473         @param arg: the arguments of the symbol.         @param args: the arguments of the symbol.
474         @type arg: C{list}         @type args: C{list}
475         @param shape: the shape         @param shape: the shape
476         @type shape: C{tuple} of C{int}         @type shape: C{tuple} of C{int}
477         @param dim: spatial dimension of the symbol. If dim=C{None} the spatial dimension is undefined.           @param dim: spatial dimension of the symbol. If dim=C{None} the spatial dimension is undefined.  
# Line 535  class Symbol(object): Line 514  class Symbol(object):
514         """         """
515         the shape of the symbol.         the shape of the symbol.
516    
517         @return : the shape of the symbol.         @return: the shape of the symbol.
518         @rtype: C{tuple} of C{int}         @rtype: C{tuple} of C{int}
519         """         """
520         return self.__shape         return self.__shape
# Line 544  class Symbol(object): Line 523  class Symbol(object):
523         """         """
524         the spatial dimension         the spatial dimension
525    
526         @return : the spatial dimension         @return: the spatial dimension
527         @rtype: C{int} if the dimension is defined. Otherwise C{None} is returned.         @rtype: C{int} if the dimension is defined. Otherwise C{None} is returned.
528         """         """
529         return self.__dim         return self.__dim
# Line 568  class Symbol(object): Line 547  class Symbol(object):
547         """         """
548         substitutes symbols in the arguments of this object and returns the result as a list.         substitutes symbols in the arguments of this object and returns the result as a list.
549    
550         @param argvals: L{Symbols} and their substitutes. The L{Symbol} u in the expression defining this object is replaced by argvals[u].         @param argvals: L{Symbol} and their substitutes. The L{Symbol} u in the expression defining this object is replaced by argvals[u].
551         @type argvals: C{dict} with keywords of type L{Symbol}.         @type argvals: C{dict} with keywords of type L{Symbol}.
552         @rtype: C{list} of objects         @rtype: C{list} of objects
553         @return: list of the object assigned to the arguments through substitution or for the arguments which are not L{Symbols} the value assigned to the argument at instantiation.         @return: list of the object assigned to the arguments through substitution or for the arguments which are not L{Symbol} the value assigned to the argument at instantiation.
554         """         """
555         out=[]         out=[]
556         for a in self.getArgument():         for a in self.getArgument():
# Line 597  class Symbol(object): Line 576  class Symbol(object):
576            else:            else:
577                s=pokeShape(s)+arg.getShape()                s=pokeShape(s)+arg.getShape()
578                if len(s)>0:                if len(s)>0:
579                   out.append(numarray.zeros(s),numarray.Float)                   out.append(numarray.zeros(s),numarray.Float64)
580                else:                else:
581                   out.append(a)                   out.append(a)
582         return out         return out
# Line 687  class Symbol(object): Line 666  class Symbol(object):
666         else:         else:
667            s=self.getShape()+arg.getShape()            s=self.getShape()+arg.getShape()
668            if len(s)>0:            if len(s)>0:
669               return numarray.zeros(s,numarray.Float)               return numarray.zeros(s,numarray.Float64)
670            else:            else:
671               return 0.               return 0.
672    
# Line 695  class Symbol(object): Line 674  class Symbol(object):
674         """         """
675         returns -self.         returns -self.
676    
677         @return:  a S{Symbol} representing the negative of the object         @return:  a L{Symbol} representing the negative of the object
678         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
679         """         """
680         return self*(-1.)         return self*(-1.)
# Line 704  class Symbol(object): Line 683  class Symbol(object):
683         """         """
684         returns +self.         returns +self.
685    
686         @return:  a S{Symbol} representing the positive of the object         @return:  a L{Symbol} representing the positive of the object
687         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
688         """         """
689         return self*(1.)         return self*(1.)
690    
691     def __abs__(self):     def __abs__(self):
692         """         """
693         returns a S{Symbol} representing the absolute value of the object.         returns a L{Symbol} representing the absolute value of the object.
694         """         """
695         return Abs_Symbol(self)         return Abs_Symbol(self)
696    
# Line 721  class Symbol(object): Line 700  class Symbol(object):
700    
701         @param other: object to be added to this object         @param other: object to be added to this object
702         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
703         @return:  a S{Symbol} representing the sum of this object and C{other}         @return:  a L{Symbol} representing the sum of this object and C{other}
704         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
705         """         """
706         return add(self,other)         return add(self,other)
# Line 732  class Symbol(object): Line 711  class Symbol(object):
711    
712         @param other: object this object is added to         @param other: object this object is added to
713         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
714         @return: a S{Symbol} representing the sum of C{other} and this object object         @return: a L{Symbol} representing the sum of C{other} and this object object
715         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
716         """         """
717         return add(other,self)         return add(other,self)
# Line 743  class Symbol(object): Line 722  class Symbol(object):
722    
723         @param other: object to be subtracted from this object         @param other: object to be subtracted from this object
724         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
725         @return: a S{Symbol} representing the difference of C{other} and this object         @return: a L{Symbol} representing the difference of C{other} and this object
726         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
727         """         """
728         return add(self,-other)         return add(self,-other)
# Line 754  class Symbol(object): Line 733  class Symbol(object):
733    
734         @param other: object this object is been subtracted from         @param other: object this object is been subtracted from
735         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
736         @return: a S{Symbol} representing the difference of this object and C{other}.         @return: a L{Symbol} representing the difference of this object and C{other}.
737         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
738         """         """
739         return add(-self,other)         return add(-self,other)
# Line 765  class Symbol(object): Line 744  class Symbol(object):
744    
745         @param other: object to be mutiplied by this object         @param other: object to be mutiplied by this object
746         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
747         @return: a S{Symbol} representing the product of the object and C{other}.         @return: a L{Symbol} representing the product of the object and C{other}.
748         @rtype: L{DependendSymbol} or 0 if other is identical to zero.         @rtype: L{DependendSymbol} or 0 if other is identical to zero.
749         """         """
750         return mult(self,other)         return mult(self,other)
# Line 776  class Symbol(object): Line 755  class Symbol(object):
755    
756         @param other: object this object is multiplied with         @param other: object this object is multiplied with
757         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
758         @return: a S{Symbol} representing the product of C{other} and the object.         @return: a L{Symbol} representing the product of C{other} and the object.
759         @rtype: L{DependendSymbol} or 0 if other is identical to zero.         @rtype: L{DependendSymbol} or 0 if other is identical to zero.
760         """         """
761         return mult(other,self)         return mult(other,self)
# Line 787  class Symbol(object): Line 766  class Symbol(object):
766    
767         @param other: object dividing this object         @param other: object dividing this object
768         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
769         @return: a S{Symbol} representing the quotient of this object and C{other}         @return: a L{Symbol} representing the quotient of this object and C{other}
770         @rtype: L{DependendSymbol}         @rtype: L{DependendSymbol}
771         """         """
772         return quotient(self,other)         return quotient(self,other)
# Line 798  class Symbol(object): Line 777  class Symbol(object):
777    
778         @param other: object dividing this object         @param other: object dividing this object
779         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
780         @return: a S{Symbol} representing the quotient of C{other} and this object         @return: a L{Symbol} representing the quotient of C{other} and this object
781         @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.         @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.
782         """         """
783         return quotient(other,self)         return quotient(other,self)
# Line 809  class Symbol(object): Line 788  class Symbol(object):
788    
789         @param other: exponent         @param other: exponent
790         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
791         @return: a S{Symbol} representing the power of this object to C{other}         @return: a L{Symbol} representing the power of this object to C{other}
792         @rtype: L{DependendSymbol} or 1 if C{other} is identical to zero.         @rtype: L{DependendSymbol} or 1 if C{other} is identical to zero.
793         """         """
794         return power(self,other)         return power(self,other)
# Line 820  class Symbol(object): Line 799  class Symbol(object):
799    
800         @param other: basis         @param other: basis
801         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
802         @return: a S{Symbol} representing the power of C{other} to this object         @return: a L{Symbol} representing the power of C{other} to this object
803         @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.         @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.
804         """         """
805         return power(other,self)         return power(other,self)
806    
807       def __getitem__(self,index):
808           """
809           returns the slice defined by index
810    
811           @param index: defines a
812           @type index: C{slice} or C{int} or a C{tuple} of them
813           @return: a L{Symbol} representing the slice defined by index
814           @rtype: L{DependendSymbol}
815           """
816           return GetSlice_Symbol(self,index)
817    
818  class DependendSymbol(Symbol):  class DependendSymbol(Symbol):
819     """     """
820     DependendSymbol extents L{Symbol} by modifying the == operator to allow two instances to be equal.     DependendSymbol extents L{Symbol} by modifying the == operator to allow two instances to be equal.
821     Two DependendSymbol are equal if they have the same shape, the same arguments and one of them has an unspecified spatial dimension or the spatial dimension is identical       Two DependendSymbol are equal if they have the same shape, the same arguments and one of them has an unspecified spatial dimension or the spatial dimension is identical  
822        
823     Example:     Example::
824        
825     u1=Symbol(shape=(3,4),dim=2,args=[4.])       u1=Symbol(shape=(3,4),dim=2,args=[4.])
826     u2=Symbol(shape=(3,4),dim=2,args=[4.])       u2=Symbol(shape=(3,4),dim=2,args=[4.])
827     print u1==u2       print u1==u2
828     False       False
829        
830        but       but::
831    
832     u1=DependendSymbol(shape=(3,4),dim=2,args=[4.])       u1=DependendSymbol(shape=(3,4),dim=2,args=[4.])
833     u2=DependendSymbol(shape=(3,4),dim=2,args=[4.])       u2=DependendSymbol(shape=(3,4),dim=2,args=[4.])
834     u3=DependendSymbol(shape=(2,),dim=2,args=[4.])         u3=DependendSymbol(shape=(2,),dim=2,args=[4.])  
835     print u1==u2, u1==u3       print u1==u2, u1==u3
836     True False       True False
837    
838     @note: DependendSymbol should be used as return value of functions with L{Symbol} arguments. This will allow the optimizer to remove redundant function calls.     @note: DependendSymbol should be used as return value of functions with L{Symbol} arguments. This will allow the optimizer to remove redundant function calls.
839     """     """
# Line 875  class DependendSymbol(Symbol): Line 865  class DependendSymbol(Symbol):
865  #=========================================================  #=========================================================
866  #  Unary operations prserving the shape  #  Unary operations prserving the shape
867  #========================================================  #========================================================
868    class GetSlice_Symbol(DependendSymbol):
869       """
870       L{Symbol} representing getting a slice for a L{Symbol}
871       """
872       def __init__(self,arg,index):
873          """
874          initialization of wherePositive L{Symbol} with argument arg
875          @param arg: argument
876          @type arg: L{Symbol}.
877          @param index: defines index
878          @type index: C{slice} or C{int} or a C{tuple} of them
879          @raises IndexError: if length of index is larger than rank of arg or a index start or stop is out of range
880          @raises ValueError: if a step is given
881          """
882          if not isinstance(index,tuple): index=(index,)
883          if len(index)>arg.getRank():
884               raise IndexError,"GetSlice_Symbol: index out of range."
885          sh=()
886          index2=()
887          for i in range(len(index)):
888             ix=index[i]
889             if isinstance(ix,int):
890                if ix<0 or ix>=arg.getShape()[i]:
891                   raise ValueError,"GetSlice_Symbol: index out of range."
892                index2=index2+(ix,)
893             else:
894               if not ix.step==None:
895                 raise ValueError,"GetSlice_Symbol: steping is not supported."
896               if ix.start==None:
897                  s=0
898               else:
899                  s=ix.start
900               if ix.stop==None:
901                  e=arg.getShape()[i]
902               else:
903                  e=ix.stop
904                  if e>arg.getShape()[i]:
905                     raise IndexError,"GetSlice_Symbol: index out of range."
906               index2=index2+(slice(s,e),)
907               if e>s:
908                   sh=sh+(e-s,)
909               elif s>e:
910                   raise IndexError,"GetSlice_Symbol: slice start must be less or equal slice end"
911          for i in range(len(index),arg.getRank()):
912              index2=index2+(slice(0,arg.getShape()[i]),)
913              sh=sh+(arg.getShape()[i],)
914          super(GetSlice_Symbol, self).__init__(args=[arg,index2],shape=sh,dim=arg.getDim())
915    
916       def getMyCode(self,argstrs,format="escript"):
917          """
918          returns a program code that can be used to evaluate the symbol.
919    
920          @param argstrs: gives for each argument a string representing the argument for the evaluation.
921          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
922          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
923          @type format: C{str}
924          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
925          @rtype: C{str}
926          @raise NotImplementedError: if the requested format is not available
927          """
928          if format=="escript" or format=="str"  or format=="text":
929             return "%s.__getitem__(%s)"%(argstrs[0],argstrs[1])
930          else:
931             raise NotImplementedError,"GetItem_Symbol does not provide program code for format %s."%format
932    
933       def substitute(self,argvals):
934          """
935          assigns new values to symbols in the definition of the symbol.
936          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
937    
938          @param argvals: new values assigned to symbols
939          @type argvals: C{dict} with keywords of type L{Symbol}.
940          @return: result of the substitution process. Operations are executed as much as possible.
941          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
942          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
943          """
944          if argvals.has_key(self):
945             arg=argvals[self]
946             if self.isAppropriateValue(arg):
947                return arg
948             else:
949                raise TypeError,"%s: new value is not appropriate."%str(self)
950          else:
951             args=self.getSubstitutedArguments(argvals)
952             arg=args[0]
953             index=args[1]
954             return arg.__getitem__(index)
955    
956  def log10(arg):  def log10(arg):
957     """     """
958     returns base-10 logarithm of argument arg     returns base-10 logarithm of argument arg
959    
960     @param arg: argument     @param arg: argument
961     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
962     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
963     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
964     """     """
965     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 903  def wherePositive(arg): Line 981  def wherePositive(arg):
981    
982     @param arg: argument     @param arg: argument
983     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
984     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
985     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
986     """     """
987     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
988        if arg.rank==0:        out=numarray.greater(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
989           if arg>0:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
990             return numarray.array(1.)        return out
          else:  
            return numarray.array(0.)  
       else:  
          return numarray.greater(arg,numarray.zeros(arg.shape,numarray.Float))  
991     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
992        return arg._wherePositive()        return arg._wherePositive()
993     elif isinstance(arg,float):     elif isinstance(arg,float):
# Line 953  class WherePositive_Symbol(DependendSymb Line 1027  class WherePositive_Symbol(DependendSymb
1027        @type format: C{str}        @type format: C{str}
1028        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1029        @rtype: C{str}        @rtype: C{str}
1030        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1031        """        """
1032        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1033            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 989  def whereNegative(arg): Line 1063  def whereNegative(arg):
1063    
1064     @param arg: argument     @param arg: argument
1065     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1066     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1067     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1068     """     """
1069     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
1070        if arg.rank==0:        out=numarray.less(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1071           if arg<0:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1072             return numarray.array(1.)        return out
          else:  
            return numarray.array(0.)  
       else:  
          return numarray.less(arg,numarray.zeros(arg.shape,numarray.Float))  
1073     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
1074        return arg._whereNegative()        return arg._whereNegative()
1075     elif isinstance(arg,float):     elif isinstance(arg,float):
# Line 1039  class WhereNegative_Symbol(DependendSymb Line 1109  class WhereNegative_Symbol(DependendSymb
1109        @type format: C{str}        @type format: C{str}
1110        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1111        @rtype: C{str}        @rtype: C{str}
1112        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1113        """        """
1114        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1115            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1075  def whereNonNegative(arg): Line 1145  def whereNonNegative(arg):
1145    
1146     @param arg: argument     @param arg: argument
1147     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1148     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1149     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1150     """     """
1151     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
1152        if arg.rank==0:        out=numarray.greater_equal(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1153           if arg<0:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1154             return numarray.array(0.)        return out
          else:  
            return numarray.array(1.)  
       else:  
          return numarray.greater_equal(arg,numarray.zeros(arg.shape,numarray.Float))  
1155     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
1156        return arg._whereNonNegative()        return arg._whereNonNegative()
1157     elif isinstance(arg,float):     elif isinstance(arg,float):
# Line 1109  def whereNonPositive(arg): Line 1175  def whereNonPositive(arg):
1175    
1176     @param arg: argument     @param arg: argument
1177     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1178     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1179     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1180     """     """
1181     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
1182        if arg.rank==0:        out=numarray.less_equal(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1183           if arg>0:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1184             return numarray.array(0.)        return out
          else:  
            return numarray.array(1.)  
       else:  
          return numarray.less_equal(arg,numarray.zeros(arg.shape,numarray.Float))*1.  
1185     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
1186        return arg._whereNonPositive()        return arg._whereNonPositive()
1187     elif isinstance(arg,float):     elif isinstance(arg,float):
# Line 1145  def whereZero(arg,tol=0.): Line 1207  def whereZero(arg,tol=0.):
1207     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1208     @param tol: tolerance. values with absolute value less then tol are accepted as zero.     @param tol: tolerance. values with absolute value less then tol are accepted as zero.
1209     @type tol: C{float}     @type tol: C{float}
1210     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1211     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1212     """     """
1213     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
1214        if arg.rank==0:        out=numarray.less_equal(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float64))*1.
1215           if abs(arg)<=tol:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1216             return numarray.array(1.)        return out
          else:  
            return numarray.array(0.)  
       else:  
          return numarray.less_equal(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float))*1.  
1217     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
1218        if tol>0.:        return arg._whereZero(tol)
          return whereNegative(abs(arg)-tol)  
       else:  
          return arg._whereZero()  
1219     elif isinstance(arg,float):     elif isinstance(arg,float):
1220        if abs(arg)<=tol:        if abs(arg)<=tol:
1221          return 1.          return 1.
# Line 1198  class WhereZero_Symbol(DependendSymbol): Line 1253  class WhereZero_Symbol(DependendSymbol):
1253        @type format: C{str}        @type format: C{str}
1254        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1255        @rtype: C{str}        @rtype: C{str}
1256        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1257        """        """
1258        if format=="escript" or format=="str"  or format=="text":        if format=="escript" or format=="str"  or format=="text":
1259           return "whereZero(%s,tol=%s)"%(argstrs[0],argstrs[1])           return "whereZero(%s,tol=%s)"%(argstrs[0],argstrs[1])
# Line 1232  def whereNonZero(arg,tol=0.): Line 1287  def whereNonZero(arg,tol=0.):
1287    
1288     @param arg: argument     @param arg: argument
1289     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1290     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1291     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1292     """     """
1293     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
1294        if arg.rank==0:        out=numarray.greater(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float64))*1.
1295          if abs(arg)>tol:        if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1296             return numarray.array(1.)        return out
         else:  
            return numarray.array(0.)  
       else:  
          return numarray.greater(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float))*1.  
1297     elif isinstance(arg,escript.Data):     elif isinstance(arg,escript.Data):
1298        if tol>0.:        return arg._whereNonZero(tol)
          return 1.-whereZero(arg,tol)  
       else:  
          return arg._whereNonZero()  
1299     elif isinstance(arg,float):     elif isinstance(arg,float):
1300        if abs(arg)>tol:        if abs(arg)>tol:
1301          return 1.          return 1.
# Line 1269  def sin(arg): Line 1317  def sin(arg):
1317    
1318     @param arg: argument     @param arg: argument
1319     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1320     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1321     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1322     """     """
1323     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1307  class Sin_Symbol(DependendSymbol): Line 1355  class Sin_Symbol(DependendSymbol):
1355        @type format: C{str}        @type format: C{str}
1356        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1357        @rtype: C{str}        @rtype: C{str}
1358        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1359        """        """
1360        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1361            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1359  def cos(arg): Line 1407  def cos(arg):
1407    
1408     @param arg: argument     @param arg: argument
1409     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1410     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1411     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1412     """     """
1413     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1397  class Cos_Symbol(DependendSymbol): Line 1445  class Cos_Symbol(DependendSymbol):
1445        @type format: C{str}        @type format: C{str}
1446        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1447        @rtype: C{str}        @rtype: C{str}
1448        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1449        """        """
1450        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1451            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1449  def tan(arg): Line 1497  def tan(arg):
1497    
1498     @param arg: argument     @param arg: argument
1499     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1500     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1501     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1502     """     """
1503     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1487  class Tan_Symbol(DependendSymbol): Line 1535  class Tan_Symbol(DependendSymbol):
1535        @type format: C{str}        @type format: C{str}
1536        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1537        @rtype: C{str}        @rtype: C{str}
1538        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1539        """        """
1540        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1541            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1539  def asin(arg): Line 1587  def asin(arg):
1587    
1588     @param arg: argument     @param arg: argument
1589     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1590     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1591     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1592     """     """
1593     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1577  class Asin_Symbol(DependendSymbol): Line 1625  class Asin_Symbol(DependendSymbol):
1625        @type format: C{str}        @type format: C{str}
1626        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1627        @rtype: C{str}        @rtype: C{str}
1628        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1629        """        """
1630        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1631            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1629  def acos(arg): Line 1677  def acos(arg):
1677    
1678     @param arg: argument     @param arg: argument
1679     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1680     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1681     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1682     """     """
1683     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1667  class Acos_Symbol(DependendSymbol): Line 1715  class Acos_Symbol(DependendSymbol):
1715        @type format: C{str}        @type format: C{str}
1716        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1717        @rtype: C{str}        @rtype: C{str}
1718        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1719        """        """
1720        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1721            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1719  def atan(arg): Line 1767  def atan(arg):
1767    
1768     @param arg: argument     @param arg: argument
1769     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1770     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1771     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1772     """     """
1773     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1757  class Atan_Symbol(DependendSymbol): Line 1805  class Atan_Symbol(DependendSymbol):
1805        @type format: C{str}        @type format: C{str}
1806        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1807        @rtype: C{str}        @rtype: C{str}
1808        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1809        """        """
1810        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1811            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1809  def sinh(arg): Line 1857  def sinh(arg):
1857    
1858     @param arg: argument     @param arg: argument
1859     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1860     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1861     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1862     """     """
1863     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1847  class Sinh_Symbol(DependendSymbol): Line 1895  class Sinh_Symbol(DependendSymbol):
1895        @type format: C{str}        @type format: C{str}
1896        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1897        @rtype: C{str}        @rtype: C{str}
1898        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1899        """        """
1900        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1901            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1899  def cosh(arg): Line 1947  def cosh(arg):
1947    
1948     @param arg: argument     @param arg: argument
1949     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1950     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1951     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
1952     """     """
1953     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 1937  class Cosh_Symbol(DependendSymbol): Line 1985  class Cosh_Symbol(DependendSymbol):
1985        @type format: C{str}        @type format: C{str}
1986        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1987        @rtype: C{str}        @rtype: C{str}
1988        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
1989        """        """
1990        if isinstance(argstrs,list):        if isinstance(argstrs,list):
1991            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 1989  def tanh(arg): Line 2037  def tanh(arg):
2037    
2038     @param arg: argument     @param arg: argument
2039     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2040     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2041     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2042     """     """
2043     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2027  class Tanh_Symbol(DependendSymbol): Line 2075  class Tanh_Symbol(DependendSymbol):
2075        @type format: C{str}        @type format: C{str}
2076        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2077        @rtype: C{str}        @rtype: C{str}
2078        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2079        """        """
2080        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2081            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2079  def asinh(arg): Line 2127  def asinh(arg):
2127    
2128     @param arg: argument     @param arg: argument
2129     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2130     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2131     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2132     """     """
2133     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2117  class Asinh_Symbol(DependendSymbol): Line 2165  class Asinh_Symbol(DependendSymbol):
2165        @type format: C{str}        @type format: C{str}
2166        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2167        @rtype: C{str}        @rtype: C{str}
2168        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2169        """        """
2170        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2171            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2169  def acosh(arg): Line 2217  def acosh(arg):
2217    
2218     @param arg: argument     @param arg: argument
2219     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2220     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2221     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2222     """     """
2223     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2207  class Acosh_Symbol(DependendSymbol): Line 2255  class Acosh_Symbol(DependendSymbol):
2255        @type format: C{str}        @type format: C{str}
2256        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2257        @rtype: C{str}        @rtype: C{str}
2258        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2259        """        """
2260        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2261            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2259  def atanh(arg): Line 2307  def atanh(arg):
2307    
2308     @param arg: argument     @param arg: argument
2309     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2310     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2311     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2312     """     """
2313     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2297  class Atanh_Symbol(DependendSymbol): Line 2345  class Atanh_Symbol(DependendSymbol):
2345        @type format: C{str}        @type format: C{str}
2346        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2347        @rtype: C{str}        @rtype: C{str}
2348        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2349        """        """
2350        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2351            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2349  def exp(arg): Line 2397  def exp(arg):
2397    
2398     @param arg: argument     @param arg: argument
2399     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2400     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2401     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2402     """     """
2403     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2387  class Exp_Symbol(DependendSymbol): Line 2435  class Exp_Symbol(DependendSymbol):
2435        @type format: C{str}        @type format: C{str}
2436        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2437        @rtype: C{str}        @rtype: C{str}
2438        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2439        """        """
2440        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2441            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2439  def sqrt(arg): Line 2487  def sqrt(arg):
2487    
2488     @param arg: argument     @param arg: argument
2489     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2490     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2491     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2492     """     """
2493     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2477  class Sqrt_Symbol(DependendSymbol): Line 2525  class Sqrt_Symbol(DependendSymbol):
2525        @type format: C{str}        @type format: C{str}
2526        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2527        @rtype: C{str}        @rtype: C{str}
2528        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2529        """        """
2530        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2531            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2529  def log(arg): Line 2577  def log(arg):
2577    
2578     @param arg: argument     @param arg: argument
2579     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2580     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2581     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2582     """     """
2583     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2567  class Log_Symbol(DependendSymbol): Line 2615  class Log_Symbol(DependendSymbol):
2615        @type format: C{str}        @type format: C{str}
2616        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2617        @rtype: C{str}        @rtype: C{str}
2618        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2619        """        """
2620        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2621            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2619  def sign(arg): Line 2667  def sign(arg):
2667    
2668     @param arg: argument     @param arg: argument
2669     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2670     @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2671     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2672     """     """
2673     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2667  class Abs_Symbol(DependendSymbol): Line 2715  class Abs_Symbol(DependendSymbol):
2715        @type format: C{str}        @type format: C{str}
2716        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2717        @rtype: C{str}        @rtype: C{str}
2718        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2719        """        """
2720        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2721            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2719  def minval(arg): Line 2767  def minval(arg):
2767    
2768     @param arg: argument     @param arg: argument
2769     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2770     @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2771     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2772     """     """
2773     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2760  class Minval_Symbol(DependendSymbol): Line 2808  class Minval_Symbol(DependendSymbol):
2808        @type format: C{str}        @type format: C{str}
2809        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2810        @rtype: C{str}        @rtype: C{str}
2811        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2812        """        """
2813        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2814            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2796  def maxval(arg): Line 2844  def maxval(arg):
2844    
2845     @param arg: argument     @param arg: argument
2846     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2847     @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2848     @raises TypeError: if the type of the argument is not expected.     @raises TypeError: if the type of the argument is not expected.
2849     """     """
2850     if isinstance(arg,numarray.NumArray):     if isinstance(arg,numarray.NumArray):
# Line 2837  class Maxval_Symbol(DependendSymbol): Line 2885  class Maxval_Symbol(DependendSymbol):
2885        @type format: C{str}        @type format: C{str}
2886        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2887        @rtype: C{str}        @rtype: C{str}
2888        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
2889        """        """
2890        if isinstance(argstrs,list):        if isinstance(argstrs,list):
2891            argstrs=argstrs[0]            argstrs=argstrs[0]
# Line 2873  def length(arg): Line 2921  def length(arg):
2921    
2922     @param arg: argument     @param arg: argument
2923     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.     @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2924     @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.     @rtype: C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2925     """     """
2926     return sqrt(inner(arg,arg))     return sqrt(inner(arg,arg))
2927    
2928    def trace(arg,axis_offset=0):
2929       """
2930       returns the trace of arg which the sum of arg[k,k] over k.
2931    
2932       @param arg: argument
2933       @type arg: L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2934       @param axis_offset: C{axis_offset} to components to sum over. C{axis_offset} must be non-negative and less than the rank of arg +1. The dimensions on component
2935                      C{axis_offset} and axis_offset+1 must be equal.
2936       @type axis_offset: C{int}
2937       @return: trace of arg. The rank of the returned object is minus 2 of the rank of arg.
2938       @rtype: L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2939       """
2940       if isinstance(arg,numarray.NumArray):
2941          sh=arg.shape
2942          if len(sh)<2:
2943            raise ValueError,"rank of argument must be greater than 1"
2944          if axis_offset<0 or axis_offset>len(sh)-2:
2945            raise ValueError,"axis_offset must be between 0 and %s"%len(sh)-2
2946          s1=1
2947          for i in range(axis_offset): s1*=sh[i]
2948          s2=1
2949          for i in range(axis_offset+2,len(sh)): s2*=sh[i]
2950          if not sh[axis_offset] == sh[axis_offset+1]:
2951            raise ValueError,"dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
2952          arg_reshaped=numarray.reshape(arg,(s1,sh[axis_offset],sh[axis_offset],s2))
2953          out=numarray.zeros([s1,s2],numarray.Float64)
2954          for i1 in range(s1):
2955            for i2 in range(s2):
2956                for j in range(sh[axis_offset]): out[i1,i2]+=arg_reshaped[i1,j,j,i2]
2957          out.resize(sh[:axis_offset]+sh[axis_offset+2:])
2958          return out
2959       elif isinstance(arg,escript.Data):
2960          if arg.getRank()<2:
2961            raise ValueError,"rank of argument must be greater than 1"
2962          if axis_offset<0 or axis_offset>arg.getRank()-2:
2963            raise ValueError,"axis_offset must be between 0 and %s"%arg.getRank()-2
2964          s=list(arg.getShape())        
2965          if not s[axis_offset] == s[axis_offset+1]:
2966            raise ValueError,"dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
2967          return arg._trace(axis_offset)
2968       elif isinstance(arg,float):
2969          raise TypeError,"illegal argument type float."
2970       elif isinstance(arg,int):
2971          raise TypeError,"illegal argument type int."
2972       elif isinstance(arg,Symbol):
2973          return Trace_Symbol(arg,axis_offset)
2974       else:
2975          raise TypeError,"Unknown argument type."
2976    
2977    class Trace_Symbol(DependendSymbol):
2978       """
2979       L{Symbol} representing the result of the trace function
2980       """
2981       def __init__(self,arg,axis_offset=0):
2982          """
2983          initialization of trace L{Symbol} with argument arg
2984          @param arg: argument of function
2985          @type arg: L{Symbol}.
2986          @param axis_offset: C{axis_offset} to components to sum over. C{axis_offset} must be non-negative and less than the rank of arg +1. The dimensions on component
2987                      C{axis_offset} and axis_offset+1 must be equal.
2988          @type axis_offset: C{int}
2989          """
2990          if arg.getRank()<2:
2991            raise ValueError,"rank of argument must be greater than 1"
2992          if axis_offset<0 or axis_offset>arg.getRank()-2:
2993            raise ValueError,"axis_offset must be between 0 and %s"%arg.getRank()-2
2994          s=list(arg.getShape())        
2995          if not s[axis_offset] == s[axis_offset+1]:
2996            raise ValueError,"dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
2997          super(Trace_Symbol,self).__init__(args=[arg,axis_offset],shape=tuple(s[0:axis_offset]+s[axis_offset+2:]),dim=arg.getDim())
2998    
2999       def getMyCode(self,argstrs,format="escript"):
3000          """
3001          returns a program code that can be used to evaluate the symbol.
3002    
3003          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3004          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3005          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3006          @type format: C{str}
3007          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3008          @rtype: C{str}
3009          @raise NotImplementedError: if the requested format is not available
3010          """
3011          if format=="escript" or format=="str"  or format=="text":
3012             return "trace(%s,axis_offset=%s)"%(argstrs[0],argstrs[1])
3013          else:
3014             raise NotImplementedError,"Trace_Symbol does not provide program code for format %s."%format
3015    
3016       def substitute(self,argvals):
3017          """
3018          assigns new values to symbols in the definition of the symbol.
3019          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3020    
3021          @param argvals: new values assigned to symbols
3022          @type argvals: C{dict} with keywords of type L{Symbol}.
3023          @return: result of the substitution process. Operations are executed as much as possible.
3024          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3025          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3026          """
3027          if argvals.has_key(self):
3028             arg=argvals[self]
3029             if self.isAppropriateValue(arg):
3030                return arg
3031             else:
3032                raise TypeError,"%s: new value is not appropriate."%str(self)
3033          else:
3034             arg=self.getSubstitutedArguments(argvals)
3035             return trace(arg[0],axis_offset=arg[1])
3036    
3037       def diff(self,arg):
3038          """
3039          differential of this object
3040    
3041          @param arg: the derivative is calculated with respect to arg
3042          @type arg: L{escript.Symbol}
3043          @return: derivative with respect to C{arg}
3044          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3045          """
3046          if arg==self:
3047             return identity(self.getShape())
3048          else:
3049             return trace(self.getDifferentiatedArguments(arg)[0],axis_offset=self.getArgument()[1])
3050    
3051    def transpose(arg,axis_offset=None):
3052       """
3053       returns the transpose of arg by swaping the first C{axis_offset} and the last rank-axis_offset components.
3054    
3055       @param arg: argument
3056       @type arg: L{escript.Data}, L{Symbol}, L{numarray.NumArray}, C{float}, C{int}
3057       @param axis_offset: the first C{axis_offset} components are swapped with rest. If C{axis_offset} must be non-negative and less or equal the rank of arg.
3058                           if C{axis_offset} is not present C{int(r/2)} where r is the rank of arg is used.
3059       @type axis_offset: C{int}
3060       @return: transpose of arg
3061       @rtype: L{escript.Data}, L{Symbol}, L{numarray.NumArray},C{float}, C{int} depending on the type of arg.
3062       """
3063       if isinstance(arg,numarray.NumArray):
3064          if axis_offset==None: axis_offset=int(arg.rank/2)
3065          return numarray.transpose(arg,axes=range(axis_offset,arg.rank)+range(0,axis_offset))
3066       elif isinstance(arg,escript.Data):
3067          r=arg.getRank()
3068          if axis_offset==None: axis_offset=int(r/2)
3069          if axis_offset<0 or axis_offset>r:
3070            raise ValueError,"axis_offset must be between 0 and %s"%r
3071          return arg._transpose(axis_offset)
3072       elif isinstance(arg,float):
3073          if not ( axis_offset==0 or axis_offset==None):
3074            raise ValueError,"axis_offset must be 0 for float argument"
3075          return arg
3076       elif isinstance(arg,int):
3077          if not ( axis_offset==0 or axis_offset==None):
3078            raise ValueError,"axis_offset must be 0 for int argument"
3079          return float(arg)
3080       elif isinstance(arg,Symbol):
3081          if axis_offset==None: axis_offset=int(arg.getRank()/2)
3082          return Transpose_Symbol(arg,axis_offset)
3083       else:
3084          raise TypeError,"Unknown argument type."
3085    
3086    class Transpose_Symbol(DependendSymbol):
3087       """
3088       L{Symbol} representing the result of the transpose function
3089       """
3090       def __init__(self,arg,axis_offset=None):
3091          """
3092          initialization of transpose L{Symbol} with argument arg
3093    
3094          @param arg: argument of function
3095          @type arg: L{Symbol}.
3096          @param axis_offset: the first C{axis_offset} components are swapped with rest. If C{axis_offset} must be non-negative and less or equal the rank of arg.
3097                           if C{axis_offset} is not present C{int(r/2)} where r is the rank of arg is used.
3098          @type axis_offset: C{int}
3099          """
3100          if axis_offset==None: axis_offset=int(arg.getRank()/2)
3101          if axis_offset<0 or axis_offset>arg.getRank():
3102            raise ValueError,"axis_offset must be between 0 and %s"%arg.getRank()
3103          s=arg.getShape()
3104          super(Transpose_Symbol,self).__init__(args=[arg,axis_offset],shape=s[axis_offset:]+s[:axis_offset],dim=arg.getDim())
3105    
3106       def getMyCode(self,argstrs,format="escript"):
3107          """
3108          returns a program code that can be used to evaluate the symbol.
3109    
3110          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3111          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3112          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3113          @type format: C{str}
3114          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3115          @rtype: C{str}
3116          @raise NotImplementedError: if the requested format is not available
3117          """
3118          if format=="escript" or format=="str"  or format=="text":
3119             return "transpose(%s,axis_offset=%s)"%(argstrs[0],argstrs[1])
3120          else:
3121             raise NotImplementedError,"Transpose_Symbol does not provide program code for format %s."%format
3122    
3123       def substitute(self,argvals):
3124          """
3125          assigns new values to symbols in the definition of the symbol.
3126          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3127    
3128          @param argvals: new values assigned to symbols
3129          @type argvals: C{dict} with keywords of type L{Symbol}.
3130          @return: result of the substitution process. Operations are executed as much as possible.
3131          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3132          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3133          """
3134          if argvals.has_key(self):
3135             arg=argvals[self]
3136             if self.isAppropriateValue(arg):
3137                return arg
3138             else:
3139                raise TypeError,"%s: new value is not appropriate."%str(self)
3140          else:
3141             arg=self.getSubstitutedArguments(argvals)
3142             return transpose(arg[0],axis_offset=arg[1])
3143    
3144       def diff(self,arg):
3145          """
3146          differential of this object
3147    
3148          @param arg: the derivative is calculated with respect to arg
3149          @type arg: L{escript.Symbol}
3150          @return: derivative with respect to C{arg}
3151          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3152          """
3153          if arg==self:
3154             return identity(self.getShape())
3155          else:
3156             return transpose(self.getDifferentiatedArguments(arg)[0],axis_offset=self.getArgument()[1])
3157    
3158    def swap_axes(arg,axis0=0,axis1=1):
3159       """
3160       returns the swap of arg by swaping the components axis0 and axis1
3161    
3162       @param arg: argument
3163       @type arg: L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
3164       @param axis0: axis. C{axis0} must be non-negative and less than the rank of arg.
3165       @type axis0: C{int}
3166       @param axis1: axis. C{axis1} must be non-negative and less than the rank of arg.
3167       @type axis1: C{int}
3168       @return: C{arg} with swaped components
3169       @rtype: L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
3170       """
3171       if axis0 > axis1:
3172          axis0,axis1=axis1,axis0
3173       if isinstance(arg,numarray.NumArray):
3174          return numarray.swapaxes(arg,axis0,axis1)
3175       elif isinstance(arg,escript.Data):
3176          return arg._swap_axes(axis0,axis1)
3177       elif isinstance(arg,float):
3178          raise TyepError,"float argument is not supported."
3179       elif isinstance(arg,int):
3180          raise TyepError,"int argument is not supported."
3181       elif isinstance(arg,Symbol):
3182          return SwapAxes_Symbol(arg,axis0,axis1)
3183       else:
3184          raise TypeError,"Unknown argument type."
3185    
3186    class SwapAxes_Symbol(DependendSymbol):
3187       """
3188       L{Symbol} representing the result of the swap function
3189       """
3190       def __init__(self,arg,axis0=0,axis1=1):
3191          """
3192          initialization of swap L{Symbol} with argument arg
3193    
3194          @param arg: argument
3195          @type arg: L{Symbol}.
3196          @param axis0: axis. C{axis0} must be non-negative and less than the rank of arg.
3197          @type axis0: C{int}
3198          @param axis1: axis. C{axis1} must be non-negative and less than the rank of arg.
3199          @type axis1: C{int}
3200          """
3201          if arg.getRank()<2:
3202             raise ValueError,"argument must have at least rank 2."
3203          if axis0<0 or axis0>arg.getRank()-1:
3204             raise ValueError,"axis0 must be between 0 and %s"%arg.getRank()-1
3205          if axis1<0 or axis1>arg.getRank()-1:
3206             raise ValueError,"axis1 must be between 0 and %s"%arg.getRank()-1
3207          if axis0 == axis1:
3208             raise ValueError,"axis indices must be different."
3209          if axis0 > axis1:
3210             axis0,axis1=axis1,axis0
3211          s=arg.getShape()
3212          s_out=[]
3213          for i in range(len(s)):
3214             if i == axis0:
3215                s_out.append(s[axis1])
3216             elif i == axis1:
3217                s_out.append(s[axis0])
3218             else:
3219                s_out.append(s[i])
3220          super(SwapAxes_Symbol,self).__init__(args=[arg,axis0,axis1],shape=tuple(s_out),dim=arg.getDim())
3221    
3222       def getMyCode(self,argstrs,format="escript"):
3223          """
3224          returns a program code that can be used to evaluate the symbol.
3225    
3226          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3227          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3228          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3229          @type format: C{str}
3230          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3231          @rtype: C{str}
3232          @raise NotImplementedError: if the requested format is not available
3233          """
3234          if format=="escript" or format=="str"  or format=="text":
3235             return "swap(%s,axis_offset=%s)"%(argstrs[0],argstrs[1])
3236          else:
3237             raise NotImplementedError,"SwapAxes_Symbol does not provide program code for format %s."%format
3238    
3239       def substitute(self,argvals):
3240          """
3241          assigns new values to symbols in the definition of the symbol.
3242          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3243    
3244          @param argvals: new values assigned to symbols
3245          @type argvals: C{dict} with keywords of type L{Symbol}.
3246          @return: result of the substitution process. Operations are executed as much as possible.
3247          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3248          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3249          """
3250          if argvals.has_key(self):
3251             arg=argvals[self]
3252             if self.isAppropriateValue(arg):
3253                return arg
3254             else:
3255                raise TypeError,"%s: new value is not appropriate."%str(self)
3256          else:
3257             arg=self.getSubstitutedArguments(argvals)
3258             return swap_axes(arg[0],axis0=arg[1],axis1=arg[2])
3259    
3260       def diff(self,arg):
3261          """
3262          differential of this object
3263    
3264          @param arg: the derivative is calculated with respect to arg
3265          @type arg: L{escript.Symbol}
3266          @return: derivative with respect to C{arg}
3267          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3268          """
3269          if arg==self:
3270             return identity(self.getShape())
3271          else:
3272             return swap_axes(self.getDifferentiatedArguments(arg)[0],axis0=self.getArgument()[1],axis1=self.getArgument()[2])
3273    
3274    def symmetric(arg):
3275        """
3276        returns the symmetric part of the square matrix arg. This is (arg+transpose(arg))/2
3277    
3278        @param arg: square matrix. Must have rank 2 or 4 and be square.
3279        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3280        @return: symmetric part of arg
3281        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
3282        """
3283        if isinstance(arg,numarray.NumArray):
3284          if arg.rank==2:
3285            if not (arg.shape[0]==arg.shape[1]):
3286               raise ValueError,"argument must be square."
3287          elif arg.rank==4:
3288            if not (arg.shape[0]==arg.shape[2] and arg.shape[1]==arg.shape[3]):
3289               raise ValueError,"argument must be square."
3290          else:
3291            raise ValueError,"rank 2 or 4 is required."
3292          return (arg+transpose(arg))/2
3293        elif isinstance(arg,escript.Data):
3294          if arg.getRank()==2:
3295            if not (arg.getShape()[0]==arg.getShape()[1]):
3296               raise ValueError,"argument must be square."
3297            return arg._symmetric()
3298          elif arg.getRank()==4:
3299            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3300               raise ValueError,"argument must be square."
3301            return arg._symmetric()
3302          else:
3303            raise ValueError,"rank 2 or 4 is required."
3304        elif isinstance(arg,float):
3305          return arg
3306        elif isinstance(arg,int):
3307          return float(arg)
3308        elif isinstance(arg,Symbol):
3309          if arg.getRank()==2:
3310            if not (arg.getShape()[0]==arg.getShape()[1]):
3311               raise ValueError,"argument must be square."
3312          elif arg.getRank()==4:
3313            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3314               raise ValueError,"argument must be square."
3315          else:
3316            raise ValueError,"rank 2 or 4 is required."
3317          return (arg+transpose(arg))/2
3318        else:
3319          raise TypeError,"symmetric: Unknown argument type."
3320    
3321    def nonsymmetric(arg):
3322        """
3323        returns the nonsymmetric part of the square matrix arg. This is (arg-transpose(arg))/2
3324    
3325        @param arg: square matrix. Must have rank 2 or 4 and be square.
3326        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3327        @return: nonsymmetric part of arg
3328        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
3329        """
3330        if isinstance(arg,numarray.NumArray):
3331          if arg.rank==2:
3332            if not (arg.shape[0]==arg.shape[1]):
3333               raise ValueError,"nonsymmetric: argument must be square."
3334          elif arg.rank==4:
3335            if not (arg.shape[0]==arg.shape[2] and arg.shape[1]==arg.shape[3]):
3336               raise ValueError,"nonsymmetric: argument must be square."
3337          else:
3338            raise ValueError,"nonsymmetric: rank 2 or 4 is required."
3339          return (arg-transpose(arg))/2
3340        elif isinstance(arg,escript.Data):
3341          if arg.getRank()==2:
3342            if not (arg.getShape()[0]==arg.getShape()[1]):
3343               raise ValueError,"argument must be square."
3344            return arg._nonsymmetric()
3345          elif arg.getRank()==4:
3346            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3347               raise ValueError,"argument must be square."
3348            return arg._nonsymmetric()
3349          else:
3350            raise ValueError,"rank 2 or 4 is required."
3351        elif isinstance(arg,float):
3352          return arg
3353        elif isinstance(arg,int):
3354          return float(arg)
3355        elif isinstance(arg,Symbol):
3356          if arg.getRank()==2:
3357            if not (arg.getShape()[0]==arg.getShape()[1]):
3358               raise ValueError,"nonsymmetric: argument must be square."
3359          elif arg.getRank()==4:
3360            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3361               raise ValueError,"nonsymmetric: argument must be square."
3362          else:
3363            raise ValueError,"nonsymmetric: rank 2 or 4 is required."
3364          return (arg-transpose(arg))/2
3365        else:
3366          raise TypeError,"nonsymmetric: Unknown argument type."
3367    
3368    def inverse(arg):
3369        """
3370        returns the inverse of the square matrix arg.
3371    
3372        @param arg: square matrix. Must have rank 2 and the first and second dimension must be equal.
3373        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3374        @return: inverse arg_inv of the argument. It will be matrix_mult(inverse(arg),arg) almost equal to kronecker(arg.getShape()[0])
3375        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
3376        @note: for L{escript.Data} objects the dimension is restricted to 3.
3377        """
3378        import numarray.linear_algebra # This statement should be after the next statement but then somehow numarray is gone.
3379        if isinstance(arg,numarray.NumArray):
3380          return numarray.linear_algebra.inverse(arg)
3381        elif isinstance(arg,escript.Data):
3382          return escript_inverse(arg)
3383        elif isinstance(arg,float):
3384          return 1./arg
3385        elif isinstance(arg,int):
3386          return 1./float(arg)
3387        elif isinstance(arg,Symbol):
3388          return Inverse_Symbol(arg)
3389        else:
3390          raise TypeError,"inverse: Unknown argument type."
3391    
3392    def escript_inverse(arg): # this should be escript._inverse and use LAPACK
3393          "arg is a Data objects!!!"
3394          if not arg.getRank()==2:
3395            raise ValueError,"escript_inverse: argument must have rank 2"
3396          s=arg.getShape()      
3397          if not s[0] == s[1]:
3398            raise ValueError,"escript_inverse: argument must be a square matrix."
3399          out=escript.Data(0.,s,arg.getFunctionSpace())
3400          if s[0]==1:
3401              if inf(abs(arg[0,0]))==0: # in c this should be done point wise as abs(arg[0,0](i))<=0.
3402                  raise ZeroDivisionError,"escript_inverse: argument not invertible"
3403              out[0,0]=1./arg[0,0]
3404          elif s[0]==2:
3405              A11=arg[0,0]
3406              A12=arg[0,1]
3407              A21=arg[1,0]
3408              A22=arg[1,1]
3409              D = A11*A22-A12*A21
3410              if inf(abs(D))==0: # in c this should be done point wise as abs(D(i))<=0.
3411                  raise ZeroDivisionError,"escript_inverse: argument not invertible"
3412              D=1./D
3413              out[0,0]= A22*D
3414              out[1,0]=-A21*D
3415              out[0,1]=-A12*D
3416              out[1,1]= A11*D
3417          elif s[0]==3:
3418              A11=arg[0,0]
3419              A21=arg[1,0]
3420              A31=arg[2,0]
3421              A12=arg[0,1]
3422              A22=arg[1,1]
3423              A32=arg[2,1]
3424              A13=arg[0,2]
3425              A23=arg[1,2]
3426              A33=arg[2,2]
3427              D  =  A11*(A22*A33-A23*A32)+ A12*(A31*A23-A21*A33)+A13*(A21*A32-A31*A22)
3428              if inf(abs(D))==0: # in c this should be done point wise as abs(D(i))<=0.
3429                  raise ZeroDivisionError,"escript_inverse: argument not invertible"
3430              D=1./D
3431              out[0,0]=(A22*A33-A23*A32)*D
3432              out[1,0]=(A31*A23-A21*A33)*D
3433              out[2,0]=(A21*A32-A31*A22)*D
3434              out[0,1]=(A13*A32-A12*A33)*D
3435              out[1,1]=(A11*A33-A31*A13)*D
3436              out[2,1]=(A12*A31-A11*A32)*D
3437              out[0,2]=(A12*A23-A13*A22)*D
3438              out[1,2]=(A13*A21-A11*A23)*D
3439              out[2,2]=(A11*A22-A12*A21)*D
3440          else:
3441             raise TypeError,"escript_inverse: only matrix dimensions 1,2,3 are supported right now."
3442          return out
3443    
3444    class Inverse_Symbol(DependendSymbol):
3445       """
3446       L{Symbol} representing the result of the inverse function
3447       """
3448       def __init__(self,arg):
3449          """
3450          initialization of inverse L{Symbol} with argument arg
3451          @param arg: argument of function
3452          @type arg: L{Symbol}.
3453          """
3454          if not arg.getRank()==2:
3455            raise ValueError,"Inverse_Symbol:: argument must have rank 2"
3456          s=arg.getShape()
3457          if not s[0] == s[1]:
3458            raise ValueError,"Inverse_Symbol:: argument must be a square matrix."
3459          super(Inverse_Symbol,self).__init__(args=[arg],shape=s,dim=arg.getDim())
3460    
3461       def getMyCode(self,argstrs,format="escript"):
3462          """
3463          returns a program code that can be used to evaluate the symbol.
3464    
3465          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3466          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3467          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3468          @type format: C{str}
3469          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3470          @rtype: C{str}
3471          @raise NotImplementedError: if the requested format is not available
3472          """
3473          if format=="escript" or format=="str"  or format=="text":
3474             return "inverse(%s)"%argstrs[0]
3475          else:
3476             raise NotImplementedError,"Inverse_Symbol does not provide program code for format %s."%format
3477    
3478       def substitute(self,argvals):
3479          """
3480          assigns new values to symbols in the definition of the symbol.
3481          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3482    
3483          @param argvals: new values assigned to symbols
3484          @type argvals: C{dict} with keywords of type L{Symbol}.
3485          @return: result of the substitution process. Operations are executed as much as possible.
3486          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3487          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3488          """
3489          if argvals.has_key(self):
3490             arg=argvals[self]
3491             if self.isAppropriateValue(arg):
3492                return arg
3493             else:
3494                raise TypeError,"%s: new value is not appropriate."%str(self)
3495          else:
3496             arg=self.getSubstitutedArguments(argvals)
3497             return inverse(arg[0])
3498    
3499       def diff(self,arg):
3500          """
3501          differential of this object
3502    
3503          @param arg: the derivative is calculated with respect to arg
3504          @type arg: L{escript.Symbol}
3505          @return: derivative with respect to C{arg}
3506          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3507          """
3508          if arg==self:
3509             return identity(self.getShape())
3510          else:
3511             return -matrix_mult(matrix_mult(self,self.getDifferentiatedArguments(arg)[0]),self)
3512    
3513    def eigenvalues(arg):
3514        """
3515        returns the eigenvalues of the square matrix arg.
3516    
3517        @param arg: square matrix. Must have rank 2 and the first and second dimension must be equal.
3518                    arg must be symmetric, ie. transpose(arg)==arg (this is not checked).
3519        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3520        @return: the eigenvalues in increasing order.
3521        @rtype: L{numarray.NumArray},L{escript.Data}, L{Symbol} depending on the input.
3522        @note: for L{escript.Data} and L{Symbol} objects the dimension is restricted to 3.
3523        """
3524        if isinstance(arg,numarray.NumArray):
3525          out=numarray.linear_algebra.eigenvalues((arg+numarray.transpose(arg))/2.)
3526          out.sort()
3527          return out
3528        elif isinstance(arg,escript.Data):
3529          return arg._eigenvalues()
3530        elif isinstance(arg,Symbol):
3531          if not arg.getRank()==2:
3532            raise ValueError,"eigenvalues: argument must have rank 2"
3533          s=arg.getShape()      
3534          if not s[0] == s[1]:
3535            raise ValueError,"eigenvalues: argument must be a square matrix."
3536          if s[0]==1:
3537              return arg[0]
3538          elif s[0]==2:
3539              arg1=symmetric(arg)
3540              A11=arg1[0,0]
3541              A12=arg1[0,1]
3542              A22=arg1[1,1]
3543              trA=(A11+A22)/2.
3544              A11-=trA
3545              A22-=trA
3546              s=sqrt(A12**2-A11*A22)
3547              return trA+s*numarray.array([-1.,1.],type=numarray.Float64)
3548          elif s[0]==3:
3549              arg1=symmetric(arg)
3550              A11=arg1[0,0]
3551              A12=arg1[0,1]
3552              A22=arg1[1,1]
3553              A13=arg1[0,2]
3554              A23=arg1[1,2]
3555              A33=arg1[2,2]
3556              trA=(A11+A22+A33)/3.
3557              A11-=trA
3558              A22-=trA
3559              A33-=trA
3560              A13_2=A13**2
3561              A23_2=A23**2
3562              A12_2=A12**2
3563              p=A13_2+A23_2+A12_2+(A11**2+A22**2+A33**2)/2.
3564              q=A13_2*A22+A23_2*A11+A12_2*A33-A11*A22*A33-2*A12*A23*A13
3565              sq_p=sqrt(p/3.)
3566              alpha_3=acos(clip(-q*(sq_p+whereZero(p,0.)*1.e-15)**(-3.)/2.,-1.,1.))/3.  # whereZero is protection against divison by zero
3567              sq_p*=2.
3568              f=cos(alpha_3)               *numarray.array([0.,0.,1.],type=numarray.Float64) \
3569               -cos(alpha_3+numarray.pi/3.)*numarray.array([0.,1.,0.],type=numarray.Float64) \
3570               -cos(alpha_3-numarray.pi/3.)*numarray.array([1.,0.,0.],type=numarray.Float64)
3571              return trA+sq_p*f
3572          else:
3573             raise TypeError,"eigenvalues: only matrix dimensions 1,2,3 are supported right now."
3574        elif isinstance(arg,float):
3575          return arg
3576        elif isinstance(arg,int):
3577          return float(arg)
3578        else:
3579          raise TypeError,"eigenvalues: Unknown argument type."
3580    
3581    def eigenvalues_and_eigenvectors(arg):
3582        """
3583        returns the eigenvalues and eigenvectors of the square matrix arg.
3584    
3585        @param arg: square matrix. Must have rank 2 and the first and second dimension must be equal.
3586                    arg must be symmetric, ie. transpose(arg)==arg (this is not checked).
3587        @type arg: L{escript.Data}
3588        @return: the eigenvalues and eigenvectors. The eigenvalues are ordered by increasing value. The
3589                 eigenvectors are orthogonal and normalized. If V are the eigenvectors than V[:,i] is
3590                 the eigenvector coresponding to the i-th eigenvalue.
3591        @rtype: L{tuple} of L{escript.Data}.
3592        @note: The dimension is restricted to 3.
3593        """
3594        if isinstance(arg,numarray.NumArray):
3595          raise TypeError,"eigenvalues_and_eigenvectors is not supporting numarray arguments"
3596        elif isinstance(arg,escript.Data):
3597          return arg._eigenvalues_and_eigenvectors()
3598        elif isinstance(arg,Symbol):
3599          raise TypeError,"eigenvalues_and_eigenvectors is not supporting Symbol arguments"
3600        elif isinstance(arg,float):
3601          return (numarray.array([[arg]],numarray.Float),numarray.ones((1,1),numarray.Float))
3602        elif isinstance(arg,int):
3603          return (numarray.array([[arg]],numarray.Float),numarray.ones((1,1),numarray.Float))
3604        else:
3605          raise TypeError,"eigenvalues: Unknown argument type."
3606  #=======================================================  #=======================================================
3607  #  Binary operations:  #  Binary operations:
3608  #=======================================================  #=======================================================
# Line 2936  class Add_Symbol(DependendSymbol): Line 3662  class Add_Symbol(DependendSymbol):
3662        @type format: C{str}        @type format: C{str}
3663        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3664        @rtype: C{str}        @rtype: C{str}
3665        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
3666        """        """
3667        if format=="str" or format=="text":        if format=="str" or format=="text":
3668           return "(%s)+(%s)"%(argstrs[0],argstrs[1])           return "(%s)+(%s)"%(argstrs[0],argstrs[1])
# Line 2995  def mult(arg0,arg1): Line 3721  def mult(arg0,arg1):
3721         """         """
3722         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
3723         if testForZero(args[0]) or testForZero(args[1]):         if testForZero(args[0]) or testForZero(args[1]):
3724            return numarray.zeros(pokeShape(args[0]),numarray.Float)            return numarray.zeros(pokeShape(args[0]),numarray.Float64)
3725         else:         else:
3726            if isinstance(args[0],Symbol) or isinstance(args[1],Symbol) :            if isinstance(args[0],Symbol) or isinstance(args[1],Symbol) :
3727                return Mult_Symbol(args[0],args[1])                return Mult_Symbol(args[0],args[1])
# Line 3035  class Mult_Symbol(DependendSymbol): Line 3761  class Mult_Symbol(DependendSymbol):
3761        @type format: C{str}        @type format: C{str}
3762        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3763        @rtype: C{str}        @rtype: C{str}
3764        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
3765        """        """
3766        if format=="str" or format=="text":        if format=="str" or format=="text":
3767           return "(%s)*(%s)"%(argstrs[0],argstrs[1])           return "(%s)*(%s)"%(argstrs[0],argstrs[1])
# Line 3095  def quotient(arg0,arg1): Line 3821  def quotient(arg0,arg1):
3821         """         """
3822         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
3823         if testForZero(args[0]):         if testForZero(args[0]):
3824            return numarray.zeros(pokeShape(args[0]),numarray.Float)            return numarray.zeros(pokeShape(args[0]),numarray.Float64)
3825         elif isinstance(args[0],Symbol):         elif isinstance(args[0],Symbol):
3826            if isinstance(args[1],Symbol):            if isinstance(args[1],Symbol):
3827               return Quotient_Symbol(args[0],args[1])               return Quotient_Symbol(args[0],args[1])
# Line 3140  class Quotient_Symbol(DependendSymbol): Line 3866  class Quotient_Symbol(DependendSymbol):
3866        @type format: C{str}        @type format: C{str}
3867        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3868        @rtype: C{str}        @rtype: C{str}
3869        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
3870        """        """
3871        if format=="str" or format=="text":        if format=="str" or format=="text":
3872           return "(%s)/(%s)"%(argstrs[0],argstrs[1])           return "(%s)/(%s)"%(argstrs[0],argstrs[1])
# Line 3201  def power(arg0,arg1): Line 3927  def power(arg0,arg1):
3927         """         """
3928         args=matchShape(arg0,arg1)         args=matchShape(arg0,arg1)
3929         if testForZero(args[0]):         if testForZero(args[0]):
3930            return numarray.zeros(args[0],numarray.Float)            return numarray.zeros(pokeShape(args[0]),numarray.Float64)
3931         elif testForZero(args[1]):         elif testForZero(args[1]):
3932            return numarray.ones(args[0],numarray.Float)            return numarray.ones(pokeShape(args[1]),numarray.Float64)
3933         elif isinstance(args[0],Symbol) or isinstance(args[1],Symbol):         elif isinstance(args[0],Symbol) or isinstance(args[1],Symbol):
3934            return Power_Symbol(args[0],args[1])            return Power_Symbol(args[0],args[1])
3935         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 3244  class Power_Symbol(DependendSymbol): Line 3970  class Power_Symbol(DependendSymbol):
3970        @type format: C{str}        @type format: C{str}
3971        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3972        @rtype: C{str}        @rtype: C{str}
3973        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
3974        """        """
3975        if format=="escript" or format=="str" or format=="text":        if format=="escript" or format=="str" or format=="text":
3976           return "(%s)**(%s)"%(argstrs[0],argstrs[1])           return "(%s)**(%s)"%(argstrs[0],argstrs[1])
# Line 3304  def maximum(*args): Line 4030  def maximum(*args):
4030         if out==None:         if out==None:
4031            out=a            out=a
4032         else:         else:
4033            m=whereNegative(out-a)            diff=add(a,-out)
4034            out=m*a+(1.-m)*out            out=add(out,mult(wherePositive(diff),diff))
4035      return out      return out
4036        
4037  def minimum(*arg):  def minimum(*args):
4038      """      """
4039      the minimum over arguments args      the minimum over arguments args
4040    
# Line 3322  def minimum(*arg): Line 4048  def minimum(*arg):
4048         if out==None:         if out==None:
4049            out=a            out=a
4050         else:         else:
4051            m=whereNegative(out-a)            diff=add(a,-out)
4052            out=m*out+(1.-m)*a            out=add(out,mult(whereNegative(diff),diff))
4053      return out      return out
4054    
4055    def clip(arg,minval=0.,maxval=1.):
4056        """
4057        cuts the values of arg between minval and maxval
4058    
4059        @param arg: argument
4060        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float}
4061        @param minval: lower range
4062        @type minval: C{float}
4063        @param maxval: upper range
4064        @type maxval: C{float}
4065        @return: is on object with all its value between minval and maxval. value of the argument that greater then minval and
4066                 less then maxval are unchanged.
4067        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{int} or C{float} depending on the input
4068        @raise ValueError: if minval>maxval
4069        """
4070        if minval>maxval:
4071           raise ValueError,"minval = %s must be less then maxval %s"%(minval,maxval)
4072        return minimum(maximum(minval,arg),maxval)
4073    
4074        
4075  def inner(arg0,arg1):  def inner(arg0,arg1):
4076      """      """
# Line 3340  def inner(arg0,arg1): Line 4086  def inner(arg0,arg1):
4086      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4087      @param arg1: second argument      @param arg1: second argument
4088      @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}      @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4089      @return : the inner product of arg0 and arg1 at each data point      @return: the inner product of arg0 and arg1 at each data point
4090      @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
4091      @raise ValueError: if the shapes of the arguments are not identical      @raise ValueError: if the shapes of the arguments are not identical
4092      """      """
# Line 3348  def inner(arg0,arg1): Line 4094  def inner(arg0,arg1):
4094      sh1=pokeShape(arg1)      sh1=pokeShape(arg1)
4095      if not sh0==sh1:      if not sh0==sh1:
4096          raise ValueError,"inner: shape of arguments does not match"          raise ValueError,"inner: shape of arguments does not match"
4097      return generalTensorProduct(arg0,arg1,offset=len(sh0))      return generalTensorProduct(arg0,arg1,axis_offset=len(sh0))
4098    
4099    def outer(arg0,arg1):
4100        """
4101        the outer product of the two argument:
4102    
4103        out[t,s]=arg0[t]*arg1[s]
4104    
4105        where
4106    
4107            - s runs through arg0.Shape
4108            - t runs through arg1.Shape
4109    
4110        @param arg0: first argument
4111        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4112        @param arg1: second argument
4113        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4114        @return: the outer product of arg0 and arg1 at each data point
4115        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4116        """
4117        return generalTensorProduct(arg0,arg1,axis_offset=0)
4118    
4119  def matrixmult(arg0,arg1):  def matrixmult(arg0,arg1):
4120      """      """
4121        see L{matrix_mult}
4122        """
4123        return matrix_mult(arg0,arg1)
4124    
4125    def matrix_mult(arg0,arg1):
4126        """
4127      matrix-matrix or matrix-vector product of the two argument:      matrix-matrix or matrix-vector product of the two argument:
4128    
4129      out[s0]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0]      out[s0]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0]
4130    
4131            or      or
4132    
4133      out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0,s1]      out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0,s1]
4134    
4135      The second dimension of arg0 and the length of arg1 must match.      The second dimension of arg0 and the first dimension of arg1 must match.
4136    
4137      @param arg0: first argument of rank 2      @param arg0: first argument of rank 2
4138      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
# Line 3376  def matrixmult(arg0,arg1): Line 4148  def matrixmult(arg0,arg1):
4148          raise ValueError,"first argument must have rank 2"          raise ValueError,"first argument must have rank 2"
4149      if not len(sh1)==2 and not len(sh1)==1:      if not len(sh1)==2 and not len(sh1)==1:
4150          raise ValueError,"second argument must have rank 1 or 2"          raise ValueError,"second argument must have rank 1 or 2"
4151      return generalTensorProduct(arg0,arg1,offset=1)      return generalTensorProduct(arg0,arg1,axis_offset=1)
4152    
4153  def outer(arg0,arg1):  def tensormult(arg0,arg1):
4154      """      """
4155      the outer product of the two argument:      use L{tensor_mult}
   
     out[t,s]=arg0[t]*arg1[s]  
   
     where s runs through arg0.Shape  
           t runs through arg1.Shape  
   
     @param arg0: first argument  
     @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}  
     @param arg1: second argument  
     @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}  
     @return: the outer product of arg0 and arg1 at each data point  
     @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input  
4156      """      """
4157      return generalTensorProduct(arg0,arg1,offset=0)      return tensor_mult(arg0,arg1)
   
4158    
4159  def tensormult(arg0,arg1):  def tensor_mult(arg0,arg1):
4160      """      """
4161      the tensor product of the two argument:      the tensor product of the two argument:
   
4162            
4163      for arg0 of rank 2 this is      for arg0 of rank 2 this is
4164    
4165      out[s0]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0]        out[s0]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0]  
4166    
4167                   or      or
4168    
4169      out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0,s1]      out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[r0,s1]
4170    
# Line 3415  def tensormult(arg0,arg1): Line 4173  def tensormult(arg0,arg1):
4173    
4174      out[s0,s1,s2,s3]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1,s2,s3]      out[s0,s1,s2,s3]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1,s2,s3]
4175                                
4176                   or      or
4177    
4178      out[s0,s1,s2]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1,s2]      out[s0,s1,s2]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1,s2]
4179    
4180                   or      or
4181    
4182      out[s0,s1]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1]      out[s0,s1]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[r0,r1]
4183    
4184      In the first case the the second dimension of arg0 and the length of arg1 must match and        In the first case the the second dimension of arg0 and the last dimension of arg1 must match and  
4185      in the second case the two last dimensions of arg0 must match the shape of arg1.      in the second case the two last dimensions of arg0 must match the two first dimensions of arg1.
4186    
4187      @param arg0: first argument of rank 2 or 4      @param arg0: first argument of rank 2 or 4
4188      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
# Line 3436  def tensormult(arg0,arg1): Line 4194  def tensormult(arg0,arg1):
4194      sh0=pokeShape(arg0)      sh0=pokeShape(arg0)
4195      sh1=pokeShape(arg1)      sh1=pokeShape(arg1)
4196      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):      if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4197         return generalTensorProduct(arg0,arg1,offset=1)         return generalTensorProduct(arg0,arg1,axis_offset=1)
4198      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):
4199         return generalTensorProduct(arg0,arg1,offset=2)         return generalTensorProduct(arg0,arg1,axis_offset=2)
4200      else:      else:
4201          raise ValueError,"tensormult: first argument must have rank 2 or 4"          raise ValueError,"tensor_mult: first argument must have rank 2 or 4"
4202    
4203  def generalTensorProduct(arg0,arg1,offset=0):  def generalTensorProduct(arg0,arg1,axis_offset=0):
4204      """      """
4205      generalized tensor product      generalized tensor product
4206    
4207      out[s,t]=S{Sigma}_r arg0[s,r]*arg1[r,t]      out[s,t]=S{Sigma}_r arg0[s,r]*arg1[r,t]
4208    
4209      where s runs through arg0.Shape[:arg0.Rank-offset]      where
           r runs trough arg0.Shape[:offset]  
           t runs through arg1.Shape[offset:]  
4210    
4211      In the first case the the second dimension of arg0 and the length of arg1 must match and            - s runs through arg0.Shape[:arg0.Rank-axis_offset]
4212      in the second case the two last dimensions of arg0 must match the shape of arg1.          - r runs trough arg0.Shape[:axis_offset]
4213            - t runs through arg1.Shape[axis_offset:]
4214    
4215      @param arg0: first argument      @param arg0: first argument
4216      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4217      @param arg1: second argument of shape greater of 1 or 2 depending on rank of arg0      @param arg1: second argument
4218      @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}      @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4219      @return: the general tensor product of arg0 and arg1 at each data point.      @return: the general tensor product of arg0 and arg1 at each data point.
4220      @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
# Line 3467  def generalTensorProduct(arg0,arg1,offse Line 4224  def generalTensorProduct(arg0,arg1,offse
4224      # at this stage arg0 and arg0 are both numarray.NumArray or escript.Data or Symbols      # at this stage arg0 and arg0 are both numarray.NumArray or escript.Data or Symbols
4225      if isinstance(arg0,numarray.NumArray):      if isinstance(arg0,numarray.NumArray):
4226         if isinstance(arg1,Symbol):         if isinstance(arg1,Symbol):
4227             return GeneralTensorProduct_Symbol(arg0,arg1,offset)             return GeneralTensorProduct_Symbol(arg0,arg1,axis_offset)
4228         else:         else:
4229             if not arg0.shape[arg0.rank-offset:]==arg1.shape[:offset]:             if not arg0.shape[arg0.rank-axis_offset:]==arg1.shape[:axis_offset]:
4230                 raise ValueError,"generalTensorProduct: dimensions of last %s components in left argument don't match the first %s components in the right argument."%(offset,offset)                 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)
4231             arg0_c=arg0.copy()             arg0_c=arg0.copy()
4232             arg1_c=arg1.copy()             arg1_c=arg1.copy()
4233             sh0,sh1=arg0.shape,arg1.shape             sh0,sh1=arg0.shape,arg1.shape
4234             d0,d1,d01=1,1,1             d0,d1,d01=1,1,1
4235             for i in sh0[:arg0.rank-offset]: d0*=i             for i in sh0[:arg0.rank-axis_offset]: d0*=i
4236             for i in sh1[offset:]: d1*=i             for i in sh1[axis_offset:]: d1*=i
4237             for i in sh1[:offset]: d01*=i             for i in sh1[:axis_offset]: d01*=i
4238             arg0_c.resize((d0,d01))             arg0_c.resize((d0,d01))
4239             arg1_c.resize((d01,d1))             arg1_c.resize((d01,d1))
4240             out=numarray.zeros((d0,d1),numarray.Float)             out=numarray.zeros((d0,d1),numarray.Float64)
4241             for i0 in range(d0):             for i0 in range(d0):
4242                      for i1 in range(d1):                      for i1 in range(d1):
4243                           out[i0,i1]=numarray.sum(arg0_c[i0,:]*arg1_c[:,i1])                           out[i0,i1]=numarray.sum(arg0_c[i0,:]*arg1_c[:,i1])
4244             out.resize(sh0[:arg0.rank-offset]+sh1[offset:])             out.resize(sh0[:arg0.rank-axis_offset]+sh1[axis_offset:])
4245             return out             return out
4246      elif isinstance(arg0,escript.Data):      elif isinstance(arg0,escript.Data):
4247         if isinstance(arg1,Symbol):         if isinstance(arg1,Symbol):
4248             return GeneralTensorProduct_Symbol(arg0,arg1,offset)             return GeneralTensorProduct_Symbol(arg0,arg1,axis_offset)
4249         else:         else:
4250             return escript_generalTensorProduct(arg0,arg1,offset) # this calls has to be replaced by escript._generalTensorProduct(arg0,arg1,offset)             return escript_generalTensorProduct(arg0,arg1,axis_offset) # this calls has to be replaced by escript._generalTensorProduct(arg0,arg1,axis_offset)
4251      else:            else:      
4252         return GeneralTensorProduct_Symbol(arg0,arg1,offset)         return GeneralTensorProduct_Symbol(arg0,arg1,axis_offset)
4253                                    
4254  class GeneralTensorProduct_Symbol(DependendSymbol):  class GeneralTensorProduct_Symbol(DependendSymbol):
4255     """     """
4256     Symbol representing the quotient of two arguments.     Symbol representing the general tensor product of two arguments
4257     """     """
4258     def __init__(self,arg0,arg1,offset=0):     def __init__(self,arg0,arg1,axis_offset=0):
4259         """         """
4260         initialization of L{Symbol} representing the quotient of two arguments         initialization of L{Symbol} representing the general tensor product of two arguments.
4261    
4262         @param arg0: numerator         @param arg0: first argument
4263         @type arg0: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type arg0: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4264         @param arg1: denominator         @param arg1: second argument
4265         @type arg1: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.         @type arg1: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4266         @raise ValueError: if both arguments do not have the same shape.         @raise ValueError: illegal dimension
4267         @note: if both arguments have a spatial dimension, they must equal.         @note: if both arguments have a spatial dimension, they must equal.
4268         """         """
4269         sh_arg0=pokeShape(arg0)         sh_arg0=pokeShape(arg0)
4270         sh_arg1=pokeShape(arg1)         sh_arg1=pokeShape(arg1)
4271         sh0=sh_arg0[:len(sh_arg0)-offset]         sh0=sh_arg0[:len(sh_arg0)-axis_offset]
4272         sh01=sh_arg0[len(sh_arg0)-offset:]         sh01=sh_arg0[len(sh_arg0)-axis_offset:]
4273         sh10=sh_arg1[:offset]         sh10=sh_arg1[:axis_offset]
4274         sh1=sh_arg1[offset:]         sh1=sh_arg1[axis_offset:]
4275         if not sh01==sh10:         if not sh01==sh10:
4276             raise ValueError,"dimensions of last %s components in left argument don't match the first %s components in the right argument."%(offset,offset)             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)
4277         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0+sh1,args=[arg0,arg1,offset])         DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0+sh1,args=[arg0,arg1,axis_offset])
4278    
4279     def getMyCode(self,argstrs,format="escript"):     def getMyCode(self,argstrs,format="escript"):
4280        """        """
# Line 3529  class GeneralTensorProduct_Symbol(Depend Line 4286  class GeneralTensorProduct_Symbol(Depend
4286        @type format: C{str}        @type format: C{str}
4287        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.        @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4288        @rtype: C{str}        @rtype: C{str}
4289        @raise: NotImplementedError: if the requested format is not available        @raise NotImplementedError: if the requested format is not available
4290        """        """
4291        if format=="escript" or format=="str" or format=="text":        if format=="escript" or format=="str" or format=="text":
4292           return "generalTensorProduct(%s,%s,offset=%s)"%(argstrs[0],argstrs[1],argstrs[2])           return "generalTensorProduct(%s,%s,axis_offset=%s)"%(argstrs[0],argstrs[1],argstrs[2])
4293        else:        else:
4294           raise NotImplementedError,"%s does not provide program code for format %s."%(str(self),format)           raise NotImplementedError,"%s does not provide program code for format %s."%(str(self),format)
4295    
# Line 3557  class GeneralTensorProduct_Symbol(Depend Line 4314  class GeneralTensorProduct_Symbol(Depend
4314           args=self.getSubstitutedArguments(argvals)           args=self.getSubstitutedArguments(argvals)
4315           return generalTensorProduct(args[0],args[1],args[2])           return generalTensorProduct(args[0],args[1],args[2])
4316    
4317  def escript_generalTensorProduct(arg0,arg1,offset): # this should be escript._generalTensorProduct  def escript_generalTensorProduct(arg0,arg1,axis_offset,transpose=0):
4318      "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!!!"
4319      # calculate the return shape:      return C_GeneralTensorProduct(arg0, arg1, axis_offset, transpose)
4320      shape0=arg0.getShape()[:arg0.getRank()-offset]  
4321      shape01=arg0.getShape()[arg0.getRank()-offset:]  def transposed_matrix_mult(arg0,arg1):
4322      shape10=arg1.getShape()[:offset]      """
4323      shape1=arg1.getShape()[offset:]      transposed(matrix)-matrix or transposed(matrix)-vector product of the two argument:
4324      if not shape01==shape10:  
4325          raise ValueError,"dimensions of last %s components in left argument don't match the first %s components in the right argument."%(offset,offset)      out[s0]=S{Sigma}_{r0} arg0[r0,s0]*arg1[r0]
4326    
4327      # create return value:      or
4328      out=escript.Data(0.,tuple(shape0+shape1),arg0.getFunctionSpace())  
4329      #      out[s0,s1]=S{Sigma}_{r0} arg0[r0,s0]*arg1[r0,s1]
4330      s0=[[]]  
4331      for k in shape0:      The function call transposed_matrix_mult(arg0,arg1) is equivalent to matrix_mult(transpose(arg0),arg1).
4332            s=[]  
4333            for j in s0:      The first dimension of arg0 and arg1 must match.
4334                  for i in range(k): s.append(j+[slice(i,i)])  
4335            s0=s      @param arg0: first argument of rank 2
4336      s1=[[]]      @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4337      for k in shape1:      @param arg1: second argument of at least rank 1
4338            s=[]      @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4339            for j in s1:      @return: the product of the transposed of arg0 and arg1 at each data point
4340                  for i in range(k): s.append(j+[slice(i,i)])      @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4341            s1=s      @raise ValueError: if the shapes of the arguments are not appropriate
4342      s01=[[]]      """
4343      for k in shape01:      sh0=pokeShape(arg0)
4344            s=[]      sh1=pokeShape(arg1)
4345            for j in s01:      if not len(sh0)==2 :
4346                  for i in range(k): s.append(j+[slice(i,i)])          raise ValueError,"first argument must have rank 2"
4347            s01=s      if not len(sh1)==2 and not len(sh1)==1:
4348            raise ValueError,"second argument must have rank 1 or 2"
4349      for i0 in s0:      return generalTransposedTensorProduct(arg0,arg1,axis_offset=1)
4350         for i1 in s1:  
4351           s=escript.Scalar(0.,arg0.getFunctionSpace())  def transposed_tensor_mult(arg0,arg1):
4352           for i01 in s01:      """
4353              s+=arg0.__getitem__(tuple(i0+i01))*arg1.__getitem__(tuple(i01+i1))      the tensor product of the transposed of the first and the second argument
4354           out.__setitem__(tuple(i0+i1),s)      
4355      return out      for arg0 of rank 2 this is
4356    
4357        out[s0]=S{Sigma}_{r0} arg0[r0,s0]*arg1[r0]  
4358    
4359        or
4360    
4361        out[s0,s1]=S{Sigma}_{r0} arg0[r0,s0]*arg1[r0,s1]
4362    
4363      
4364        and for arg0 of rank 4 this is
4365    
4366        out[s0,s1,s2,s3]=S{Sigma}_{r0,r1} arg0[r0,r1,s0,s1]*arg1[r0,r1,s2,s3]
4367                  
4368        or
4369    
4370        out[s0,s1,s2]=S{Sigma}_{r0,r1} arg0[r0,r1,s0,s1]*arg1[r0,r1,s2]
4371    
4372        or
4373    
4374        out[s0,s1]=S{Sigma}_{r0,r1} arg0[r0,r1,s0,s1]*arg1[r0,r1]
4375    
4376        In the first case the the first dimension of arg0 and the first dimension of arg1 must match and  
4377        in the second case the two first dimensions of arg0 must match the two first dimension of arg1.
4378    
4379        The function call transposed_tensor_mult(arg0,arg1) is equivalent to tensor_mult(transpose(arg0),arg1).
4380    
4381        @param arg0: first argument of rank 2 or 4
4382        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4383        @param arg1: second argument of shape greater of 1 or 2 depending on rank of arg0
4384        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4385        @return: the tensor product of tarnsposed of arg0 and arg1 at each data point
4386        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4387        """
4388        sh0=pokeShape(arg0)
4389        sh1=pokeShape(arg1)
4390        if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4391           return generalTransposedTensorProduct(arg0,arg1,axis_offset=1)
4392        elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):
4393           return generalTransposedTensorProduct(arg0,arg1,axis_offset=2)
4394        else:
4395            raise ValueError,"first argument must have rank 2 or 4"
4396    
4397    def generalTransposedTensorProduct(arg0,arg1,axis_offset=0):
4398        """
4399        generalized tensor product of transposed of arg0 and arg1:
4400    
4401        out[s,t]=S{Sigma}_r arg0[r,s]*arg1[r,t]
4402    
4403        where
4404    
4405            - s runs through arg0.Shape[axis_offset:]
4406            - r runs trough arg0.Shape[:axis_offset]
4407            - t runs through arg1.Shape[axis_offset:]
4408    
4409        The function call generalTransposedTensorProduct(arg0,arg1,axis_offset) is equivalent
4410        to generalTensorProduct(transpose(arg0,arg0.rank-axis_offset),arg1,axis_offset).
4411    
4412        @param arg0: first argument
4413        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4414        @param arg1: second argument
4415        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4416        @return: the general tensor product of transposed(arg0) and arg1 at each data point.
4417        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4418        """
4419        if isinstance(arg0,float) and isinstance(arg1,float): return arg1*arg0
4420        arg0,arg1=matchType(arg0,arg1)
4421        # at this stage arg0 and arg0 are both numarray.NumArray or escript.Data or Symbols
4422        if isinstance(arg0,numarray.NumArray):
4423           if isinstance(arg1,Symbol):
4424               return GeneralTransposedTensorProduct_Symbol(arg0,arg1,axis_offset)
4425           else:
4426               if not arg0.shape[:axis_offset]==arg1.shape[:axis_offset]:
4427                   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)
4428               arg0_c=arg0.copy()
4429               arg1_c=arg1.copy()
4430               sh0,sh1=arg0.shape,arg1.shape
4431               d0,d1,d01=1,1,1
4432               for i in sh0[axis_offset:]: d0*=i
4433               for i in sh1[axis_offset:]: d1*=i
4434               for i in sh0[:axis_offset]: d01*=i
4435               arg0_c.resize((d01,d0))
4436               arg1_c.resize((d01,d1))
4437               out=numarray.zeros((d0,d1),numarray.Float64)
4438               for i0 in range(d0):
4439                        for i1 in range(d1):
4440                             out[i0,i1]=numarray.sum(arg0_c[:,i0]*arg1_c[:,i1])
4441               out.resize(sh0[axis_offset:]+sh1[axis_offset:])
4442               return out
4443        elif isinstance(arg0,escript.Data):
4444           if isinstance(arg1,Symbol):
4445               return GeneralTransposedTensorProduct_Symbol(arg0,arg1,axis_offset)
4446           else:
4447               return escript_generalTransposedTensorProduct(arg0,arg1,axis_offset) # this calls has to be replaced by escript._generalTensorProduct(arg0,arg1,axis_offset)
4448        else:      
4449           return GeneralTransposedTensorProduct_Symbol(arg0,arg1,axis_offset)
4450                    
4451    class GeneralTransposedTensorProduct_Symbol(DependendSymbol):
4452       """
4453       Symbol representing the general tensor product of the transposed of arg0 and arg1
4454       """
4455       def __init__(self,arg0,arg1,axis_offset=0):
4456           """
4457           initialization of L{Symbol} representing tensor product of the transposed of arg0 and arg1
4458    
4459           @param arg0: first argument
4460           @type arg0: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4461           @param arg1: second argument
4462           @type arg1: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4463           @raise ValueError: inconsistent dimensions of arguments.
4464           @note: if both arguments have a spatial dimension, they must equal.
4465           """
4466           sh_arg0=pokeShape(arg0)
4467           sh_arg1=pokeShape(arg1)
4468           sh01=sh_arg0[:axis_offset]
4469           sh10=sh_arg1[:axis_offset]
4470           sh0=sh_arg0[axis_offset:]
4471           sh1=sh_arg1[axis_offset:]
4472           if not sh01==sh10:
4473               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)
4474           DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0+sh1,args=[arg0,arg1,axis_offset])
4475    
4476       def getMyCode(self,argstrs,format="escript"):
4477          """
4478          returns a program code that can be used to evaluate the symbol.
4479    
4480          @param argstrs: gives for each argument a string representing the argument for the evaluation.
4481          @type argstrs: C{list} of length 2 of C{str}.
4482          @param format: specifies the format to be used. At the moment only "escript", "str" and "text" are supported.
4483          @type format: C{str}
4484          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4485          @rtype: C{str}
4486          @raise NotImplementedError: if the requested format is not available
4487          """
4488          if format=="escript" or format=="str" or format=="text":
4489             return "generalTransposedTensorProduct(%s,%s,axis_offset=%s)"%(argstrs[0],argstrs[1],argstrs[2])
4490          else:
4491             raise NotImplementedError,"%s does not provide program code for format %s."%(str(self),format)
4492    
4493       def substitute(self,argvals):
4494          """
4495          assigns new values to symbols in the definition of the symbol.
4496          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
4497    
4498          @param argvals: new values assigned to symbols
4499          @type argvals: C{dict} with keywords of type L{Symbol}.
4500          @return: result of the substitution process. Operations are executed as much as possible.
4501          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
4502          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
4503          """
4504          if argvals.has_key(self):
4505             arg=argvals[self]
4506             if self.isAppropriateValue(arg):
4507                return arg
4508             else:
4509                raise TypeError,"%s: new value is not appropriate."%str(self)
4510          else:
4511             args=self.getSubstitutedArguments(argvals)
4512             return generalTransposedTensorProduct(args[0],args[1],args[2])
4513    
4514    def escript_generalTransposedTensorProduct(arg0,arg1,axis_offset): # this should be escript._generalTransposedTensorProduct
4515        "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"
4516        return C_GeneralTensorProduct(arg0, arg1, axis_offset, 1)
4517    
4518    def matrix_transposed_mult(arg0,arg1):
4519        """
4520        matrix-transposed(matrix) product of the two argument:
4521    
4522        out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[s1,r0]
4523    
4524        The function call matrix_transposed_mult(arg0,arg1) is equivalent to matrix_mult(arg0,transpose(arg1)).
4525    
4526        The last dimensions of arg0 and arg1 must match.
4527    
4528        @param arg0: first argument of rank 2
4529        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4530        @param arg1: second argument of rank 2
4531        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4532        @return: the product of arg0 and the transposed of arg1 at each data point
4533        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4534        @raise ValueError: if the shapes of the arguments are not appropriate
4535        """
4536        sh0=pokeShape(arg0)
4537        sh1=pokeShape(arg1)
4538        if not len(sh0)==2 :
4539            raise ValueError,"first argument must have rank 2"
4540        if not len(sh1)==2 and not len(sh1)==1:
4541            raise ValueError,"second argument must have rank 1 or 2"
4542        return generalTensorTransposedProduct(arg0,arg1,axis_offset=1)
4543    
4544    def tensor_transposed_mult(arg0,arg1):
4545        """
4546        the tensor product of the first and the transpose of the second argument
4547        
4548        for arg0 of rank 2 this is
4549    
4550        out[s0,s1]=S{Sigma}_{r0} arg0[s0,r0]*arg1[s1,r0]
4551    
4552        and for arg0 of rank 4 this is
4553    
4554        out[s0,s1,s2,s3]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[s2,s3,r0,r1]
4555                  
4556        or
4557    
4558        out[s0,s1,s2]=S{Sigma}_{r0,r1} arg0[s0,s1,r0,r1]*arg1[s2,r0,r1]
4559    
4560        In the first case the the second dimension of arg0 and arg1 must match and  
4561        in the second case the two last dimensions of arg0 must match the two last dimension of arg1.
4562    
4563        The function call tensor_transpose_mult(arg0,arg1) is equivalent to tensor_mult(arg0,transpose(arg1)).
4564    
4565        @param arg0: first argument of rank 2 or 4
4566        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4567        @param arg1: second argument of shape greater of 1 or 2 depending on rank of arg0
4568        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
4569        @return: the tensor product of tarnsposed of arg0 and arg1 at each data point
4570        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4571        """
4572        sh0=pokeShape(arg0)
4573        sh1=pokeShape(arg1)
4574        if len(sh0)==2 and ( len(sh1)==2 or len(sh1)==1 ):
4575           return generalTensorTransposedProduct(arg0,arg1,axis_offset=1)
4576        elif len(sh0)==4 and (len(sh1)==2 or len(sh1)==3 or len(sh1)==4):
4577           return generalTensorTransposedProduct(arg0,arg1,axis_offset=2)
4578        else:
4579            raise ValueError,"first argument must have rank 2 or 4"
4580    
4581    def generalTensorTransposedProduct(arg0,arg1,axis_offset=0):
4582        """
4583        generalized tensor product of transposed of arg0 and arg1:
4584    
4585        out[s,t]=S{Sigma}_r arg0[s,r]*arg1[t,r]
4586    
4587        where
4588    
4589            - s runs through arg0.Shape[:arg0.Rank-axis_offset]
4590            - r runs trough arg0.Shape[arg1.Rank-axis_offset:]
4591            - t runs through arg1.Shape[arg1.Rank-axis_offset:]
4592    
4593        The function call generalTensorTransposedProduct(arg0,arg1,axis_offset) is equivalent
4594        to generalTensorProduct(arg0,transpose(arg1,arg1.Rank-axis_offset),axis_offset).
4595    
4596        @param arg0: first argument
4597        @type arg0: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4598        @param arg1: second argument
4599        @type arg1: L{numarray.NumArray}, L{escript.Data}, L{Symbol}, C{float}, C{int}
4600        @return: the general tensor product of transposed(arg0) and arg1 at each data point.
4601        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
4602        """
4603        if isinstance(arg0,float) and isinstance(arg1,float): return arg1*arg0
4604        arg0,arg1=matchType(arg0,arg1)
4605        # at this stage arg0 and arg0 are both numarray.NumArray or escript.Data or Symbols
4606        if isinstance(arg0,numarray.NumArray):
4607           if isinstance(arg1,Symbol):
4608               return GeneralTensorTransposedProduct_Symbol(arg0,arg1,axis_offset)
4609           else:
4610               if not arg0.shape[arg0.rank-axis_offset:]==arg1.shape[arg1.rank-axis_offset:]:
4611                   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)
4612               arg0_c=arg0.copy()
4613               arg1_c=arg1.copy()
4614               sh0,sh1=arg0.shape,arg1.shape
4615               d0,d1,d01=1,1,1
4616               for i in sh0[:arg0.rank-axis_offset]: d0*=i
4617               for i in sh1[:arg1.rank-axis_offset]: d1*=i
4618               for i in sh1[arg1.rank-axis_offset:]: d01*=i
4619               arg0_c.resize((d0,d01))
4620               arg1_c.resize((d1,d01))
4621               out=numarray.zeros((d0,d1),numarray.Float64)
4622               for i0 in range(d0):
4623                        for i1 in range(d1):
4624                             out[i0,i1]=numarray.sum(arg0_c[i0,:]*arg1_c[i1,:])
4625               out.resize(sh0[:arg0.rank-axis_offset]+sh1[:arg1.rank-axis_offset])
4626               return out
4627        elif isinstance(arg0,escript.Data):
4628           if isinstance(arg1,Symbol):
4629               return GeneralTensorTransposedProduct_Symbol(arg0,arg1,axis_offset)
4630           else:
4631               return escript_generalTensorTransposedProduct(arg0,arg1,axis_offset) # this calls has to be replaced by escript._generalTensorProduct(arg0,arg1,axis_offset)
4632        else:      
4633           return GeneralTensorTransposedProduct_Symbol(arg0,arg1,axis_offset)
4634                    
4635    class GeneralTensorTransposedProduct_Symbol(DependendSymbol):
4636       """
4637       Symbol representing the general tensor product of arg0 and the transpose of arg1
4638       """
4639       def __init__(self,arg0,arg1,axis_offset=0):
4640           """
4641           initialization of L{Symbol} representing the general tensor product of arg0 and the transpose of arg1
4642    
4643           @param arg0: first argument
4644           @type arg0: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4645           @param arg1: second argument
4646           @type arg1: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
4647           @raise ValueError: inconsistent dimensions of arguments.
4648           @note: if both arguments have a spatial dimension, they must equal.
4649           """
4650           sh_arg0=pokeShape(arg0)
4651           sh_arg1=pokeShape(arg1)
4652           sh0=sh_arg0[:len(sh_arg0)-axis_offset]
4653           sh01=sh_arg0[len(sh_arg0)-axis_offset:]
4654           sh10=sh_arg1[len(sh_arg1)-axis_offset:]
4655           sh1=sh_arg1[:len(sh_arg1)-axis_offset]
4656           if not sh01==sh10:
4657               raise ValueError,"dimensions of last %s components in left argument don't match the last %s components in the right argument."%(axis_offset,axis_offset)
4658           DependendSymbol.__init__(self,dim=commonDim(arg0,arg1),shape=sh0+sh1,args=[arg0,arg1,axis_offset])
4659    
4660       def getMyCode(self,argstrs,format="escript"):
4661          """
4662          returns a program code that can be used to evaluate the symbol.
4663    
4664          @param argstrs: gives for each argument a string representing the argument for the evaluation.
4665          @type argstrs: C{list} of length 2 of C{str}.
4666          @param format: specifies the format to be used. At the moment only "escript", "str" and "text" are supported.
4667          @type format: C{str}
4668          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4669          @rtype: C{str}
4670          @raise NotImplementedError: if the requested format is not available
4671          """
4672          if format=="escript" or format=="str" or format=="text":
4673             return "generalTensorTransposedProduct(%s,%s,axis_offset=%s)"%(argstrs[0],argstrs[1],argstrs[2])
4674          else:
4675             raise NotImplementedError,"%s does not provide program code for format %s."%(str(self),format)
4676    
4677       def substitute(self,argvals):
4678          """
4679          assigns new values to symbols in the definition of the symbol.
4680          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
4681    
4682          @param argvals: new values assigned to symbols
4683          @type argvals: C{dict} with keywords of type L{Symbol}.
4684          @return: result of the substitution process. Operations are executed as much as possible.
4685          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
4686          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
4687          """
4688          if argvals.has_key(self):
4689             arg=argvals[self]
4690             if self.isAppropriateValue(arg):
4691                return arg
4692             else:
4693                raise TypeError,"%s: new value is not appropriate."%str(self)
4694          else:
4695             args=self.getSubstitutedArguments(argvals)
4696             return generalTensorTransposedProduct(args[0],args[1],args[2])
4697    
4698    def escript_generalTensorTransposedProduct(arg0,arg1,axis_offset): # this should be escript._generalTensorTransposedProduct
4699        "arg0 and arg1 are both Data objects but not neccesrily on the same function space. they could be identical!!!"
4700        return C_GeneralTensorProduct(arg0, arg1, axis_offset, 2)
4701    
4702  #=========================================================  #=========================================================
4703  #   some little helpers  #  functions dealing with spatial dependency
4704  #=========================================================  #=========================================================
4705  def grad(arg,where=None):  def grad(arg,where=None):
4706      """      """
4707      Returns the spatial gradient of arg at where.      Returns the spatial gradient of arg at where.
4708    
4709        If C{g} is the returned object, then
4710    
4711          - if C{arg} is rank 0 C{g[s]} is the derivative of C{arg} with respect to the C{s}-th spatial dimension.
4712          - if C{arg} is rank 1 C{g[i,s]} is the derivative of C{arg[i]} with respect to the C{s}-th spatial dimension.
4713          - if C{arg} is rank 2 C{g[i,j,s]} is the derivative of C{arg[i,j]} with respect to the C{s}-th spatial dimension.
4714          - if C{arg} is rank 3 C{g[i,j,k,s]} is the derivative of C{arg[i,j,k]} with respect to the C{s}-th spatial dimension.
4715    
4716      @param arg:   Data object representing the function which gradient      @param arg: function which gradient to be calculated. Its rank has to be less than 3.
4717                    to be calculated.      @type arg: L{escript.Data} or L{Symbol}
4718      @param where: FunctionSpace in which the gradient will be calculated.      @param where: FunctionSpace in which the gradient will be calculated.
4719                    If not present or C{None} an appropriate default is used.                    If not present or C{None} an appropriate default is used.
4720        @type where: C{None} or L{escript.FunctionSpace}
4721        @return: gradient of arg.
4722        @rtype: L{escript.Data} or L{Symbol}
4723      """      """
4724      if isinstance(arg,Symbol):      if isinstance(arg,Symbol):
4725         return Grad_Symbol(arg,where)         return Grad_Symbol(arg,where)
# Line 3617  def grad(arg,where=None): Line 4729  def grad(arg,where=None):
4729         else:         else:
4730            return arg._grad(where)            return arg._grad(where)
4731      else:      else:
4732        raise TypeError,"grad: Unknown argument type."         raise TypeError,"grad: Unknown argument type."
4733    
4734    class Grad_Symbol(DependendSymbol):
4735       """
4736       L{Symbol} representing the result of the gradient operator
4737       """
4738       def __init__(self,arg,where=None):
4739          """
4740          initialization of gradient L{Symbol} with argument arg
4741          @param arg: argument of function
4742          @type arg: L{Symbol}.
4743          @param where: FunctionSpace in which the gradient will be calculated.
4744                      If not present or C{None} an appropriate default is used.
4745          @type where: C{None} or L{escript.FunctionSpace}
4746          """
4747          d=arg.getDim()
4748          if d==None:
4749             raise ValueError,"argument must have a spatial dimension"
4750          super(Grad_Symbol,self).__init__(args=[arg,where],shape=arg.getShape()+(d,),dim=d)
4751    
4752       def getMyCode(self,argstrs,format="escript"):
4753          """
4754          returns a program code that can be used to evaluate the symbol.
4755    
4756          @param argstrs: gives for each argument a string representing the argument for the evaluation.
4757          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
4758          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
4759          @type format: C{str}
4760          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4761          @rtype: C{str}
4762          @raise NotImplementedError: if the requested format is not available
4763          """
4764          if format=="escript" or format=="str"  or format=="text":
4765             return "grad(%s,where=%s)"%(argstrs[0],argstrs[1])
4766          else:
4767             raise NotImplementedError,"Trace_Symbol does not provide program code for format %s."%format
4768    
4769       def substitute(self,argvals):
4770          """
4771          assigns new values to symbols in the definition of the symbol.
4772          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
4773    
4774          @param argvals: new values assigned to symbols
4775          @type argvals: C{dict} with keywords of type L{Symbol}.
4776          @return: result of the substitution process. Operations are executed as much as possible.
4777          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
4778          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
4779          """
4780          if argvals.has_key(self):
4781             arg=argvals[self]
4782             if self.isAppropriateValue(arg):
4783                return arg
4784             else:
4785                raise TypeError,"%s: new value is not appropriate."%str(self)
4786          else:
4787             arg=self.getSubstitutedArguments(argvals)
4788             return grad(arg[0],where=arg[1])
4789    
4790       def diff(self,arg):
4791          """
4792          differential of this object
4793    
4794          @param arg: the derivative is calculated with respect to arg
4795          @type arg: L{escript.Symbol}
4796          @return: derivative with respect to C{arg}
4797          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
4798          """
4799          if arg==self:
4800             return identity(self.getShape())
4801          else:
4802             return grad(self.getDifferentiatedArguments(arg)[0],where=self.getArgument()[1])
4803    
4804  def integrate(arg,where=None):  def integrate(arg,where=None):
4805      """      """
4806      Return the integral if the function represented by Data object arg over      Return the integral of the function C{arg} over its domain. If C{where} is present C{arg} is interpolated to C{where}
4807      its domain.      before integration.
4808    
4809      @param arg:   Data object representing the function which is integrated.      @param arg:   the function which is integrated.
4810        @type arg: L{escript.Data} or L{Symbol}
4811      @param where: FunctionSpace in which the integral is calculated.      @param where: FunctionSpace in which the integral is calculated.
4812                    If not present or C{None} an appropriate default is used.                    If not present or C{None} an appropriate default is used.
4813        @type where: C{None} or L{escript.FunctionSpace}
4814        @return: integral of arg.
4815        @rtype: C{float}, C{numarray.NumArray} or L{Symbol}
4816      """      """
4817      if isinstance(arg,numarray.NumArray):      if isinstance(arg,Symbol):
         if checkForZero(arg):  
            return arg  
         else:  
            raise TypeError,"integrate: cannot intergrate argument"  
     elif isinstance(arg,float):  
         if checkForZero(arg):  
            return arg  
         else:  
            raise TypeError,"integrate: cannot intergrate argument"  
     elif isinstance(arg,int):  
         if checkForZero(arg):  
            return float(arg)  
         else:  
            raise TypeError,"integrate: cannot intergrate argument"  
     elif isinstance(arg,Symbol):  
4818         return Integrate_Symbol(arg,where)         return Integrate_Symbol(arg,where)
4819      elif isinstance(arg,escript.Data):      elif isinstance(arg,escript.Data):
4820         if not where==None: arg=escript.Data(arg,where)         if not where==None: arg=escript.Data(arg,where)
# Line 3654  def integrate(arg,where=None): Line 4825  def integrate(arg,where=None):
4825      else:      else:
4826        raise TypeError,"integrate: Unknown argument type."        raise TypeError,"integrate: Unknown argument type."
4827    
4828    class Integrate_Symbol(DependendSymbol):
4829       """
4830       L{Symbol} representing the result of the spatial integration operator
4831       """
4832       def __init__(self,arg,where=None):
4833          """
4834          initialization of integration L{Symbol} with argument arg
4835          @param arg: argument of the integration
4836          @type arg: L{Symbol}.
4837          @param where: FunctionSpace in which the integration will be calculated.
4838                      If not present or C{None} an appropriate default is used.
4839          @type where: C{None} or L{escript.FunctionSpace}
4840          """
4841          super(Integrate_Symbol,self).__init__(args=[arg,where],shape=arg.getShape(),dim=arg.getDim())
4842    
4843       def getMyCode(self,argstrs,format="escript"):
4844          """
4845          returns a program code that can be used to evaluate the symbol.
4846    
4847          @param argstrs: gives for each argument a string representing the argument for the evaluation.
4848          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
4849          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
4850          @type format: C{str}
4851          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4852          @rtype: C{str}
4853          @raise NotImplementedError: if the requested format is not available
4854          """
4855          if format=="escript" or format=="str"  or format=="text":
4856             return "integrate(%s,where=%s)"%(argstrs[0],argstrs[1])
4857          else:
4858             raise NotImplementedError,"Trace_Symbol does not provide program code for format %s."%format
4859    
4860       def substitute(self,argvals):
4861          """
4862          assigns new values to symbols in the definition of the symbol.
4863          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
4864    
4865          @param argvals: new values assigned to symbols
4866          @type argvals: C{dict} with keywords of type L{Symbol}.
4867          @return: result of the substitution process. Operations are executed as much as possible.
4868          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
4869          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
4870          """
4871          if argvals.has_key(self):
4872             arg=argvals[self]
4873             if self.isAppropriateValue(arg):
4874                return arg
4875             else:
4876                raise TypeError,"%s: new value is not appropriate."%str(self)
4877          else:
4878             arg=self.getSubstitutedArguments(argvals)
4879             return integrate(arg[0],where=arg[1])
4880    
4881       def diff(self,arg):
4882          """
4883          differential of this object
4884    
4885          @param arg: the derivative is calculated with respect to arg
4886          @type arg: L{escript.Symbol}
4887          @return: derivative with respect to C{arg}
4888          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
4889          """
4890          if arg==self:
4891             return identity(self.getShape())
4892          else:
4893             return integrate(self.getDifferentiatedArguments(arg)[0],where=self.getArgument()[1])
4894    
4895    
4896  def interpolate(arg,where):  def interpolate(arg,where):
4897      """      """
4898      Interpolates the function into the FunctionSpace where.      interpolates the function into the FunctionSpace where.
4899    
4900      @param arg:    interpolant      @param arg: interpolant
4901      @param where:  FunctionSpace to interpolate to      @type arg: L{escript.Data} or L{Symbol}
4902        @param where: FunctionSpace to be interpolated to
4903        @type where: L{escript.FunctionSpace}
4904        @return: interpolated argument
4905        @rtype: C{escript.Data} or L{Symbol}
4906      """      """
4907      if testForZero(arg):      if isinstance(arg,Symbol):
4908        return 0         return Interpolate_Symbol(arg,where)
     elif isinstance(arg,Symbol):  
        return Interpolated_Symbol(arg,where)  
4909      else:      else:
4910         return escript.Data(arg,where)         return escript.Data(arg,where)
4911    
4912  def div(arg,where=None):  class Interpolate_Symbol(DependendSymbol):
4913      """     """
4914      Returns the divergence of arg at where.     L{Symbol} representing the result of the interpolation operator
4915       """
4916       def __init__(self,arg,where):
4917          """
4918          initialization of interpolation L{Symbol} with argument arg
4919          @param arg: argument of the interpolation
4920          @type arg: L{Symbol}.
4921          @param where: FunctionSpace into which the argument is interpolated.
4922          @type where: L{escript.FunctionSpace}
4923          """
4924          super(Interpolate_Symbol,self).__init__(args=[arg,where],shape=arg.getShape(),dim=arg.getDim())
4925    
4926      @param arg:   Data object representing the function which gradient to     def getMyCode(self,argstrs,format="escript"):
4927                    be calculated.        """
4928      @param where: FunctionSpace in which the gradient will be calculated.        returns a program code that can be used to evaluate the symbol.
                   If not present or C{None} an appropriate default is used.  
     """  
     g=grad(arg,where)  
     return trace(g,axis0=g.getRank()-2,axis1=g.getRank()-1)  
4929    
4930  def jump(arg):        @param argstrs: gives for each argument a string representing the argument for the evaluation.
4931      """        @type argstrs: C{str} or a C{list} of length 1 of C{str}.
4932      Returns the jump of arg across a continuity.        @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
4933          @type format: C{str}
4934          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
4935          @rtype: C{str}
4936          @raise NotImplementedError: if the requested format is not available
4937          """
4938          if format=="escript" or format=="str"  or format=="text":
4939             return "interpolate(%s,where=%s)"%(argstrs[0],argstrs[1])
4940          else:
4941             raise NotImplementedError,"Trace_Symbol does not provide program code for format %s."%format
4942    
4943      @param arg:   Data object representing the function which gradient     def substitute(self,argvals):
4944                    to be calculated.        """
4945      """        assigns new values to symbols in the definition of the symbol.
4946      d=arg.getDomain()        The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
     return arg.interpolate(escript.FunctionOnContactOne())-arg.interpolate(escript.FunctionOnContactZero())  
4947    
4948  #=============================        @param argvals: new values assigned to symbols
4949  #        @type argvals: C{dict} with keywords of type L{Symbol}.
4950  # wrapper for various functions: if the argument has attribute the function name        @return: result of the substitution process. Operations are executed as much as possible.
4951  # as an argument it calls the corresponding methods. Otherwise the corresponding        @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
4952  # numarray function is called.        @raise TypeError: if a value for a L{Symbol} cannot be substituted.
4953          """
4954          if argvals.has_key(self):
4955             arg=argvals[self]
4956             if self.isAppropriateValue(arg):
4957                return arg
4958             else:
4959                raise TypeError,"%s: new value is not appropriate."%str(self)
4960          else:
4961             arg=self.getSubstitutedArguments(argvals)
4962             return interpolate(arg[0],where=arg[1])
4963    
4964       def diff(self,arg):
4965          """
4966          differential of this object
4967    
4968          @param arg: the derivative is calculated with respect to arg
4969          @type arg: L{escript.Symbol}
4970          @return: derivative with respect to C{arg}
4971          @rtype: L{Symbol} but other types such as L{escript.Data}, L{numarray.NumArray}  are possible.
4972          """
4973          if arg==self:
4974             return identity(self.getShape())
4975          else:
4976             return interpolate(self.getDifferentiatedArguments(arg)[0],where=self.getArgument()[1])
4977    
 # functions involving the underlying Domain:  
4978    
4979  def transpose(arg,axis=None):  def div(arg,where=None):
4980      """      """
4981      Returns the transpose of the Data object arg.      returns the divergence of arg at where.
4982    
4983      @param arg:      @param arg: function which divergence to be calculated. Its shape has to be (d,) where d is the spatial dimension.
4984        @type arg: L{escript.Data} or L{Symbol}
4985        @param where: FunctionSpace in which the divergence will be calculated.
4986                      If not present or C{None} an appropriate default is used.
4987        @type where: C{None} or L{escript.FunctionSpace}
4988        @return: divergence of arg.
4989        @rtype: L{escript.Data} or L{Symbol}
4990      """      """
     if axis==None:  
        r=0  
        if hasattr(arg,"getRank"): r=arg.getRank()  
        if hasattr(arg,"rank"): r=arg.rank  
        axis=r/2  
4991      if isinstance(arg,Symbol):      if isinstance(arg,Symbol):
4992         return Transpose_Symbol(arg,axis=r)          dim=arg.getDim()
4993      if isinstance(arg,escript.Data):      elif isinstance(arg,escript.Data):
4994         # hack for transpose          dim=arg.getDomain().getDim()
        r=arg.getRank()  
        if r!=2: raise ValueError,"Tranpose only avalaible for rank 2 objects"  
        s=arg.getShape()  
        out=escript.Data(0.,(s[1],s[0]),arg.getFunctionSpace())  
        for i in range(s[0]):  
           for j in range(s[1]):  
              out[j,i]=arg[i,j]  
        return out  
        # end hack for transpose  
        return arg.transpose(axis)  
4995      else:      else:
4996         return numarray.transpose(arg,axis=axis)          raise TypeError,"div: argument type not supported"
4997        if not arg.getShape()==(dim,):
4998          raise ValueError,"div: expected shape is (%s,)"%dim
4999        return trace(grad(arg,where))
5000    
5001  def trace(arg,axis0=0,axis1=1):  def jump(arg,domain=None):
5002      """      """
5003      Return      returns the jump of arg across the continuity of the domain
5004    
5005      @param arg:      @param arg: argument
5006        @type arg: L{escript.Data} or L{Symbol}
5007        @param domain: the domain where the discontinuity is located. If domain is not present or equal to C{None}
5008                       the domain of arg is used. If arg is a L{Symbol} the domain must be present.
5009        @type domain: C{None} or L{escript.Domain}
5010        @return: jump of arg
5011        @rtype: L{escript.Data} or L{Symbol}
5012      """      """
5013      if isinstance(arg,Symbol):      if domain==None: domain=arg.getDomain()
5014         s=list(arg.getShape())              return interpolate(arg,escript.FunctionOnContactOne(domain))-interpolate(arg,escript.FunctionOnContactZero(domain))
        s=tuple(s[0:axis0]+s[axis0+1:axis1]+s[axis1+1:])  
        return Trace_Symbol(arg,axis0=axis0,axis1=axis1)  
     elif isinstance(arg,escript.Data):  
        # hack for trace  
        s=arg.getShape()  
        if s[axis0]!=s[axis1]:  
            raise ValueError,"illegal axis in trace"  
        out=escript.Scalar(0.,arg.getFunctionSpace())  
        for i in range(s[axis0]):  
           out+=arg[i,i]  
        return out  
        # end hack for trace  
     else:  
        return numarray.trace(arg,axis0=axis0,axis1=axis1)  
5015    
5016    def L2(arg):
5017        """
5018        returns the L2 norm of arg at where
5019        
5020        @param arg: function which L2 to be calculated.
5021        @type arg: L{escript.Data} or L{Symbol}
5022        @return: L2 norm of arg.
5023        @rtype: L{float} or L{Symbol}
5024        @note: L2(arg) is equivalent to sqrt(integrate(inner(arg,arg)))
5025        """
5026        return sqrt(integrate(inner(arg,arg)))
5027    #=============================
5028    #
5029    
5030  def reorderComponents(arg,index):  def reorderComponents(arg,index):
5031      """      """
5032      resorts the component of arg according to index      resorts the component of arg according to index
5033    
5034      """      """
5035      pass      raise NotImplementedError
5036  #  #
5037  # $Log: util.py,v $  # $Log: util.py,v $
5038  # Revision 1.14.2.16  2005/10/19 06:09:57  gross  # Revision 1.14.2.16  2005/10/19 06:09:57  gross

Legend:
Removed from v.341  
changed lines
  Added in v.813

  ViewVC Help
Powered by ViewVC 1.1.26