Mailman 3 June 2009 - NumPy-Discussion

Tutorial topics for SciPy'09 Conference
by Fernando Perez 01 Jul '09

01 Jul '09

Hi all, The time for the Scipy'09 conference is rapidly approaching, and we would like to both announce the plan for tutorials and solicit feedback from everyone on topics of interest. Broadly speaking, the plan is something along the lines of what we had last year: one continuous 2-day tutorial aimed at introductory users, starting from the very basics, and in parallel a set of 'advanced' tutorials, consisting of a series of 2-hour sessions on specific topics. We will request that the presenters for the advanced tutorials keep the 'tutorial' word very much in mind, so that the sessions really contain hands-on learning work and not simply a 2-hour long slide presentation. We will thus require that all the tutorials will be based on tools that the attendees can install at least 2 weeks in advance on all platforms (no "I released it last night" software). With that in mind, we'd like feedback from all of you on possible topics for the advanced tutorials. We have space for 8 slots total, and here are in no particular order some possible topics. At this point there are no guarantees yet that we can get presentations for these, but we'd like to establish a first list of preferred topics to try and secure the presentations as soon as possible. This is simply a list of candiate topics that various people have informally suggested so far: - Mayavi/TVTK - Advanced topics in matplotlib - Statistics with Scipy - The TimeSeries scikit - Designing scientific interfaces with Traits - Advanced numpy - Sparse Linear Algebra with Scipy - Structured and record arrays in numpy - Cython - Sage - general tutorial - Sage - specific topics, suggestions welcome - Using GPUs with PyCUDA - Testing strategies for scientific codes - Parallel processing and mpi4py - Graph theory with Networkx - Design patterns for efficient iterator-based scientific codes. - Symbolic computing with sympy We'd like to hear from any ideas on other possible topics of interest, and we'll then run a doodle poll to gather quantitative feedback with the final list of candidates. Many thanks, f

1 2

Advice on converting Numarray C extension?
by Russell E. Owen 01 Jul '09

01 Jul '09

I have an old Numarray C extension (or, rather, a Python package containing a C extension) that I would like to convert to numpy (in a way that is likely to be supported long-term). Options I have found include: - Use the new numpy extension. This seems likely to be fast and future-proof. But I am finding the documentation slow going. Does anyone know of a simple example (e.g. read in an array, create a new array)? - Use the Numarray compatible C API. Simple (and takes advantage of the nice Numarray tutorial example for documentation), but will this be supported in the long term? - Switch to ctypes. Simple in concept. But I'm wondering if I can get distutils to build the resulting package. - Use SWIG. I have some experience with it, but not with numpy arrays. - Use Cython to replace the C code. No idea if this is a long-term supported package. Another option is to try to rewrite in pure python. Perhaps the numpy indexing is sophisticated enough to allow an efficient solution. The C extension computes a radial profile from a 2-d masked array: radProf(r)= sum of all unmasked pixels at radius r about some specified center index I can easily generate (and cache) a 2-d array of radius index, but is it possible to use that to efficiently generate the desired sum? Any opinions? -- Russell

3 9

Re: [Numpy-discussion] permutation symbol
by David Goldsmith 30 Jun '09

30 Jun '09

No, but my guess is that it might be useful in numerical General Relativity, e.g., as well as pedagogically in teaching differential geometry and tensor algebra. DG --- On Tue, 6/30/09, Nils Wagner <nwagner(a)iam.uni-stuttgart.de> wrote: > From: Nils Wagner <nwagner(a)iam.uni-stuttgart.de> > Subject: Re: [Numpy-discussion] permutation symbol > To: "Discussion of Numerical Python" <numpy-discussion(a)scipy.org> > Date: Tuesday, June 30, 2009, 11:50 AM > On Tue, 30 Jun 2009 11:30:51 -0700 > (PDT) > David Goldsmith <d_l_goldsmith(a)yahoo.com> > wrote: > > > > Great, Nils! Now, can you generalize it to N-D > for us? > >;-) > > > > DG > > Just curious - Do you have any application for N-D case in > > mind ? > > Nils > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion >

2 1

Re: [Numpy-discussion] permutation symbol
by David Goldsmith 30 Jun '09

30 Jun '09

Great, Nils! Now, can you generalize it to N-D for us? ;-) DG --- On Tue, 6/30/09, Nils Wagner <nwagner(a)iam.uni-stuttgart.de> wrote: > From: Nils Wagner <nwagner(a)iam.uni-stuttgart.de> > Subject: Re: [Numpy-discussion] permutation symbol > To: "Discussion of Numerical Python" <numpy-discussion(a)scipy.org> > Date: Tuesday, June 30, 2009, 11:26 AM > On Tue, 30 Jun 2009 11:10:39 -0600 > Charles R Harris <charlesr.harris(a)gmail.com> > wrote: > > On Tue, Jun 30, 2009 at 10:56 AM, Charles R Harris > < > > charlesr.harris(a)gmail.com> > wrote: > > > >> > >> > >> On Tue, Jun 30, 2009 at 10:40 AM, Nils Wagner > < > >> nwagner(a)iam.uni-stuttgart.de> > wrote: > >> > >>> On Tue, 30 Jun 2009 10:27:05 -0600 > >>> Charles R Harris <charlesr.harris(a)gmail.com> > wrote: > >>> > On Tue, Jun 30, 2009 at 5:11 AM, Nils > Wagner > >>> > <nwagner(a)iam.uni-stuttgart.de>wrote: > >>> > > >>> >> On Tue, 30 Jun 2009 11:22:34 +0200 > >>> >> "Nils Wagner" <nwagner(a)iam.uni-stuttgart.de> > wrote: > >>> >> > >>> >>> Hi all, > >>> >>> > >>> >>> How can I build the following > product with numpy > >>> >>> > >>> >>> q_i = \varepsilon_{ijk} q_{kj} > >>> >>> > >>> >>> where \varepsilon_{ijk} > denotes the permutation > >>>symbol. > >>> >>> > >>> >>> Nils > >>> >>> > >>> >> Sorry for replying to myself. > >>> >> The permutation symbol is also known > as the > >>>Levi-Civita > >>> >>symbol. > >>> >> I found an explicit expression at > >>> >> http://en.wikipedia.org/wiki/Levi-Civita_symbol > >>> >> > >>> >> How do I build the product of the > Levi-Civita symbol > >>> >>\varepsilon_{ijk} and > >>> >> the two dimensional array > >>> >> q_{kj}, i,j,k = 1,2,3 ? > >>> >> > >>> > > >>> > Write it out explicitly. It essentially > >>>antisymmetrizes > >>> >q and the three off > >>> > diagonal elements can then be treated as > a vector. > >>> >Depending on how q is > >>> > formed and the resulting vector is used > there may be > >>> >other things you can do > >>> > when you use it in a more general > expression. If this > >>>is > >>> >part of a general > >>> > calculation there might be other ways of > expressing > >>>it. > >>> > > >>> > Chuck > >>> > >>> Hi Chuck, > >>> > >>> Thank you for your response. > >>> The problem at hand is described in a paper by > Angeles > >>> namely equation (17c) in > >>> "Automatic computation of the screw parameters > of > >>> rigid-body motions. > >>> Part I: Finitely-separated positions" > >>> Journal of Dynamic systems, Measurement and > Control, > >>>Vol. > >>> 108 (1986) pp. 32-38 > >>> > >> > >> You can solve this problem using quaternions also, > in > >>which case it reduces > >> to an eigenvalue problem. You will note that such > things > >>as PCA are used in > >> the papers that reference the cited work so you > can't > >>really get around that > >> bit of inefficiency. > >> > > > > Here's a reference to the quaternion approach: > > http://people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf. > > >You can > > get the translation part from the motion of the > >centroid. > > > > If you are into abstractions you will note that the > >problem reduces to > > minimising a quadratic form in the quaternion > >components. The rest is just > > algebra ;) > > > > Chuck > > It turns out that the product is simply an invariant of a > 3 \times 3 matrix. > > from numpy import array, zeros, identity > from numpy.linalg import norm > > > def vect(A): > """ linear invariant of a 3 x 3 > matrix """ > tmp = zeros(3,float) > tmp[0] = 0.5*(A[2,1]-A[1,2]) > tmp[1] = 0.5*(A[0,2]-A[2,0]) > tmp[2] = 0.5*(A[1,0]-A[0,1]) > > return tmp > > Q = array([[0,0,-1],[-1,0,0],[0,1,0]]) > > q = vect(Q) > print q > > > Nils > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion >

3 2

Re: [Numpy-discussion] permutation symbol
by David Goldsmith 30 Jun '09

30 Jun '09

Hi! I can't help but wonder if - out there somewhere - there's an optimized algorithm for implementing the general N-D Permutation Symbol: "The symbol can be generalized to an arbitrary number of elements, in which case the permutation symbol is (-1)^(i(p)), where i(p) is the number of transpositions of pairs of elements (i.e., permutation inversions) that must be composed to build up the permutation p (Skiena 1990)." (perhaps in the referenced paper - anyone have "free" access to such reprints online) which we could add to NumPy, including it's specializations to the 3-D and 2-D cases as convenience functions. (I'd write it given the algorithm, but combinatorics was never my strong suit, i.e., I'm sure someone else on this list could figure out at least a brute force algorithm more efficiently than yours truly.) Alternatively, isn't this "just" a tensor product with one factor a constant tensor, and thus we could implement it that way using an existing resource? DG --- On Tue, 6/30/09, Charles R Harris <charlesr.harris(a)gmail.com> wrote: > From: Charles R Harris <charlesr.harris(a)gmail.com> > Subject: Re: [Numpy-discussion] permutation symbol > To: "Discussion of Numerical Python" <numpy-discussion(a)scipy.org> > Date: Tuesday, June 30, 2009, 10:10 AM > > > On Tue, Jun 30, 2009 at 10:56 AM, > Charles R Harris <charlesr.harris(a)gmail.com> > wrote: > > > > On > Tue, Jun 30, 2009 at 10:40 AM, Nils Wagner <nwagner(a)iam.uni-stuttgart.de> > wrote: > > > On Tue, 30 Jun 2009 10:27:05 -0600 > > Charles R Harris <charlesr.harris(a)gmail.com> > wrote: > > > On Tue, Jun 30, 2009 at 5:11 AM, Nils Wagner > > > <nwagner(a)iam.uni-stuttgart.de>wrote: > > > > > >> On Tue, 30 Jun 2009 11:22:34 +0200 > > >> "Nils Wagner" <nwagner(a)iam.uni-stuttgart.de> > wrote: > > >> > > >>> Hi all, > > >>> > > >>> How can I build the following product with > numpy > > >>> > > >>> q_i = \varepsilon_{ijk} q_{kj} > > >>> > > >>> where \varepsilon_{ijk} denotes the > permutation symbol. > > >>> > > >>> Nils > > >>> > > >> Sorry for replying to myself. > > >> The permutation symbol is also known as the > Levi-Civita > > >>symbol. > > >> I found an explicit expression at > > >> http://en.wikipedia.org/wiki/Levi-Civita_symbol > > >> > > >> How do I build the product of the Levi-Civita > symbol > > >>\varepsilon_{ijk} and > > >> the two dimensional array > > >> q_{kj}, i,j,k = 1,2,3 ? > > >> > > > > > > Write it out explicitly. It essentially > antisymmetrizes > > >q and the three off > > > diagonal elements can then be treated as a vector. > > >Depending on how q is > > > formed and the resulting vector is used there may be > > >other things you can do > > > when you use it in a more general expression. If this > is > > >part of a general > > > calculation there might be other ways of expressing > it. > > > > > > Chuck > > > > Hi Chuck, > > > > Thank you for your response. > > The problem at hand is described in a paper by Angeles > > namely equation (17c) in > > "Automatic computation of the screw parameters of > > rigid-body motions. > > Part I: Finitely-separated positions" > > Journal of Dynamic systems, Measurement and Control, Vol. > > 108 (1986) pp. 32-38 > > > You can solve this problem using quaternions also, in which > case it reduces to an eigenvalue problem. You will note that > such things as PCA are used in the papers that reference the > cited work so you can't really get around that bit of > inefficiency. > > > > Here's a reference to the quaternion approach: http://people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf. > You can get the translation part from the motion of the > centroid. > > > If you are into abstractions you will note that the problem > reduces to minimising a quadratic form in the quaternion > components. The rest is just algebra ;) > > Chuck > > > > -----Inline Attachment Follows----- > > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion >

1 0

Code for compatibility with gfortran, license issues
by David Cournapeau 30 Jun '09

30 Jun '09

Hi, I started working on a new approach for windows 64 bits support, to be able to combine gfortran and visual studio. Basically, I am reimplementing the needed functions from libgfortran so that it can be built with MS compiler, but I cannot hope to do that without using gfortran sources (under the GPL 3), if only for the signature/API of the functions. Is it ok to include those sources in numpy/scipy, given that they will only be used in those cases where libgfortran would have to be used anyway ? I got pretty far very quickly using this method (full numpy works except for a few unit tests, scipy.linalg, scipy.sparse, scipy.fftpack are already working), so it would be nice if it was possible :) cheers, David

5 6

Deprecating function_base.unique()
by Neil Crighton 29 Jun '09

29 Jun '09

Hi, There was some discussion (e.g. http://article.gmane.org/gmane.comp.python.numeric.general/30629) about changes to the arraysetops module to consolidate the separate unique/non-unique functions and rename setmember1d to in1d. There's a patch that makes these changes in ticket 1133 (http://projects.scipy.org/numpy/ticket/1133) The patch also deprecates function_base.unique and renames arraysetops.unique1d to arraysetops.unique as per this discussion: http://article.gmane.org/gmane.comp.python.numeric.general/31068 Is deprecating function_base.unique the right thing to do? Alternatively, we could deprecate unique1d and transfer its functionality to function_base.unique, but since many of the setops functions use unique, arraysetops seems a more natural place for it. Does anyone have any thoughts about this? Neil

1 0

numpy/numexpr performance (particle simulation)
by Prashant Saxena 29 Jun '09

29 Jun '09

Hi, I am doing a little test using numpy and numexpr to do a particle simulation. I never used either of them much and this is the first time I have to go deeper. Here is the code: import numpy as np import numexpr as nexpr class Particle( object ): def __init__( self, id ): self.position = [0.0, 0.0, 0.0] self.color = [0.5, 0.5, 0.5] self.direction = [0.0, 0.0, 0.0] self..id = id self.lifeSpan = 0 self.type = 0 class Emitter( object ): def __init__( self ): self.particles = np.empty([], dtype=Particle()) self.currentParticleId = 0 self.numParticles = 1000 self.emissionRate = 10 self.position = [0.0, 0.0, 0.0] self.rotation = [0.0, 0.0, 0.0] # Add a single particle def addParticle( self ): """ Add a single particle in the emitter. """ if self.currentParticleId < self.numParticles: self.particles = np.append( self.particles, Particle( self..currentParticleId ) ) self.currentParticleId+=1 ####################################################### Problem 1: self.particles = np.empty([], dtype=Particle()) In "Emitter" class, how do I initialize a numpy array of "Particle" type? Problem 2: If problem 1 can be solved, is it possible to use numexpr to alter the position of each particle's position in emitter.particles array. For example, multiply x,y and z of each particle by using three random values. In other words self.particles[0].position*= random(1.0) self.particles[1].position*= random(2.0) self.particles[2].position*= random(3.0) If problem 1 cannot be solved then may be I could use something like this using simple list: a = [Particle(0), Particle(1)] How do I modify the position of each particle in array"a" using numexpr? I would be using 1-5 million particles where each particle's attribute postion, rotation and color will be modified using some complex equations. The above code is just the begining hence any tips for performance would be greatly apperciated. Regards Prashant Yahoo! recommends that you upgrade to the new and safer Internet Explorer 8. http://downloads.yahoo.com/in/internetexplorer/

3 2

ANN: SciPy 2009 student sponsorship
by Jarrod Millman 29 Jun '09

29 Jun '09

I am pleased to announce that the Python Software Foundation is sponsoring 10 students' travel, registration, and accommodation for the SciPy 2009 conference (Aug. 18-23). The focus of the conference is both on scientific libraries and tools developed with Python and on scientific or engineering achievements using Python. If you're in college or a graduate program, please check out the details here: http://conference.scipy.org/student-funding About the conference -------------------- SciPy 2009, the 8th Python in Science conference, will be held from August 18-23, 2009 at Caltech in Pasadena, CA, USA. The conference starts with two days of tutorials to the scientific Python tools. There will be two tracks, one for introduction of the basic tools to beginners, and one for more advanced tools. The tutorials will be followed by two days of talks. Both days of talks will begin with a keynote address. The first day’s keynote will be given by Peter Norvig, the Director of Research at Google; while, the second keynote will be delivered by Jon Guyer, a Materials Scientist in the Thermodynamics and Kinetics Group at NIST. The program committee will select the remaining talks from submissions to our call for papers. All selected talks will be included in our conference proceedings edited by the program committee. After the talks each day we will provide several rooms for impromptu birds of a feather discussions. Finally, the last two days of the conference will be used for a number of coding sprints on the major software projects in our community. For the 8th consecutive year, the conference will bring together the developers and users of the open source software stack for scientific computing with Python. Attendees have the opportunity to review the available tools and how they apply to specific problems. By providing a forum for developers to share their Python expertise with the wider commercial, academic, and research communities, this conference fosters collaboration and facilitates the sharing of software components, techniques, and a vision for high level language use in scientific computing. For further information, please visit the conference homepage: http://conference.scipy.org. Important Dates --------------- * Friday, July 3: Abstracts Due * Friday, July 10: Announce accepted talks, post schedule * Friday, July 10: Early Registration ends * Tuesday-Wednesday, August 18-19: Tutorials * Thursday-Friday, August 20-21: Conference * Saturday-Sunday, August 22-23: Sprints * Friday, September 4: Papers for proceedings due Executive Committee ------------------- * Jarrod Millman, UC Berkeley, USA (Conference Chair) * Gaël Varoquaux, INRIA Saclay, France (Program Co-Chair) * Stéfan van der Walt, University of Stellenbosch, South Africa (Program Co-Chair) * Fernando Pérez, UC Berkeley, USA (Tutorial Chair)

1 0

stderr
by Charles R Harris 28 Jun '09

28 Jun '09

There are a few spots where messages are printed to stderr. Some of these look almost like debugging stuff, for instance NPY_NO_EXPORT void format_@name@(char *buf, size_t buflen, @name@ val, unsigned int prec) { /* XXX: Find a correct size here for format string */ char format[64], *res; size_t i, cnt; PyOS_snprintf(format, sizeof(format), _FMT1, prec); res = NumPyOS_ascii_format@type@(buf, buflen, format, val, 0); if (res == NULL) { fprintf(stderr, "Error while formatting\n"); return; } /* If nothing but digits after sign, append ".0" */ cnt = strlen(buf); for (i = (val < 0) ? 1 : 0; i < cnt; ++i) { if (!isdigit(Py_CHARMASK(buf[i]))) { break; } } if (i == cnt && buflen >= cnt + 3) { strcpy(&buf[cnt],".0"); } } Do we want to raise an error here? Alternatively, we could use an assert. Chuck

3 3