Mailman 3 August 2012 - NumPy-Discussion

__array_wrap__ and dot
by Aronne Merrelli 22 Aug '12

22 Aug '12

Hello list, I posted an issue many months ago [1] about confusing __array_wrap__ behavior in a subclass of np.matrix. Since that time I wasn't really using the code that used the matrix subclass, but lately I have starting using the code so I've run into the issue again. Looking into a little more detail, I think the root of the problem is how __array_wrap__ interacts with the dot functions. Operations like M*2 will be processed differently depending on whether M is a ndarray or matrix subclass; it appears that a matrix will always call np.dot for that operation, while an ndarray will go through np.multiply. Since dot acts differently than multiply here, you get different results between the two subclasses. I updated the gist to show the problem [2]. The output I get on my system is pasted at the bottom. It looks like __array_wrap__ gets called sometimes, depending on the type of one of the arguments. I'm not seeing any logical pattern to it, which makes it confusing. My expectation (perhaps incorrect) is that __array_wrap__ should be invoked any time dot is called. I guess that depends on whether dot should be considered a ufunc - I don't know the answer to that question. If someone with more knowledge of the internals could take a look at this, that would be great - I think this is a bug, but I'm not really sure what is the expected behavior here. Thanks, Aronne [1] http://mail.scipy.org/pipermail/numpy-discussion/2011-December/059665.html [2] https://gist.github.com/1511354 Output from run_test (see the gist code) using ipython: In [1]: import matrix_array_wrap_test as mtest; import numpy as np In [2]: np.__version__ Out[2]: '1.6.1' In [3]: np.dot Out[3]: <function numpy.core._dotblas.dot> In [4]: mtest.run_test() aw called? | creation | use Yes | o=ArrSubClass(2) | o2 = o + 2 Yes | o=ArrSubClass(2) | o2 = o + 2.0 Yes | o=ArrSubClass(2) | o2 = o * 2 Yes | o=ArrSubClass(2) | o2 = o * 2.0 Yes | o=ArrSubClass(2) | o2 = o * o.T No | o=ArrSubClass(2) | o2 = np.dot(o,2) No | o=ArrSubClass(2) | o2 = np.dot(o,2.0) No | o=ArrSubClass(2) | o2 = np.dot(o,o.T) Yes | o=ArrSubClass(2) | o2 = np.core.multiarray.dot(o,2) Yes | o=ArrSubClass(2) | o2 = np.core.multiarray.dot(o,2.0) No | o=ArrSubClass(2) | o2 = np.core.multiarray.dot(o,o.T) Yes | o=MatSubClass([1, 1]) | o2 = o + 2 Yes | o=MatSubClass([1, 1]) | o2 = o + 2.0 Yes | o=MatSubClass([1, 1]) | o2 = o * 2 No | o=MatSubClass([1, 1]) | o2 = o * 2.0 No | o=MatSubClass([1, 1]) | o2 = o * o.T Yes | o=MatSubClass([1, 1]) | o2 = np.dot(o,2) No | o=MatSubClass([1, 1]) | o2 = np.dot(o,2.0) No | o=MatSubClass([1, 1]) | o2 = np.dot(o,o.T) Yes | o=MatSubClass([1, 1]) | o2 = np.core.multiarray.dot(o,2) Yes | o=MatSubClass([1, 1]) | o2 = np.core.multiarray.dot(o,2.0) No | o=MatSubClass([1, 1]) | o2 = np.core.multiarray.dot(o,o.T) Yes | o=MatSubClass([1.0, 1.0]) | o2 = o + 2 Yes | o=MatSubClass([1.0, 1.0]) | o2 = o + 2.0 No | o=MatSubClass([1.0, 1.0]) | o2 = o * 2 No | o=MatSubClass([1.0, 1.0]) | o2 = o * 2.0 No | o=MatSubClass([1.0, 1.0]) | o2 = o * o.T No | o=MatSubClass([1.0, 1.0]) | o2 = np.dot(o,2) No | o=MatSubClass([1.0, 1.0]) | o2 = np.dot(o,2.0) No | o=MatSubClass([1.0, 1.0]) | o2 = np.dot(o,o.T) Yes | o=MatSubClass([1.0, 1.0]) | o2 = np.core.multiarray.dot(o,2) Yes | o=MatSubClass([1.0, 1.0]) | o2 = np.core.multiarray.dot(o,2.0) No | o=MatSubClass([1.0, 1.0]) | o2 = np.core.multiarray.dot(o,o.T)

1 0

[ANN] carray 0.5 released
by Francesc Alted 21 Aug '12

21 Aug '12

Announcing carray 0.5 ===================== What's new ---------- carray 0.5 supports completely transparent storage on-disk in addition to memory. That means that everything that can be done with an in-memory container can be done using the disk instead. The advantages of a disk-based container is that your addressable space is much larger than just your available memory. Also, as carray is based on a chunked and compressed data layout based on the super-fast Blosc compression library, and the different cache levels existing in both modern operating systems and the internal carray machinery, the data access speed is very good. The format chosen for the persistence layer is based on the 'bloscpack' library (thanks to Valentin Haenel for his inspiration) and described in 'persistence.rst', although not everything has been implemented yet. You may want to contribute by proposing enhancements to it. See: https://github.com/FrancescAlted/carray/wiki/PersistenceProposal CAVEAT: The bloscpack format is still evolving, so don't trust on forward compatibility of the format, at least until 1.0, where the internal format will be declared frozen. For more detailed info, see the release notes in: https://github.com/FrancescAlted/carray/wiki/Release-0.5 What it is ---------- carray is a chunked container for numerical data. Chunking allows for efficient enlarging/shrinking of data container. In addition, it can also be compressed for reducing memory/disk needs. The compression process is carried out internally by Blosc, a high-performance compressor that is optimized for binary data. carray can use numexpr internally so as to accelerate many vector and query operations (although it can use pure NumPy for doing so too). numexpr can use optimize the memory usage and use several cores for doing the computations, so it is blazing fast. Moreover, with the introduction of a carray/ctable disk-based container (in version 0.5), it can be used for seamlessly performing out-of-core computations. carray comes with an exhaustive test suite and fully supports both 32-bit and 64-bit platforms. Also, it is typically tested on both UNIX and Windows operating systems. Resources --------- Visit the main carray site repository at: http://github.com/FrancescAlted/carray You can download a source package from: http://carray.pytables.org/download Manual: http://carray.pytables.org/docs/manual Home of Blosc compressor: http://blosc.pytables.org User's mail list: carray(a)googlegroups.com http://groups.google.com/group/carray ---- Enjoy! -- Francesc Alted

1 0

A step toward merging odeint and ode
by Fabrice Silva 21 Aug '12

21 Aug '12

I made a pull request [1] to integrate the LSODA solver that is used in odeint into the modular scipy.integrate.ode generic class. In a similar way as for vode, it just wraps the already present lsoda.f file (see .pyf file) and exposes it within an IntegratorBase subclass adjusting the coefficients before calling lsoda. Note that lsoda provide automatic switching between stiff and non-stiff methods, feature that is not present in the available vode integrator. Final note: tests are ok! Regards, [1] https://github.com/scipy/scipy/pull/273 -- Fabrice Silva

2 6

Preventing lossy cast for new float dtypes ?
by David Cournapeau 20 Aug '12

20 Aug '12

Hi, I have started toying with implementing a quad precision dtype for numpy on supported platforms, using the __float128 + quadmath lib from gcc. I have noticed invalid (and unexpected) downcast to long double in some cases, especially for ufuncs (e.g. when I don't define my own ufunc for a given operation). Looking down in numpy ufunc machinery, I can see that the issue is coming from the assumption that long double is the highest precision possible for a float type, and the only way I can 'fix' this is to define kind to a value != 'f' in my dtype definition (in which case I get an expected invalid cast exception). Is there a way to still avoid those casts while keeping the 'f' kind ? thanks, David

2 1

Difference between np.loadtxt depending on whether you supply a file object or a filename
by Andrew Nelson 20 Aug '12

20 Aug '12

Dear list, I observe a difference when I try to load a 2D numpy array from a file object compared to if I supply a filename viz: >>> np.version.version '1.5.1' >>> f=open('fit_theoretical.txt') >>> a=np.loadtxt(f) >>> a.shape (1000,) >>> a=np.loadtxt('fit_theoretical.txt') >>> a.shape (500, 2) This strikes me as unexpected, it's not a documented behaviour. Any ideas? cheers, Andrew. -- _____________________________________ Dr. Andrew Nelson _____________________________________

2 1

Functions for stride manipulations
by Sebastian Berg 19 Aug '12

19 Aug '12

Hey, Inspired by an existing PR into numpy, I created two functions based on stride_tricks which I thought might be useful inside numpy. So if I get any feedback/pointers, I would add some tests and create a PR. The first function rolling_window is to create for every point of the original array, a view of its neighbourhood in arbitrary dimensions with some options to make it easy to use. For example for a 3 dim array: rolling_window(a, (3,3,3), asteps=(3,3,3)).max((3,4,5)) Gives the maximum for all non-overlapping 3x3x3 subarrays and: rolling_window(a, 3, axes=0).max(-1) would create the rolling maximum over all 3 successive values along the 0th axis. Plus a function permute_axes which allows to give a tuple for switching axes and can also merge them. A (10, 4, 3) shaped array would give: permute_axes(a, (2, 1, 0)).shape -> (3, 4, 10) permute_axes(a, (1, (0, 2))).shape -> (4, 30) A copy is created when necessary, this basically allows to do multiple swap_axes and a reshape in a single call. I hope this might be useful... Regards, Sebastian

1 0

Re: [Numpy-discussion] 64bit infrastructure
by Sayth Renshaw 17 Aug '12

17 Aug '12

> > Hi Sayth, > > On Thu, Aug 16, 2012 at 8:51 PM, Sayth Renshaw <flebber.crue(a)gmail.com> wrote: >> Hi >> >> I was just wondering if the current absence of 64bit builds was as a result >> of an infrastructure or funding concern. >> >> I like many now have only 64 bit systems and am using the unofficial builds. >> >> If it was something we could raise funds to resolve, how much? Or is it more >> complex than that? > > Which platform, windows? > > Ondrej > > Yes I have a Windows 7 64Bit install and a ubuntu 64 Bit install. Sayth

2 1

64bit infrastructure
by Sayth Renshaw 17 Aug '12

17 Aug '12

Hi I was just wondering if the current absence of 64bit builds was as a result of an infrastructure or funding concern. I like many now have only 64 bit systems and am using the unofficial builds. If it was something we could raise funds to resolve, how much? Or is it more complex than that? Sayth

3 2

Correlated distributions (?)
by Andrea Gavana 17 Aug '12

17 Aug '12

Hi All, once again, my apologies for a (possibly) very ignorant question, my google-fu is failing me... also because I am not sure of what exactly I should look for. My problem is relatively simple. Let's assume I have two Python objects, A and B, and one of their attributes can assume a value of "True" or "False" depending on the results of a uniform random distribution sample, i.e.: probability_A = 0.95 probability_B = 0.86 A.has_failed = False B.has_failed = False if numpy.random.random() < probability_A: A.has_failed = True if numpy.random.random() < probability_B: B.has_failed = True Now, I know that there is a correlation factor between the failing/not failing of A and the failing/not failing of B. Specifically, If A fails, then B should have 80% more chance of failing, but I have been banging my head to find out how I should modify the "probability_B" number (or the extremes of the uniform distribution, if that makes sense) in order to reflect that correlation. I have been looking at correlated distributions, but it appears that most of the results I have found relate to normal distributions, there is very little about non-normal (and especially uniform) distributions. It's also very likely that I am not looking in the right direction, so I would appreciate any suggestion you may share. Thank you in advance. Andrea. "Imagination Is The Only Weapon In The War Against Reality." http://xoomer.alice.it/infinity77/

4 3

Enhanced array destruction needed for C flexible numeric type
by Lee McCuller 17 Aug '12

17 Aug '12

I've been implementing a native C++ code with all of the functionality of the Python Uncertanties package (http://packages.python.org/uncertainties/). Automatic error propagation (and differentiation) is extremely useful, but native speed is needed to use in things like function fitting. I have created the type and successfully bound it to numpy (though still missing a lot of UFuncs), but I have a hangup in that I have no clean way to call a destructor on the numeric type (it is variable length via a held pointer). I could keep it in the array as a python object and get xdecref to do it, but I would much rather A: not have the overhead and B: I really want to keep the GIL free. All of the other PyArray_ArrFuncs which manage the type appear to allow pairing a destructor when memory is overwritten, and constructors are all well supported (like copyswap, setitem and fillwithscalar). I'd like to propose adding a flag for PyArray_Descr.flags NPY_ITEM_USES_DESTRUCTOR and a function PyArray_ArrFuncs.destructor with void destruct(void *data, npy_intp dstride, npy_intp len, void *arr) where data points into the data array to be affected, dstride is the data stride, len is number of affected elements, and arr is the array object containing data. It may be possible to change the call structure of another element when this flag is set (to preserve ABI compatibility for existing types). This flag would by inherited when NPY_FROM_FIELDS is set. Subclassing PyArray_Descr may also allow preservation of the API. The current implementation of NPY_ITEM_REFCOUNT could potentially use this more-general destructor implementation to call xdecref. I could potentialy implement this in a branch, but I'd like to know what other's think of it as well as how best to do it to minimize its effect on the rest of the code and external libraries. For anyone interested in the details of this uncertainties code - My error propagation class are a C++ numeric type, which carry around a list of partial derivatives between new values and all of the values used to compute them. This list is a C++ vector which needs its destructor called. My bindings are hand-written, but should provide a good template for binding to other C++ numeric types. For instance, it would not be hard to extend them also add numpy support for arbitrary precision math types, which would also need their C++ destructors called.

1 0