Mailman 3 November 2008 - NumPy-Discussion

PyArray_ContiguousFromObject & PyDECREF
by Michiel Jan Laurens de Hoon 17 May '11

17 May '11

Dear pythoneers, In Chapter 12 ("Writing a C extension") in the Numerical Python manual, in the first example PyArray_ContiguousFromObject is used together with PyDECREF, while in the last example PyArray_ContiguousFromObject is used but the arrays are not PyDECREF'ed before the function returns. From Chapter 13 ("C API Reference") it seems that a PyDECREF is needed in both cases. Or am I missing something here? Is PyDECREF not needed for some reason in the last example? Thanks in advance, Michiel de Hoon U Tokyo First example: ============== static PyObject * trace(PyObject *self, PyObject *args) { PyObject *input; PyArrayObject *array; double sum; int i, n; if (!PyArg_ParseTuple(args, "O", &input)) return NULL; array = (PyArrayObject *) PyArray_ContiguousFromObject(input, PyArray_DOUBLE, 2, 2); if (array == NULL) return NULL; n = array->dimensions[0]; if (n > array->dimensions[1]) n = array->dimensions[1]; sum = 0.; for (i = 0; i < n; i++) sum += *(double *)(array->data + i*array->strides[0] + i*array->strides[1]); Py_DECREF(array); return PyFloat_FromDouble(sum); } Last example: ============= static PyObject * matrix_vector(PyObject *self, PyObject *args) { PyObject *input1, *input2; PyArrayObject *matrix, *vector, *result; int dimensions[1]; double factor[1]; double real_zero[1] = {0.}; long int_one[1] = {1}; long dim0[1], dim1[1]; extern dgemv_(char *trans, long *m, long *n, double *alpha, double *a, long *lda, double *x, long *incx, double *beta, double *Y, long *incy); if (!PyArg_ParseTuple(args, "dOO", factor, &input1, &input2)) return NULL; matrix = (PyArrayObject *) PyArray_ContiguousFromObject(input1, PyArray_DOUBLE, 2, 2); if (matrix == NULL) return NULL; vector = (PyArrayObject *) PyArray_ContiguousFromObject(input2, PyArray_DOUBLE, 1, 1); if (vector == NULL) return NULL; if (matrix->dimensions[1] != vector->dimensions[0]) { PyErr_SetString(PyExc_ValueError, "array dimensions are not compatible"); return NULL; } dimensions[0] = matrix->dimensions[0]; result = (PyArrayObject *)PyArray_FromDims(1, dimensions, PyArray_DOUBLE); if (result == NULL) return NULL; dim0[0] = (long)matrix->dimensions[0]; dim1[0] = (long)matrix->dimensions[1]; dgemv_("T", dim1, dim0, factor, (double *)matrix->data, dim1, (double *)vector->data, int_one, real_zero, (double *)result->data, int_one); return PyArray_Return(result); } -- Michiel de Hoon, Assistant Professor University of Tokyo, Institute of Medical Science Human Genome Center 4-6-1 Shirokane-dai, Minato-ku Tokyo 108-8639 Japan http://bonsai.ims.u-tokyo.ac.jp/~mdehoon

3 2

Test failures on 2.6
by T J 18 Nov '10

18 Nov '10

Hi, I'm getting a couple of test failures with Python 2.6, Numpy 1.2.0, Nose 0.10.4: nose version 0.10.4 ..........................................................................................................................................................................................................................................................................................................................................................................................................................................F................K........................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................./share/home/me/usr/lib/python2.6/site-packages/numpy/lib/tests/test_io.py:68: SyntaxWarning: assertion is always true, perhaps remove parentheses? assert(c.readlines(), ........................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................./share/home/me/usr/lib/python2.6/site-packages/numpy/ma/tests/test_core.py:1315: SyntaxWarning: assertion is always true, perhaps remove parentheses? assert(store._mask, True) /home/me/usr/lib/python2.6/site-packages/numpy/ma/tests/test_core.py:1322: SyntaxWarning: assertion is always true, perhaps remove parentheses? assert(store._mask, True) /home/me/usr/lib/python2.6/site-packages/numpy/ma/tests/test_core.py:1989: SyntaxWarning: assertion is always true, perhaps remove parentheses? assert(test.mask, [0,1,0,0,0,0,0,0,0,0]) ...............................................E................................................................................................................................................................................ ====================================================================== ERROR: Tests the min/max functions with explicit outputs ---------------------------------------------------------------------- Traceback (most recent call last): File "/home/me/usr/lib/python2.6/site-packages/numpy/ma/tests/test_core.py", line 653, in test_minmax_funcs_with_output result = npfunc(xm,axis=0,out=nout) File "/home/me/usr/lib/python2.6/site-packages/numpy/core/fromnumeric.py", line 1525, in amin return amin(axis, out) File "/home/me/usr/lib/python2.6/site-packages/numpy/ma/core.py", line 2978, in min np.putmask(out, newmask, np.nan) ValueError: cannot convert float NaN to integer ====================================================================== FAIL: test_umath.TestComplexFunctions.test_against_cmath ---------------------------------------------------------------------- Traceback (most recent call last): File "/home/me/usr/lib/python2.6/site-packages/nose-0.10.4-py2.6.egg/nose/case.py", line 182, in runTest self.test(*self.arg) File "/home/me/usr/lib/python2.6/site-packages/numpy/core/tests/test_umath.py", line 268, in test_against_cmath assert abs(a - b) < atol, "%s %s: %s; cmath: %s"%(fname,p,a,b) AssertionError: arcsin 2: (1.57079632679-1.31695789692j); cmath: (1.57079632679+1.31695789692j) ---------------------------------------------------------------------- Ran 1726 tests in 8.856s FAILED (KNOWNFAIL=1, errors=1, failures=1) <nose.result.TextTestResult run=1726 errors=1 failures=1>

2 1

Unpleasant behavior with poly1d and numpy scalar multiplication
by Fernando Perez 15 Feb '10

15 Feb '10

Hi all, consider this little script: from numpy import poly1d, float, float32 p=poly1d([1.,2.]) three=float(3) three32=float32(3) print 'three*p:',three*p print 'three32*p:',three32*p print 'p*three32:',p*three32 which produces when run: In [3]: run pol1d.py three*p: 3 x + 6 three32*p: [ 3. 6.] p*three32: 3 x + 6 The fact that multiplication between poly1d objects and numbers is: - non-commutative when the numbers are numpy scalars - different for the same number if it is a python float vs a numpy scalar is rather unpleasant, and I can see this causing hard to find bugs, depending on whether your code gets a parameter that came as a python float or a numpy one. This was found today by a colleague on numpy 1.0.4.dev3937. It feels like a bug to me, do others agree? Or is it consistent with a part of the zen of numpy I've missed thus far? Thanks, f

6 17

Compact way of performing array math with specified result type?
by Russell E. Owen 08 Feb '10

08 Feb '10

I often find myself doing simple math on sequences of numbers (which might or might not be numpy arrays) where I want the result (and thus the inputs) coerced to a particular data type. I'd like to be able to say: numpy.divide(seq1, seq2, dtype=float) but ufuncs don't allow on to specify a result type. So I do this instead: numpy.array(seq1, dtype=float) / numpy.array(seq2, dtype=float) Is there a more compact solution (without having to create the result array first and supply it as an argument)? -- Russell

4 3

Is anyone knowledgeable about dll deployment on windows ?
by David Cournapeau 01 Dec '09

01 Dec '09

Hi, My head hurt trying to understand dll management with windows. I wanted to find a sane way to use dll for numscons, but I can't see how to do this, so I was wondering if anyone on this ML had any deep knowledge on how to install dll, and reuse them with python extensions ? The problem is the following: - you have some dll (say MKL) in a path (say C:\program files\MKL). - how do you tell another dll (more exactly a python extension, a .pyd) to look for mkl dll into the MKL path ? It seems that there is no rpath-like capability on windows, and the only way to load dll is to put them into some windows directory, or to put the dll path into PATH. Is it really the only way ? cheers, David

8 25

Slicing/selection in multiple dimensions simultaneously
by Mike Ressler 17 Mar '09

17 Mar '09

The following seems to be a wart: is it expected? Set up a 10x10 array and some indexing arrays: a=arange(100) a.shape=(10,10) q=array([0,2,4,6,8]) r=array([0,5]) Suppose I want to extract only the "even" numbered rows from a - then print a[q,:] <works - output deleted> Every fifth column: print a[:,r] <works - output deleted> Only the even rows of every fifth column: print a[q,r] --------------------------------------------------------------------------- <type 'exceptions.ValueError'> Traceback (most recent call last) /.../.../.../<ipython console> in <module>() <type 'exceptions.ValueError'>: shape mismatch: objects cannot be broadcast to a single shape But, this works: print a[q,:][:,r] [[ 0 5] [20 25] [40 45] [60 65] [80 85]] So why does the a[q,r] form have problems? Thanks for your insights. Mike -- mike.ressler(a)alum.mit.edu

8 16

Automatic differentiation (was Re: second-order gradient)
by Rob Clewley 11 Mar '09

11 Mar '09

> Maybe we should focus on writing a decent 'deriv' function then. I > know Konrad Hinsen's Scientific had a derivatives package > (Scientific.Functions.Derivatives) that implemented automatic > differentiation: > > http://en.wikipedia.org/wiki/Automatic_differentiation That would be great, but wouldn't that be best suited as a utility requiring Sympy? You'll want to take advantage of all sorts of symbolic classes, especially for any source code transformation approach. IMO Hinsen's implementation isn't a very efficient or attractive solution to AD given the great existing C/C++ codes out there. Maybe we should be looking to provide a python interface to an existing open source package such as ADOL-C, but I'm all in favour of a new pure python approach too. What would be perfect is to have a single interface to a python AD package that would support a faster implementation if the user wished to install a C/C++ package, otherwise would default to a pure python equivalent. -Rob

4 4

linalg.norm missing an 'axis' kwarg?!
by Hans Meine 10 Feb '09

10 Feb '09

Hi, I have a 2D matrix comprising a sequence of vectors, and I want to compute the norm of each vector. np.linalg.norm seems to be the best bet, but it does not support axis. Wouldn't this be a nice feature? Greetings, Hans

4 8

Re: [Numpy-discussion] Complete LAPACK needed (Frank Lagor)
by Frank Lagor 14 Jan '09

14 Jan '09

Thanks so much for your help, David. I'm sorry I did not receive your posts previously -- I have the digest mode on and there is a bit of a delay. I'll try to change my options next time I post a request. Thanks so much again, Frank

2 1

ANN: numpy.i - added managed deallocation to ARGOUTVIEW_ARRAY1 (ARGOUTVIEWM_ARRAY1)
by Egor Zindy 08 Jan '09

08 Jan '09

Dear List, after I tried to write a simple ARGOUTVIEW_ARRAY1 example (see http://code.google.com/p/ezwidgets/wiki/NumpySWIGMinGW#A_simple_ARGOUTVIEW_… ), I started wondering about memory deallocation. Turns out a few clever people already did all the thinking (see http://blog.enthought.com/?p=62 ) and a few more clever people use this in a swig file (see http://niftilib.sourceforge.net/pynifti, file nifticlib.i). All this concentrated knowledge helped me implement a single ARGOUTVIEWM_ARRAY1 fragment in numpy.i to do some testing (see attached numpy.i). How to use it? In yourfile.i, the %init function has a few lines added to it: %init %{ import_array(); /* initialize the new Python type for memory deallocation */ _MyDeallocType.tp_new = PyType_GenericNew; if (PyType_Ready(&_MyDeallocType) < 0) return; %} ... and that's it! then just use ARGOUTVIEWM_ARRAY1 instead of ARGOUTVIEW_ARRAY1 and python does the deallocation for you when the python array is destroyed (see the examples attached). Everything compiles, but it would be nice to get rid of the warnings... ezalloc_wrap.c:2763: warning: initialization from incompatible pointer type ezalloc_wrap.c: In function `_wrap_alloc_managed': ezalloc_wrap.c:2844: warning: assignment from incompatible pointer type writing build\temp.win32-2.5\Release\_ezalloc.def Compilation: python setup_alloc.py build Testing: The attached test_alloc.py does 2048 allocations of 1MB each for managed and unmanaged arrays. Output on my XP laptop with 1GB RAM as follows: ARGOUTVIEWM_ARRAY1 (managed arrays) - 2048 allocations (1048576 bytes each) Done! ARGOUTVIEW_ARRAY1 (unmanaged, leaking) - 2048 allocations (1048576 bytes each) Step 482 failed TODO: Find a better name for the methods (not sure I like ARGOUTVIEWM_ARRAY), then do the missing fragments (2D arrays), clear the warnings and verify the method for a possible inclusion in the official numpy.i file. Thank you for reading! Regards, Egor

2 6

NumPy-Discussion November 2008