Re: [Numpy-discussion] failure to register ufunc loops for user defined types
Hi Geoffrey, On Mon, Dec 5, 2011 at 12:37 AM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 6:45 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 6:59 PM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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 (164)/((164)+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 irv...@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() Thanks. Could you put this work into a separate branch, say fixuserloops, and enter a
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Mon, Dec 5, 2011 at 6:59 AM, Charles R Harris charlesr.har...@gmail.com wrote: Hi Geoffrey, On Mon, Dec 5, 2011 at 12:37 AM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 6:45 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 6:59 PM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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 (164)/((164)+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 irv...@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
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 9:45 PM, Charles R Harris charlesr.har...@gmail.com wrote: 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. I agree that those will be useful, but I am worried about adding more stuff in multiarray. User-types should really be separated from multiarray. Ideally, they should be plugins but separated from multiarray would be a good first step. I realize it is a bit unfair to have this ready for Geoffray's code changes, but depending on the timelines for the 2.0.0 milestone, I think this would be a useful thing to have. Otherwise, if some ABI/API changes are needed after 2.0, we will be dragged down with this for years. I am willing to spend time on this. Geoffray, does this sound acceptable to you ? David ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 11:37 PM, Geoffrey Irving irv...@naml.us wrote: snip 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 irv...@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() To support this properly, I think we would need to convert needs_api into an enum with this hybrid mode as another case. While it isn't done currently, I was imagining using a thread pool to multithread the trivially data-parallel operations when needs_api is false, and I suspect the PyGILState_Ensure/Release would trigger undefined behavior in a thread created entirely outside of the Python system. For comparison, I created a special mechanism for simplified multi-threaded exceptions in the nditer in the 'errmsg' parameter: http://docs.scipy.org/doc/numpy/reference/c-api.iterator.html#NpyIter_GetIterNext Worth considering is also the fact that the PyGILState API is incompatible with multiple embedded interpreters. Maybe that's not something anyone does with NumPy, though. -Mark Geoffrey ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On 5 December 2011 17:25, Mark Wiebe mwwi...@gmail.com wrote: On Sun, Dec 4, 2011 at 11:37 PM, Geoffrey Irving irv...@naml.us wrote: snip 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 irv...@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() To support this properly, I think we would need to convert needs_api into an enum with this hybrid mode as another case. While it isn't done currently, I was imagining using a thread pool to multithread the trivially data-parallel operations when needs_api is false, and I suspect the PyGILState_Ensure/Release would trigger undefined behavior in a thread created entirely outside of the Python system. PyGILState_Ensure/Release can be safely used by non-python threads with the only requirement that the GIL has been initialized previously in the main thread (PyEval_InitThreads). For comparison, I created a special mechanism for simplified multi-threaded exceptions in the nditer in the 'errmsg' parameter: http://docs.scipy.org/doc/numpy/reference/c-api.iterator.html#NpyIter_GetIterNext Worth considering is also the fact that the PyGILState API is incompatible with multiple embedded interpreters. Maybe that's not something anyone does with NumPy, though. -Mark Geoffrey ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Mon, Dec 5, 2011 at 8:58 AM, David Cournapeau courn...@gmail.com wrote: On Sun, Dec 4, 2011 at 9:45 PM, Charles R Harris charlesr.har...@gmail.com wrote: 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. I agree that those will be useful, but I am worried about adding more stuff in multiarray. User-types should really be separated from multiarray. Ideally, they should be plugins but separated from multiarray would be a good first step. I think the object and datetime dtypes should also be moved out of the core multiarray module at some point. The user-type mechanism could be improved a lot based on Martin's feedback after he did the quaternion implementation, and needs further expansion to be able to support object and datetime arrays as currently implemented. I realize it is a bit unfair to have this ready for Geoffray's code changes, but depending on the timelines for the 2.0.0 milestone, I think this would be a useful thing to have. Otherwise, if some ABI/API changes are needed after 2.0, we will be dragged down with this for years. I am willing to spend time on this. Geoffray, does this sound acceptable to you ? A rational type could be added without breaking the ABI, in the same way it was done for datetime and half in 1.6. I think the revamp of the user-type mechanism needs its own NEP design document, because changing it will be a very delicate operation in dealing with how it interacts with the NumPy core, and making it much more programmer-friendly will take a fair number of design iterations. -Mark David ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Mon, Dec 5, 2011 at 9:37 AM, mark florisson markflorisso...@gmail.comwrote: On 5 December 2011 17:25, Mark Wiebe mwwi...@gmail.com wrote: On Sun, Dec 4, 2011 at 11:37 PM, Geoffrey Irving irv...@naml.us wrote: snip 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 irv...@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() To support this properly, I think we would need to convert needs_api into an enum with this hybrid mode as another case. While it isn't done currently, I was imagining using a thread pool to multithread the trivially data-parallel operations when needs_api is false, and I suspect the PyGILState_Ensure/Release would trigger undefined behavior in a thread created entirely outside of the Python system. PyGILState_Ensure/Release can be safely used by non-python threads with the only requirement that the GIL has been initialized previously in the main thread (PyEval_InitThreads). Is there a way this could efficiently be used to propagate any errors back to the main thread, for example using TBB as the thread pool? The innermost task code which calls the inner loop can't call PyErr_Occurred() without first calling PyGILState_Ensure itself, which would kill utilization. Maybe this is an ABI problem in NumPy that needs to be fixed, to mandate that inner loops always return an error code and disallow them from setting the Python exception state without returning failure. -Mark For comparison, I created a special mechanism for simplified multi-threaded exceptions in the nditer in the 'errmsg' parameter: http://docs.scipy.org/doc/numpy/reference/c-api.iterator.html#NpyIter_GetIterNext Worth considering is also the fact that the PyGILState API is incompatible with multiple embedded interpreters. Maybe that's not something anyone does with NumPy, though. -Mark Geoffrey ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On 5 December 2011 17:48, Mark Wiebe mwwi...@gmail.com wrote: On Mon, Dec 5, 2011 at 9:37 AM, mark florisson markflorisso...@gmail.com wrote: On 5 December 2011 17:25, Mark Wiebe mwwi...@gmail.com wrote: On Sun, Dec 4, 2011 at 11:37 PM, Geoffrey Irving irv...@naml.us wrote: snip 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 irv...@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() To support this properly, I think we would need to convert needs_api into an enum with this hybrid mode as another case. While it isn't done currently, I was imagining using a thread pool to multithread the trivially data-parallel operations when needs_api is false, and I suspect the PyGILState_Ensure/Release would trigger undefined behavior in a thread created entirely outside of the Python system. PyGILState_Ensure/Release can be safely used by non-python threads with the only requirement that the GIL has been initialized previously in the main thread (PyEval_InitThreads). Is there a way this could efficiently be used to propagate any errors back to the main thread, for example using TBB as the thread pool? The innermost task code which calls the inner loop can't call PyErr_Occurred() without first calling PyGILState_Ensure itself, which would kill utilization. No, there is no way these things can be efficient, as the GIL is likely contented anyway (I wasn't making a point for these functions, just wanted to clarify). There is in fact the additional problem that PyGILState_Ensure would initialize a threadstate, you set an exception, and when you call PyGILState_Release the threadstate gets deleted along with the exception, before you will even have a chance to check with PyErr_Occurred(). For cython.parallel I worked around this by calling PyGILState_Ensure (to initialize the thread state), followed immediately by Py_BEGIN_ALLOW_THREADS before starting any work. You then have to fetch the exception and restore it in another thread when you want to propagate it. It's a total mess, it's inefficient and if you can avoid it you should. Maybe this is an ABI problem in NumPy that needs to be fixed, to mandate that inner loops always return an error code and disallow them from setting the Python exception state without returning failure. That would likely be the best thing. -Mark For comparison, I created a special mechanism for simplified multi-threaded exceptions in the nditer in the 'errmsg' parameter: http://docs.scipy.org/doc/numpy/reference/c-api.iterator.html#NpyIter_GetIterNext Worth considering is also the fact that the PyGILState API is incompatible with multiple embedded interpreters. Maybe that's not something anyone does with NumPy, though. -Mark Geoffrey ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On 5 December 2011 17:57, mark florisson markflorisso...@gmail.com wrote: On 5 December 2011 17:48, Mark Wiebe mwwi...@gmail.com wrote: On Mon, Dec 5, 2011 at 9:37 AM, mark florisson markflorisso...@gmail.com wrote: On 5 December 2011 17:25, Mark Wiebe mwwi...@gmail.com wrote: On Sun, Dec 4, 2011 at 11:37 PM, Geoffrey Irving irv...@naml.us wrote: snip 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 irv...@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() To support this properly, I think we would need to convert needs_api into an enum with this hybrid mode as another case. While it isn't done currently, I was imagining using a thread pool to multithread the trivially data-parallel operations when needs_api is false, and I suspect the PyGILState_Ensure/Release would trigger undefined behavior in a thread created entirely outside of the Python system. PyGILState_Ensure/Release can be safely used by non-python threads with the only requirement that the GIL has been initialized previously in the main thread (PyEval_InitThreads). Is there a way this could efficiently be used to propagate any errors back to the main thread, for example using TBB as the thread pool? The innermost task code which calls the inner loop can't call PyErr_Occurred() without first calling PyGILState_Ensure itself, which would kill utilization. No, there is no way these things can be efficient, as the GIL is likely contented anyway (I wasn't making a point for these functions, just wanted to clarify). There is in fact the additional problem that PyGILState_Ensure would initialize a threadstate, you set an exception, and when you call PyGILState_Release the threadstate gets deleted along with the exception, before you will even have a chance to check with PyErr_Occurred(). To clarify, this case will only happen if you're doing this from a non-Python thread that doesn't have a threadstate to begin with. For cython.parallel I worked around this by calling PyGILState_Ensure (to initialize the thread state), followed immediately by Py_BEGIN_ALLOW_THREADS before starting any work. You then have to fetch the exception and restore it in another thread when you want to propagate it. It's a total mess, it's inefficient and if you can avoid it you should. Maybe this is an ABI problem in NumPy that needs to be fixed, to mandate that inner loops always return an error code and disallow them from setting the Python exception state without returning failure. That would likely be the best thing. -Mark For comparison, I created a special mechanism for simplified multi-threaded exceptions in the nditer in the 'errmsg' parameter: http://docs.scipy.org/doc/numpy/reference/c-api.iterator.html#NpyIter_GetIterNext Worth considering is also the fact that the PyGILState API is incompatible with multiple embedded interpreters. Maybe that's not something anyone does with NumPy, though. -Mark Geoffrey ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sat, Dec 3, 2011 at 8:14 PM, Geoffrey Irving irv...@naml.us wrote:
Hello,
I'm trying to add a fixed precision rational number dtype to numpy,
and am running into an issue trying to register ufunc loops. The code
in question looks like
int npy_rational = PyArray_RegisterDataType(rational_descr);
PyObject* equal = ... // extract equal object from the imported numpy
module
int types[3] = {npy_rational,npy_rational,NPY_BOOL};
if
(PyUFunc_RegisterLoopForType((PyUFuncObject*)ufunc,npy_rational,rational_ufunc_##name,_types,0)0)
return 0;
In Python 2.6.7 with the latest numpy from git, I get
from rational import *
i = array([rational(5,3)])
i
array([5/3], dtype=rational)
equal(i,i)
Traceback (most recent call last):
File stdin, line 1, in module
TypeError: ufunc 'equal' not supported for the input types, and
the inputs could not be safely coerced to any supported types
according to the casting rule ''safe''
The same thing happens with (rational,rational)-rational ufuncs like
multiply.
The full extension module code is here:
https://github.com/girving/poker/blob/rational/rational.cpp
I realize this isn't much information to go on, but let me know if
anything comes to mind in terms of possible reasons or further tests
to run. Unfortunately it looks like the ufunc ntypes and types
properties aren't updated based on user-defined loops, so I'm not yet
sure if the problem is in registry or resolution.
It's also possible someone else hit this before:
http://projects.scipy.org/numpy/ticket/1913.
I haven't tried adding a new type and can't offer any suggestions. But
there was a recent implementation of a quaternion type that might be worth
looking at for comparison. You can find it here http://tinyurl.com/83kesqd.
Chuck
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
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;
}
}
On Sun, Dec 4, 2011 at 9:29 AM, Charles R Harris
charlesr.har...@gmail.com wrote:
On Sat, Dec 3, 2011 at 8:14 PM, Geoffrey Irving irv...@naml.us wrote:
Hello,
I'm trying to add a fixed precision rational number dtype to numpy,
and am running into an issue trying to register ufunc loops. The code
in question looks like
int npy_rational = PyArray_RegisterDataType(rational_descr);
PyObject* equal = ... // extract equal object from the imported numpy
module
int types[3] = {npy_rational,npy_rational,NPY_BOOL};
if
(PyUFunc_RegisterLoopForType((PyUFuncObject*)ufunc,npy_rational,rational_ufunc_##name,_types,0)0)
return 0;
In Python 2.6.7 with the latest numpy from git, I get
from rational import *
i = array([rational(5,3)])
i
array([5/3], dtype=rational)
equal(i,i)
Traceback (most recent call last):
File stdin, line 1, in module
TypeError: ufunc 'equal' not supported for the input types, and
the inputs could not be safely coerced to any supported types
according to the casting rule ''safe''
The same thing happens with (rational,rational)-rational ufuncs like
multiply.
The full extension module code is here:
https://github.com/girving/poker/blob/rational/rational.cpp
I realize this isn't much information to go on, but let me know if
anything comes to mind in terms of possible reasons or further tests
to run. Unfortunately it looks like the ufunc ntypes and types
properties aren't updated based on user-defined loops, so I'm not yet
sure if the problem is in registry or resolution.
It's also possible someone else hit this before:
http://projects.scipy.org/numpy/ticket/1913.
I haven't tried adding a new type and can't offer any suggestions. But there
was a recent implementation of a quaternion type that might be worth looking
at for comparison. You can find it here.
Chuck
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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. Speaking of tests... I was wondering if you could be talked into putting together a simple user type for including in the tests? snip Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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 (164)/((164)+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. As for the patch below, I wouldn't bother looking at it until I get the rest of the bugs out of the way (whether they're in my code or numpy). Geoffrey - diff --git a/numpy/core/src/umath/ufunc_type_resolution.c b/numpy/core/src/umath/ufunc_type_resolution.c index 0d6cf19..4e81e92 100644 --- a/numpy/core/src/umath/ufunc_type_resolution.c +++ b/numpy/core/src/umath/ufunc_type_resolution.c @@ -1656,7 +1656,7 @@ linear_search_userloop_type_resolver(PyUFuncObject *self, /* Found a match */ case 1: set_ufunc_loop_data_types(self, op, out_dtype, types); -return 0; +return 1; } funcdata = funcdata-next; ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 6:59 PM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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 (164)/((164)+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? As for the patch below, I wouldn't bother looking at it until I get the rest of the bugs out of the way (whether they're in my code or numpy). snip Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] failure to register ufunc loops for user defined types
On Sun, Dec 4, 2011 at 6:45 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 6:59 PM, Geoffrey Irving irv...@naml.us wrote: On Sun, Dec 4, 2011 at 5:18 PM, Charles R Harris charlesr.har...@gmail.com wrote: On Sun, Dec 4, 2011 at 5:41 PM, Geoffrey Irving irv...@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 (164)/((164)+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 irv...@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 ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
