[PYTHON MATRIX-SIG] Final matrix object renaming and packaging

[James Hugunin]

   Date: Tue, 16 Jan 1996 12:10:28 -0500
   From: hinsenk@ERE.UMontreal.CA (Hinsen Konrad)

      1) Konrad's numeric patches to the python core (incorporated into the
      working version of python by Guido) will be required and included with
      the distribution.  These will be the only patches to the python core
      required.

   There has been a slight modification in the meantime (there will be
   only 'j' for imaginary constants, not 'i') and two corrections that
   perhaps not everybody has received. I'll prepare a new set of patch
   files for the beta release.

Great, just upload the patch files to "the site" when you're ready.


      This PyArray C type will not implement automatic type-coercion (unlike the
      current implementation).  The reason for this is that I have decided
      type-coercion is a major pain for arrays of floats (as opposed to
      doubles) and I happen to use arrays of floats for most of my work.  If
      somebody can give me a good suggestion for how to keep automatic
      type-coercion and have Matrix_f(1,2,3)*1.2 == Matrix_f(1.2,2.4,3.6)
      then I might consider changing this decision.  See later note on
      Array.py for an alternative.

   One way to solve your "float" problem would be to write
      Matrix_f(1,2,3)*Matrix_f(1.2),
   i.e. cast 1.2 explicitly to float. I don't see how you can avoid
   that cast anyway; without type coercion, Matrix_f(1,2,3)*1.2 would
   be an error.

I'm willing to be a little bit sloppy with some of the details of
types, and force the type of a scalar argument to correspond to the
type of the matrix.  This is essentially what I'm doing anyway when I
say Matrix_f(1.2) [Note, I'll have another post to come about
constructor function].

   I'd prefer to have type coercion everywhere for consistency.  If the
   "high-level" classes Array and Matrix are really fast enough (about
   which I have doubts, see below), then it doesn't matter that much, but
   in general I'd give top priority to consistency and put the burden of
   more complicated coding on those who want top speed (which of course
   must be possible).

We seem to just have a difference of opinion here, I'll wait and see
if more opinions come in.


      3) Two dynamically linkable modules called "umathmodule.c", and "ieee_umathmodule.c"

   I don't think that using "ieee" in the name is a good idea. After all,
   there is no guarantee that this will really use IEEE arithmetic (a
   Vax, for example, doesn't have IEEE arithmetic hardware). Why not
   just "fast_umath"?

Good point, actually I was originally considering "slow_umath" and
"fast_umath", but I knew you wouldn't like that.  So, the new proposed
name is "umath" and "fast_umath".


      4) Two python objects, "Array.py" and "Matrix.py"

      Array is essentially a python binding around the underlying C type,
      and this will also provide for automatic type-coercion and will
      generally assume that it is only working with arrays of type long,
      double, and complex double (the three types of arrays that are
      equivalent to python objects).  In my initial tests of this object on
      LLNL's simple benchmark, I found that the performance was only 5%
      slower than using the C object directly.

   Nevertheless the Python wrapper could cause a significant overhead
   for small arrays, which will also be used heavily. For example,
   I am planning to change my vector class (3d vectors for geometry
   etc.) to use an array of length 3 as its internal representation.
   Do you have any data on the speed penalty for such small arrays?

No data yet on that, if you want to send me a nice simple "real"
benchmark code for 3d vectors I'd be happy to add it to my bag of
tricks.  Right now everything is being optimized for LLNL style
physics code because they sent me such a nice benchmark to run.

      Matrix will inherit almost everything from Array, however it will be
      limited to 2 or fewer dimensions, and m1*m2 where m1 and m2 are
      matrices will perform matrix style multiplication.  If the linear
      algebra people would like, other changes can be made (ie. ~m1 ==
      m1.transpose(), ...).  Based on the experiments with Array, the
      performance penalty for this approach should be minimal.

   I suppose some discussion will be necessary about the details of
   this class, which until now didn't exist. Given that many people
   seem to prefer a distinction between row and column vectors, it
   would be a good idea to include separate constructors for them.

Now that I've brought things out into python classes, I think I'll
probably just let somebody who's into linear algebra code worry about
those sorts of refinements.

      5) A standard library "Numeric.py" which will be the standard way of
      importing multiarray, Array, Matrix, umath, etc.  It will include the
      inverted trig functions ("sec", "csc", "arcsec", ...) as well as most
      of the standard functions currently in Matrix.py

   I am not sure it is a good idea to have such an all-encompassing
   module. It would contain a rather arbitrary selection of things,
   i.e. those which exist at the time it is introduced. I expect
   that other numerical modules will appear later. So it would be
   better to group things according to functions or applications.

I disagree here.  I find that there is a certain core of
functions that are very nice to be able to grab with a single import
statement.  I would hate to have to start my typical script with
import umath, UMathMore, Matrix, Array, ArrayConstructors,
ArrayUtilities, ...



      9?) A "numericmodule.c" which contains reasonably efficient fft,
      matrix inversion, convolution, random numbers, eigenvalues and
      filtering functions stolen from existing C code so that the package
      can be viewed as a "complete" numerical computing system?

   Again this should be several modules: linear algebra, random
   numbers, fft/convolution. We won't achieve anything "complete"
   anyway, so let's not pretend we do.

You're obviously right here.  Part of what I was hoping to do with
this was to create a simple module so that people (even those without
a FORTRAN compiler) could take the inverse of a matrix in a reasonable
amount of time, and possibly to set up some sort of API, so that when
somebody does write the ultimate linear algebra module they will be
likely include a function called "inverse" so that I could then
trivially plug in this new more powerful system.


=================
MATRIX-SIG  - SIG on Matrix Math for Python

send messages to: matrix-sig@python.org
administrivia to: matrix-sig-request@python.org
=================