Mailman 3 February 2008 - NumPy-Discussion

resizing arrays
by Neal Becker 15 Feb '08

15 Feb '08

I have a situation where I'm going to save some data sets to plot, but I don't know a-priori how many sets there will be. I'm using this code: try: shape = list(phase_plots.shape) shape[0] += 1 phase_plots.resize (shape, refcheck=0) except NameError: phase_plots = empty ((1, 2*iterations+1, l)) This works, I'm just wondering if this is a reasonable approach or if maybe something else would be better (or more efficient).

2 1

Re: [Numpy-discussion] String sort
by Francesc Altet 15 Feb '08

15 Feb '08

A Wednesday 13 February 2008, Charles R Harris escrigué: > On Feb 13, 2008 10:56 AM, Francesc Altet <faltet(a)carabos.com> wrote: > > Be warned, I'd like to stress out that these are my figures for my > > _own laptop_. It would be nice if you can verify all of this with > > other achitectures (your Core2 machine seems different enough). I > > can run the benchmarks on Windows (installed in the same laptop) > > too. Tell me if you are interested on me doing this. > > Its easy enough to test if you compile from svn, just add your new > copy function and change the name in this line: > > #copy=copy_string, copy_ucs4# > > to use your function instead of copy_string. I've spoken too fast. I've never compiled NumPy on Windows before, and had problems when trying to compile it using MSVC 7.1 and a recent copy of the repository. Well, in any case, I've exercised the Opteron platform, using gcc 4.1.3 (i.e. the one that can optimize newqsort correctly), and this brings new light to our study. From the plot (attached), it can be drawn the next conclusions: 1) copy_string2 (the combination of manual copy and memcpy) is not better than memcpy for *any* value of the string length in our Opteron platform. Also, the improvements with Pentium4 was not that big (<20%). In consequence, I'd choose to always use memcpy and discard copy_string2. 2) Curiously enough, the indirect sort in Opterons is *always* faster than newqsort+memcpy. For large values of string lengths (> 256), the speed-up can be up to 3x, which is a lot. And I'd say that this keeps true also for most modern processors (read Core2, Barcelona). For older processors (Pentium4), the indirect method can be slower than direct plot for small lengths, but by a very few extent (<10%). Conclusion 2 makes me wonder if it wouldn't be useful the introduction of a new flag in sort, like: """ `indirect` - Use the indirect method for sorting. This requires more memory (4/8 bytes per array element), but for sorting arrays of strings it is almost always faster than the direct approach (default). Beware: this is not the case when using numerical values, where the use of this method for sorting is not recommendable. """ Agreed, that could introduce some confusion, but as this flag would be 'False' by default, not many people should bother about its existence, and can definitely help to people who cares about sorting string performance. Cheers, -- >0,0< Francesc Altet http://www.carabos.com/ V V Cárabos Coop. V. Enjoy Data "-" ------------------------------------------------------- -- >0,0< Francesc Altet http://www.carabos.com/ V V Cárabos Coop. V. Enjoy Data "-"

2 11

MLPY - Machine Learning Py - Python/NumPy based package for machine learning
by Davide Albanese 15 Feb '08

15 Feb '08

*Machine Learning Py* (MLPY) is a *Python/NumPy* based package for machine learning. The package now includes: * *Support Vector Machines* (linear, gaussian, polinomial, terminated ramps) for 2-class problems * *Fisher Discriminant Analysis* for 2-class problems * *Iterative Relief* for feature weighting for 2-class problems * *Feature Ranking* methods based on Recursive Feature Elimination (rfe, onerfe, erfe, bisrfe, sqrtrfe) and Recursive Forward Selection (rfs) * *Input Data* functions * *Confidence Interval* functions Requires Python <http://www.python.org/> >= 2.4 and NumPy <http://www.scipy.org/> >= 1.0.3.* MLPY* is a project of MPBA Group <http://mpa.fbk.eu/> (mpa.fbk.eu) at Fondazione Bruno Kessler (www.fbk.eu). <http://www.fbk.eu/>* MLPY* is free software. It is licensed under the GNU General Public License (GPL) version 3 <http://www.gnu.org/licenses/gpl-3.0.html>. HomePage: mlpy.fbk.eu

3 10

f2py: sharing F90 module data between modules
by Garry Willgoose 14 Feb '08

14 Feb '08

I have a suite of fortran modules that I want to wrap with f2py independently (so they appear to python as seperate imports) but where each module has access to another fortran module (which contains global data that is shared between the suite of fortran modules). I currently compile all the fortran modules, including the global data, in a single f2py statement so all the fortran code gets imported in a single statement import tsim where tsim is the python name of the aggregated modules. What I would like to do is as below, where fm1, fm2 are two example modules that will be called from python, and fm3 has the global data shared between fm1, fm2. <fm1.f90> module fm1 contains subroutine modify(myvalue) use fm3 real :: myvalue value=myvalue end subroutine modify real function result() use fm3 result=value end function result end module fm1 <fm2.f90> module fm2 contains subroutine modify(myvalue) use fm3 real, intent(in) :: myvalue value=myvalue end subroutine modify real function result() use fm3 result=value end function result end module fm2 <fm3.f90> module fm3 real :: value=3.0 end module fm3 The result of running the following script is import fm1,fm2 result1=fm1.fm1.result() result2=fm2.fm2.result() print 'original',result1,result2 x=10.0 fm1.fm1.modify(x) result1=fm1.fm1.result() result2=fm2.fm2.result() print 'fm1 allocate',result1,result2 x=20.0 fm2.fm2.modify(x) result1=fm1.fm1.result() result2=fm2.fm2.result() print 'fm2 allocate',result1,result2 is original 3.0 3.0 fm1 allocate 10.0 3.0 fm2 allocate 10.0 20.0 which clearly shows that fm1, and fm2 have independent images of fm3 (I vaguely recall some time back finding docs that said that f2py/ python loaded each imported into seperate memory and common/modules couldn't be shared between python modules). As I mentioned if fm1,fm2,fm3 are all compiled in a single f2py statement (so they are in a single python module) then this works fine. I guess the question is there any way that I can get fm3 to be shared between fm1 and fm2? The reasons for wanting to do this are because I'm developing a plug-in like architecture for environmental modelling where the user can develop new fortran modules (suitably f2py'ed) that can just be dropped into the module search path but still have access to the global data (subject to fortran module interfaces, etc). ==================================================================== Prof Garry Willgoose, Australian Professorial Fellow in Environmental Engineering, Director, Centre for Climate Impact Management (C2IM), School of Engineering, The University of Newcastle, Callaghan, 2308 Australia. Centre webpage: www.c2im.org.au Phone: (International) +61 2 4921 6050 (Tues-Fri AM); +61 2 6545 9574 (Fri PM-Mon) FAX: (International) +61 2 4921 6991 (Uni); +61 2 6545 9574 (personal and Telluric) Env. Engg. Secretary: (International) +61 2 4921 6042 email: garry.willgoose(a)newcastle.edu.au; g.willgoose(a)telluricresearch.com email-for-life: garry.willgoose(a)alum.mit.edu personal webpage: www.telluricresearch.com/garry ==================================================================== "Do not go where the path may lead, go instead where there is no path and leave a trail" Ralph Waldo Emerson ====================================================================

4 5

test failures
by Charles R Harris 14 Feb '08

14 Feb '08

Hi Stefan, I believe these come from your latest commit. File "/usr/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 25, in test_ufunclike Failed example: nx.sign(a) Expected: array([ 1., -1., 0., 0., 1., -1.]) Got: array([ 1., -1., -1., 0., 1., -1.]) ********************************************************************** File "/usr/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 40, in test_ufunclike Failed example: nx.sign(a, y) Expected: array([ True, True, False, False, True, True], dtype=bool) Got: array([ True, True, True, False, True, True], dtype=bool) ********************************************************************** File "/usr/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 43, in test_ufunclike Failed example: y Expected: array([ True, True, False, False, True, True], dtype=bool) Got: array([ True, True, True, False, True, True], dtype=bool) ............................................................................................................................................... ====================================================================== FAIL: Ticket #588 ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/lib/python2.5/site-packages/numpy/core/tests/test_regression.py", line 734, in check_dot_negative_stride assert_equal(np.dot(x,z),np.dot(x,y2)) File "/usr/lib/python2.5/site-packages/numpy/testing/utils.py", line 143, in assert_equal return assert_array_equal(actual, desired, err_msg) File "/usr/lib/python2.5/site-packages/numpy/testing/utils.py", line 225, in assert_array_equal verbose=verbose, header='Arrays are not equal') File "/usr/lib/python2.5/site-packages/numpy/testing/utils.py", line 217, in assert_array_compare assert cond, msg AssertionError: Arrays are not equal (mismatch 100.0%) x: array([[ 55924.]]) y: array([[ 640000.]]) ====================================================================== FAIL: Ticket #396 ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/lib/python2.5/site-packages/numpy/core/tests/test_regression.py", line 600, in check_poly1d_nan_roots self.failUnlessRaises(np.linalg.LinAlgError,getattr,p,"r") AssertionError: LinAlgError not raised ---------------------------------------------------------------------- Ran 805 tests in 2.983s FAILED (failures=2) Out[1]: <unittest._TextTestResult run=805 errors=0 failures=2> Chuck

2 2

Re: [Numpy-discussion] String sort
by Francesc Altet 14 Feb '08

14 Feb '08

A Thursday 14 February 2008, escriguéreu: > > In any case, if anybody have access to an Opteron machine and gcc > > 4.2.3, it would be great if he can run the benchmark and contribute > > his feedback. > > Here it is with gcc 4.2.3 on an Opteron 246 (2.0 GHz): > > uller:~$ ./sort423_O2 # with -O2 > Benchmark with 1000000 strings of size 15 > C qsort with C style compare: 0.770000 > C qsort with Python style compare: 0.740000 > NumPy newqsort: 0.630000 > > uller:~$ ./sort423_O3 # with -O3 > Benchmark with 1000000 strings of size 15 > C qsort with C style compare: 0.640000 > C qsort with Python style compare: 0.660000 > NumPy newqsort: 0.400000 And here are my timings with gcc 4.1.3 and using a similar Opteron than yours (270 @ 2.0 GHz): With -O2: Benchmark with 1000000 strings of size 15 C qsort with C style compare: 0.750000 C qsort with Python style compare: 0.700000 NumPy newqsort: 0.690000 With -O3: Benchmark with 1000000 strings of size 15 C qsort with C style compare: 0.670000 C qsort with Python style compare: 0.620000 NumPy newqsort: 0.380000 So, it seems clear that the GCC people has fixed in 4.2.3 the problem with the optimizer introduced in 4.2.1. Very good! By the way, it's nice to see the wide range of platforms that this list allows to test out :-) Cheers, -- >0,0< Francesc Altet http://www.carabos.com/ V V Cárabos Coop. V. Enjoy Data "-"

3 3

Proxy array class and units
by Dan Goodman 14 Feb '08

14 Feb '08

Hi all, I'm looking for some advice on how to implement a 'unit checking proxy array'. Background: I'm writing a package to run simulations which make extensive use of linear algebra, for which I'm using numpy. However - it is important to my package that quantities can have dimesions, so I've written a class Quantity subclassed from float which includes the physical dimensions and raises exceptions etc. if you try to add volts to kilograms or whatever. This seems to work fine, but I also wanted arrays containing physical dimensions. I wrote a class qarray subclassed from numpy.ndarray which mostly just copies the functionality of an array with dtype=object. Here is the issue. The code for my package internally doesn't use these qarray objects, it uses numpy arrays because it would obviously be terribly slow to be repeatedly checking if the dimensions were consistent all the time. However, I want to present a (somewhat) consistent strict unit-checking front to the user of this package. As one part of this, I would like a class that appears to act exactly like a qarray (which itself is supposed to act very much like a numpy array) but maintains a connection with its numpy array counterpart in my package. Suppose I have an array X in my package, then I want an object Y which the user is free to manipulate like an array, but if they do something like Y[0]=2*volt then the object Y does (1) check that Y[0] has volt units, (2) if so, update Y[0] and also X[0]. My question is about what is the best way to implement something like this. So far, what I've done is written a proxy array class (call it pqarray say) derived from qarray (which is pretty much just array with dtype=object but not exactly as it does some other things too). When my package creates a pqarray for the user, it does the following: take the underlying numpy array X; create Y a qarray based on X with the correct units (this necessarily copies the data in X because the dtype is different); create Z from X, Y, a pqarray that knows that Y is based on X. The class pqarray simply overrides the __setitem__ method to update itself and X if the units are correct (or raise an exception). Now, this does actually work, but it's not terribly elegant and there is a problem. The standard way that this proxy array will be accessed is through another object, let's just call it P. P has a property val (it has to be a property in general, because the internal numpy array could change at any moment). If you write P.val it returns a pqarray corresponding to the internal representation of the data. Let's suppose P.val has length 1000. Suppose the user writes: for i in range(1000): P.val[i] = i*volt Now what happens is that this creates a new pqarray 1000 times, and this involves copying 1000 elements each time for something of the order of 1M operations total when it should really only be of the order of 1000 operations total. Now, obviously given the way it's written, the user 'ought' to write: val = P.val for i in range(1000): val[i] = i*volt or even better (and what you would do in practice): P.val = qarray(range(1000))*volt and then it will be O(1000) operations, but (a) the user shouldn't need to know this, and (b) my pqarray class seems horrifically ugly and I want something better. It seems to me that my options are: (1) Use my current solution, and instruct the user about the potential problem (which probably won't come up too often anyway). (2) Write pqarray to be a subclass of numpy.ndarray so that no data is copied in creating the pqarray, and give pqarray __getitem__ and __setitem__ methods that return a Quantity with the right units and do unit checking, respectively. I have actually implemented this, and it sort of works but I have some worries about it. For example, if Y is a pqarray of the numpy array X then it seems print Y will give the same result as print X (i.e. it forgets the units). I could presumably look this up and override this behaviour but that's a game which has no end (or at least, a very distant end that might end up with me having reimplemented the whole of numpy.ndarray). (3) Have some sort of caching mechanism, so that P.val returns the cached value unless the underlying numpy array is changed. The problem with this is that you have to keep track of whether or not the underlying numpy array has changed which could be fiddly and rather unsafe. (4) Give up hoping that P.val can be a subclass of an array, make it a new class that keeps a copy of the necessary units and the underlying array, and can do some basic array operations. Since it's not an array type object itself, the user can be told not to expect it to behave like one. This would be quite a sad solution though, because I want the user to be able to write things like mean(P.val) and have them just work as they'd expect (which works at the moment). (5) Something I haven't thought of yet...? Any advice of any sort would be very much appreciated. I did search for relevant information but nothing much came up. Plenty about how to do subclassing of numpy array's, but that's not really my problem here - it's more specific than that. Many thanks in advance, -- Dan Goodman http://thesamovar.net/contact

2 2

String sort
by Francesc Altet 14 Feb '08

14 Feb '08

Hi, I'm a bit confused that the sort method of a string character doesn't allow a mergesort: >>> s = numpy.empty(10, "S10") >>> s.sort(kind="merge") TypeError: desired sort not supported for this type However, by looking at the numpy sources, it seems that the only implemented method for sorting array strings is "merge" (I presume because it is stable). So, perhaps the message above should be fixed. Also, in the context of my work in indexing, and because of the slowness of the current implementation in NumPy, I've ended with an implementation of the quicksort method for 1-D array strings. For moderately large arrays, it is about 2.5x-3x faster than the (supposedly) mergesort version in NumPy, not only due to the quicksort, but also because I've implemented a couple of macros for efficient string swapping and copy. If this is of interest for NumPy developers, tell me and I will provide the code. Cheers, -- >0,0< Francesc Altet http://www.carabos.com/ V V Cárabos Coop. V. Enjoy Data "-"

12 67

Building a static libnumpy
by Jussi Enkovaara 14 Feb '08

14 Feb '08

Hi, I need to use numpy in an environment which does not support shared libraries. Previously, I have used the old Numeric where I created a Makefile to build a static Numeric library which was later on linked to the python interpreter. As far as I understood, numpy uses sort of extended distutils. I was wondering if this extended distutils enables building of static libraries or do I have to go the cumbersome Makefile-way again? Regards, Jussi Enkovaara

5 6

Question for Travis
by Charles R Harris 13 Feb '08

13 Feb '08

Travis, I notice that you used PyDataMem_NEW, PyDimMem_NEW, and friends to allocate memory in the sort routines. Is there a good reason to use these rather than malloc? Chuck

2 1