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

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

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

trunk/esys2/escript/py_src/util.py revision 102 by jgs, Wed Dec 15 07:08:39 2004 UTC trunk/escript/py_src/util.py revision 637 by gross, Thu Mar 23 10:55:31 2006 UTC
# Line 1  Line 1 
1  # $Id$  # $Id$
2    
 ## @file util.py  
   
3  """  """
4  @brief Utility functions for escript  Utility functions for escript
5    
6    @var __author__: name of author
7    @var __copyright__: copyrights
8    @var __license__: licence agreement
9    @var __url__: url entry point on documentation
10    @var __version__: version
11    @var __date__: date of the version
12  """  """
13                                                                                                                                                                                                        
14    __author__="Lutz Gross, l.gross@uq.edu.au"
15    __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"
21    __version__="$Revision$"
22    __date__="$Date$"
23    
24    
25    import math
26  import numarray  import numarray
27  import escript  import escript
28  #  import os
29  #   escript constants (have to be consistent witj utilC.h  
30  #  #=========================================================
31  UNKNOWN=-1  #   some helpers:
32  EPSILON=1.e-15  #=========================================================
33  Pi=numarray.pi  def saveVTK(filename,domain=None,**data):
 # some solver options:  
 NO_REORDERING=0  
 MINIMUM_FILL_IN=1  
 NESTED_DISSECTION=2  
 # solver methods  
 DEFAULT_METHOD=0  
 DIRECT=1  
 CHOLEVSKY=2  
 PCG=3  
 CR=4  
 CGS=5  
 BICGSTAB=6  
 SSOR=7  
 ILU0=8  
 ILUT=9  
 JACOBI=10  
 GMRES=11  
 PRES20=12  
   
 METHOD_KEY="method"  
 SYMMETRY_KEY="symmetric"  
 TOLERANCE_KEY="tolerance"  
   
 # supported file formats:  
 VRML=1  
 PNG=2  
 JPEG=3  
 JPG=3  
 PS=4  
 OOGL=5  
 BMP=6  
 TIFF=7  
 OPENINVENTOR=8  
 RENDERMAN=9  
 PNM=10  
 #  
 # wrapper for various functions: if the argument has attribute the function name  
 # as an argument it calls the correspong methods. Otherwise the coresponsing numarray  
 # function is called.  
 #  
 # functions involving the underlying Domain:  
 #  
 def grad(arg,where=None):  
34      """      """
35      @brief returns the spatial gradient of the Data object arg      writes a L{Data} objects into a files using the the VTK XML file format.
36    
37      @param arg: Data object representing the function which gradient to be calculated.      Example:
38      @param where: FunctionSpace in which the gradient will be. If None Function(dom) where dom is the  
39                    domain of the Data object arg.         tmp=Scalar(..)
40           v=Vector(..)
41           saveVTK("solution.xml",temperature=tmp,velovity=v)
42    
43        tmp and v are written into "solution.xml" where tmp is named "temperature" and v is named "velovity"
44    
45        @param filename: file name of the output file
46        @type filename: C{str}
47        @param domain: domain of the L{Data} object. If not specified, the domain of the given L{Data} objects is used.
48        @type domain: L{escript.Domain}
49        @keyword <name>: writes the assigned value to the VTK file using <name> as identifier.
50        @type <name>: L{Data} object.
51        @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.
52      """      """
53      if where==None:      if domain==None:
54         return arg.grad()         for i in data.keys():
55              if not data[i].isEmpty(): domain=data[i].getFunctionSpace().getDomain()
56        if domain==None:
57            raise ValueError,"no domain detected."
58      else:      else:
59         return arg.grad(where)          domain.saveVTK(filename,data)
60    
61  def integrate(arg):  def saveDX(filename,domain=None,**data):
62      """      """
63      @brief return the integral if the function represented by Data object arg over its domain.      writes a L{Data} objects into a files using the the DX file format.
64    
65        Example:
66    
67      @param arg         tmp=Scalar(..)
68           v=Vector(..)
69           saveDX("solution.dx",temperature=tmp,velovity=v)
70    
71        tmp and v are written into "solution.dx" where tmp is named "temperature" and v is named "velovity".
72    
73        @param filename: file name of the output file
74        @type filename: C{str}
75        @param domain: domain of the L{Data} object. If not specified, the domain of the given L{Data} objects is used.
76        @type domain: L{escript.Domain}
77        @keyword <name>: writes the assigned value to the DX file using <name> as identifier. The identifier can be used to select the data set when data are imported into DX.
78        @type <name>: L{Data} object.
79        @note: The data objects have to be defined on the same domain. They may not be in the same L{FunctionSpace} but one cannot expect that all L{FunctionSpace} can be mixed. Typically, data on the boundary and data on the interior cannot be mixed.
80      """      """
81      return arg.integrate()      if domain==None:
82           for i in data.keys():
83              if not data[i].isEmpty(): domain=data[i].getFunctionSpace().getDomain()
84        if domain==None:
85            raise ValueError,"no domain detected."
86        else:
87            domain.saveDX(filename,data)
88    
89  def interpolate(arg,where):  def kronecker(d=3):
90       """
91       return the kronecker S{delta}-symbol
92    
93       @param d: dimension or an object that has the C{getDim} method defining the dimension
94       @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
95       @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise
96       @rtype d: L{numarray.NumArray} or L{escript.Data} of rank 2.
97       """
98       return identityTensor(d)
99    
100    def identity(shape=()):
101       """
102       return the shape x shape identity tensor
103    
104       @param shape: input shape for the identity tensor
105       @type shape: C{tuple} of C{int}
106       @return: array of shape shape x shape with M{u[i,k]=1} for M{i=k} and M{u[i,k]=0} otherwise for len(shape)=1 and, for len(shape)=2: M{u[i,j,k,l]=1} for M{i=k and j=l} and M{u[i,j,k,l]=0} otherwise.
107       @rtype: L{numarray.NumArray} of rank 1, rankk 2 or rank 4.
108       @raise ValueError: if len(shape)>2.
109       """
110       if len(shape)>0:
111          out=numarray.zeros(shape+shape,numarray.Float64)
112          if len(shape)==1:
113              for i0 in range(shape[0]):
114                 out[i0,i0]=1.
115          elif len(shape)==2:
116              for i0 in range(shape[0]):
117                 for i1 in range(shape[1]):
118                    out[i0,i1,i0,i1]=1.
119          else:
120              raise ValueError,"identity: length of shape is restricted to 2."
121       else:
122          out=1.
123       return out
124    
125    def identityTensor(d=3):
126       """
127       return the dxd identity matrix
128    
129       @param d: dimension or an object that has the C{getDim} method defining the dimension
130       @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
131       @return: the object u of rank 2 with M{u[i,j]=1} for M{i=j} and M{u[i,j]=0} otherwise
132       @rtype d: L{numarray.NumArray} or L{escript.Data} of rank 2
133       """
134       if isinstance(d,escript.FunctionSpace):
135           return escript.Data(identity((d.getDim(),)),d)
136       elif isinstance(d,escript.Domain):
137           return identity((d.getDim(),))
138       else:
139           return identity((d,))
140    
141    def identityTensor4(d=3):
142       """
143       return the dxdxdxd identity tensor
144    
145       @param d: dimension or an object that has the C{getDim} method defining the dimension
146       @type d: C{int} or any object with a C{getDim} method
147       @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
148       @rtype d: L{numarray.NumArray} or L{escript.Data} of rank 4.
149       """
150       if isinstance(d,escript.FunctionSpace):
151           return escript.Data(identity((d.getDim(),d.getDim())),d)
152       elif isinstance(d,escript.Domain):
153           return identity((d.getDim(),d.getDim()))
154       else:
155           return identity((d,d))
156    
157    def unitVector(i=0,d=3):
158       """
159       return a unit vector u of dimension d with nonzero index i:
160    
161       @param i: index
162       @type i: C{int}
163       @param d: dimension or an object that has the C{getDim} method defining the dimension
164       @type d: C{int}, L{escript.Domain} or L{escript.FunctionSpace}
165       @return: the object u of rank 1 with M{u[j]=1} for M{j=i} and M{u[i]=0} otherwise
166       @rtype d: L{numarray.NumArray} or L{escript.Data} of rank 1
167       """
168       return kronecker(d)[i]
169    
170    #=========================================================================
171    #   global reduction operations (these functions have no symbolic version)
172    #=========================================================================
173    def Lsup(arg):
174      """      """
175      @brief interpolates the function represented by Data object arg into the FunctionSpace where.      returns the Lsup-norm of argument arg. This is the maximum absolute value over all data points.
176        This function is equivalent to sup(abs(arg)).
177    
178      @param arg      @param arg: argument
179      @param where      @type arg: C{float}, C{int}, L{escript.Data}, L{numarray.NumArray}.
180        @return: maximum value of the absolute value of arg over all components and all data points
181        @rtype: C{float}
182        @raise TypeError: if type of arg cannot be processed
183      """      """
184      return arg.interpolate(where)      if isinstance(arg,numarray.NumArray):
185            return sup(abs(arg))
186        elif isinstance(arg,escript.Data):
187            return arg._Lsup()
188        elif isinstance(arg,float):
189            return abs(arg)
190        elif isinstance(arg,int):
191            return abs(float(arg))
192        else:
193          raise TypeError,"Lsup: Unknown argument type."
194    
195  # functions returning Data objects:  def sup(arg):
196        """
197        returns the maximum value over all data points.
198    
199  def transpose(arg,axis=None):      @param arg: argument
200        @type arg: C{float}, C{int}, L{escript.Data}, L{numarray.NumArray}.
201        @return: maximum value of arg over all components and all data points
202        @rtype: C{float}
203        @raise TypeError: if type of arg cannot be processed
204      """      """
205      @brief returns the transpose of the Data object arg.      if isinstance(arg,numarray.NumArray):
206            return arg.max()
207        elif isinstance(arg,escript.Data):
208            return arg._sup()
209        elif isinstance(arg,float):
210            return arg
211        elif isinstance(arg,int):
212            return float(arg)
213        else:
214          raise TypeError,"sup: Unknown argument type."
215    
216      @param arg  def inf(arg):
217      """      """
218      if isinstance(arg,escript.Data):      returns the maximum value over all data points.
219         if axis==None: axis=arg.getRank()/2  
220         return arg.transpose(axis)      @param arg: argument
221        @type arg: C{float}, C{int}, L{escript.Data}, L{numarray.NumArray}.
222        @return : minimum value of arg over all components and all data points
223        @rtype: C{float}
224        @raise TypeError: if type of arg cannot be processed
225        """
226        if isinstance(arg,numarray.NumArray):
227            return arg.min()
228        elif isinstance(arg,escript.Data):
229            return arg._inf()
230        elif isinstance(arg,float):
231            return arg
232        elif isinstance(arg,int):
233            return float(arg)
234      else:      else:
235         if axis==None: axis=arg.rank/2        raise TypeError,"inf: Unknown argument type."
236         return numarray.transpose(arg,axis=axis)  
237    
238  def trace(arg):  #=========================================================================
239    #   some little helpers
240    #=========================================================================
241    def pokeShape(arg):
242      """      """
243      @brief      identifies the shape of its argument
244    
245      @param arg      @param arg: a given object
246        @type arg: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, C{Symbol}
247        @return: the shape of the argument
248        @rtype: C{tuple} of C{int}
249        @raise TypeError: if type of arg cannot be processed
250      """      """
251      if isinstance(arg,escript.Data):  
252         return arg.trace()      if isinstance(arg,numarray.NumArray):
253            return arg.shape
254        elif isinstance(arg,escript.Data):
255            return arg.getShape()
256        elif isinstance(arg,float):
257            return ()
258        elif isinstance(arg,int):
259            return ()
260        elif isinstance(arg,Symbol):
261            return arg.getShape()
262      else:      else:
263         return numarray.trace(arg)        raise TypeError,"pokeShape: cannot identify shape"
264    
265  def exp(arg):  def pokeDim(arg):
266      """      """
267      @brief      identifies the spatial dimension of its argument
268    
269      @param arg      @param arg: a given object
270        @type arg: any
271        @return: the spatial dimension of the argument, if available, or C{None}
272        @rtype: C{int} or C{None}
273      """      """
274    
275      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
276         return arg.exp()          return arg.getFunctionSpace().getDim()
277        elif isinstance(arg,Symbol):
278            return arg.getDim()
279      else:      else:
280         return numarray.exp(arg)          return None
281    
282  def sqrt(arg):  def commonShape(arg0,arg1):
283      """      """
284      @brief      returns a shape to which arg0 can be extendent from the right and arg1 can be extended from the left.
285    
286      @param arg      @param arg0: an object with a shape (see L{pokeShape})
287        @param arg1: an object with a shape (see L{pokeShape})
288        @return: the shape of arg0 or arg1 such that the left port equals the shape of arg0 and the right end equals the shape of arg1.
289        @rtype: C{tuple} of C{int}
290        @raise ValueError: if no shape can be found.
291      """      """
292      if isinstance(arg,escript.Data):      sh0=pokeShape(arg0)
293         return arg.sqrt()      sh1=pokeShape(arg1)
294        if len(sh0)<len(sh1):
295           if not sh0==sh1[:len(sh0)]:
296                 raise ValueError,"argument 0 cannot be extended to the shape of argument 1"
297           return sh1
298        elif len(sh0)>len(sh1):
299           if not sh1==sh0[:len(sh1)]:
300                 raise ValueError,"argument 1 cannot be extended to the shape of argument 0"
301           return sh0
302      else:      else:
303         return numarray.sqrt(arg)         if not sh0==sh1:
304                 raise ValueError,"argument 1 and argument 0 have not the same shape."
305           return sh0
306    
307  def sin(arg):  def commonDim(*args):
308      """      """
309      @brief      identifies, if possible, the spatial dimension across a set of objects which may or my not have a spatial dimension.
310    
311      @param arg      @param *args: given objects
312        @return: the spatial dimension of the objects with identifiable dimension (see L{pokeDim}). If none the objects has
313                 a spatial dimension C{None} is returned.
314        @rtype: C{int} or C{None}
315        @raise ValueError: if the objects with identifiable dimension don't have the same spatial dimension.
316      """      """
317      if isinstance(arg,escript.Data):      out=None
318         return arg.sin()      for a in args:
319           d=pokeDim(a)
320           if not out==None:
321              if not (d==None or out==d):
322                 raise ValueError,"dimension of arguments don't match"
323           else:
324              out=d
325        return out
326    
327    def testForZero(arg):
328        """
329        test the argument for being identical to Zero
330    
331        @param arg: a given object
332        @type arg: typically L{numarray.NumArray},L{escript.Data},C{float}, C{int}
333        @return : True if the argument is identical to zero.
334        @rtype : C{bool}
335        """
336        if isinstance(arg,numarray.NumArray):
337           return not Lsup(arg)>0.
338        elif isinstance(arg,escript.Data):
339           return False
340        elif isinstance(arg,float):
341           return not Lsup(arg)>0.
342        elif isinstance(arg,int):
343           return not Lsup(arg)>0.
344        elif isinstance(arg,Symbol):
345           return False
346      else:      else:
347         return numarray.sin(arg)         return False
348    
349  def tan(arg):  def matchType(arg0=0.,arg1=0.):
350      """      """
351      @brief      converting arg0 and arg1 both to the same type L{numarray.NumArray} or L{escript.Data} or, if one of arg0 or arg1 is of type L{Symbol}, the other one to be of type L{numarray.NumArray} or L{escript.Data}.
352    
353      @param arg      @param arg0: first argument
354        @type arg0: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, C{Symbol}
355        @param arg1: second argument
356        @type arg1: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, C{Symbol}
357        @return: a tuple representing arg0 and arg1 with the same type or with one of them being a L{Symbol}
358        @rtype: C{tuple} of two L{numarray.NumArray}, two L{escript.Data}, a C{Symbol} and one of the types L{numarray.NumArray} or L{escript.Data}.  
359        @raise TypeError: if type of arg0 or arg1 cannot be processed
360      """      """
361      if isinstance(arg,escript.Data):      if isinstance(arg0,numarray.NumArray):
362         return arg.tan()         if isinstance(arg1,numarray.NumArray):
363              pass
364           elif isinstance(arg1,escript.Data):
365              arg0=escript.Data(arg0,arg1.getFunctionSpace())
366           elif isinstance(arg1,float):
367              arg1=numarray.array(arg1,type=numarray.Float64)
368           elif isinstance(arg1,int):
369              arg1=numarray.array(float(arg1),type=numarray.Float64)
370           elif isinstance(arg1,Symbol):
371              pass
372           else:
373              raise TypeError,"function: Unknown type of second argument."    
374        elif isinstance(arg0,escript.Data):
375           if isinstance(arg1,numarray.NumArray):
376              arg1=escript.Data(arg1,arg0.getFunctionSpace())
377           elif isinstance(arg1,escript.Data):
378              pass
379           elif isinstance(arg1,float):
380              arg1=escript.Data(arg1,(),arg0.getFunctionSpace())
381           elif isinstance(arg1,int):
382              arg1=escript.Data(float(arg1),(),arg0.getFunctionSpace())
383           elif isinstance(arg1,Symbol):
384              pass
385           else:
386              raise TypeError,"function: Unknown type of second argument."    
387        elif isinstance(arg0,Symbol):
388           if isinstance(arg1,numarray.NumArray):
389              pass
390           elif isinstance(arg1,escript.Data):
391              pass
392           elif isinstance(arg1,float):
393              arg1=numarray.array(arg1,type=numarray.Float64)
394           elif isinstance(arg1,int):
395              arg1=numarray.array(float(arg1),type=numarray.Float64)
396           elif isinstance(arg1,Symbol):
397              pass
398           else:
399              raise TypeError,"function: Unknown type of second argument."    
400        elif isinstance(arg0,float):
401           if isinstance(arg1,numarray.NumArray):
402              arg0=numarray.array(arg0,type=numarray.Float64)
403           elif isinstance(arg1,escript.Data):
404              arg0=escript.Data(arg0,arg1.getFunctionSpace())
405           elif isinstance(arg1,float):
406              arg0=numarray.array(arg0,type=numarray.Float64)
407              arg1=numarray.array(arg1,type=numarray.Float64)
408           elif isinstance(arg1,int):
409              arg0=numarray.array(arg0,type=numarray.Float64)
410              arg1=numarray.array(float(arg1),type=numarray.Float64)
411           elif isinstance(arg1,Symbol):
412              arg0=numarray.array(arg0,type=numarray.Float64)
413           else:
414              raise TypeError,"function: Unknown type of second argument."    
415        elif isinstance(arg0,int):
416           if isinstance(arg1,numarray.NumArray):
417              arg0=numarray.array(float(arg0),type=numarray.Float64)
418           elif isinstance(arg1,escript.Data):
419              arg0=escript.Data(float(arg0),arg1.getFunctionSpace())
420           elif isinstance(arg1,float):
421              arg0=numarray.array(float(arg0),type=numarray.Float64)
422              arg1=numarray.array(arg1,type=numarray.Float64)
423           elif isinstance(arg1,int):
424              arg0=numarray.array(float(arg0),type=numarray.Float64)
425              arg1=numarray.array(float(arg1),type=numarray.Float64)
426           elif isinstance(arg1,Symbol):
427              arg0=numarray.array(float(arg0),type=numarray.Float64)
428           else:
429              raise TypeError,"function: Unknown type of second argument."    
430      else:      else:
431         return numarray.tan(arg)        raise TypeError,"function: Unknown type of first argument."    
432    
433        return arg0,arg1
434    
435    def matchShape(arg0,arg1):
436        """
437        
438    
439        If shape is not given the shape "largest" shape of args is used.
440    
441        @param args: a given ob
442        @type arg: typically L{numarray.NumArray},L{escript.Data},C{float}, C{int}
443        @return: True if the argument is identical to zero.
444        @rtype: C{list} of C{int}
445        """
446        sh=commonShape(arg0,arg1)
447        sh0=pokeShape(arg0)
448        sh1=pokeShape(arg1)
449        if len(sh0)<len(sh):
450           return outer(arg0,numarray.ones(sh[len(sh0):],numarray.Float64)),arg1
451        elif len(sh1)<len(sh):
452           return arg0,outer(arg1,numarray.ones(sh[len(sh1):],numarray.Float64))
453        else:
454           return arg0,arg1
455    #=========================================================
456    #   symbolic tool box starts here:
457    #=========================================================
458    class Symbol(object):
459       """
460       Symbol class.
461    
462       Symbol class objects provide the same functionality as L{numarray.NumArray} and L{escript.Data} objects
463       but they do not have a value and therefore cannot be plotted or visualize. The main purpose is the possibilty
464       calculate derivatives with respect to other Symbols used to define a Symbol.
465    
466       """
467       def __init__(self,shape=(),args=[],dim=None):
468           """
469           Creates an instance of a symbol of a given shape. The symbol may depending on a list of arguments args which may be
470           symbols or any other object.
471    
472           @param arg: the arguments of the symbol.
473           @type arg: C{list}
474           @param shape: the shape
475           @type shape: C{tuple} of C{int}
476           @param dim: spatial dimension of the symbol. If dim=C{None} the spatial dimension is undefined.  
477           @type dim: C{None} or C{int}  
478    
479           """
480           if len(shape)>4:
481               raise ValueError,"Symbol supports only tensors up to order 4"
482           self.__args=args
483           self.__shape=shape
484           self.__dim=dim
485    
486       def getArgument(self,i=None):
487           """
488           returns the i-th argument of the symbol
489    
490           @param i: index of the argument requested.
491           @type i: C{int} or C{None}
492           @raise IndexError: if the requested index does not exist
493           @return: the vlaue of the i-th argument or i is not specified the list of all arguments.
494           @rtype: a single object or a list of objects
495           """
496           if i==None:
497              return self.__args
498           else:
499              if i<0 or i>=len(self.__args):
500                 raise IndexError,"there are only %s arguments"%len(self.__args)
501              return self.__args[i]
502    
503       def getRank(self):
504           """
505           the rank of the symbol
506    
507           @return: the rank of the symbol. This is length of the shape
508           @rtype: C{int}
509           """
510           return len(self.getShape())
511    
512       def getShape(self):
513           """
514           the shape of the symbol.
515    
516           @return : the shape of the symbol.
517           @rtype: C{tuple} of C{int}
518           """
519           return self.__shape
520    
521       def getDim(self):
522           """
523           the spatial dimension
524    
525           @return : the spatial dimension
526           @rtype: C{int} if the dimension is defined. Otherwise C{None} is returned.
527           """
528           return self.__dim
529      
530       def __str__(self):
531           """
532           a string representation of the symbol.
533           @return: a string representation of the object
534           @rtype: C{str}
535           """
536           args=[]
537           for arg in self.getArgument():
538              args.append(str(arg))
539           try:
540               out=self.getMyCode(args,format="str")
541           except NotImplementedError:
542               out="<Symbol %s>"%id(self)
543           return out
544          
545       def getSubstitutedArguments(self,argvals):
546           """
547           substitutes symbols in the arguments of this object and returns the result as a list.
548    
549           @param argvals: L{Symbols} and their substitutes. The L{Symbol} u in the expression defining this object is replaced by argvals[u].
550           @type argvals: C{dict} with keywords of type L{Symbol}.
551           @rtype: C{list} of objects
552           @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.
553           """
554           out=[]
555           for a in self.getArgument():
556              if isinstance(a,Symbol):
557                 out.append(a.substitute(argvals))
558              else:
559                 out.append(a)
560           return out
561    
562       def getDifferentiatedArguments(self,arg):
563           """
564           applifies differentials to the arguments of this object and returns the result as a list.
565    
566           @param arg: the derivative is calculated with respect to arg
567           @type arg: typically L{escript.Symbol} but can also be C{float}, L{escript.Data}, L{numarray.NumArray} depending the involved functions and data.
568           @rtype: C{list} of objects
569           @return: list of object obtained by calculating the derivatives of the argumenst with respct to arg
570           """
571           out=[]
572           for a in self.getArgument():
573              if isinstance(a,Symbol):
574                 out.append(a.substitute(argvals))
575              else:
576                  s=pokeShape(s)+arg.getShape()
577                  if len(s)>0:
578                     out.append(numarray.zeros(s),numarray.Float64)
579                  else:
580                     out.append(a)
581           return out
582    
583       def isAppropriateValue(self,arg):
584          """
585          checks if the given argument arg can be used as a substitution of this object. The method checks
586          the shape of arg and, if the spatial dimension is defined, the spatial dimension of arg.    
587    
588          @param arg: a given object
589          @type arg: L{numarray.NumArray},L{escript.Data},C{float}, C{int}, C{Symbol}
590          @return: True if arg is a suitbale object to be used for substitution. Otherwise False is returned.
591          @rtype: C{bool}
592          """
593          if isinstance(arg,numarray.NumArray):
594              return arg.shape==self.getShape()
595          elif isinstance(arg,escript.Data):
596              if self.getDim()==None:
597                  return arg.getShape()==self.getShape()
598              elif self.getDim()==arg.getFunctionSpace().getDim():
599                  return arg.getShape()==self.getShape()
600              else:
601                  return False
602          elif isinstance(arg,Symbol):
603              if self.getDim()==None:
604                  return arg.getShape()==self.getShape()
605              elif self.getDim()==arg.getDim():
606                  return arg.getShape()==self.getShape()
607              else:
608                  return False
609          elif isinstance(arg,float):
610              return ()==self.getShape()
611          elif isinstance(arg,int):
612              return ()==self.getShape()
613          else:
614             return False
615    
616       def getMyCode(self,argstrs,format="escript"):
617           """
618           returns a program code that can be used to evaluate the symbol.
619    
620           @param argstrs: gives for each argument a string representing the argument for the evaluation.
621           @type argstrs: C{list} of C{str}.
622           @param format: specifies the format to be used. At the moment only "escript", "str" and "text" are supported.
623           @type format: C{str}
624           @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
625           @rtype: C{str}
626           @raise NotImplementedError: if no implementation for the given format is available
627           @note: This method has to be overwritten by subclasses.
628           """
629           raise NotImplementedError,"no code for %s representation available"%format
630    
631       def substitute(self,argvals):
632          """  
633          assigns new values to symbols in the definition of the symbol.
634          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
635    
636          @param argvals: new values assigned to symbols
637          @type argvals: C{dict} with keywords of type L{Symbol}.
638          @return: result of the substitution process. Operations are executed as much as possible.
639          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
640          @note: this method has to be overwritten by a particular L{Symbol}
641          @raise NotImplementedError: if no implementation for the given format is available
642          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
643          """
644          if argvals.has_key(self):
645             arg=argvals[self]
646             if self.isAppropriateValue(arg):
647                return arg
648             else:
649                raise TypeError,"Symbol: new value is not appropriate."
650          else:
651             raise NotImplementedError,"no substitution in %s avialable"%str(self)
652    
653       def diff(self,arg):
654           """
655           returns the derivative of the symbol with respect to L{Symbol} arg
656    
657           @param arg: the derivative is calculated with respect to arg
658           @type arg: typically L{escript.Symbol} but can also be C{float}, L{escript.Data}, L{numarray.NumArray} depending the involved functions and data.
659           @return: derivative with respect to C{arg}
660           @rtype: typically L{escript.Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
661           @note: this method is overwritten by a particular L{Symbol}
662           """
663           if arg==self:
664              return identity(self.getShape())
665           else:
666              s=self.getShape()+arg.getShape()
667              if len(s)>0:
668                 return numarray.zeros(s,numarray.Float64)
669              else:
670                 return 0.
671    
672       def __neg__(self):
673           """
674           returns -self.
675    
676           @return:  a S{Symbol} representing the negative of the object
677           @rtype: L{DependendSymbol}
678           """
679           return self*(-1.)
680    
681       def __pos__(self):
682           """
683           returns +self.
684    
685           @return:  a S{Symbol} representing the positive of the object
686           @rtype: L{DependendSymbol}
687           """
688           return self*(1.)
689    
690       def __abs__(self):
691           """
692           returns a S{Symbol} representing the absolute value of the object.
693           """
694           return Abs_Symbol(self)
695    
696       def __add__(self,other):
697           """
698           add another object to this object
699    
700           @param other: object to be added to this object
701           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
702           @return:  a S{Symbol} representing the sum of this object and C{other}
703           @rtype: L{DependendSymbol}
704           """
705           return add(self,other)
706    
707       def __radd__(self,other):
708           """
709           add this object to another object
710    
711           @param other: object this object is added to
712           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
713           @return: a S{Symbol} representing the sum of C{other} and this object object
714           @rtype: L{DependendSymbol}
715           """
716           return add(other,self)
717    
718       def __sub__(self,other):
719           """
720           subtracts another object from this object
721    
722           @param other: object to be subtracted from this object
723           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
724           @return: a S{Symbol} representing the difference of C{other} and this object
725           @rtype: L{DependendSymbol}
726           """
727           return add(self,-other)
728    
729       def __rsub__(self,other):
730           """
731           subtracts this object from another object
732    
733           @param other: object this object is been subtracted from
734           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
735           @return: a S{Symbol} representing the difference of this object and C{other}.
736           @rtype: L{DependendSymbol}
737           """
738           return add(-self,other)
739    
740       def __mul__(self,other):
741           """
742           multiplies this object with other object
743    
744           @param other: object to be mutiplied by this object
745           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
746           @return: a S{Symbol} representing the product of the object and C{other}.
747           @rtype: L{DependendSymbol} or 0 if other is identical to zero.
748           """
749           return mult(self,other)
750    
751       def __rmul__(self,other):
752           """
753           multiplies this object with other object
754    
755           @param other: object this object is multiplied with
756           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
757           @return: a S{Symbol} representing the product of C{other} and the object.
758           @rtype: L{DependendSymbol} or 0 if other is identical to zero.
759           """
760           return mult(other,self)
761    
762       def __div__(self,other):
763           """
764           divides this object by other object
765    
766           @param other: object dividing this object
767           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
768           @return: a S{Symbol} representing the quotient of this object and C{other}
769           @rtype: L{DependendSymbol}
770           """
771           return quotient(self,other)
772    
773       def __rdiv__(self,other):
774           """
775           divides this object by other object
776    
777           @param other: object dividing this object
778           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
779           @return: a S{Symbol} representing the quotient of C{other} and this object
780           @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.
781           """
782           return quotient(other,self)
783    
784       def __pow__(self,other):
785           """
786           raises this object to the power of other
787    
788           @param other: exponent
789           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
790           @return: a S{Symbol} representing the power of this object to C{other}
791           @rtype: L{DependendSymbol} or 1 if C{other} is identical to zero.
792           """
793           return power(self,other)
794    
795       def __rpow__(self,other):
796           """
797           raises an object to the power of this object
798    
799           @param other: basis
800           @type other: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray}.
801           @return: a S{Symbol} representing the power of C{other} to this object
802           @rtype: L{DependendSymbol} or 0 if C{other} is identical to zero.
803           """
804           return power(other,self)
805    
806       def __getitem__(self,index):
807           """
808           returns the slice defined by index
809    
810           @param index: defines a
811           @type index: C{slice} or C{int} or a C{tuple} of them
812           @return: a S{Symbol} representing the slice defined by index
813           @rtype: L{DependendSymbol}
814           """
815           return GetSlice_Symbol(self,index)
816    
817    class DependendSymbol(Symbol):
818       """
819       DependendSymbol extents L{Symbol} by modifying the == operator to allow two instances to be equal.
820       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  
821      
822       Example:
823      
824       u1=Symbol(shape=(3,4),dim=2,args=[4.])
825       u2=Symbol(shape=(3,4),dim=2,args=[4.])
826       print u1==u2
827       False
828      
829          but  
830    
831       u1=DependendSymbol(shape=(3,4),dim=2,args=[4.])
832       u2=DependendSymbol(shape=(3,4),dim=2,args=[4.])
833       u3=DependendSymbol(shape=(2,),dim=2,args=[4.])  
834       print u1==u2, u1==u3
835       True False
836    
837       @note: DependendSymbol should be used as return value of functions with L{Symbol} arguments. This will allow the optimizer to remove redundant function calls.
838       """
839       def __eq__(self,other):
840          """
841          checks if other equals self
842    
843          @param other: any object
844          @return: True if other has the same class like self, and the shape, the spatial diemsnion and the arguments are equal.
845          @rtype: C{bool}
846          """
847          if isinstance(other,DependendSymbol):
848             if self.__class__==other.__class__:  
849                if self.getShape()==other.getShape():
850                  if self.getArgument()==other.getArgument():
851                     if self.getDim()==None or other.getDim()==None or self.getDim()==other.getDim():
852                        return True
853          return False
854    
855       def __ne__(self,other):
856          """
857          checks if other equals self
858    
859          @param other: any object
860          @return: Flase if other has the same class like self, and the shape, the spatial diemsnion and the arguments are equal.
861          @rtype: C{bool}
862          """
863          return not self==other
864    #=========================================================
865    #  Unary operations prserving the shape
866    #========================================================
867    class GetSlice_Symbol(DependendSymbol):
868       """
869       L{Symbol} representing getting a slice for a L{Symbol}
870       """
871       def __init__(self,arg,index):
872          """
873          initialization of wherePositive L{Symbol} with argument arg
874          @param arg: argument
875          @type arg: L{Symbol}.
876          @param index: defines index
877          @type index: C{slice} or C{int} or a C{tuple} of them
878          @raises IndexError: if length of index is larger than rank of arg or a index start or stop is out of range
879          @raises ValueError: if a step is given
880          """
881          if not isinstance(index,tuple): index=(index,)
882          if len(index)>arg.getRank():
883               raise IndexError,"GetSlice_Symbol: index out of range."
884          sh=()
885          index2=()
886          for i in range(len(index)):
887             ix=index[i]
888             if isinstance(ix,int):
889                if ix<0 or ix>=arg.getShape()[i]:
890                   raise ValueError,"GetSlice_Symbol: index out of range."
891                index2=index2+(ix,)
892             else:
893               if not ix.step==None:
894                 raise ValueError,"GetSlice_Symbol: steping is not supported."
895               if ix.start==None:
896                  s=0
897               else:
898                  s=ix.start
899               if ix.stop==None:
900                  e=arg.getShape()[i]
901               else:
902                  e=ix.stop
903                  if e>arg.getShape()[i]:
904                     raise IndexError,"GetSlice_Symbol: index out of range."
905               index2=index2+(slice(s,e),)
906               if e>s:
907                   sh=sh+(e-s,)
908               elif s>e:
909                   raise IndexError,"GetSlice_Symbol: slice start must be less or equal slice end"
910          for i in range(len(index),arg.getRank()):
911              index2=index2+(slice(0,arg.getShape()[i]),)
912              sh=sh+(arg.getShape()[i],)
913          super(GetSlice_Symbol, self).__init__(args=[arg,index2],shape=sh,dim=arg.getDim())
914    
915       def getMyCode(self,argstrs,format="escript"):
916          """
917          returns a program code that can be used to evaluate the symbol.
918    
919          @param argstrs: gives for each argument a string representing the argument for the evaluation.
920          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
921          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
922          @type format: C{str}
923          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
924          @rtype: C{str}
925          @raise: NotImplementedError: if the requested format is not available
926          """
927          if format=="escript" or format=="str"  or format=="text":
928             return "%s.__getitem__(%s)"%(argstrs[0],argstrs[1])
929          else:
930             raise NotImplementedError,"GetItem_Symbol does not provide program code for format %s."%format
931    
932       def substitute(self,argvals):
933          """
934          assigns new values to symbols in the definition of the symbol.
935          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
936    
937          @param argvals: new values assigned to symbols
938          @type argvals: C{dict} with keywords of type L{Symbol}.
939          @return: result of the substitution process. Operations are executed as much as possible.
940          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
941          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
942          """
943          if argvals.has_key(self):
944             arg=argvals[self]
945             if self.isAppropriateValue(arg):
946                return arg
947             else:
948                raise TypeError,"%s: new value is not appropriate."%str(self)
949          else:
950             args=self.getSubstitutedArguments(argvals)
951             arg=args[0]
952             index=args[1]
953             return arg.__getitem__(index)
954    
955    def log10(arg):
956       """
957       returns base-10 logarithm of argument arg
958    
959       @param arg: argument
960       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
961       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
962       @raises TypeError: if the type of the argument is not expected.
963       """
964       if isinstance(arg,numarray.NumArray):
965          return numarray.log10(arg)
966       elif isinstance(arg,escript.Data):
967          return arg._log10()
968       elif isinstance(arg,float):
969          return math.log10(arg)
970       elif isinstance(arg,int):
971          return math.log10(float(arg))
972       elif isinstance(arg,Symbol):
973          return log(arg)/log(10.)
974       else:
975          raise TypeError,"log10: Unknown argument type."
976    
977    def wherePositive(arg):
978       """
979       returns mask of positive values of argument arg
980    
981       @param arg: argument
982       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
983       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
984       @raises TypeError: if the type of the argument is not expected.
985       """
986       if isinstance(arg,numarray.NumArray):
987          out=numarray.greater(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
988          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
989          return out
990       elif isinstance(arg,escript.Data):
991          return arg._wherePositive()
992       elif isinstance(arg,float):
993          if arg>0:
994            return 1.
995          else:
996            return 0.
997       elif isinstance(arg,int):
998          if arg>0:
999            return 1.
1000          else:
1001            return 0.
1002       elif isinstance(arg,Symbol):
1003          return WherePositive_Symbol(arg)
1004       else:
1005          raise TypeError,"wherePositive: Unknown argument type."
1006    
1007    class WherePositive_Symbol(DependendSymbol):
1008       """
1009       L{Symbol} representing the result of the mask of positive values function
1010       """
1011       def __init__(self,arg):
1012          """
1013          initialization of wherePositive L{Symbol} with argument arg
1014          @param arg: argument of function
1015          @type arg: typically L{Symbol}.
1016          """
1017          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1018    
1019       def getMyCode(self,argstrs,format="escript"):
1020          """
1021          returns a program code that can be used to evaluate the symbol.
1022    
1023          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1024          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1025          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1026          @type format: C{str}
1027          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1028          @rtype: C{str}
1029          @raise: NotImplementedError: if the requested format is not available
1030          """
1031          if isinstance(argstrs,list):
1032              argstrs=argstrs[0]
1033          if format=="escript" or format=="str"  or format=="text":
1034             return "wherePositive(%s)"%argstrs
1035          else:
1036             raise NotImplementedError,"WherePositive_Symbol does not provide program code for format %s."%format
1037    
1038       def substitute(self,argvals):
1039          """
1040          assigns new values to symbols in the definition of the symbol.
1041          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1042    
1043          @param argvals: new values assigned to symbols
1044          @type argvals: C{dict} with keywords of type L{Symbol}.
1045          @return: result of the substitution process. Operations are executed as much as possible.
1046          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1047          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1048          """
1049          if argvals.has_key(self):
1050             arg=argvals[self]
1051             if self.isAppropriateValue(arg):
1052                return arg
1053             else:
1054                raise TypeError,"%s: new value is not appropriate."%str(self)
1055          else:
1056             arg=self.getSubstitutedArguments(argvals)[0]
1057             return wherePositive(arg)
1058    
1059    def whereNegative(arg):
1060       """
1061       returns mask of positive values of argument arg
1062    
1063       @param arg: argument
1064       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1065       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1066       @raises TypeError: if the type of the argument is not expected.
1067       """
1068       if isinstance(arg,numarray.NumArray):
1069          out=numarray.less(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1070          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1071          return out
1072       elif isinstance(arg,escript.Data):
1073          return arg._whereNegative()
1074       elif isinstance(arg,float):
1075          if arg<0:
1076            return 1.
1077          else:
1078            return 0.
1079       elif isinstance(arg,int):
1080          if arg<0:
1081            return 1.
1082          else:
1083            return 0.
1084       elif isinstance(arg,Symbol):
1085          return WhereNegative_Symbol(arg)
1086       else:
1087          raise TypeError,"whereNegative: Unknown argument type."
1088    
1089    class WhereNegative_Symbol(DependendSymbol):
1090       """
1091       L{Symbol} representing the result of the mask of positive values function
1092       """
1093       def __init__(self,arg):
1094          """
1095          initialization of whereNegative L{Symbol} with argument arg
1096          @param arg: argument of function
1097          @type arg: typically L{Symbol}.
1098          """
1099          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1100    
1101       def getMyCode(self,argstrs,format="escript"):
1102          """
1103          returns a program code that can be used to evaluate the symbol.
1104    
1105          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1106          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1107          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1108          @type format: C{str}
1109          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1110          @rtype: C{str}
1111          @raise: NotImplementedError: if the requested format is not available
1112          """
1113          if isinstance(argstrs,list):
1114              argstrs=argstrs[0]
1115          if format=="escript" or format=="str"  or format=="text":
1116             return "whereNegative(%s)"%argstrs
1117          else:
1118             raise NotImplementedError,"WhereNegative_Symbol does not provide program code for format %s."%format
1119    
1120       def substitute(self,argvals):
1121          """
1122          assigns new values to symbols in the definition of the symbol.
1123          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1124    
1125          @param argvals: new values assigned to symbols
1126          @type argvals: C{dict} with keywords of type L{Symbol}.
1127          @return: result of the substitution process. Operations are executed as much as possible.
1128          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1129          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1130          """
1131          if argvals.has_key(self):
1132             arg=argvals[self]
1133             if self.isAppropriateValue(arg):
1134                return arg
1135             else:
1136                raise TypeError,"%s: new value is not appropriate."%str(self)
1137          else:
1138             arg=self.getSubstitutedArguments(argvals)[0]
1139             return whereNegative(arg)
1140    
1141    def whereNonNegative(arg):
1142       """
1143       returns mask of non-negative values of argument arg
1144    
1145       @param arg: argument
1146       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1147       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1148       @raises TypeError: if the type of the argument is not expected.
1149       """
1150       if isinstance(arg,numarray.NumArray):
1151          out=numarray.greater_equal(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1152          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1153          return out
1154       elif isinstance(arg,escript.Data):
1155          return arg._whereNonNegative()
1156       elif isinstance(arg,float):
1157          if arg<0:
1158            return 0.
1159          else:
1160            return 1.
1161       elif isinstance(arg,int):
1162          if arg<0:
1163            return 0.
1164          else:
1165            return 1.
1166       elif isinstance(arg,Symbol):
1167          return 1.-whereNegative(arg)
1168       else:
1169          raise TypeError,"whereNonNegative: Unknown argument type."
1170    
1171    def whereNonPositive(arg):
1172       """
1173       returns mask of non-positive values of argument arg
1174    
1175       @param arg: argument
1176       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1177       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1178       @raises TypeError: if the type of the argument is not expected.
1179       """
1180       if isinstance(arg,numarray.NumArray):
1181          out=numarray.less_equal(arg,numarray.zeros(arg.shape,numarray.Float64))*1.
1182          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1183          return out
1184       elif isinstance(arg,escript.Data):
1185          return arg._whereNonPositive()
1186       elif isinstance(arg,float):
1187          if arg>0:
1188            return 0.
1189          else:
1190            return 1.
1191       elif isinstance(arg,int):
1192          if arg>0:
1193            return 0.
1194          else:
1195            return 1.
1196       elif isinstance(arg,Symbol):
1197          return 1.-wherePositive(arg)
1198       else:
1199          raise TypeError,"whereNonPositive: Unknown argument type."
1200    
1201    def whereZero(arg,tol=0.):
1202       """
1203       returns mask of zero entries of argument arg
1204    
1205       @param arg: argument
1206       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1207       @param tol: tolerance. values with absolute value less then tol are accepted as zero.
1208       @type tol: C{float}
1209       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1210       @raises TypeError: if the type of the argument is not expected.
1211       """
1212       if isinstance(arg,numarray.NumArray):
1213          out=numarray.less_equal(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float64))*1.
1214          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1215          return out
1216       elif isinstance(arg,escript.Data):
1217          if tol>0.:
1218             return whereNegative(abs(arg)-tol)
1219          else:
1220             return arg._whereZero()
1221       elif isinstance(arg,float):
1222          if abs(arg)<=tol:
1223            return 1.
1224          else:
1225            return 0.
1226       elif isinstance(arg,int):
1227          if abs(float(arg))<=tol:
1228            return 1.
1229          else:
1230            return 0.
1231       elif isinstance(arg,Symbol):
1232          return WhereZero_Symbol(arg,tol)
1233       else:
1234          raise TypeError,"whereZero: Unknown argument type."
1235    
1236    class WhereZero_Symbol(DependendSymbol):
1237       """
1238       L{Symbol} representing the result of the mask of zero entries function
1239       """
1240       def __init__(self,arg,tol=0.):
1241          """
1242          initialization of whereZero L{Symbol} with argument arg
1243          @param arg: argument of function
1244          @type arg: typically L{Symbol}.
1245          """
1246          DependendSymbol.__init__(self,args=[arg,tol],shape=arg.getShape(),dim=arg.getDim())
1247    
1248       def getMyCode(self,argstrs,format="escript"):
1249          """
1250          returns a program code that can be used to evaluate the symbol.
1251    
1252          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1253          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1254          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1255          @type format: C{str}
1256          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1257          @rtype: C{str}
1258          @raise: NotImplementedError: if the requested format is not available
1259          """
1260          if format=="escript" or format=="str"  or format=="text":
1261             return "whereZero(%s,tol=%s)"%(argstrs[0],argstrs[1])
1262          else:
1263             raise NotImplementedError,"WhereZero_Symbol does not provide program code for format %s."%format
1264    
1265       def substitute(self,argvals):
1266          """
1267          assigns new values to symbols in the definition of the symbol.
1268          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1269    
1270          @param argvals: new values assigned to symbols
1271          @type argvals: C{dict} with keywords of type L{Symbol}.
1272          @return: result of the substitution process. Operations are executed as much as possible.
1273          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1274          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1275          """
1276          if argvals.has_key(self):
1277             arg=argvals[self]
1278             if self.isAppropriateValue(arg):
1279                return arg
1280             else:
1281                raise TypeError,"%s: new value is not appropriate."%str(self)
1282          else:
1283             arg=self.getSubstitutedArguments(argvals)
1284             return whereZero(arg[0],arg[1])
1285    
1286    def whereNonZero(arg,tol=0.):
1287       """
1288       returns mask of values different from zero of argument arg
1289    
1290       @param arg: argument
1291       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1292       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1293       @raises TypeError: if the type of the argument is not expected.
1294       """
1295       if isinstance(arg,numarray.NumArray):
1296          out=numarray.greater(abs(arg)-tol,numarray.zeros(arg.shape,numarray.Float64))*1.
1297          if isinstance(out,float): out=numarray.array(out,type=numarray.Float64)
1298          return out
1299       elif isinstance(arg,escript.Data):
1300          if tol>0.:
1301             return 1.-whereZero(arg,tol)
1302          else:
1303             return arg._whereNonZero()
1304       elif isinstance(arg,float):
1305          if abs(arg)>tol:
1306            return 1.
1307          else:
1308            return 0.
1309       elif isinstance(arg,int):
1310          if abs(float(arg))>tol:
1311            return 1.
1312          else:
1313            return 0.
1314       elif isinstance(arg,Symbol):
1315          return 1.-whereZero(arg,tol)
1316       else:
1317          raise TypeError,"whereNonZero: Unknown argument type."
1318    
1319    def sin(arg):
1320       """
1321       returns sine of argument arg
1322    
1323       @param arg: argument
1324       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1325       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1326       @raises TypeError: if the type of the argument is not expected.
1327       """
1328       if isinstance(arg,numarray.NumArray):
1329          return numarray.sin(arg)
1330       elif isinstance(arg,escript.Data):
1331          return arg._sin()
1332       elif isinstance(arg,float):
1333          return math.sin(arg)
1334       elif isinstance(arg,int):
1335          return math.sin(arg)
1336       elif isinstance(arg,Symbol):
1337          return Sin_Symbol(arg)
1338       else:
1339          raise TypeError,"sin: Unknown argument type."
1340    
1341    class Sin_Symbol(DependendSymbol):
1342       """
1343       L{Symbol} representing the result of the sine function
1344       """
1345       def __init__(self,arg):
1346          """
1347          initialization of sin L{Symbol} with argument arg
1348          @param arg: argument of function
1349          @type arg: typically L{Symbol}.
1350          """
1351          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1352    
1353       def getMyCode(self,argstrs,format="escript"):
1354          """
1355          returns a program code that can be used to evaluate the symbol.
1356    
1357          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1358          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1359          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1360          @type format: C{str}
1361          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1362          @rtype: C{str}
1363          @raise: NotImplementedError: if the requested format is not available
1364          """
1365          if isinstance(argstrs,list):
1366              argstrs=argstrs[0]
1367          if format=="escript" or format=="str"  or format=="text":
1368             return "sin(%s)"%argstrs
1369          else:
1370             raise NotImplementedError,"Sin_Symbol does not provide program code for format %s."%format
1371    
1372       def substitute(self,argvals):
1373          """
1374          assigns new values to symbols in the definition of the symbol.
1375          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1376    
1377          @param argvals: new values assigned to symbols
1378          @type argvals: C{dict} with keywords of type L{Symbol}.
1379          @return: result of the substitution process. Operations are executed as much as possible.
1380          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1381          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1382          """
1383          if argvals.has_key(self):
1384             arg=argvals[self]
1385             if self.isAppropriateValue(arg):
1386                return arg
1387             else:
1388                raise TypeError,"%s: new value is not appropriate."%str(self)
1389          else:
1390             arg=self.getSubstitutedArguments(argvals)[0]
1391             return sin(arg)
1392    
1393       def diff(self,arg):
1394          """
1395          differential of this object
1396    
1397          @param arg: the derivative is calculated with respect to arg
1398          @type arg: L{escript.Symbol}
1399          @return: derivative with respect to C{arg}
1400          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1401          """
1402          if arg==self:
1403             return identity(self.getShape())
1404          else:
1405             myarg=self.getArgument()[0]
1406             val=matchShape(cos(myarg),self.getDifferentiatedArguments(arg)[0])
1407             return val[0]*val[1]
1408    
1409  def cos(arg):  def cos(arg):
1410       """
1411       returns cosine of argument arg
1412    
1413       @param arg: argument
1414       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1415       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1416       @raises TypeError: if the type of the argument is not expected.
1417       """
1418       if isinstance(arg,numarray.NumArray):
1419          return numarray.cos(arg)
1420       elif isinstance(arg,escript.Data):
1421          return arg._cos()
1422       elif isinstance(arg,float):
1423          return math.cos(arg)
1424       elif isinstance(arg,int):
1425          return math.cos(arg)
1426       elif isinstance(arg,Symbol):
1427          return Cos_Symbol(arg)
1428       else:
1429          raise TypeError,"cos: Unknown argument type."
1430    
1431    class Cos_Symbol(DependendSymbol):
1432       """
1433       L{Symbol} representing the result of the cosine function
1434       """
1435       def __init__(self,arg):
1436          """
1437          initialization of cos L{Symbol} with argument arg
1438          @param arg: argument of function
1439          @type arg: typically L{Symbol}.
1440          """
1441          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1442    
1443       def getMyCode(self,argstrs,format="escript"):
1444          """
1445          returns a program code that can be used to evaluate the symbol.
1446    
1447          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1448          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1449          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1450          @type format: C{str}
1451          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1452          @rtype: C{str}
1453          @raise: NotImplementedError: if the requested format is not available
1454          """
1455          if isinstance(argstrs,list):
1456              argstrs=argstrs[0]
1457          if format=="escript" or format=="str"  or format=="text":
1458             return "cos(%s)"%argstrs
1459          else:
1460             raise NotImplementedError,"Cos_Symbol does not provide program code for format %s."%format
1461    
1462       def substitute(self,argvals):
1463          """
1464          assigns new values to symbols in the definition of the symbol.
1465          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1466    
1467          @param argvals: new values assigned to symbols
1468          @type argvals: C{dict} with keywords of type L{Symbol}.
1469          @return: result of the substitution process. Operations are executed as much as possible.
1470          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1471          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1472          """
1473          if argvals.has_key(self):
1474             arg=argvals[self]
1475             if self.isAppropriateValue(arg):
1476                return arg
1477             else:
1478                raise TypeError,"%s: new value is not appropriate."%str(self)
1479          else:
1480             arg=self.getSubstitutedArguments(argvals)[0]
1481             return cos(arg)
1482    
1483       def diff(self,arg):
1484          """
1485          differential of this object
1486    
1487          @param arg: the derivative is calculated with respect to arg
1488          @type arg: L{escript.Symbol}
1489          @return: derivative with respect to C{arg}
1490          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1491          """
1492          if arg==self:
1493             return identity(self.getShape())
1494          else:
1495             myarg=self.getArgument()[0]
1496             val=matchShape(-sin(myarg),self.getDifferentiatedArguments(arg)[0])
1497             return val[0]*val[1]
1498    
1499    def tan(arg):
1500       """
1501       returns tangent of argument arg
1502    
1503       @param arg: argument
1504       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1505       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1506       @raises TypeError: if the type of the argument is not expected.
1507       """
1508       if isinstance(arg,numarray.NumArray):
1509          return numarray.tan(arg)
1510       elif isinstance(arg,escript.Data):
1511          return arg._tan()
1512       elif isinstance(arg,float):
1513          return math.tan(arg)
1514       elif isinstance(arg,int):
1515          return math.tan(arg)
1516       elif isinstance(arg,Symbol):
1517          return Tan_Symbol(arg)
1518       else:
1519          raise TypeError,"tan: Unknown argument type."
1520    
1521    class Tan_Symbol(DependendSymbol):
1522       """
1523       L{Symbol} representing the result of the tangent function
1524       """
1525       def __init__(self,arg):
1526          """
1527          initialization of tan L{Symbol} with argument arg
1528          @param arg: argument of function
1529          @type arg: typically L{Symbol}.
1530          """
1531          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1532    
1533       def getMyCode(self,argstrs,format="escript"):
1534          """
1535          returns a program code that can be used to evaluate the symbol.
1536    
1537          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1538          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1539          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1540          @type format: C{str}
1541          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1542          @rtype: C{str}
1543          @raise: NotImplementedError: if the requested format is not available
1544          """
1545          if isinstance(argstrs,list):
1546              argstrs=argstrs[0]
1547          if format=="escript" or format=="str"  or format=="text":
1548             return "tan(%s)"%argstrs
1549          else:
1550             raise NotImplementedError,"Tan_Symbol does not provide program code for format %s."%format
1551    
1552       def substitute(self,argvals):
1553          """
1554          assigns new values to symbols in the definition of the symbol.
1555          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1556    
1557          @param argvals: new values assigned to symbols
1558          @type argvals: C{dict} with keywords of type L{Symbol}.
1559          @return: result of the substitution process. Operations are executed as much as possible.
1560          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1561          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1562          """
1563          if argvals.has_key(self):
1564             arg=argvals[self]
1565             if self.isAppropriateValue(arg):
1566                return arg
1567             else:
1568                raise TypeError,"%s: new value is not appropriate."%str(self)
1569          else:
1570             arg=self.getSubstitutedArguments(argvals)[0]
1571             return tan(arg)
1572    
1573       def diff(self,arg):
1574          """
1575          differential of this object
1576    
1577          @param arg: the derivative is calculated with respect to arg
1578          @type arg: L{escript.Symbol}
1579          @return: derivative with respect to C{arg}
1580          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1581          """
1582          if arg==self:
1583             return identity(self.getShape())
1584          else:
1585             myarg=self.getArgument()[0]
1586             val=matchShape(1./cos(myarg)**2,self.getDifferentiatedArguments(arg)[0])
1587             return val[0]*val[1]
1588    
1589    def asin(arg):
1590       """
1591       returns inverse sine of argument arg
1592    
1593       @param arg: argument
1594       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1595       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1596       @raises TypeError: if the type of the argument is not expected.
1597       """
1598       if isinstance(arg,numarray.NumArray):
1599          return numarray.arcsin(arg)
1600       elif isinstance(arg,escript.Data):
1601          return arg._asin()
1602       elif isinstance(arg,float):
1603          return math.asin(arg)
1604       elif isinstance(arg,int):
1605          return math.asin(arg)
1606       elif isinstance(arg,Symbol):
1607          return Asin_Symbol(arg)
1608       else:
1609          raise TypeError,"asin: Unknown argument type."
1610    
1611    class Asin_Symbol(DependendSymbol):
1612       """
1613       L{Symbol} representing the result of the inverse sine function
1614       """
1615       def __init__(self,arg):
1616          """
1617          initialization of asin L{Symbol} with argument arg
1618          @param arg: argument of function
1619          @type arg: typically L{Symbol}.
1620          """
1621          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1622    
1623       def getMyCode(self,argstrs,format="escript"):
1624          """
1625          returns a program code that can be used to evaluate the symbol.
1626    
1627          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1628          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1629          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1630          @type format: C{str}
1631          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1632          @rtype: C{str}
1633          @raise: NotImplementedError: if the requested format is not available
1634          """
1635          if isinstance(argstrs,list):
1636              argstrs=argstrs[0]
1637          if format=="escript" or format=="str"  or format=="text":
1638             return "asin(%s)"%argstrs
1639          else:
1640             raise NotImplementedError,"Asin_Symbol does not provide program code for format %s."%format
1641    
1642       def substitute(self,argvals):
1643          """
1644          assigns new values to symbols in the definition of the symbol.
1645          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1646    
1647          @param argvals: new values assigned to symbols
1648          @type argvals: C{dict} with keywords of type L{Symbol}.
1649          @return: result of the substitution process. Operations are executed as much as possible.
1650          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1651          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1652          """
1653          if argvals.has_key(self):
1654             arg=argvals[self]
1655             if self.isAppropriateValue(arg):
1656                return arg
1657             else:
1658                raise TypeError,"%s: new value is not appropriate."%str(self)
1659          else:
1660             arg=self.getSubstitutedArguments(argvals)[0]
1661             return asin(arg)
1662    
1663       def diff(self,arg):
1664          """
1665          differential of this object
1666    
1667          @param arg: the derivative is calculated with respect to arg
1668          @type arg: L{escript.Symbol}
1669          @return: derivative with respect to C{arg}
1670          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1671          """
1672          if arg==self:
1673             return identity(self.getShape())
1674          else:
1675             myarg=self.getArgument()[0]
1676             val=matchShape(1./sqrt(1.-myarg**2),self.getDifferentiatedArguments(arg)[0])
1677             return val[0]*val[1]
1678    
1679    def acos(arg):
1680       """
1681       returns inverse cosine of argument arg
1682    
1683       @param arg: argument
1684       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1685       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1686       @raises TypeError: if the type of the argument is not expected.
1687       """
1688       if isinstance(arg,numarray.NumArray):
1689          return numarray.arccos(arg)
1690       elif isinstance(arg,escript.Data):
1691          return arg._acos()
1692       elif isinstance(arg,float):
1693          return math.acos(arg)
1694       elif isinstance(arg,int):
1695          return math.acos(arg)
1696       elif isinstance(arg,Symbol):
1697          return Acos_Symbol(arg)
1698       else:
1699          raise TypeError,"acos: Unknown argument type."
1700    
1701    class Acos_Symbol(DependendSymbol):
1702       """
1703       L{Symbol} representing the result of the inverse cosine function
1704       """
1705       def __init__(self,arg):
1706          """
1707          initialization of acos L{Symbol} with argument arg
1708          @param arg: argument of function
1709          @type arg: typically L{Symbol}.
1710          """
1711          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1712    
1713       def getMyCode(self,argstrs,format="escript"):
1714          """
1715          returns a program code that can be used to evaluate the symbol.
1716    
1717          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1718          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1719          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1720          @type format: C{str}
1721          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1722          @rtype: C{str}
1723          @raise: NotImplementedError: if the requested format is not available
1724          """
1725          if isinstance(argstrs,list):
1726              argstrs=argstrs[0]
1727          if format=="escript" or format=="str"  or format=="text":
1728             return "acos(%s)"%argstrs
1729          else:
1730             raise NotImplementedError,"Acos_Symbol does not provide program code for format %s."%format
1731    
1732       def substitute(self,argvals):
1733          """
1734          assigns new values to symbols in the definition of the symbol.
1735          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1736    
1737          @param argvals: new values assigned to symbols
1738          @type argvals: C{dict} with keywords of type L{Symbol}.
1739          @return: result of the substitution process. Operations are executed as much as possible.
1740          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1741          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1742          """
1743          if argvals.has_key(self):
1744             arg=argvals[self]
1745             if self.isAppropriateValue(arg):
1746                return arg
1747             else:
1748                raise TypeError,"%s: new value is not appropriate."%str(self)
1749          else:
1750             arg=self.getSubstitutedArguments(argvals)[0]
1751             return acos(arg)
1752    
1753       def diff(self,arg):
1754          """
1755          differential of this object
1756    
1757          @param arg: the derivative is calculated with respect to arg
1758          @type arg: L{escript.Symbol}
1759          @return: derivative with respect to C{arg}
1760          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1761          """
1762          if arg==self:
1763             return identity(self.getShape())
1764          else:
1765             myarg=self.getArgument()[0]
1766             val=matchShape(-1./sqrt(1.-myarg**2),self.getDifferentiatedArguments(arg)[0])
1767             return val[0]*val[1]
1768    
1769    def atan(arg):
1770       """
1771       returns inverse tangent of argument arg
1772    
1773       @param arg: argument
1774       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1775       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1776       @raises TypeError: if the type of the argument is not expected.
1777       """
1778       if isinstance(arg,numarray.NumArray):
1779          return numarray.arctan(arg)
1780       elif isinstance(arg,escript.Data):
1781          return arg._atan()
1782       elif isinstance(arg,float):
1783          return math.atan(arg)
1784       elif isinstance(arg,int):
1785          return math.atan(arg)
1786       elif isinstance(arg,Symbol):
1787          return Atan_Symbol(arg)
1788       else:
1789          raise TypeError,"atan: Unknown argument type."
1790    
1791    class Atan_Symbol(DependendSymbol):
1792       """
1793       L{Symbol} representing the result of the inverse tangent function
1794       """
1795       def __init__(self,arg):
1796          """
1797          initialization of atan L{Symbol} with argument arg
1798          @param arg: argument of function
1799          @type arg: typically L{Symbol}.
1800          """
1801          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1802    
1803       def getMyCode(self,argstrs,format="escript"):
1804          """
1805          returns a program code that can be used to evaluate the symbol.
1806    
1807          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1808          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1809          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1810          @type format: C{str}
1811          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1812          @rtype: C{str}
1813          @raise: NotImplementedError: if the requested format is not available
1814          """
1815          if isinstance(argstrs,list):
1816              argstrs=argstrs[0]
1817          if format=="escript" or format=="str"  or format=="text":
1818             return "atan(%s)"%argstrs
1819          else:
1820             raise NotImplementedError,"Atan_Symbol does not provide program code for format %s."%format
1821    
1822       def substitute(self,argvals):
1823          """
1824          assigns new values to symbols in the definition of the symbol.
1825          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1826    
1827          @param argvals: new values assigned to symbols
1828          @type argvals: C{dict} with keywords of type L{Symbol}.
1829          @return: result of the substitution process. Operations are executed as much as possible.
1830          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1831          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1832          """
1833          if argvals.has_key(self):
1834             arg=argvals[self]
1835             if self.isAppropriateValue(arg):
1836                return arg
1837             else:
1838                raise TypeError,"%s: new value is not appropriate."%str(self)
1839          else:
1840             arg=self.getSubstitutedArguments(argvals)[0]
1841             return atan(arg)
1842    
1843       def diff(self,arg):
1844          """
1845          differential of this object
1846    
1847          @param arg: the derivative is calculated with respect to arg
1848          @type arg: L{escript.Symbol}
1849          @return: derivative with respect to C{arg}
1850          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1851          """
1852          if arg==self:
1853             return identity(self.getShape())
1854          else:
1855             myarg=self.getArgument()[0]
1856             val=matchShape(1./(1+myarg**2),self.getDifferentiatedArguments(arg)[0])
1857             return val[0]*val[1]
1858    
1859    def sinh(arg):
1860       """
1861       returns hyperbolic sine of argument arg
1862    
1863       @param arg: argument
1864       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1865       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1866       @raises TypeError: if the type of the argument is not expected.
1867       """
1868       if isinstance(arg,numarray.NumArray):
1869          return numarray.sinh(arg)
1870       elif isinstance(arg,escript.Data):
1871          return arg._sinh()
1872       elif isinstance(arg,float):
1873          return math.sinh(arg)
1874       elif isinstance(arg,int):
1875          return math.sinh(arg)
1876       elif isinstance(arg,Symbol):
1877          return Sinh_Symbol(arg)
1878       else:
1879          raise TypeError,"sinh: Unknown argument type."
1880    
1881    class Sinh_Symbol(DependendSymbol):
1882       """
1883       L{Symbol} representing the result of the hyperbolic sine function
1884       """
1885       def __init__(self,arg):
1886          """
1887          initialization of sinh L{Symbol} with argument arg
1888          @param arg: argument of function
1889          @type arg: typically L{Symbol}.
1890          """
1891          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1892    
1893       def getMyCode(self,argstrs,format="escript"):
1894          """
1895          returns a program code that can be used to evaluate the symbol.
1896    
1897          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1898          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1899          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1900          @type format: C{str}
1901          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1902          @rtype: C{str}
1903          @raise: NotImplementedError: if the requested format is not available
1904          """
1905          if isinstance(argstrs,list):
1906              argstrs=argstrs[0]
1907          if format=="escript" or format=="str"  or format=="text":
1908             return "sinh(%s)"%argstrs
1909          else:
1910             raise NotImplementedError,"Sinh_Symbol does not provide program code for format %s."%format
1911    
1912       def substitute(self,argvals):
1913          """
1914          assigns new values to symbols in the definition of the symbol.
1915          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
1916    
1917          @param argvals: new values assigned to symbols
1918          @type argvals: C{dict} with keywords of type L{Symbol}.
1919          @return: result of the substitution process. Operations are executed as much as possible.
1920          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
1921          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
1922          """
1923          if argvals.has_key(self):
1924             arg=argvals[self]
1925             if self.isAppropriateValue(arg):
1926                return arg
1927             else:
1928                raise TypeError,"%s: new value is not appropriate."%str(self)
1929          else:
1930             arg=self.getSubstitutedArguments(argvals)[0]
1931             return sinh(arg)
1932    
1933       def diff(self,arg):
1934          """
1935          differential of this object
1936    
1937          @param arg: the derivative is calculated with respect to arg
1938          @type arg: L{escript.Symbol}
1939          @return: derivative with respect to C{arg}
1940          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
1941          """
1942          if arg==self:
1943             return identity(self.getShape())
1944          else:
1945             myarg=self.getArgument()[0]
1946             val=matchShape(cosh(myarg),self.getDifferentiatedArguments(arg)[0])
1947             return val[0]*val[1]
1948    
1949    def cosh(arg):
1950       """
1951       returns hyperbolic cosine of argument arg
1952    
1953       @param arg: argument
1954       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
1955       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
1956       @raises TypeError: if the type of the argument is not expected.
1957       """
1958       if isinstance(arg,numarray.NumArray):
1959          return numarray.cosh(arg)
1960       elif isinstance(arg,escript.Data):
1961          return arg._cosh()
1962       elif isinstance(arg,float):
1963          return math.cosh(arg)
1964       elif isinstance(arg,int):
1965          return math.cosh(arg)
1966       elif isinstance(arg,Symbol):
1967          return Cosh_Symbol(arg)
1968       else:
1969          raise TypeError,"cosh: Unknown argument type."
1970    
1971    class Cosh_Symbol(DependendSymbol):
1972       """
1973       L{Symbol} representing the result of the hyperbolic cosine function
1974       """
1975       def __init__(self,arg):
1976          """
1977          initialization of cosh L{Symbol} with argument arg
1978          @param arg: argument of function
1979          @type arg: typically L{Symbol}.
1980          """
1981          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
1982    
1983       def getMyCode(self,argstrs,format="escript"):
1984          """
1985          returns a program code that can be used to evaluate the symbol.
1986    
1987          @param argstrs: gives for each argument a string representing the argument for the evaluation.
1988          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
1989          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
1990          @type format: C{str}
1991          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
1992          @rtype: C{str}
1993          @raise: NotImplementedError: if the requested format is not available
1994          """
1995          if isinstance(argstrs,list):
1996              argstrs=argstrs[0]
1997          if format=="escript" or format=="str"  or format=="text":
1998             return "cosh(%s)"%argstrs
1999          else:
2000             raise NotImplementedError,"Cosh_Symbol does not provide program code for format %s."%format
2001    
2002       def substitute(self,argvals):
2003          """
2004          assigns new values to symbols in the definition of the symbol.
2005          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2006    
2007          @param argvals: new values assigned to symbols
2008          @type argvals: C{dict} with keywords of type L{Symbol}.
2009          @return: result of the substitution process. Operations are executed as much as possible.
2010          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2011          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2012          """
2013          if argvals.has_key(self):
2014             arg=argvals[self]
2015             if self.isAppropriateValue(arg):
2016                return arg
2017             else:
2018                raise TypeError,"%s: new value is not appropriate."%str(self)
2019          else:
2020             arg=self.getSubstitutedArguments(argvals)[0]
2021             return cosh(arg)
2022    
2023       def diff(self,arg):
2024          """
2025          differential of this object
2026    
2027          @param arg: the derivative is calculated with respect to arg
2028          @type arg: L{escript.Symbol}
2029          @return: derivative with respect to C{arg}
2030          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2031          """
2032          if arg==self:
2033             return identity(self.getShape())
2034          else:
2035             myarg=self.getArgument()[0]
2036             val=matchShape(sinh(myarg),self.getDifferentiatedArguments(arg)[0])
2037             return val[0]*val[1]
2038    
2039    def tanh(arg):
2040       """
2041       returns hyperbolic tangent of argument arg
2042    
2043       @param arg: argument
2044       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2045       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2046       @raises TypeError: if the type of the argument is not expected.
2047       """
2048       if isinstance(arg,numarray.NumArray):
2049          return numarray.tanh(arg)
2050       elif isinstance(arg,escript.Data):
2051          return arg._tanh()
2052       elif isinstance(arg,float):
2053          return math.tanh(arg)
2054       elif isinstance(arg,int):
2055          return math.tanh(arg)
2056       elif isinstance(arg,Symbol):
2057          return Tanh_Symbol(arg)
2058       else:
2059          raise TypeError,"tanh: Unknown argument type."
2060    
2061    class Tanh_Symbol(DependendSymbol):
2062       """
2063       L{Symbol} representing the result of the hyperbolic tangent function
2064       """
2065       def __init__(self,arg):
2066          """
2067          initialization of tanh L{Symbol} with argument arg
2068          @param arg: argument of function
2069          @type arg: typically L{Symbol}.
2070          """
2071          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2072    
2073       def getMyCode(self,argstrs,format="escript"):
2074          """
2075          returns a program code that can be used to evaluate the symbol.
2076    
2077          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2078          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2079          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2080          @type format: C{str}
2081          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2082          @rtype: C{str}
2083          @raise: NotImplementedError: if the requested format is not available
2084          """
2085          if isinstance(argstrs,list):
2086              argstrs=argstrs[0]
2087          if format=="escript" or format=="str"  or format=="text":
2088             return "tanh(%s)"%argstrs
2089          else:
2090             raise NotImplementedError,"Tanh_Symbol does not provide program code for format %s."%format
2091    
2092       def substitute(self,argvals):
2093          """
2094          assigns new values to symbols in the definition of the symbol.
2095          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2096    
2097          @param argvals: new values assigned to symbols
2098          @type argvals: C{dict} with keywords of type L{Symbol}.
2099          @return: result of the substitution process. Operations are executed as much as possible.
2100          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2101          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2102          """
2103          if argvals.has_key(self):
2104             arg=argvals[self]
2105             if self.isAppropriateValue(arg):
2106                return arg
2107             else:
2108                raise TypeError,"%s: new value is not appropriate."%str(self)
2109          else:
2110             arg=self.getSubstitutedArguments(argvals)[0]
2111             return tanh(arg)
2112    
2113       def diff(self,arg):
2114          """
2115          differential of this object
2116    
2117          @param arg: the derivative is calculated with respect to arg
2118          @type arg: L{escript.Symbol}
2119          @return: derivative with respect to C{arg}
2120          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2121          """
2122          if arg==self:
2123             return identity(self.getShape())
2124          else:
2125             myarg=self.getArgument()[0]
2126             val=matchShape(1./cosh(myarg)**2,self.getDifferentiatedArguments(arg)[0])
2127             return val[0]*val[1]
2128    
2129    def asinh(arg):
2130       """
2131       returns inverse hyperbolic sine of argument arg
2132    
2133       @param arg: argument
2134       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2135       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2136       @raises TypeError: if the type of the argument is not expected.
2137       """
2138       if isinstance(arg,numarray.NumArray):
2139          return numarray.arcsinh(arg)
2140       elif isinstance(arg,escript.Data):
2141          return arg._asinh()
2142       elif isinstance(arg,float):
2143          return numarray.arcsinh(arg)
2144       elif isinstance(arg,int):
2145          return numarray.arcsinh(float(arg))
2146       elif isinstance(arg,Symbol):
2147          return Asinh_Symbol(arg)
2148       else:
2149          raise TypeError,"asinh: Unknown argument type."
2150    
2151    class Asinh_Symbol(DependendSymbol):
2152       """
2153       L{Symbol} representing the result of the inverse hyperbolic sine function
2154       """
2155       def __init__(self,arg):
2156          """
2157          initialization of asinh L{Symbol} with argument arg
2158          @param arg: argument of function
2159          @type arg: typically L{Symbol}.
2160          """
2161          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2162    
2163       def getMyCode(self,argstrs,format="escript"):
2164          """
2165          returns a program code that can be used to evaluate the symbol.
2166    
2167          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2168          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2169          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2170          @type format: C{str}
2171          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2172          @rtype: C{str}
2173          @raise: NotImplementedError: if the requested format is not available
2174          """
2175          if isinstance(argstrs,list):
2176              argstrs=argstrs[0]
2177          if format=="escript" or format=="str"  or format=="text":
2178             return "asinh(%s)"%argstrs
2179          else:
2180             raise NotImplementedError,"Asinh_Symbol does not provide program code for format %s."%format
2181    
2182       def substitute(self,argvals):
2183          """
2184          assigns new values to symbols in the definition of the symbol.
2185          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2186    
2187          @param argvals: new values assigned to symbols
2188          @type argvals: C{dict} with keywords of type L{Symbol}.
2189          @return: result of the substitution process. Operations are executed as much as possible.
2190          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2191          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2192          """
2193          if argvals.has_key(self):
2194             arg=argvals[self]
2195             if self.isAppropriateValue(arg):
2196                return arg
2197             else:
2198                raise TypeError,"%s: new value is not appropriate."%str(self)
2199          else:
2200             arg=self.getSubstitutedArguments(argvals)[0]
2201             return asinh(arg)
2202    
2203       def diff(self,arg):
2204          """
2205          differential of this object
2206    
2207          @param arg: the derivative is calculated with respect to arg
2208          @type arg: L{escript.Symbol}
2209          @return: derivative with respect to C{arg}
2210          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2211          """
2212          if arg==self:
2213             return identity(self.getShape())
2214          else:
2215             myarg=self.getArgument()[0]
2216             val=matchShape(1./sqrt(myarg**2+1),self.getDifferentiatedArguments(arg)[0])
2217             return val[0]*val[1]
2218    
2219    def acosh(arg):
2220       """
2221       returns inverse hyperolic cosine of argument arg
2222    
2223       @param arg: argument
2224       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2225       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2226       @raises TypeError: if the type of the argument is not expected.
2227       """
2228       if isinstance(arg,numarray.NumArray):
2229          return numarray.arccosh(arg)
2230       elif isinstance(arg,escript.Data):
2231          return arg._acosh()
2232       elif isinstance(arg,float):
2233          return numarray.arccosh(arg)
2234       elif isinstance(arg,int):
2235          return numarray.arccosh(float(arg))
2236       elif isinstance(arg,Symbol):
2237          return Acosh_Symbol(arg)
2238       else:
2239          raise TypeError,"acosh: Unknown argument type."
2240    
2241    class Acosh_Symbol(DependendSymbol):
2242       """
2243       L{Symbol} representing the result of the inverse hyperolic cosine function
2244       """
2245       def __init__(self,arg):
2246          """
2247          initialization of acosh L{Symbol} with argument arg
2248          @param arg: argument of function
2249          @type arg: typically L{Symbol}.
2250          """
2251          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2252    
2253       def getMyCode(self,argstrs,format="escript"):
2254          """
2255          returns a program code that can be used to evaluate the symbol.
2256    
2257          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2258          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2259          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2260          @type format: C{str}
2261          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2262          @rtype: C{str}
2263          @raise: NotImplementedError: if the requested format is not available
2264          """
2265          if isinstance(argstrs,list):
2266              argstrs=argstrs[0]
2267          if format=="escript" or format=="str"  or format=="text":
2268             return "acosh(%s)"%argstrs
2269          else:
2270             raise NotImplementedError,"Acosh_Symbol does not provide program code for format %s."%format
2271    
2272       def substitute(self,argvals):
2273          """
2274          assigns new values to symbols in the definition of the symbol.
2275          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2276    
2277          @param argvals: new values assigned to symbols
2278          @type argvals: C{dict} with keywords of type L{Symbol}.
2279          @return: result of the substitution process. Operations are executed as much as possible.
2280          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2281          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2282          """
2283          if argvals.has_key(self):
2284             arg=argvals[self]
2285             if self.isAppropriateValue(arg):
2286                return arg
2287             else:
2288                raise TypeError,"%s: new value is not appropriate."%str(self)
2289          else:
2290             arg=self.getSubstitutedArguments(argvals)[0]
2291             return acosh(arg)
2292    
2293       def diff(self,arg):
2294          """
2295          differential of this object
2296    
2297          @param arg: the derivative is calculated with respect to arg
2298          @type arg: L{escript.Symbol}
2299          @return: derivative with respect to C{arg}
2300          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2301          """
2302          if arg==self:
2303             return identity(self.getShape())
2304          else:
2305             myarg=self.getArgument()[0]
2306             val=matchShape(1./sqrt(myarg**2-1),self.getDifferentiatedArguments(arg)[0])
2307             return val[0]*val[1]
2308    
2309    def atanh(arg):
2310       """
2311       returns inverse hyperbolic tangent of argument arg
2312    
2313       @param arg: argument
2314       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2315       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2316       @raises TypeError: if the type of the argument is not expected.
2317       """
2318       if isinstance(arg,numarray.NumArray):
2319          return numarray.arctanh(arg)
2320       elif isinstance(arg,escript.Data):
2321          return arg._atanh()
2322       elif isinstance(arg,float):
2323          return numarray.arctanh(arg)
2324       elif isinstance(arg,int):
2325          return numarray.arctanh(float(arg))
2326       elif isinstance(arg,Symbol):
2327          return Atanh_Symbol(arg)
2328       else:
2329          raise TypeError,"atanh: Unknown argument type."
2330    
2331    class Atanh_Symbol(DependendSymbol):
2332       """
2333       L{Symbol} representing the result of the inverse hyperbolic tangent function
2334       """
2335       def __init__(self,arg):
2336          """
2337          initialization of atanh L{Symbol} with argument arg
2338          @param arg: argument of function
2339          @type arg: typically L{Symbol}.
2340          """
2341          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2342    
2343       def getMyCode(self,argstrs,format="escript"):
2344          """
2345          returns a program code that can be used to evaluate the symbol.
2346    
2347          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2348          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2349          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2350          @type format: C{str}
2351          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2352          @rtype: C{str}
2353          @raise: NotImplementedError: if the requested format is not available
2354          """
2355          if isinstance(argstrs,list):
2356              argstrs=argstrs[0]
2357          if format=="escript" or format=="str"  or format=="text":
2358             return "atanh(%s)"%argstrs
2359          else:
2360             raise NotImplementedError,"Atanh_Symbol does not provide program code for format %s."%format
2361    
2362       def substitute(self,argvals):
2363          """
2364          assigns new values to symbols in the definition of the symbol.
2365          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2366    
2367          @param argvals: new values assigned to symbols
2368          @type argvals: C{dict} with keywords of type L{Symbol}.
2369          @return: result of the substitution process. Operations are executed as much as possible.
2370          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2371          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2372          """
2373          if argvals.has_key(self):
2374             arg=argvals[self]
2375             if self.isAppropriateValue(arg):
2376                return arg
2377             else:
2378                raise TypeError,"%s: new value is not appropriate."%str(self)
2379          else:
2380             arg=self.getSubstitutedArguments(argvals)[0]
2381             return atanh(arg)
2382    
2383       def diff(self,arg):
2384          """
2385          differential of this object
2386    
2387          @param arg: the derivative is calculated with respect to arg
2388          @type arg: L{escript.Symbol}
2389          @return: derivative with respect to C{arg}
2390          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2391          """
2392          if arg==self:
2393             return identity(self.getShape())
2394          else:
2395             myarg=self.getArgument()[0]
2396             val=matchShape(1./(1.-myarg**2),self.getDifferentiatedArguments(arg)[0])
2397             return val[0]*val[1]
2398    
2399    def exp(arg):
2400       """
2401       returns exponential of argument arg
2402    
2403       @param arg: argument
2404       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2405       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2406       @raises TypeError: if the type of the argument is not expected.
2407       """
2408       if isinstance(arg,numarray.NumArray):
2409          return numarray.exp(arg)
2410       elif isinstance(arg,escript.Data):
2411          return arg._exp()
2412       elif isinstance(arg,float):
2413          return math.exp(arg)
2414       elif isinstance(arg,int):
2415          return math.exp(arg)
2416       elif isinstance(arg,Symbol):
2417          return Exp_Symbol(arg)
2418       else:
2419          raise TypeError,"exp: Unknown argument type."
2420    
2421    class Exp_Symbol(DependendSymbol):
2422       """
2423       L{Symbol} representing the result of the exponential function
2424       """
2425       def __init__(self,arg):
2426          """
2427          initialization of exp L{Symbol} with argument arg
2428          @param arg: argument of function
2429          @type arg: typically L{Symbol}.
2430          """
2431          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2432    
2433       def getMyCode(self,argstrs,format="escript"):
2434          """
2435          returns a program code that can be used to evaluate the symbol.
2436    
2437          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2438          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2439          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2440          @type format: C{str}
2441          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2442          @rtype: C{str}
2443          @raise: NotImplementedError: if the requested format is not available
2444          """
2445          if isinstance(argstrs,list):
2446              argstrs=argstrs[0]
2447          if format=="escript" or format=="str"  or format=="text":
2448             return "exp(%s)"%argstrs
2449          else:
2450             raise NotImplementedError,"Exp_Symbol does not provide program code for format %s."%format
2451    
2452       def substitute(self,argvals):
2453          """
2454          assigns new values to symbols in the definition of the symbol.
2455          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2456    
2457          @param argvals: new values assigned to symbols
2458          @type argvals: C{dict} with keywords of type L{Symbol}.
2459          @return: result of the substitution process. Operations are executed as much as possible.
2460          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2461          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2462          """
2463          if argvals.has_key(self):
2464             arg=argvals[self]
2465             if self.isAppropriateValue(arg):
2466                return arg
2467             else:
2468                raise TypeError,"%s: new value is not appropriate."%str(self)
2469          else:
2470             arg=self.getSubstitutedArguments(argvals)[0]
2471             return exp(arg)
2472    
2473       def diff(self,arg):
2474          """
2475          differential of this object
2476    
2477          @param arg: the derivative is calculated with respect to arg
2478          @type arg: L{escript.Symbol}
2479          @return: derivative with respect to C{arg}
2480          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2481          """
2482          if arg==self:
2483             return identity(self.getShape())
2484          else:
2485             myarg=self.getArgument()[0]
2486             val=matchShape(self,self.getDifferentiatedArguments(arg)[0])
2487             return val[0]*val[1]
2488    
2489    def sqrt(arg):
2490       """
2491       returns square root of argument arg
2492    
2493       @param arg: argument
2494       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2495       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2496       @raises TypeError: if the type of the argument is not expected.
2497       """
2498       if isinstance(arg,numarray.NumArray):
2499          return numarray.sqrt(arg)
2500       elif isinstance(arg,escript.Data):
2501          return arg._sqrt()
2502       elif isinstance(arg,float):
2503          return math.sqrt(arg)
2504       elif isinstance(arg,int):
2505          return math.sqrt(arg)
2506       elif isinstance(arg,Symbol):
2507          return Sqrt_Symbol(arg)
2508       else:
2509          raise TypeError,"sqrt: Unknown argument type."
2510    
2511    class Sqrt_Symbol(DependendSymbol):
2512       """
2513       L{Symbol} representing the result of the square root function
2514       """
2515       def __init__(self,arg):
2516          """
2517          initialization of sqrt L{Symbol} with argument arg
2518          @param arg: argument of function
2519          @type arg: typically L{Symbol}.
2520          """
2521          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2522    
2523       def getMyCode(self,argstrs,format="escript"):
2524          """
2525          returns a program code that can be used to evaluate the symbol.
2526    
2527          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2528          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2529          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2530          @type format: C{str}
2531          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2532          @rtype: C{str}
2533          @raise: NotImplementedError: if the requested format is not available
2534          """
2535          if isinstance(argstrs,list):
2536              argstrs=argstrs[0]
2537          if format=="escript" or format=="str"  or format=="text":
2538             return "sqrt(%s)"%argstrs
2539          else:
2540             raise NotImplementedError,"Sqrt_Symbol does not provide program code for format %s."%format
2541    
2542       def substitute(self,argvals):
2543          """
2544          assigns new values to symbols in the definition of the symbol.
2545          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2546    
2547          @param argvals: new values assigned to symbols
2548          @type argvals: C{dict} with keywords of type L{Symbol}.
2549          @return: result of the substitution process. Operations are executed as much as possible.
2550          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2551          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2552          """
2553          if argvals.has_key(self):
2554             arg=argvals[self]
2555             if self.isAppropriateValue(arg):
2556                return arg
2557             else:
2558                raise TypeError,"%s: new value is not appropriate."%str(self)
2559          else:
2560             arg=self.getSubstitutedArguments(argvals)[0]
2561             return sqrt(arg)
2562    
2563       def diff(self,arg):
2564          """
2565          differential of this object
2566    
2567          @param arg: the derivative is calculated with respect to arg
2568          @type arg: L{escript.Symbol}
2569          @return: derivative with respect to C{arg}
2570          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2571          """
2572          if arg==self:
2573             return identity(self.getShape())
2574          else:
2575             myarg=self.getArgument()[0]
2576             val=matchShape(0.5/self,self.getDifferentiatedArguments(arg)[0])
2577             return val[0]*val[1]
2578    
2579    def log(arg):
2580       """
2581       returns natural logarithm of argument arg
2582    
2583       @param arg: argument
2584       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2585       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2586       @raises TypeError: if the type of the argument is not expected.
2587       """
2588       if isinstance(arg,numarray.NumArray):
2589          return numarray.log(arg)
2590       elif isinstance(arg,escript.Data):
2591          return arg._log()
2592       elif isinstance(arg,float):
2593          return math.log(arg)
2594       elif isinstance(arg,int):
2595          return math.log(arg)
2596       elif isinstance(arg,Symbol):
2597          return Log_Symbol(arg)
2598       else:
2599          raise TypeError,"log: Unknown argument type."
2600    
2601    class Log_Symbol(DependendSymbol):
2602       """
2603       L{Symbol} representing the result of the natural logarithm function
2604       """
2605       def __init__(self,arg):
2606          """
2607          initialization of log L{Symbol} with argument arg
2608          @param arg: argument of function
2609          @type arg: typically L{Symbol}.
2610          """
2611          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2612    
2613       def getMyCode(self,argstrs,format="escript"):
2614          """
2615          returns a program code that can be used to evaluate the symbol.
2616    
2617          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2618          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2619          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2620          @type format: C{str}
2621          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2622          @rtype: C{str}
2623          @raise: NotImplementedError: if the requested format is not available
2624          """
2625          if isinstance(argstrs,list):
2626              argstrs=argstrs[0]
2627          if format=="escript" or format=="str"  or format=="text":
2628             return "log(%s)"%argstrs
2629          else:
2630             raise NotImplementedError,"Log_Symbol does not provide program code for format %s."%format
2631    
2632       def substitute(self,argvals):
2633          """
2634          assigns new values to symbols in the definition of the symbol.
2635          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2636    
2637          @param argvals: new values assigned to symbols
2638          @type argvals: C{dict} with keywords of type L{Symbol}.
2639          @return: result of the substitution process. Operations are executed as much as possible.
2640          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2641          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2642          """
2643          if argvals.has_key(self):
2644             arg=argvals[self]
2645             if self.isAppropriateValue(arg):
2646                return arg
2647             else:
2648                raise TypeError,"%s: new value is not appropriate."%str(self)
2649          else:
2650             arg=self.getSubstitutedArguments(argvals)[0]
2651             return log(arg)
2652    
2653       def diff(self,arg):
2654          """
2655          differential of this object
2656    
2657          @param arg: the derivative is calculated with respect to arg
2658          @type arg: L{escript.Symbol}
2659          @return: derivative with respect to C{arg}
2660          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2661          """
2662          if arg==self:
2663             return identity(self.getShape())
2664          else:
2665             myarg=self.getArgument()[0]
2666             val=matchShape(1./arg,self.getDifferentiatedArguments(arg)[0])
2667             return val[0]*val[1]
2668    
2669    def sign(arg):
2670       """
2671       returns sign of argument arg
2672    
2673       @param arg: argument
2674       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2675       @rtype:C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2676       @raises TypeError: if the type of the argument is not expected.
2677       """
2678       if isinstance(arg,numarray.NumArray):
2679          return wherePositive(arg)-whereNegative(arg)
2680       elif isinstance(arg,escript.Data):
2681          return arg._sign()
2682       elif isinstance(arg,float):
2683          if arg>0:
2684            return 1.
2685          elif arg<0:
2686            return -1.
2687          else:
2688            return 0.
2689       elif isinstance(arg,int):
2690          if float(arg)>0:
2691            return 1.
2692          elif float(arg)<0:
2693            return -1.
2694          else:
2695            return 0.
2696       elif isinstance(arg,Symbol):
2697          return wherePositive(arg)-whereNegative(arg)
2698       else:
2699          raise TypeError,"sign: Unknown argument type."
2700    
2701    class Abs_Symbol(DependendSymbol):
2702       """
2703       L{Symbol} representing the result of the absolute value function
2704       """
2705       def __init__(self,arg):
2706          """
2707          initialization of abs L{Symbol} with argument arg
2708          @param arg: argument of function
2709          @type arg: typically L{Symbol}.
2710          """
2711          DependendSymbol.__init__(self,args=[arg],shape=arg.getShape(),dim=arg.getDim())
2712    
2713       def getMyCode(self,argstrs,format="escript"):
2714          """
2715          returns a program code that can be used to evaluate the symbol.
2716    
2717          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2718          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2719          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2720          @type format: C{str}
2721          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2722          @rtype: C{str}
2723          @raise: NotImplementedError: if the requested format is not available
2724          """
2725          if isinstance(argstrs,list):
2726              argstrs=argstrs[0]
2727          if format=="escript" or format=="str"  or format=="text":
2728             return "abs(%s)"%argstrs
2729          else:
2730             raise NotImplementedError,"Abs_Symbol does not provide program code for format %s."%format
2731    
2732       def substitute(self,argvals):
2733          """
2734          assigns new values to symbols in the definition of the symbol.
2735          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2736    
2737          @param argvals: new values assigned to symbols
2738          @type argvals: C{dict} with keywords of type L{Symbol}.
2739          @return: result of the substitution process. Operations are executed as much as possible.
2740          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2741          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2742          """
2743          if argvals.has_key(self):
2744             arg=argvals[self]
2745             if self.isAppropriateValue(arg):
2746                return arg
2747             else:
2748                raise TypeError,"%s: new value is not appropriate."%str(self)
2749          else:
2750             arg=self.getSubstitutedArguments(argvals)[0]
2751             return abs(arg)
2752    
2753       def diff(self,arg):
2754          """
2755          differential of this object
2756    
2757          @param arg: the derivative is calculated with respect to arg
2758          @type arg: L{escript.Symbol}
2759          @return: derivative with respect to C{arg}
2760          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
2761          """
2762          if arg==self:
2763             return identity(self.getShape())
2764          else:
2765             myarg=self.getArgument()[0]
2766             val=matchShape(sign(myarg),self.getDifferentiatedArguments(arg)[0])
2767             return val[0]*val[1]
2768    
2769    def minval(arg):
2770       """
2771       returns minimum value over all components of arg at each data point
2772    
2773       @param arg: argument
2774       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2775       @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2776       @raises TypeError: if the type of the argument is not expected.
2777       """
2778       if isinstance(arg,numarray.NumArray):
2779          if arg.rank==0:
2780             return float(arg)
2781          else:
2782             return arg.min()
2783       elif isinstance(arg,escript.Data):
2784          return arg._minval()
2785       elif isinstance(arg,float):
2786          return arg
2787       elif isinstance(arg,int):
2788          return float(arg)
2789       elif isinstance(arg,Symbol):
2790          return Minval_Symbol(arg)
2791       else:
2792          raise TypeError,"minval: Unknown argument type."
2793    
2794    class Minval_Symbol(DependendSymbol):
2795       """
2796       L{Symbol} representing the result of the minimum value function
2797       """
2798       def __init__(self,arg):
2799          """
2800          initialization of minimum value L{Symbol} with argument arg
2801          @param arg: argument of function
2802          @type arg: typically L{Symbol}.
2803          """
2804          DependendSymbol.__init__(self,args=[arg],shape=(),dim=arg.getDim())
2805    
2806       def getMyCode(self,argstrs,format="escript"):
2807          """
2808          returns a program code that can be used to evaluate the symbol.
2809    
2810          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2811          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2812          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2813          @type format: C{str}
2814          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2815          @rtype: C{str}
2816          @raise: NotImplementedError: if the requested format is not available
2817          """
2818          if isinstance(argstrs,list):
2819              argstrs=argstrs[0]
2820          if format=="escript" or format=="str"  or format=="text":
2821             return "minval(%s)"%argstrs
2822          else:
2823             raise NotImplementedError,"Minval_Symbol does not provide program code for format %s."%format
2824    
2825       def substitute(self,argvals):
2826          """
2827          assigns new values to symbols in the definition of the symbol.
2828          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2829    
2830          @param argvals: new values assigned to symbols
2831          @type argvals: C{dict} with keywords of type L{Symbol}.
2832          @return: result of the substitution process. Operations are executed as much as possible.
2833          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2834          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2835          """
2836          if argvals.has_key(self):
2837             arg=argvals[self]
2838             if self.isAppropriateValue(arg):
2839                return arg
2840             else:
2841                raise TypeError,"%s: new value is not appropriate."%str(self)
2842          else:
2843             arg=self.getSubstitutedArguments(argvals)[0]
2844             return minval(arg)
2845    
2846    def maxval(arg):
2847       """
2848       returns maximum value over all components of arg at each data point
2849    
2850       @param arg: argument
2851       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2852       @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2853       @raises TypeError: if the type of the argument is not expected.
2854       """
2855       if isinstance(arg,numarray.NumArray):
2856          if arg.rank==0:
2857             return float(arg)
2858          else:
2859             return arg.max()
2860       elif isinstance(arg,escript.Data):
2861          return arg._maxval()
2862       elif isinstance(arg,float):
2863          return arg
2864       elif isinstance(arg,int):
2865          return float(arg)
2866       elif isinstance(arg,Symbol):
2867          return Maxval_Symbol(arg)
2868       else:
2869          raise TypeError,"maxval: Unknown argument type."
2870    
2871    class Maxval_Symbol(DependendSymbol):
2872       """
2873       L{Symbol} representing the result of the maximum value function
2874       """
2875       def __init__(self,arg):
2876          """
2877          initialization of maximum value L{Symbol} with argument arg
2878          @param arg: argument of function
2879          @type arg: typically L{Symbol}.
2880          """
2881          DependendSymbol.__init__(self,args=[arg],shape=(),dim=arg.getDim())
2882    
2883       def getMyCode(self,argstrs,format="escript"):
2884          """
2885          returns a program code that can be used to evaluate the symbol.
2886    
2887          @param argstrs: gives for each argument a string representing the argument for the evaluation.
2888          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
2889          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
2890          @type format: C{str}
2891          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
2892          @rtype: C{str}
2893          @raise: NotImplementedError: if the requested format is not available
2894          """
2895          if isinstance(argstrs,list):
2896              argstrs=argstrs[0]
2897          if format=="escript" or format=="str"  or format=="text":
2898             return "maxval(%s)"%argstrs
2899          else:
2900             raise NotImplementedError,"Maxval_Symbol does not provide program code for format %s."%format
2901    
2902       def substitute(self,argvals):
2903          """
2904          assigns new values to symbols in the definition of the symbol.
2905          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
2906    
2907          @param argvals: new values assigned to symbols
2908          @type argvals: C{dict} with keywords of type L{Symbol}.
2909          @return: result of the substitution process. Operations are executed as much as possible.
2910          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
2911          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
2912          """
2913          if argvals.has_key(self):
2914             arg=argvals[self]
2915             if self.isAppropriateValue(arg):
2916                return arg
2917             else:
2918                raise TypeError,"%s: new value is not appropriate."%str(self)
2919          else:
2920             arg=self.getSubstitutedArguments(argvals)[0]
2921             return maxval(arg)
2922    
2923    def length(arg):
2924       """
2925       returns length/Euclidean norm of argument arg at each data point
2926    
2927       @param arg: argument
2928       @type arg: C{float}, L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2929       @rtype:C{float}, L{escript.Data}, L{Symbol} depending on the type of arg.
2930       """
2931       return sqrt(inner(arg,arg))
2932    
2933    def trace(arg,axis_offset=0):
2934       """
2935       returns the trace of arg which the sum of arg[k,k] over k.
2936    
2937       @param arg: argument
2938       @type arg: L{escript.Data}, L{Symbol}, L{numarray.NumArray}.
2939       @param axis_offset: 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
2940                      axis_offset and axis_offset+1 must be equal.
2941       @type axis_offset: C{int}
2942       @return: trace of arg. The rank of the returned object is minus 2 of the rank of arg.
2943       @rtype: L{escript.Data}, L{Symbol}, L{numarray.NumArray} depending on the type of arg.
2944       """
2945       if isinstance(arg,numarray.NumArray):
2946          sh=arg.shape
2947          if len(sh)<2:
2948            raise ValueError,"trace: rank of argument must be greater than 1"
2949          if axis_offset<0 or axis_offset>len(sh)-2:
2950            raise ValueError,"trace: axis_offset must be between 0 and %s"%len(sh)-2
2951          s1=1
2952          for i in range(axis_offset): s1*=sh[i]
2953          s2=1
2954          for i in range(axis_offset+2,len(sh)): s2*=sh[i]
2955          if not sh[axis_offset] == sh[axis_offset+1]:
2956            raise ValueError,"trace: dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
2957          arg_reshaped=numarray.reshape(arg,(s1,sh[axis_offset],sh[axis_offset],s2))
2958          out=numarray.zeros([s1,s2],numarray.Float64)
2959          for i1 in range(s1):
2960            for i2 in range(s2):
2961                for j in range(sh[axis_offset]): out[i1,i2]+=arg_reshaped[i1,j,j,i2]
2962          out.resize(sh[:axis_offset]+sh[axis_offset+2:])
2963          return out
2964       elif isinstance(arg,escript.Data):
2965          return escript_trace(arg,axis_offset)
2966       elif isinstance(arg,float):
2967          raise TypeError,"trace: illegal argument type float."
2968       elif isinstance(arg,int):
2969          raise TypeError,"trace: illegal argument type int."
2970       elif isinstance(arg,Symbol):
2971          return Trace_Symbol(arg,axis_offset)
2972       else:
2973          raise TypeError,"trace: Unknown argument type."
2974    
2975    def escript_trace(arg,axis_offset): # this should be escript._trace
2976          "arg si a Data objects!!!"
2977          if arg.getRank()<2:
2978            raise ValueError,"escript_trace: rank of argument must be greater than 1"
2979          if axis_offset<0 or axis_offset>arg.getRank()-2:
2980            raise ValueError,"escript_trace: axis_offset must be between 0 and %s"%arg.getRank()-2
2981          s=list(arg.getShape())        
2982          if not s[axis_offset] == s[axis_offset+1]:
2983            raise ValueError,"escript_trace: dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
2984          out=escript.Data(0.,tuple(s[0:axis_offset]+s[axis_offset+2:]),arg.getFunctionSpace())
2985          if arg.getRank()==2:
2986             for i0 in range(s[0]):
2987                out+=arg[i0,i0]
2988          elif arg.getRank()==3:
2989             if axis_offset==0:
2990                for i0 in range(s[0]):
2991                      for i2 in range(s[2]):
2992                             out[i2]+=arg[i0,i0,i2]
2993             elif axis_offset==1:
2994                for i0 in range(s[0]):
2995                   for i1 in range(s[1]):
2996                             out[i0]+=arg[i0,i1,i1]
2997          elif arg.getRank()==4:
2998             if axis_offset==0:
2999                for i0 in range(s[0]):
3000                      for i2 in range(s[2]):
3001                         for i3 in range(s[3]):
3002                             out[i2,i3]+=arg[i0,i0,i2,i3]
3003             elif axis_offset==1:
3004                for i0 in range(s[0]):
3005                   for i1 in range(s[1]):
3006                         for i3 in range(s[3]):
3007                             out[i0,i3]+=arg[i0,i1,i1,i3]
3008             elif axis_offset==2:
3009                for i0 in range(s[0]):
3010                   for i1 in range(s[1]):
3011                      for i2 in range(s[2]):
3012                             out[i0,i1]+=arg[i0,i1,i2,i2]
3013          return out
3014    class Trace_Symbol(DependendSymbol):
3015       """
3016       L{Symbol} representing the result of the trace function
3017       """
3018       def __init__(self,arg,axis_offset=0):
3019          """
3020          initialization of trace L{Symbol} with argument arg
3021          @param arg: argument of function
3022          @type arg: L{Symbol}.
3023          @param axis_offset: 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
3024                      axis_offset and axis_offset+1 must be equal.
3025          @type axis_offset: C{int}
3026          """
3027          if arg.getRank()<2:
3028            raise ValueError,"Trace_Symbol: rank of argument must be greater than 1"
3029          if axis_offset<0 or axis_offset>arg.getRank()-2:
3030            raise ValueError,"Trace_Symbol: axis_offset must be between 0 and %s"%arg.getRank()-2
3031          s=list(arg.getShape())        
3032          if not s[axis_offset] == s[axis_offset+1]:
3033            raise ValueError,"Trace_Symbol: dimensions of component %s and %s must match."%(axis_offset.axis_offset+1)
3034          super(Trace_Symbol,self).__init__(args=[arg,axis_offset],shape=tuple(s[0:axis_offset]+s[axis_offset+2:]),dim=arg.getDim())
3035    
3036       def getMyCode(self,argstrs,format="escript"):
3037          """
3038          returns a program code that can be used to evaluate the symbol.
3039    
3040          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3041          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3042          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3043          @type format: C{str}
3044          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3045          @rtype: C{str}
3046          @raise: NotImplementedError: if the requested format is not available
3047          """
3048          if format=="escript" or format=="str"  or format=="text":
3049             return "trace(%s,axis_offset=%s)"%(argstrs[0],argstrs[1])
3050          else:
3051             raise NotImplementedError,"Trace_Symbol does not provide program code for format %s."%format
3052    
3053       def substitute(self,argvals):
3054          """
3055          assigns new values to symbols in the definition of the symbol.
3056          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3057    
3058          @param argvals: new values assigned to symbols
3059          @type argvals: C{dict} with keywords of type L{Symbol}.
3060          @return: result of the substitution process. Operations are executed as much as possible.
3061          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3062          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3063          """
3064          if argvals.has_key(self):
3065             arg=argvals[self]
3066             if self.isAppropriateValue(arg):
3067                return arg
3068             else:
3069                raise TypeError,"%s: new value is not appropriate."%str(self)
3070          else:
3071             arg=self.getSubstitutedArguments(argvals)
3072             return trace(arg[0],axis_offset=arg[1])
3073    
3074       def diff(self,arg):
3075          """
3076          differential of this object
3077    
3078          @param arg: the derivative is calculated with respect to arg
3079          @type arg: L{escript.Symbol}
3080          @return: derivative with respect to C{arg}
3081          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3082          """
3083          if arg==self:
3084             return identity(self.getShape())
3085          else:
3086             return trace(self.getDifferentiatedArguments(arg)[0],axis_offset=self.getArgument()[1])
3087    
3088    def transpose(arg,axis_offset=None):
3089       """
3090       returns the transpose of arg by swaping the first axis_offset and the last rank-axis_offset components.
3091    
3092       @param arg: argument
3093       @type arg: L{escript.Data}, L{Symbol}, L{numarray.NumArray}, C{float}, C{int}
3094       @param axis_offset: the first axis_offset components are swapped with rest. If C{axis_offset} must be non-negative and less or equal the rank of arg.
3095                           if axis_offset is not present C{int(r/2)} where r is the rank of arg is used.
3096       @type axis_offset: C{int}
3097       @return: transpose of arg
3098       @rtype: L{escript.Data}, L{Symbol}, L{numarray.NumArray},C{float}, C{int} depending on the type of arg.
3099       """
3100       if isinstance(arg,numarray.NumArray):
3101          if axis_offset==None: axis_offset=int(arg.rank/2)
3102          return numarray.transpose(arg,axes=range(axis_offset,arg.rank)+range(0,axis_offset))
3103       elif isinstance(arg,escript.Data):
3104          if axis_offset==None: axis_offset=int(arg.getRank()/2)
3105          return escript_transpose(arg,axis_offset)
3106       elif isinstance(arg,float):
3107          if not ( axis_offset==0 or axis_offset==None):
3108            raise ValueError,"transpose: axis_offset must be 0 for float argument"
3109          return arg
3110       elif isinstance(arg,int):
3111          if not ( axis_offset==0 or axis_offset==None):
3112            raise ValueError,"transpose: axis_offset must be 0 for int argument"
3113          return float(arg)
3114       elif isinstance(arg,Symbol):
3115          if axis_offset==None: axis_offset=int(arg.getRank()/2)
3116          return Transpose_Symbol(arg,axis_offset)
3117       else:
3118          raise TypeError,"transpose: Unknown argument type."
3119    
3120    def escript_transpose(arg,axis_offset): # this should be escript._transpose
3121          "arg si a Data objects!!!"
3122          r=arg.getRank()
3123          if axis_offset<0 or axis_offset>r:
3124            raise ValueError,"escript_transpose: axis_offset must be between 0 and %s"%r
3125          s=arg.getShape()
3126          s_out=s[axis_offset:]+s[:axis_offset]
3127          out=escript.Data(0.,s_out,arg.getFunctionSpace())
3128          if r==4:
3129             if axis_offset==1:
3130                for i0 in range(s_out[0]):
3131                   for i1 in range(s_out[1]):
3132                      for i2 in range(s_out[2]):
3133                         for i3 in range(s_out[3]):
3134                             out[i0,i1,i2,i3]=arg[i3,i0,i1,i2]
3135             elif axis_offset==2:
3136                for i0 in range(s_out[0]):
3137                   for i1 in range(s_out[1]):
3138                      for i2 in range(s_out[2]):
3139                         for i3 in range(s_out[3]):
3140                             out[i0,i1,i2,i3]=arg[i2,i3,i0,i1]
3141             elif axis_offset==3:
3142                for i0 in range(s_out[0]):
3143                   for i1 in range(s_out[1]):
3144                      for i2 in range(s_out[2]):
3145                         for i3 in range(s_out[3]):
3146                             out[i0,i1,i2,i3]=arg[i1,i2,i3,i0]
3147             else:
3148                for i0 in range(s_out[0]):
3149                   for i1 in range(s_out[1]):
3150                      for i2 in range(s_out[2]):
3151                         for i3 in range(s_out[3]):
3152                             out[i0,i1,i2,i3]=arg[i0,i1,i2,i3]
3153          elif r==3:
3154             if axis_offset==1:
3155                for i0 in range(s_out[0]):
3156                   for i1 in range(s_out[1]):
3157                      for i2 in range(s_out[2]):
3158                             out[i0,i1,i2]=arg[i2,i0,i1]
3159             elif axis_offset==2:
3160                for i0 in range(s_out[0]):
3161                   for i1 in range(s_out[1]):
3162                      for i2 in range(s_out[2]):
3163                             out[i0,i1,i2]=arg[i1,i2,i0]
3164             else:
3165                for i0 in range(s_out[0]):
3166                   for i1 in range(s_out[1]):
3167                      for i2 in range(s_out[2]):
3168                             out[i0,i1,i2]=arg[i0,i1,i2]
3169          elif r==2:
3170             if axis_offset==1:
3171                for i0 in range(s_out[0]):
3172                   for i1 in range(s_out[1]):
3173                             out[i0,i1]=arg[i1,i0]
3174             else:
3175                for i0 in range(s_out[0]):
3176                   for i1 in range(s_out[1]):
3177                             out[i0,i1]=arg[i0,i1]
3178          elif r==1:
3179              for i0 in range(s_out[0]):
3180                   out[i0]=arg[i0]
3181          elif r==0:
3182                 out=arg+0.
3183          return out
3184    class Transpose_Symbol(DependendSymbol):
3185       """
3186       L{Symbol} representing the result of the transpose function
3187       """
3188       def __init__(self,arg,axis_offset=None):
3189          """
3190          initialization of transpose L{Symbol} with argument arg
3191    
3192          @param arg: argument of function
3193          @type arg: L{Symbol}.
3194           @param axis_offset: the first axis_offset components are swapped with rest. If C{axis_offset} must be non-negative and less or equal the rank of arg.
3195                           if axis_offset is not present C{int(r/2)} where r is the rank of arg is used.
3196          @type axis_offset: C{int}
3197          """
3198          if axis_offset==None: axis_offset=int(arg.getRank()/2)
3199          if axis_offset<0 or axis_offset>arg.getRank():
3200            raise ValueError,"escript_transpose: axis_offset must be between 0 and %s"%r
3201          s=arg.getShape()
3202          super(Transpose_Symbol,self).__init__(args=[arg,axis_offset],shape=s[axis_offset:]+s[:axis_offset],dim=arg.getDim())
3203    
3204       def getMyCode(self,argstrs,format="escript"):
3205          """
3206          returns a program code that can be used to evaluate the symbol.
3207    
3208          @param argstrs: gives for each argument a string representing the argument for the evaluation.
3209          @type argstrs: C{str} or a C{list} of length 1 of C{str}.
3210          @param format: specifies the format to be used. At the moment only "escript" ,"text" and "str" are supported.
3211          @type format: C{str}
3212          @return: a piece of program code which can be used to evaluate the expression assuming the values for the arguments are available.
3213          @rtype: C{str}
3214          @raise: NotImplementedError: if the requested format is not available
3215          """
3216          if format=="escript" or format=="str"  or format=="text":
3217             return "transpose(%s,axis_offset=%s)"%(argstrs[0],argstrs[1])
3218          else:
3219             raise NotImplementedError,"Transpose_Symbol does not provide program code for format %s."%format
3220    
3221       def substitute(self,argvals):
3222          """
3223          assigns new values to symbols in the definition of the symbol.
3224          The method replaces the L{Symbol} u by argvals[u] in the expression defining this object.
3225    
3226          @param argvals: new values assigned to symbols
3227          @type argvals: C{dict} with keywords of type L{Symbol}.
3228          @return: result of the substitution process. Operations are executed as much as possible.
3229          @rtype: L{escript.Symbol}, C{float}, L{escript.Data}, L{numarray.NumArray} depending on the degree of substitution
3230          @raise TypeError: if a value for a L{Symbol} cannot be substituted.
3231          """
3232          if argvals.has_key(self):
3233             arg=argvals[self]
3234             if self.isAppropriateValue(arg):
3235                return arg
3236             else:
3237                raise TypeError,"%s: new value is not appropriate."%str(self)
3238          else:
3239             arg=self.getSubstitutedArguments(argvals)
3240             return transpose(arg[0],axis_offset=arg[1])
3241    
3242       def diff(self,arg):
3243          """
3244          differential of this object
3245    
3246          @param arg: the derivative is calculated with respect to arg
3247          @type arg: L{escript.Symbol}
3248          @return: derivative with respect to C{arg}
3249          @rtype: typically L{Symbol} but other types such as C{float}, L{escript.Data}, L{numarray.NumArray}  are possible.
3250          """
3251          if arg==self:
3252             return identity(self.getShape())
3253          else:
3254             return transpose(self.getDifferentiatedArguments(arg)[0],axis_offset=self.getArgument()[1])
3255    def symmetric(arg):
3256      """      """
3257      @brief      returns the symmetric part of the square matrix arg. This is (arg+transpose(arg))/2
3258    
3259      @param arg      @param arg: square matrix. Must have rank 2 or 4 and be square.
3260        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3261        @return: symmetric part of arg
3262        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
3263      """      """
3264      if isinstance(arg,escript.Data):      if isinstance(arg,numarray.NumArray):
3265         return arg.cos()        if arg.rank==2:
3266            if not (arg.shape[0]==arg.shape[1]):
3267               raise ValueError,"symmetric: argument must be square."
3268          elif arg.rank==4:
3269            if not (arg.shape[0]==arg.shape[2] and arg.shape[1]==arg.shape[3]):
3270               raise ValueError,"symmetric: argument must be square."
3271          else:
3272            raise ValueError,"symmetric: rank 2 or 4 is required."
3273          return (arg+transpose(arg))/2
3274        elif isinstance(arg,escript.Data):
3275          return escript_symmetric(arg)
3276        elif isinstance(arg,float):
3277          return arg
3278        elif isinstance(arg,int):
3279          return float(arg)
3280        elif isinstance(arg,Symbol):
3281          if arg.getRank()==2:
3282            if not (arg.getShape()[0]==arg.getShape()[1]):
3283               raise ValueError,"symmetric: argument must be square."
3284          elif arg.getRank()==4:
3285            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3286               raise ValueError,"symmetric: argument must be square."
3287          else:
3288            raise ValueError,"symmetric: rank 2 or 4 is required."
3289          return (arg+transpose(arg))/2
3290      else:      else:
3291         return numarray.cos(arg)        raise TypeError,"symmetric: Unknown argument type."
3292    
3293  def maxval(arg):  def escript_symmetric(arg): # this should be implemented in c++
3294          if arg.getRank()==2:
3295            if not (arg.getShape()[0]==arg.getShape()[1]):
3296               raise ValueError,"escript_symmetric: argument must be square."
3297            out=escript.Data(0.,arg.getShape(),arg.getFunctionSpace())
3298            for i0 in range(arg.getShape()[0]):
3299               for i1 in range(arg.getShape()[1]):
3300                  out[i0,i1]=(arg[i0,i1]+arg[i1,i0])/2.
3301          elif arg.getRank()==4:
3302            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3303               raise ValueError,"escript_symmetric: argument must be square."
3304            out=escript.Data(0.,arg.getShape(),arg.getFunctionSpace())
3305            for i0 in range(arg.getShape()[0]):
3306               for i1 in range(arg.getShape()[1]):
3307                  for i2 in range(arg.getShape()[2]):
3308                     for i3 in range(arg.getShape()[3]):
3309                         out[i0,i1,i2,i3]=(arg[i0,i1,i2,i3]+arg[i2,i3,i0,i1])/2.
3310          else:
3311            raise ValueError,"escript_symmetric: rank 2 or 4 is required."
3312          return out
3313    
3314    def nonsymmetric(arg):
3315      """      """
3316      @brief      returns the nonsymmetric part of the square matrix arg. This is (arg-transpose(arg))/2
3317    
3318      @param arg      @param arg: square matrix. Must have rank 2 or 4 and be square.
3319        @type arg: L{numarray.NumArray}, L{escript.Data}, L{Symbol}
3320        @return: nonsymmetric part of arg
3321        @rtype: L{numarray.NumArray}, L{escript.Data}, L{Symbol} depending on the input
3322      """      """
3323      if isinstance(arg,escript.Data):      if isinstance(arg,numarray.NumArray):
3324         return arg.maxval()        if arg.rank==2:
3325            if not (arg.shape[0]==arg.shape[1]):
3326               raise ValueError,"nonsymmetric: argument must be square."
3327          elif arg.rank==4:
3328            if not (arg.shape[0]==arg.shape[2] and arg.shape[1]==arg.shape[3]):
3329               raise ValueError,"nonsymmetric: argument must be square."
3330          else:
3331            raise ValueError,"nonsymmetric: rank 2 or 4 is required."
3332          return (arg-transpose(arg))/2
3333        elif isinstance(arg,escript.Data):
3334          return escript_nonsymmetric(arg)
3335        elif isinstance(arg,float):
3336          return arg
3337        elif isinstance(arg,int):
3338          return float(arg)
3339        elif isinstance(arg,Symbol):
3340          if arg.getRank()==2:
3341            if not (arg.getShape()[0]==arg.getShape()[1]):
3342               raise ValueError,"nonsymmetric: argument must be square."
3343          elif arg.getRank()==4:
3344            if not (arg.getShape()[0]==arg.getShape()[2] and arg.getShape()[1]==arg.getShape()[3]):
3345               raise ValueError,"nonsymmetric: argument must be square."
3346          else:
3347            raise ValueError,"nonsymmetric: rank 2 or 4 is required."
3348          return (arg-transpose(arg))/2
3349      else:      else:
3350         return arg.max()        raise TypeError,"nonsymmetric: Unknown argument type."
3351    
3352  def minval(arg):  def escript_nonsymmetric(arg): # this should be implemented in c++
3353          if arg.getRank()==2:
3354            if not (arg.getShape()[0]==arg.getShape()[1]):
3355               raise ValueError,"escript_nonsymmetric: argument must be square."
3356            out=escript.Data(0.,arg.getShape(),arg.getFunctionSpace())
3357            for i0 in range(arg.getShape()[0]):
3358               for i1 in range(arg.getShape()[1]):
3359                  out[i0,i1]=(arg[i0,i1]-arg[i1,i0])/2.
3360          elif arg.getRank()==4: