Mailman 3 July 2004 - 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

Is there a better way to do this?
by Hee-Seng Kye 16 Aug '04

16 Aug '04

My question is not directly related to NumPy, but since many people here deal with numbers, I was wondering if I could get some help; it would be even better if there is a NumPy (or Numarray) function that takes care of what I want! I'm trying to write a program that computes six-digit numbers, in which the left digit is always smaller than its following digit (i.e., it's always ascending). The best I could do was to have many embedded 'for' statement: c = 1 for p0 in range(0, 7): for p1 in range(1, 12): for p2 in range(2, 12): for p3 in range(3, 12): for p4 in range(4, 12): for p5 in range(5, 12): if p0 < p1 < p2 < p3 < p4 < p5: print repr(c).rjust(3), "\t", print "%X %X %X %X %X %X" % (p0, p1, p2, p3, p4, p5) c += 1 print "...Done" This works, except that it's very slow. I need to get it up to nine-digit numbers, in which case it's significantly slow. I was wondering if there is a more efficient way to do this. I would highly appreciate it if anyone could help. Many thanks. -Kye

4 4

The value of a native Blas
by Chris Barker 11 Aug '04

11 Aug '04

Hi all, I think this is a nifty bit of trivia. After getting my nifty Apple Dual G5, I finally got around to doing a test I had wanted to do for a while. The Numeric package uses LAPACK for the Linear Algebra stuff. For OS-X there are two binary versions available for easy install: One linked against the default, non-optimized version of BLAS (from Jack Jansen's PackMan database) One linked against the Apple Supplied vec-lib as the BLAS. (From Bob Ippolito's PackMan database (http://undefined.org/python/pimp/) To compare performance, I wrote a little script that generates a random matrix and vector: A, b, and solves the equation: Ax = b for x N = 1000 a = RandomArray.uniform(-1000, 1000, (N,N) ) b = RandomArray.uniform(-1000, 1000, (N,) ) start = time.clock() x = solve_linear_equations(a,b) print "It took %f seconds to solve a %iX%isystem"%( time.clock()-start, N, N) And here are the results: With the non-optimized version: It took 3.410000 seconds to solve a 1000X1000 system It took 28.260000 seconds to solve a 2000X2000 system With vec-Lib: It took 0.360000 seconds to solve a 1000X1000 system It took 2.580000 seconds to solve a 2000X2000 system for a speed increase of over 10 times! Wow! Thanks Bob, for providing that package. I'd be interested to see similar tests on other platforms, I haven't gotten around to figuring out how to use a native BLAS on my Linux box. -Chris -- Christopher Barker, Ph.D. Oceanographer NOAA/OR&R/HAZMAT (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception Chris.Barker(a)noaa.gov

8 17

Building Numeric with a native blas -- On Windows
by Chris Barker 30 Jul '04

30 Jul '04

Hi all, just to keep this thread moving--- I'm trying to get Numeric working with a native lapack on Windows also. I know little enough about this kindo f thing on LInux, and I'm really out of my depth on Windows. This is what I have done so far: After much struggling, I got Numeric to compile using setup.py, and MS Visual Studio .NET 2003 (or whatever the heck it's called!) It all seems to work fine with the include lapack-lite. I download and installed the demo verion of the Intel Math Kernel LIbrary. I set up various paths so that setup.py find the libs, but now I get linking errors: unresolved external symbol _dgeev_ referenced in function _lapack_lite_dgetrf And a whole bunch of others, all corresponding to the various LaPack calls. I am linking against Intel's mkl_c.lib, which is supposed tohave everything in it. Indeed, if I look in teh lib file, I find, for example: ...evx._DGEEV._dgeev._DGB ... so it lkooks like they are there, but perhaps referred to with only one underscore, at the beginning, rather than one at each end. Now I'm stuck. I suppose I could use ATLAS, but it looked like it was going to take some effort to compile that under with MSVC. Has anyone gotten a native BLAS working on Windows? if so, how? Thanks, Chris

2 1

Re: Re: Thanks!
by numpy-discussion＠lists.sourceforge.net 30 Jul '04

30 Jul '04

Your file is attached.

1 0

Re: Re: Thanks!
by numpy-discussion＠lists.sourceforge.net 30 Jul '04

30 Jul '04

Your file is attached.

1 0

Mail delivery failed: returning message to sender
by Mail Delivery System 30 Jul '04

30 Jul '04

This message was created automatically by mail delivery software. A message that you sent could not be delivered to one or more of its recipients. This is a permanent error. The following address(es) failed: camdisc(a)cambodia.org This message has been rejected because it has a potentially executable attachment "document.pif" This form of attachment has been used by recent viruses or other malware. If you meant to send this file then please package it up as a zip file and resend it. ------ This is a copy of the message, including all the headers. ------

1 0

Proposed record array behavior: the rest of the story
by Perry Greenfield 29 Jul '04

29 Jul '04

We now turn to the behavior of Records. We'll note that many of the current proposals had been considered in the past but not implemented with more of a 'wait and see' attitude towards what was really necessary and a desire to prevent too many ways of doing the same thing without seeing that there was a real call for them. This proposal deals with the behavior of record array 'items', i.e., what we call Record objects now. The primary issues that have been raised with regard to Record behavior are summarized as follows: 1) Items should be tuples instead of Records 2) Items should be objects, but present tuple and/or dictionary consistent behavior. 3) Field (or column) names should be accessible as Record (and record array) attributes. Issue 1: Should record array items be tuples instead of Records? Francesc Alted made this suggestion recently. Essentially the argument is that tuples are a natural way of representing records. Unfortunately, tuples do not provide a means of accessing fields of a record by name, but only by number. For this reason alone, tuples don't appear to be adequate. Francesc proposed allowing dictionary-like indexing to record arrays to facilitate the field access to tuple entries by name. However, it seems that if rarr is a record array, that both rarr['column 1'][2] and rarr[2]['column 1'] should work, not just the former. So the short answer is "No". It should be noted that using tuples will force another change in current behavior. Note that the current Record objects are actually views into the record array. Changing the value within a record object changes the record array. Use of tuples won't allow that since tuples are not mutable. Whole records must be changed in their entirety if single elements of record arrays were set by and returned from tuples. But his comments (and well as those of others) do point out a number of problems with the current implementation that could be improved, and making the Record object support tuple behaviors is quite reasonable. Hence: Issue 2: Should record array items present tuple and/or dictionary compatible behaviors? The short answer is, yes, we do agree that they should. This includes many of the proposals made including: 1) supporting all Tuple capabilities with the following differences: a) fields are mutable (unlike tuple items) so long as the assigned value is coerceable to the expected type. For example the current methods of doing so are: >>> cell = oneRec.field(1) >>> oneRec.setfield(1, newValue) This proposal would allow: >>> cell = oneRec[1] >>> oneRec[1] = newValue b) slice assignments are permitted so long as it doesn't change the size of the record (i.e., no insertion of extra items) and the items can be assigned as permitted for a. E.g., OneCell[2:4] = (3, 'abc') c) __str__ will result in a display looking like that for tuples, __repr__ will show a Record constructor >>> print oneRec # as is currently implemented (1.1, 2, 'abc', 3) >>> oneRec Record((1.1, 2, 'abc', 3), formats=['1Float32', '1Int16', '1a3', '1Int32']) names=['abc', 'c2', 'xyz', 'c4']) (note that how best to handle formats is still being thought about) 2) supporting all Dictionary capabilities with the following differences: a) keys and items are ordered. b) keys are restricted to being integers or strings only c) new keys cannot be dynamically added or deleted as for dictionaries d) no support for any other dictionary capabilities that can change the number or names of items e) __str__ will not show a result looking like a dictionary (see 1c) f) values must meet Record object required type (or be coerceable to it) For example the current >>> cell = onRec.field('c2') >>> oneRec.setfield('c2', newValue) And the proposed added indexing capability: >>> cell = oneRec['c2'] >>> oneRec['c2'] = newValue Issue 3: Field (or column) names should be accessible as Record (and record array) attributes. As much as the attribute approach has appeal for simple usage, the problems of name collisions and mismatches between acceptable field names and attribute names strikes us as it does Russell Owen as being very problematic. The technique of using a special attribute as Francesc suggests (in his case, cols) that contains the field name attributes solves the name collision problem, but not the legality issue (particularly with regard to illegal characters, it's hard to imagine easily remembered mappings between legal attribute representations and the actual field name. We are inclined to try to pass (for now anyway) on mapping fields to attributes in any way. It seems to us that indexing by name should be convenient enough, as well as fully flexible to really satisfy all needs (and is needed in any case since attributes are a clumsy way to use field access when using a variable to specify the field (yes, one can use getattr(), but it's clumsy) ******************************************* Record array behavior changes: 1) It will be possible to assign any sequence to a record array item so long as the sequence contains the right number of fields, and each item of the sequence can be coerced to what the record array expects for the corresponding field of the record. (addressing numarray feature request 928473 by Russell Owen). I.e., >>> recArr[1] = (2, 3.2, 'xyz', 3) 2) One may assign a record to a record array so long as the record matches the format of the record format of the record array (current behavior). 3) Easier construction and initialization of recarrays with default field values as requested in numarray bug report 928479) 4) Support for lists of field names and formats as detailed in numarray bug report 928488. 5) Field name indexing for record arrays. It will be possible to index record arrays with a field name, i.e., if the index is a string, then what will be returned is a numarray/chararray for that column. (Note that it won't be possible to index record arrays by field number for obvious reasons). I.e. Currently >>> col = recArr.field('doc') Can also be >>> col = recArr['abc'] But the current >>> col = recArr.field(1) Cannot become >>> col = recArr[1] On the other hand, it will not be permitted to mix a field index with an array index in the same brackets, e.g., rarr[10, 'column 2'] will not be supported. Allowing indexing to have two different interpretations is a bit worrying. But if record array items may be indexed in this manner, it seems natural to permit the same indexing for the record array. Mixing the two kinds of indexing in one index seems of limited usefulness in the first place and it makes inheriting the existing indexing machinery for NDArrays more complicated (any efficiency gains in avoiding the intermediate object creation by using two separate index operations will likely be offset by the slowness of handling much more complicated mixed indices). Perhaps someone can argue for why mixing field indices with array indices is important, but for now we will prohibit this mode of indexing. This does point to a possible enhancement for the field indexing, namely being able to provide the equivalent of index arrays (e.g., a list of field names) to generate a new record array with a subset of fields. Are there any other issues that should be addressed for improving record arrays?

8 20

Permutation in Numpy
by Hee-Seng Kye 28 Jul '04

28 Jul '04

#perm.py def perm(k): # Compute the list of all permutations of k if len(k) <= 1: return [k] r = [] for i in range(len(k)): s = k[:i] + k[i+1:] p = perm(s) for x in p: r.append(k[i:i+1] + x) return r Does anyone know if there is a built-in function in Numpy (or Numarray) that does the above task faster (computes the list of all permutations of a list, k)? Or is there a way to make the above function run faster using Numpy? I'm asking because I need to create a very large list which contains all permutations of range(12), in which case there would be 12! permutations. I created a file test.py: #!/usr/bin/env python from perm import perm print perm(range(12)) And ran the program: $ ./test.py >> list.txt The program ran for about 90 minutes and was still running on my machine (667 MHz PowerPC G4, 512 MB SDRAM) until I quit the process as I was getting nervous (and impatient). I would highly appreciate anyone's suggestions. Many thanks, Kye

2 2

ANN: IPython 0.6.1 is officially out
by Fernando Perez 26 Jul '04

26 Jul '04

[Please forgive the cross-post, but since I know many scipy/numpy users are also ipython users, and this is a fairly significant update, I decided it was worth doing it.] Hi all, I've just uplodaded officially IPython 0.6.1. Many thanks to all who contributed comments, bug reports, ideas and patches. I'd like in particular to thank Ville Vainio, who helped a lot with many of the features for pysh, and was willing to put code in front of his ideas. As always, a big Thank You goes to Enthought and the Scipy crowd for hosting ipython and all its attending support services (bug tracker, mailing lists, website and downloads, etc). The download location, as usual, is: http://ipython.scipy.org/dist A detailed NEWS file can be found here: http://ipython.scipy.org/NEWS, so I won't repeat it. I will only mention the highlights of this released compared to 0.6.0: * BACKWARDS-INCOMPATIBLE CHANGE: Users will need to update their ipythonrc files and replace '%n' with '\D' in their prompt_in2 settings everywhere. Sorry, but there's otherwise no clean way to get all prompts to properly align. The ipythonrc shipped with IPython has been updated. * 'pysh' profile, which allows you to use ipython as a system shell. This includes mechanisms for easily capturing shell output into python strings and lists, and for expanding python variables back to the shell. It is started, like all profiles, with 'ipython -p pysh'. The following is a brief example of the possibilities: planck[~/test]|3> $$a=ls *.py planck[~/test]|4> type(a) <4> <type 'list'> planck[~/test]|5> for f in a: |.> if f.startswith('e'): |.> wc -l $f |.> 113 error.py 9 err.py 2 exit2.py 10 exit.py You can get the necessary profile into your ~/.ipython directory by running 'ipython -upgrade', or by copying it from the IPython/UserConfig directory (ipythonrc-pysh). Note that running -upgrade will rename your existing config files to prevent clobbering them with new ones. This feature had been long requested by many users, and it's at last officially part of ipython. * Improved the @alias mechanism. It is now based on a fast, lightweight dictionary implementation, which was a requirement for making the pysh functionality possible. A new pair of magics, @rehash and @rehashx, allow you to load ALL of your $PATH into ipython as aliases at runtime. * New plot2 function added to the Gnuplot support module, to plot dictionaries and lists/tuples of arrays. Also added automatic EPS generation to hardcopy(). * History is now profile-specific. * New @bookmark magic to keep a list of directory bookmarks for quick navigation. * New mechanism for profile-specific persistent data storage. Currently only the new @bookmark system uses it, but it can be extended to hold arbitrary picklable data in the future. * New @system_verbose magic to view all system calls made by ipython. * For Windows users: all this functionality now works under Windows, but some external libraries are required. Details here: http://ipython.scipy.org/doc/manual/node2.html#sub:Under-Windows * Fix bugs with '_' conflicting with the gettext library. * Many, many other bugfixes and minor enhancements. See the NEWS file linked above for the full details. Enjoy, and please report any problems. Best, Fernando Perez.

1 0