Re: [Numpy-discussion] Generalized ufuncs?
Robert Kern [EMAIL PROTECTED] wrote in message
news:[EMAIL PROTECTED]...
On Sun, Aug 17, 2008 at 21:55, Anne Archibald
[EMAIL PROTECTED] wrote:
2008/8/17 Robert Kern [EMAIL PROTECTED]:
I suggested that we move it to a branch for the time being so we
can
play with it and come up with examples of its use. If you have
examples that you have already written, I would love to see them.
I,
for one, am amenable to seeing this in 1.2.0, but only if we push
back
the release by at least a few weeks.
This is not a worked example, but this is exactly what is needed to
make possible the arrays of matrices functions that were discussed
some time ago. For example, users frequently want to do something
like
multiply an array of matrices by an array of matrices; that is not
currently feasible without a very large temporary array (or a loop).
But I think you were looking for examples of code using the
interface,
to see whether it worked well.
I'll take what I can get. In order of increasing utility:
1. Descriptions of use cases (as above).
2. Mockups of the Python code demonstrating the use case (e.g.
nonfunctional Python code showing a potential generalized ufunc
with inputs and outputs).
3. The C code implementing the actual functionality for the use
case.
--
Robert Kern
I have come to believe that the whole world is an enigma, a harmless
enigma that is made terrible by our own mad attempt to interpret it as
though it had an underlying truth.
-- Umberto Eco
Please find an example on inner1d in the following.
1. One use case for an inner function is described on
http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions and
http://thread.gmane.org/gmane.comp.python.numeric.general/17694.
(Another one would be the array of matrices usage mentioned above; we
have called this dot2d with signature (m,n),(n,p)-(m,p) on
http://scipy.org/scipy/numpy/ticket/887: here the matrix multiplication
would occur with respect to the last two dimensions.)
2. The mockup python code would be:
from numpy import *
N = 10
a = random.randn(3, 5, N)
b = random.randn(5, N)
# standard inner function
... inner(a, b).shape
(3, 5, 5)
# new ufunc inner1d with signature (i),(i)-(), satisfying
GeneralLoopingFunctions use case
... inner1d(a, b).shape
(3, 5)
3. Concrete implementation of inner1d in C:
/*
* This implements the function out = inner1d(in1, in2) with
* out[K] = sum_i { in1[K, i] * in2[K, i] }
* and multi-index K, as described on
* http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions
* and on http://projects.scipy.org/scipy/numpy/ticket/887.
*/
static void
DOUBLE_inner1d(char **args, intp *dimensions, intp *steps, void *func)
{
/* Standard ufunc loop length and strides. */
intp dn = dimensions[0];
intp s0 = steps[0];
intp s1 = steps[1];
intp s2 = steps[2];
intp n;
/* Additional loop length and strides for dimension i in
elementary function. */
intp di = dimensions[1];
intp i_s1 = steps[3];
intp i_s2 = steps[4];
intp i;
/* Outer loop: equivalent to standard ufuncs */
for (n = 0; n dn; n++, args[0] += s0, args[1] += s1, args[2] +=
s2) {
char *ip1 = args[0], *ip2 = args[1], *op = args[2];
/* Implement elementary function: out = sum_i { in1[i] * in2[i]
} */
npy_double sum = 0;
for (i = 0; i di; i++) {
sum += (*(npy_double *)ip1) * (*(npy_double *)ip2);
ip1 += i_s1;
ip2 += i_s2;
}
*(npy_double *)op = sum;
}
}
/* Actually create the ufunc object */
static PyUFuncGenericFunction inner1d_functions[] = { DOUBLE_inner1d };
static void *inner1d_data[] = { (void *)NULL };
static char inner1d_sigs[] = { PyArray_DOUBLE, PyArray_DOUBLE,
PyArray_DOUBLE };
static void
addInner1d(PyObject *dictionary) {
PyObject *f = PyUFunc_FromFuncAndDataAndSignature(inner1d_functions,
inner1d_data, inner1d_sigs, 1,
2, 1, PyUFunc_None, inner1d,
inner on the last dimension and
broadcast on the other dimensions,
0,
(i),(i)-());
PyDict_SetItemString(dictionary, inner1d, f);
}
___
Numpy-discussion mailing list
Numpy-discussion@scipy.org
http://projects.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Generalized ufuncs?
Charles R Harris [EMAIL PROTECTED] wrote in message news:[EMAIL PROTECTED]... On Sun, Aug 17, 2008 at 7:56 PM, Stéfan van der Walt [EMAIL PROTECTED]wrote: 2008/8/17 Robert Kern [EMAIL PROTECTED]: I suggested that we move it to a branch for the time being so we can play with it and come up with examples of its use. That branch is here: http://[EMAIL PROTECTED]/svn/numpy/branches/gen_ufuncs For an earlier thread about using vector valued ufuncs for sorts and such -- and negative reactions to the suggestion -- go herehttp://thread.gmane.org/gmane.comp.python.numeric.general/20552/focus=20560. One of the major objections was how to call such functions with the ufunc machinery and the needed machinery for type promotions, sub classes, and all that nonsense. Are these dealt with in the patch? The current numpy code for all that is a bit of a mess anyway, and it would be nice to figure out some unified interface to call through and clean up the current code in the process. In fact, I've been making some preliminary efforts in that direction by formatting the current code and working through it. Also, do we also want reduce methods and such? I suspect there is still a lot of work to do to get this whole thing up and running. Chuck -- The good news is that the patch just uses of the existing code to deal with all the tricky issues (this is why the patch is so short). By the way, sort could be implemented with the proposed specifications, its signature would be (i)-(i). I agree that it would be nice if that code could be made somewhat clearer; however, I think that this task is orthogonal to the generalized ufuncs patch, because there is no code overlap. The way the suggested implementation basically works is to remove the core dimensions from the input/output arrays, and then have the existing code handle all the intricacies over the loop dimensions. Reduce methods are currently not supported (an error is raised). Therefore, the current patch does not forestall anything and the desired functionality can be added whenever it is clear what would be best. I do not think that it would makes sense to specify/implement all possible extensions, optimizations, concrete ufuncs, morphing of existing numpy functions to ufuncs, etc. at once; presumably it is much better to start with a small but extremely flexible specification of generalized ufuncs first. Best, Hansres ___ Numpy-discussion mailing list Numpy-discussion@scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Generalized ufuncs?
I am sorry that our submission http://projects.scipy.org/scipy/numpy/ticket/887 has created some annoyance; presumably we have taken the Make contributions (e.g. code patches), (...) by submitting a 'ticket' on the Trac pages linked below on http://scipy.org/Developer_Zone somewhat too literally. It was great to receive very positive responses to our patch and to receive a very timely review; this is encouraging for submitting code to the numpy repository. I would like to add that the proposed change is not that arbitrary; it is a well-established concept -- it is outlined on the GeneralLoopingFunctions wiki page, and it is a prominent concept of perl's PDL vector library. Of course, there is still room for arguing about details. The fact that no explicit generalized ufuncs are provided, should in my opinion not be an argument why not to include the change in the 1.2.0 release. Writing extension libraries that implement such generalized ufuncs, while being able to use the standard numpy distribution, would certainly be very valuable. Furthermore, the risk for including the proposed patch in 1.2.0 is very low: existing functionality is not touched. (Except the glitch we had by declaring variables in a gcc-way.) For standard ufuncs, it should be straightforward to see that the patch does not modify the behavior at all. I wish everyone a great SciPy'08 conference and very much hope that you can include the proposed functionality in the forthcoming NumPy release. Best regards, Hansres ___ Numpy-discussion mailing list Numpy-discussion@scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Problem with the mailing list?
Charles R Harris [EMAIL PROTECTED] wrote in message news:[EMAIL PROTECTED]... Hi All, I received an email from Hans Andreas -- the gen-ufuncs guy -- and he is unable to post to the list even though subscribed. Anyone know what might be the problem? Please cc [EMAIL PROTECTED] if you reply to this post. Chuck Thanks Chuck for following up on this -- it works now! (Initially I signed up a different address than the one I used for posting.) ___ Numpy-discussion mailing list Numpy-discussion@scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Generalised ufuncs branch
Stéfan van der Walt [EMAIL PROTECTED] wrote in message news:[EMAIL
PROTECTED]...
Hi all,
I have moved the generalised ufuncs functionality off to a branch:
http://svn.scipy.org/svn/numpy/branches/gen_ufuncs
Please try it out and give us your feedback. We shall also pound on
it at the sprint during SciPy'08, and thereafter decide how and when
to best integrate it into NumPy.
Thanks to Wenjie Fu and Hans-Andreas Engel for taking the time to
think this issue through and to submit such a high-quality patch.
Regards
Stéfan
--
Motivated by Stefan's and Chuck's off-line comments, I have attached an example
implementation of a generalized ufunc:
/*
* This implements the function out = inner1d(in1, in2) with
* out[K] = sum_i { in1[K, i] * in2[K, i] }
* and multi-index K, as described on
* http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions
* and on http://projects.scipy.org/scipy/numpy/ticket/887.
*/
static void
DOUBLE_inner1d(char **args, intp *dimensions, intp *steps, void *func)
{
/* Standard ufunc loop length and strides. */
intp dn = dimensions[0];
intp s0 = steps[0];
intp s1 = steps[1];
intp s2 = steps[2];
intp n;
/* Additional loop length and strides for dimension i in elementary
function. */
intp di = dimensions[1];
intp i_s1 = steps[3];
intp i_s2 = steps[4];
intp i;
/* Outer loop: equivalent to standard ufuncs */
for (n = 0; n dn; n++, args[0] += s0, args[1] += s1, args[2] += s2) {
char *ip1 = args[0], *ip2 = args[1], *op = args[2];
/* Implement elementary function: out = sum_i { in1[i] * in2[i] } */
npy_double sum = 0;
for (i = 0; i di; i++) {
sum += (*(npy_double *)ip1) * (*(npy_double *)ip2);
ip1 += i_s1;
ip2 += i_s2;
}
*(npy_double *)op = sum;
}
}
/* Actually create the ufunc object */
static PyUFuncGenericFunction inner1d_functions[] = { DOUBLE_inner1d };
static void *inner1d_data[] = { (void *)NULL };
static char inner1d_sigs[] = { PyArray_DOUBLE, PyArray_DOUBLE, PyArray_DOUBLE };
static void
addInner1d(PyObject *dictionary) {
PyObject *f = PyUFunc_FromFuncAndDataAndSignature(inner1d_functions,
inner1d_data, inner1d_sigs, 1,
2, 1, PyUFunc_None, inner1d,
inner on the last dimension and broadcast
on the other dimensions,
0,
(i),(i)-());
PyDict_SetItemString(dictionary, inner1d, f);
}
___
Numpy-discussion mailing list
Numpy-discussion@scipy.org
http://projects.scipy.org/mailman/listinfo/numpy-discussion
