[PYTHON MATRIX-SIG] Coercion problems solved?

[Alan Watson]

I'm evaluating Numerical Python for astronomical imaging applications.
These tend to involve largish arrays. For example, state-of-the-art
visible light detectors are 2k by 4k, so even with 32-bit floats a
single image is 32 MB.

I read the discussion in the mailing list archive on the coercion
rules between arrays and scalars with some disappointment. I
appreciate that having (scalar * vector)[i] == scalar * vector[i] is a
nice property, but I saw my options as either eating the cache and
page faults, buying a computer with a bigger cache and yet more
memory, or having to wrap every scalar in a function to convert it to
a rank-0 32-bit float array. None of these prospects filled me with
unbounded joy.

Then I remembered that Python is not an assembler with structs but is
an object-oriented language. A few minutes yielded:

    import UserArray
    class array32(UserArray.UserArray):
        def __init__(self, sequence): 
            UserArray.UserArray.__init__(self, sequence, typecode='f')
        def __coerce__(self,other): 
            return (self, UserArray.UserArray(other, typecode='f'))

(I've been using Python for only two days, so this may not be quite
perfect. Also, the UserArray classes don't appear to support a copy
argument.) This allows me to implement my own `principle of least
surprise' for coercion, namely that array elements are always 32-bit
floats. I think this idea can be extended to other, more complicated,
schemes.

I've not seen this mentioned in the discussion; maybe it was too
obvious or maybe I missed it while scanning the archive. I think it
just about eliminates the problem of the current rules for those of us
who prefer to work with the smaller types.

Regards,

Alan Watson

_______________
MATRIX-SIG  - SIG on Matrix Math for Python

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