[Numpy-discussion] failure to register ufunc loops for user defined types

Geoffrey Irving irving at naml.us
Mon Dec 5 02:37:29 EST 2011


On Sun, Dec 4, 2011 at 6:45 PM, Charles R Harris
<charlesr.harris at gmail.com> wrote:
>
>
> On Sun, Dec 4, 2011 at 6:59 PM, Geoffrey Irving <irving at naml.us> wrote:
>>
>> On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris
>> <charlesr.harris at gmail.com> wrote:
>> >
>> >
>> > On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving <irving at naml.us> wrote:
>> >>
>> >> This may be the problem.  Simple diffs are pleasant.  I'm guessing
>> >> this code doesn't get a lot of testing.  Glad it's there, though!
>> >>
>> >> Geoffrey
>> >>
>> >> diff --git a/numpy/core/src/umath/ufunc_type_resolution.c
>> >> b/numpy/core/src/umath/ufunc_type_resolution.c
>> >> index 0d6cf19..a93eda1 100644
>> >> --- a/numpy/core/src/umath/ufunc_type_resolution.c
>> >> +++ b/numpy/core/src/umath/ufunc_type_resolution.c
>> >> @@ -1866,7 +1866,7 @@ linear_search_type_resolver(PyUFuncObject *self,
>> >>             case -1:
>> >>                 return -1;
>> >>             /* A loop was found */
>> >> -            case 1:
>> >> +            case 0:
>> >>                 return 0;
>> >>         }
>> >>     }
>> >>
>> >
>> > Heh. Can you verify that this fixes the problem? That function is only
>> > called once  and its return value is passed up the chain, but the
>> > documented
>> > return values of that calling function are -1, 0. So the documentation
>> > needs
>> > to be changed if this is the right thing to do.
>>
>> Actually, that patch was wrong, since
>> linear_search_userloop_type_resolver needs to return three values
>> (error, not-found, success).  A better patch follows.  I can confirm
>> that this gets me further, but I get other failures down the line, so
>> more fixes may follow.  I'll push the branch with all my fixes for
>> convenience once I have everything working.
>>
>> > Speaking of tests... I was wondering if you could be talked into putting
>> > together a simple user type for including in the tests?
>>
>> Yep, though likely not for a couple weeks.  If there's interest, I
>> could also be convinced to sanitize my entire rational class so you
>> could include that directly.  Currently it's both C++ and uses some
>> gcc specific features like __int128_t.  Basically it's
>> numerator/denominator, where both are 64 bit integers, and an
>> OverflowError is thrown if anything can't be represented as such
>> (possibly a different exception would be better in cases like
>> (1<<64)/((1<<64)+1)).  It would be easy to generalize it to rational32
>> vs. rational64 as well.
>>
>> If you want tests but not rational, it would be straightforward to
>> strip what I have down to a bare bones test case.
>>
>
> We'll see how much interest there is. If it becomes official you may get
> more feedback on features. There are some advantages to having some user
> types in numpy. One is that otherwise they tend to get lost, another is that
> having a working example or two provides a templates for others to work
> from, and finally they provide test material. Because official user types
> aren't assigned anywhere there might also be some conflicts. Maybe something
> like an extension types module would be a way around that. In any case, I
> think both rational numbers and quaternions would be useful to have and I
> hope there is some discussion of how to do that. Rationals may be a bit
> trickier than quaternions though, as usually they are used to provide exact
> arithmetic without concern for precision. I don't know how restrictive the
> 64 bit limitation will be in practice. What are you using them for?

I'm using them for frivolous analysis of poker Nash equilibria.  I'll
let others decide if it has any non-toy uses.  64 bits seems to be
enough for me, though it's possible that I'll run in trouble with
other examples.  It still exact, though, in the sense that it throws
an exception rather than doing anything weird if it overflows.  And it
has the key advantage of being orders of magnitude faster than object
arrays of Fractions.

Back to the bugs: here's a branch with all the changes I needed to get
rational arithmetic to work:

    https://github.com/girving/numpy

I discovered two more after the last email.  One is another simple 0
vs. 1 bug, and another is somewhat optional:

commit 730b05a892371d6f18d9317e5ae6dc306c0211b0
Author: Geoffrey Irving <irving at naml.us>
Date:   Sun Dec 4 20:03:46 2011 -0800

    After loops, check for PyErr_Occurred() even if needs_api is 0

    For certain types of user defined classes, casting and ufunc loops
    normally run without the Python API, but occasionally need to throw
    an error.  Currently we assume that !needs_api means no error occur.
    However, the fastest way to implement such loops is to run without
    the GIL normally and use PyGILState_Ensure/Release if an error occurs.

    In order to support this usage pattern, change all post-loop checks from

        needs_api && PyErr_Occurred()

    to simply

        PyErr_Occurred()

Geoffrey



More information about the NumPy-Discussion mailing list