[Numpy-discussion] Giving numpy the ability to multi-iterate excluding an axis

John Salvatier jsalvati@u.washington....
Wed Dec 22 02:44:28 CST 2010


This now makes sense to me, and I think it should work :D. This is all very
cool. This is going to do big things for cython and numpy.

Some hopefully constructive criticism:

When first reading through the API description, the way oa_ndim and oa_axes
work is not clear. I think your description would be clearer if you explain
what oa_ndim means (I gather something like "the number of axes over which
you wish to iterate"), currently it just says "These parameters let you
control in detail how the axes of the operand arrays get matched together
and iterated."

It's also not totally clear to me how offsetting works. What are the offsets
measured from? It seems like they are measured from another iterator, but
I'm not sure and I don't see how it gets that information.

John

On Tue, Dec 21, 2010 at 5:12 PM, Mark Wiebe <mwwiebe@gmail.com> wrote:

> On Mon, Dec 20, 2010 at 1:42 PM, John Salvatier <jsalvati@u.washington.edu
> > wrote:
>
>> A while ago, I asked a whether it was possible to multi-iterate over
>> several ndarrays but exclude a certain axis(
>> http://www.mail-archive.com/numpy-discussion@scipy.org/msg29204.html),
>> sort of a combination of PyArray_IterAllButAxis and PyArray_MultiIterNew. My
>> goal was to allow creation of relatively complex ufuncs that can allow
>> reduction or directionally dependent computation and still use broadcasting
>> (for example a moving averaging ufunc that can have changing averaging
>> parameters). I didn't get any solutions, which I take to mean that no one
>> knew how to do this.
>>
>> I am thinking about trying to make a numpy patch with this functionality,
>> and I have some questions: 1) How difficult would this kind of task be for
>> someone with non-expert C knowledge and good numpy knowledge? 2) Does anyone
>> have advice on how to do this kind of thing?
>>
>
> You may be able to do what you would like with the new iterator I've
> written.  In particular, it supports nesting multiple iterators by providing
> either pointers or offsets, and allowing you to specify any subset of the
> axes to iterate.  Here's how the code to do this in a simple 3D case might
> look, for making axis 1 the inner loop:
>
> PyArrayObject *op[2] = {a,b};
> npy_intp axes_outer[2] = {0,2}};
> npy_intp *op_axes[2];
> npy_intp axis_inner = 1;
> npy_int32 flags[2] = {NPY_ITER_READONLY, NPY_ITER_READONLY};
> NpyIter *outer, *inner;
> NpyIter_IterNext_Fn oiternext, iiternext;
> npy_intp *ooffsets;
> char **idataptrs;
>
> op_axes[0] = op_axes[1] = axes_outer;
> outer = NpyIter_MultiNew(2, op, NPY_ITER_OFFSETS,
>                            NPY_KEEPORDER, NPY_NO_CASTING, flags, NULL, 2,
> op_axes, 0);
> op_axes[0] = op_axes[1] = &axis_inner;
> inner = NpyIter_MultiNew(2, op, 0, NPY_KEEPORDER, NPY_NO_CASTING, flags,
> NULL, 1, op_axes, 0);
>
> oiternext = NpyIter_GetIterNext(outer);
> iiternext = NpyIter_GetIterNext(inner);
>
> ooffsets = (npy_intp *)NpyIter_GetDataPtrArray(outer);
> idataptrs = NpyIter_GetDataPtrArray(inner);
>
> do {
>    do {
>       char *a_data = idataptrs[0] + ooffsets[0], *b_data = idataptrs[0] +
> ooffsets[0];
>       /* Do stuff with the data */
>    } while(iiternext());
>    NpyIter_Reset(inner);
> } while(oiternext());
>
> NpyIter_Deallocate(outer);
> NpyIter_Deallocate(inner);
>
> Extending to more dimensions, or making both the inner and outer loops have
> multiple dimensions, isn't too crazy.  Is this along the lines of what you
> need?
>
> If you check out my code, note that it currently isn't exposed as NumPy API
> yet, but you can try a lot of things with the Python exposure.
>
> Cheers,
> Mark
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/numpy-discussion/attachments/20101222/3e401b2b/attachment.html 


More information about the NumPy-Discussion mailing list