[PYTHON MATRIX-SIG] Updated function list
Konrad HINSEN
hinsenk@ere.umontreal.ca
Sun, 25 Feb 1996 14:22:17 -0500
Thanks to a rainy weekend, I have completed the proposed function list.
It now covers all operations that I can think of (which doesn't mean
that nothing is missing!). The descriptions of array() and reshape()
have changed with respect to the first draft. Enjoy and ciriticize!
---------------------------------------------------------------------------
This list is an attempt to put some order into the array functions. I
have made a list of functions that should be available and proposed a
name for each function. Lines that describe new or modified operations
are marked with a '+'. Comments are at the end of the file.
As I have suggested before, it is useful to have an additional module
defining abbreviations and/or more meaningful names for often-used
combinations. These should not be part of the standard module to
reduce name space pollution and learning overhead, but they should
nevertheless be part of the standard distribution to ensure
uniformity. I'll include recommended abbreviations in the respective
sections.
The following principles have been applied:
1) Every operation that is not closely related to the implementation
details of the array object is implemented as a function, not
as a method. This is necessary because methods don't work on
scalars, which also serve as rank-0 arrays. It is also necessary
to allow the implementation of every operation in C or Python
without changing the interface. And this decision also ensures
consistency with additional modules (linear algebra etc.); such
modules can only add functions, not methods.
2) There are no functions with a variable number of array arguments.
Such functions are difficult to apply to a number of arrays that
is calculated in the code rather than constant. Instead, a sequence
of arrays is passed. This has the added advantage that a higher-dimensional
array can be passed, implying an iteration over the first axis.
Constructors
============
1) Construct an array from an arbitrary sequence object or a function:
+ array(object, shape=None, type=None)
Creates an array from the given object with a given
shape and type. If the shape and/or the type are not
specified, they are inferred from the data.
The object providing the data can be of the following
types:
1) An object with a method "toArray". This method
is called, and the resulting array is converted to
the given shape and type (if not None).
2) A nested list or tuple, or an array. The shape is inferred
from the nesting, the type from the data, unless
space and/or type are specified.
3) A file object. The file is read as a rank-2
array. If a shape and/or type is specified, the result is
changed to that shape/type.
4) A callable object. This object is called once for
each array element with the indices as arguments.
In reshaping to an explicitly given shape, the input
is reused from the beginning if necessary. It is also
allowed not to use all the input data.
This constructor replaces the current constructors
array(), copyAndReshape(), and fromFunction(), and
makes special cases like zeros() so easy that they
don't have to be standard functions any more.
Note: the current fromFunction() behaves in a subtly
different way: it calls the function only once with
index arrays as arguments. This is equivalent if
the function consists only of operators and ofuncs,
and of course much more efficient, but for most
interesting functions it doesn't work at all. In
this case I prefer generality to efficiency.
2) Construct a rank-1 array of equally spaced points:
+ arrayRange(x1, x2=None, x3=None)
Currently: function arange() does exactly the same.
Abbreviations:
zeros(n1,n2,...) equals array(0,(n1,n2,...))
ones(n1,n2,...) equals array(1,(n1,n2,...))
unitMatrix(n) equals array([1]+n*[0], (n,n))
arange equals arrayRange
Structural array operations
===========================
1) Selecting subarrays
1.1) Selecting on a "raster" along each coordinate:
Done by indexing with integers and slices
1.2) Selecting arbitrary elements specified by an array of indices i:
+ take(a, i, axis=0)
Conditions: i must be of an integer array type,
minimum.reduce(ravel(i)) >= 0,
maximum.reduce(ravel(i)) < a.shape[axis]
Currently: method a.take(i, axis=0), i can be of any type, but
must have rank 1.
1.3) Selecting the diagonal, i.e. those values where the indices
along two axes are equal (see comment 1):
+ diagonal(a, axis1=0, axis2=1)
Conditions: axis1 < len(a.shape), axis2 < len(a.shape),
a.shape[axis1] == a.shape[axis2]
2) Rearranging array elements
2.1) Changing the shape
2.1.1) General reshaping:
+ reshape(a,shape, copy=1)
shape can be of any non-nested sequence type
Condition: multiply.reduce(shape) == multiply.reduce(a.shape)
If copy==1, the returned array is a new copy, otherwise
a reference to a.
Currently: method a.reshape(length1, length2, ...)
This syntax does not provide an easy way to specify
a non-constant shape. Also, the current behaviour is
like copy=0, which can lead to unpleasant surprises
(since otherwise only indexing returns references)
and limits reshaping to contiguous arrays.
2.1.2) Reshaping to a rank-1 array:
+ ravel(a, copy=1)
equivalent to reshape(a,(multiply.reduce(a.shape),), copy)
2.1.3) Combining a group of axes into one axis (see comment 2):
+ a[i1, i2, ......, i3, i4] (see comment 3)
The double ellipsis works similar to the single one, but
contracts all the axes covered into a single axis.
2.1.4) Adding an axis of length 1:
a[..., NewAxis, ...]
2.2) Transposition:
+ transpose(a, axes)
axes is a non-nested sequence of non-negative integers
with maximum.reduce(axes) < len(a.shape)
Currently: method a.transpose(axes) does the same, but
there is a bug that sometimes gives wrong results if
an axis is used more than once in the list. It also
insists that len(axes) == len(a.shape), although there
is no need for this restriction.
3) Replicating and combining array elements
3.1) Replicating elements along an axis:
+ repeat(a, n, axis=0)
n is a non-nested integer sequence with len(n) == a.shape[axis]
Each element in a is repeated as often as indicated by
the corresponding element in n. The length of the result
along the specified axis is add.reduce(n).
Currently: function compress(n, a, axis=0) is a special case
limited to boolean (0/1) elements in n.
3.2) Concatenation of arrays along an axis:
+ concatenate((a1,a2,a3,...), axis=0)
Condition: all arrays must have the same shape for the
remaining axes.
Currently: method a.concat(a1,a2,a3,...) works only along
first axis and is difficult to apply to a variable number
of arguments.
3.3) Concatenation of arrays along a new axis:
Can be done by combining concatenate() and indexing with NewAxis.
(see comment 4)
Abbreviation:
reverse(a) = a[::-1]
Arithmetic and logical operations
=================================
1) Binary elementwise arithmetic operators
+ - * / % **
as functions: add, subtract, multiply, divide, remainder, power
2) Standard "mathematical" functions
arccos, arccosh, arcsin, arcsinh, arctan, arctanh, cos, cosh, exp,
log, log10, sin, sinh, sqrt, tan, tanh
3) Other elementwise arithmetic functions
maximum, minimum, clip
+ conjugate
Currently, conjugate() is implemented as a method both on
complex objects and array objects. It should be available
as a function, because then it can be applied to other
number objects as well. There are many algorithms that work
for both real and complex numbers providing that conjugates
are used here and there. It ought to be possible to implement
these algorithms and make them work for float and complex
objects as well as array objects. The method we have now
doesn't work on real and integer objects.
4) Comparison functions
+ equal(a,b), notEqual(a,b), greater(a,b), greaterEqual(a,b),
+ less(a,b), lessEqual(a,b)
Currently all these are methods. Apart from introducing
an arbitrary asymmetry, this makes it impossible to write
code that works both for scalars and arrays.
5) Logical and boolean operators
+ logicalAnd(a,b), logicalOr(a,b), logicalXor(a,b), logicalNot(a),
+ booleanAnd(a,b), booleanOr(a,b), booleanXor(a,b), booleanNot(a,b)
Currently there are equivalent methods with "inversed" names.
I think these names are more intuitive. About the problem of
methods, see above. For some strange reason there is currently
no logical xor operation.
6) Non-elementwise operations
6.1) Attribute functions of binary operator functions (see 1 and 5)
reduce, accumulate,
+ outer
+ Note: reduce and accumulate should return the neutral
+ element of their base operation, reshaped to agree with the
+ remaining axes, if used along an axis of length zero.
+ Example:
+ add.reduce(array(42, (2,0,1)), 1)
+ should be equal to
+ array(0, (2,1))
6.2) Dot/matrix product
+ dot(a, b, axis_a = -1, axis_b = 0)
equivalent to:
if axis_a < 0: axis_a = len(a.shape)+axis_a
if axis_b >= 0: axis_b = axis_b + len(a.shape)
add.reduce(diagonal(multiply.outer(a,b), axis_a, axis_b))
(but of course done more efficiently)
Currently: method matrixMultiply, restricted to rank-2
arrays.
Redundant method: nonzero
I can't see what this operation is good for. It seems
to be equivalent to ravel(notEqual(a,0)), which is not
an extremely frequent combination.
Abbreviations:
sum = add.reduce
cumsum = add.accumulate
product = multiply.reduce
cumproduct = multiply.accumulate
allTrue = logicalAnd.reduce
someTrue = logicalOr.reduce
ptp(a) = maximum.reduce(a)-minimum.reduce(a)
Looping, mapping, filtering, etc.
=================================
1) Looping over array elements
for element in a: ...
Looping over anything else than the first axis can be achieved
by indexing.
2) Mapping a function on one or more arrays
+ arrayMap(function, (a1, a2, a3, ..., an), ranks = None)
function must be a callable object. If ranks is specified, it
must be a sequence of integers of length n; if ranks is None,
all ranks are assumed to be one.
The frames of all arrays with respect to the rank specification
must match.
3) Filtering
+ arrayFilter(function, a, axis=0)
collects all elements along the given axis for which the
supplied function returns non-0.
See also repeat() in section "structural operations".
4) Conditional selection
+ choose(a, (b0, b2, b3, ..., bn))
a must be an array of integers between 0 and n. The result
has the same shape as a with elements selected from b0...bn
according to the value of the corresponding element of a.
Currently: method a.choose(b0, ..., bn)
Passing the arrays b0...bn as a tuple has the advantage that
this tuple can be constructed arbitrarily in the code.
5) Sorting
5.1) Sort an array
+ sortArray(a, function=lambda x:x, axis=0)
Sorts the elements along the specified axis according to the return
value of the supplied function.
Note: the default function makes sense only for rank-1 arrays.
5.2) Sort index
+ sortIndex(a, function=lambda x:x, axis=0)
returns an index array i such that take(a,i) == sort(a, function, axis)
Redundant function:
where(condition, x, y) equals choose(condition, (x,y))
Other operations
================
1) Type casts
1.1) General cast
+ a.asType(typecode)
1.2) Special cast functions that also work for scalars
+ integerArray(a)
+ floatArray(a)
+ complexArray(a)
2) Implementation dependent functions:
2.1) Obtain the amount of bytes occupied by each element
a.itemsize()
2.2) Return a byte-swapped copy
+ a.byteswapped()
2.3) Return the element type code
a.typecode()
2.4) Return the element type
+ a.elementType()
returns the type object of the Python scalar object that results
if an element is extracted. Returns None for a general object
(typecode 'O').
2.5) Checks if an array is contiguous
+ a.isContiguous()
2.1) Return data as a Python string
a.toString()
Redundant function:
a.copy()
Copying should be done with copy.copy(), just as for any other object.
Comments
========
1) APL uses a special case of transposition for selecting a
diagonal, but this is often confusing.
2) In J this is done by using ravel() with reduced rank, but
as long as we don't have functions with bounded ranks,
we need a special function.
3) I tried to find a construction that does not require
yet another syntax and this is the best I could think
of. But I am not particularly attached to it, so feel
free to think of something better.
4) APL does this with fractional indices, which is probably
the weirdest feature in APL.
-------------------------------------------------------------------------------
Konrad Hinsen | E-Mail: hinsenk@ere.umontreal.ca
Departement de chimie | Tel.: +1-514-343-6111 ext. 3953
Universite de Montreal | Fax: +1-514-343-7586
C.P. 6128, succ. Centre-Ville | Deutsch/Esperanto/English/Nederlands/
Montreal (QC) H3C 3J7 | Francais (phase experimentale)
-------------------------------------------------------------------------------
=================
MATRIX-SIG - SIG on Matrix Math for Python
send messages to: matrix-sig@python.org
administrivia to: matrix-sig-request@python.org
=================