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

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

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

revision 102 by jgs, Wed Dec 15 07:08:39 2004 UTC revision 122 by jgs, Thu Jun 9 05:38:05 2005 UTC
# Line 49  TIFF=7 Line 49  TIFF=7
49  OPENINVENTOR=8  OPENINVENTOR=8
50  RENDERMAN=9  RENDERMAN=9
51  PNM=10  PNM=10
52    
53  #  #
54  # wrapper for various functions: if the argument has attribute the function name  # wrapper for various functions: if the argument has attribute the function name
55  # as an argument it calls the correspong methods. Otherwise the coresponsing numarray  # as an argument it calls the correspong methods. Otherwise the coresponsing numarray
# Line 64  def grad(arg,where=None): Line 65  def grad(arg,where=None):
65      @param where: FunctionSpace in which the gradient will be. If None Function(dom) where dom is the      @param where: FunctionSpace in which the gradient will be. If None Function(dom) where dom is the
66                    domain of the Data object arg.                    domain of the Data object arg.
67      """      """
68      if where==None:      if isinstance(arg,escript.Data):
69         return arg.grad()         if where==None:
70              return arg.grad()
71           else:
72              return arg.grad(where)
73      else:      else:
74         return arg.grad(where)         return arg*0.
75    
76  def integrate(arg):  def integrate(arg,what=None):
77      """      """
78      @brief return the integral if the function represented by Data object arg over its domain.      @brief return the integral if the function represented by Data object arg over its domain.
79    
80      @param arg      @param arg
81      """      """
82      return arg.integrate()      if not what==None:
83           arg2=escript.Data(arg,what)
84        else:
85           arg2=arg
86        if arg2.getRank()==0:
87            return arg2.integrate()[0]
88        else:
89            return arg2.integrate()
90    
91  def interpolate(arg,where):  def interpolate(arg,where):
92      """      """
# Line 84  def interpolate(arg,where): Line 95  def interpolate(arg,where):
95      @param arg      @param arg
96      @param where      @param where
97      """      """
98      return arg.interpolate(where)      if isinstance(arg,escript.Data):
99           return arg.interpolate(where)
100        else:
101           return arg
102    
103  # functions returning Data objects:  # functions returning Data objects:
104    
# Line 95  def transpose(arg,axis=None): Line 109  def transpose(arg,axis=None):
109      @param arg      @param arg
110      """      """
111      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
112           # hack for transpose
113           r=arg.getRank()
114           if r!=2: raise ValueError,"Tranpose only avalaible for rank 2 objects"
115           s=arg.getShape()
116           out=escript.Data(0.,(s[1],s[0]),arg.getFunctionSpace())
117           for i in range(s[0]):
118              for j in range(s[1]):
119                 out[j,i]=arg[i,j]
120           return out
121           # end hack for transpose
122         if axis==None: axis=arg.getRank()/2         if axis==None: axis=arg.getRank()/2
123         return arg.transpose(axis)         return arg.transpose(axis)
124      else:      else:
# Line 108  def trace(arg): Line 132  def trace(arg):
132      @param arg      @param arg
133      """      """
134      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
135           # hack for trace
136           r=arg.getRank()
137           if r!=2: raise ValueError,"trace only avalaible for rank 2 objects"
138           s=arg.getShape()
139           out=escript.Scalar(0,arg.getFunctionSpace())
140           for i in range(min(s)):
141                 out+=arg[i,i]
142           return out
143           # end hack for trace
144         return arg.trace()         return arg.trace()
145      else:      else:
146         return numarray.trace(arg)         return numarray.trace(arg)
# Line 175  def maxval(arg): Line 208  def maxval(arg):
208      """      """
209      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
210         return arg.maxval()         return arg.maxval()
211        elif isinstance(arg,float) or isinstance(arg,int):
212           return arg
213      else:      else:
214         return arg.max()         return arg.max()
215    
# Line 186  def minval(arg): Line 221  def minval(arg):
221      """      """
222      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
223         return arg.minval()         return arg.minval()
224        elif isinstance(arg,float) or isinstance(arg,int):
225           return arg
226      else:      else:
227         return arg.max()         return arg.min()
228    
229  def length(arg):  def length(arg):
230      """      """
# Line 196  def length(arg): Line 233  def length(arg):
233      @param arg      @param arg
234      """      """
235      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
236         return arg.length()         if arg.isEmpty(): return escript.Data()
237           if arg.getRank()==0:
238              return abs(arg)
239           elif arg.getRank()==1:
240              sum=escript.Scalar(0,arg.getFunctionSpace())
241              for i in range(arg.getShape()[0]):
242                 sum+=arg[i]**2
243              return sqrt(sum)
244           elif arg.getRank()==2:
245              sum=escript.Scalar(0,arg.getFunctionSpace())
246              for i in range(arg.getShape()[0]):
247                 for j in range(arg.getShape()[1]):
248                    sum+=arg[i,j]**2
249              return sqrt(sum)
250           elif arg.getRank()==3:
251              sum=escript.Scalar(0,arg.getFunctionSpace())
252              for i in range(arg.getShape()[0]):
253                 for j in range(arg.getShape()[1]):
254                    for k in range(arg.getShape()[2]):
255                       sum+=arg[i,j,k]**2
256              return sqrt(sum)
257           elif arg.getRank()==4:
258              sum=escript.Scalar(0,arg.getFunctionSpace())
259              for i in range(arg.getShape()[0]):
260                 for j in range(arg.getShape()[1]):
261                    for k in range(arg.getShape()[2]):
262                       for l in range(arg.getShape()[3]):
263                          sum+=arg[i,j,k,l]**2
264              return sqrt(sum)
265           else:
266              raise SystemError,"length is not been implemented yet"
267           # return arg.length()
268      else:      else:
269         return sqrt((arg**2).sum())         return sqrt((arg**2).sum())
270    
271    def deviator(arg):
272        """
273        @brief
274    
275        @param arg1
276        """
277        if isinstance(arg,escript.Data):
278            shape=arg.getShape()
279        else:
280            shape=arg.shape
281        if len(shape)!=2:
282              raise ValueError,"Deviator requires rank 2 object"
283        if shape[0]!=shape[1]:
284              raise ValueError,"Deviator requires a square matrix"
285        return arg-1./(shape[0]*1.)*trace(arg)*kronecker(shape[0])
286    
287    def inner(arg1,arg2):
288        """
289        @brief
290    
291        @param arg1, arg2
292        """
293        sum=escript.Scalar(0,arg1.getFunctionSpace())
294        if arg.getRank()==0:
295              return arg1*arg2
296        elif arg.getRank()==1:
297             sum=escript.Scalar(0,arg.getFunctionSpace())
298             for i in range(arg.getShape()[0]):
299                sum+=arg1[i]*arg2[i]
300        elif arg.getRank()==2:
301            sum=escript.Scalar(0,arg.getFunctionSpace())
302            for i in range(arg.getShape()[0]):
303               for j in range(arg.getShape()[1]):
304                  sum+=arg1[i,j]*arg2[i,j]
305        elif arg.getRank()==3:
306            sum=escript.Scalar(0,arg.getFunctionSpace())
307            for i in range(arg.getShape()[0]):
308                for j in range(arg.getShape()[1]):
309                   for k in range(arg.getShape()[2]):
310                      sum+=arg1[i,j,k]*arg2[i,j,k]
311        elif arg.getRank()==4:
312            sum=escript.Scalar(0,arg.getFunctionSpace())
313            for i in range(arg.getShape()[0]):
314               for j in range(arg.getShape()[1]):
315                  for k in range(arg.getShape()[2]):
316                     for l in range(arg.getShape()[3]):
317                        sum+=arg1[i,j,k,l]*arg2[i,j,k,l]
318        else:
319              raise SystemError,"inner is not been implemented yet"
320        return sum
321    
322  def sign(arg):  def sign(arg):
323      """      """
324      @brief      @brief
# Line 229  def sup(arg): Line 348  def sup(arg):
348      """      """
349      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
350         return arg.sup()         return arg.sup()
351        elif isinstance(arg,float) or isinstance(arg,int):
352           return arg
353      else:      else:
354         return arg.max()         return arg.max()
355    
# Line 240  def inf(arg): Line 361  def inf(arg):
361      """      """
362      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
363         return arg.inf()         return arg.inf()
364        elif isinstance(arg,float) or isinstance(arg,int):
365           return arg
366      else:      else:
367         return arg.min()         return arg.min()
368    
# Line 249  def L2(arg): Line 372  def L2(arg):
372    
373      @param arg      @param arg
374      """      """
375      return arg.L2()      if isinstance(arg,escript.Data):
376           return arg.L2()
377        elif isinstance(arg,float) or isinstance(arg,int):
378           return abs(arg)
379        else:
380           return numarry.sqrt(dot(arg,arg))
381    
382  def Lsup(arg):  def Lsup(arg):
383      """      """
# Line 259  def Lsup(arg): Line 387  def Lsup(arg):
387      """      """
388      if isinstance(arg,escript.Data):      if isinstance(arg,escript.Data):
389         return arg.Lsup()         return arg.Lsup()
390        elif isinstance(arg,float) or isinstance(arg,int):
391           return abs(arg)
392      else:      else:
393         return arg.max(numarray.abs(arg))         return max(numarray.abs(arg))
394    
395    def Linf(arg):
396        """
397        @brief
398    
399        @param arg
400        """
401        if isinstance(arg,escript.Data):
402           return arg.Linf()
403        elif isinstance(arg,float) or isinstance(arg,int):
404           return abs(arg)
405        else:
406           return min(numarray.abs(arg))
407    
408  def dot(arg1,arg2):  def dot(arg1,arg2):
409      """      """
# Line 274  def dot(arg1,arg2): Line 417  def dot(arg1,arg2):
417         return arg2.dot(arg1)         return arg2.dot(arg1)
418      else:      else:
419         return numarray.dot(arg1,arg2)         return numarray.dot(arg1,arg2)
420    
421    def kronecker(d):
422       if hasattr(d,"getDim"):
423          return numarray.identity(d.getDim())
424       else:
425          return numarray.identity(d)
426    
427    def unit(i,d):
428       """
429       @brief return a unit vector of dimension d with nonzero index i
430       @param d dimension
431       @param i index
432       """
433       e = numarray.zeros((d,))
434       e[i] = 1.0
435       return e
436    
437    #
438    # $Log$
439    # Revision 1.10  2005/06/09 05:37:59  jgs
440    # Merge of development branch back to main trunk on 2005-06-09
441    #
442    # Revision 1.2.2.14  2005/05/20 04:05:23  gross
443    # some work on a darcy flow started
444    #
445    # Revision 1.2.2.13  2005/03/16 05:17:58  matt
446    # Implemented unit(idx, dim) to create cartesian unit basis vectors to
447    # complement kronecker(dim) function.
448    #
449    # Revision 1.2.2.12  2005/03/10 08:14:37  matt
450    # Added non-member Linf utility function to complement Data::Linf().
451    #
452    # Revision 1.2.2.11  2005/02/17 05:53:25  gross
453    # some bug in saveDX fixed: in fact the bug was in
454    # DataC/getDataPointShape
455    #
456    # Revision 1.2.2.10  2005/01/11 04:59:36  gross
457    # automatic interpolation in integrate switched off
458    #
459    # Revision 1.2.2.9  2005/01/11 03:38:13  gross
460    # Bug in Data.integrate() fixed for the case of rank 0. The problem is not totallly resolved as the method should return a scalar rather than a numarray object in the case of rank 0. This problem is fixed by the util.integrate wrapper.
461    #
462    # Revision 1.2.2.8  2005/01/05 04:21:41  gross
463    # FunctionSpace checking/matchig in slicing added
464    #
465    # Revision 1.2.2.7  2004/12/29 05:29:59  gross
466    # AdvectivePDE successfully tested for Peclet number 1000000. there is still a problem with setValue and Data()
467    #
468    # Revision 1.2.2.6  2004/12/24 06:05:41  gross
469    # some changes in linearPDEs to add AdevectivePDE
470    #
471    # Revision 1.2.2.5  2004/12/17 00:06:53  gross
472    # mk sets ESYS_ROOT is undefined
473    #
474    # Revision 1.2.2.4  2004/12/07 03:19:51  gross
475    # options for GMRES and PRES20 added
476    #
477    # Revision 1.2.2.3  2004/12/06 04:55:18  gross
478    # function wraper extended
479    #
480    # Revision 1.2.2.2  2004/11/22 05:44:07  gross
481    # a few more unitary functions have been added but not implemented in Data yet
482    #
483    # Revision 1.2.2.1  2004/11/12 06:58:15  gross
484    # a lot of changes to get the linearPDE class running: most important change is that there is no matrix format exposed to the user anymore. the format is chosen by the Domain according to the solver and symmetry
485    #
486    # Revision 1.2  2004/10/27 00:23:36  jgs
487    # fixed minor syntax error
488    #
489    # Revision 1.1.1.1  2004/10/26 06:53:56  jgs
490    # initial import of project esys2
491    #
492    # Revision 1.1.2.3  2004/10/26 06:43:48  jgs
493    # committing Lutz's and Paul's changes to brach jgs
494    #
495    # Revision 1.1.4.1  2004/10/20 05:32:51  cochrane
496    # Added incomplete Doxygen comments to files, or merely put the docstrings that already exist into Doxygen form.
497    #
498    # Revision 1.1  2004/08/05 03:58:27  gross
499    # Bug in Assemble_NodeCoordinates fixed
500    #
501    #

Legend:
Removed from v.102  
changed lines
  Added in v.122

  ViewVC Help
Powered by ViewVC 1.1.26