Mailman 3 August 2011 - NumPy-Discussion

3D array -> 2D array issue...
by Fred 06 Sep '11

06 Sep '11

Hi all, Let's say C1 is a 3D array, and q0 and k are 2D array. dim C1 = nx*ny*nz dim q0 = nx*ny = dim k I have to do the following: q0[0, 0] = C1[0, 0, k[0, 0]] q0[1, 1] = C1[1, 1, k[1, 1]] ... q0[i, j] = C1[i, j, k[i, j]] ... I tried q0 = C1[:, :, k] but this obviously does not work. How could I do this ala NumPy? TIA Cheers, -- Fred

5 11

Efficient way to load a 1Gb file?
by Russell E. Owen 06 Sep '11

06 Sep '11

A coworker is trying to load a 1Gb text data file into a numpy array using numpy.loadtxt, but he says it is using up all of his machine's 6Gb of RAM. Is there a more efficient way to read such text data files? -- Russell

7 13

Eigenvalues did not converge
by Rick Muller 04 Sep '11

04 Sep '11

I'm bumping into the old "Eigenvalues did not converge" error using numpy.linalg.eigh() on several different linux builds of numpy (1.4.1). The matrix is 166x166. I can compute the eigenvalues on a Macintosh build of numpy, and I can confirm that there aren't degenerate eigenvalues, and that the matrix appears to be negative definite. I've seen this before (though not for several years), and what I normally do is to build lapack with -O0. This trick did not work in the current instance. Does anyone have any tricks to getting eigh to work? Other weird things that I've noticed about this case: I can compute the eigenvalues using eigvals and eigvalsh, and can compute the eigenvals/vecs using eig(). The matrix is real symmetric, and I've tested that it's symmetric enough by forcibly symmetrizing it. Thanks in advance for any help you can offer. -- Rick Muller rpmuller(a)gmail.com 505-750-7557

3 5

converting standard array to record array
by Bryce Ready 04 Sep '11

04 Sep '11

Hello all, So i'm using numpy 1.6.0, and trying to convert a (4,4) numpy array of dtype 'f8' into a record array of this dtype: dt = np.dtype([('mat','(4,4)f8')]) > Here is the code snippet: In [21]: a = np.random.randn(4,4) > > In [22]: a.view(dt) > and the resulting error: ValueError: new type not compatible with array. > Can anyone shed some light for me on why this conversion is not possible? It is certainly technically possible, since the memory layout of the two arrays should be the same. Can anyone recommend a better way to do this conversion? Thanks in advance! -Bryce Ready

2 2

PyBUF_SIMPLE/PyBUF_FORMAT: casts to unsigned bytes
by Stefan Krah 02 Sep '11

02 Sep '11

Hello, PEP-3118 presumably intended that a PyBUF_SIMPLE request should cast the original buffer's data type to 'B' (unsigned bytes). Here is a one-dimensional example that currently occurs in Lib/test/test_multiprocessing: >>> import array, io >>> a = array.array('i', [1,2,3,4,5]) >>> m = memoryview(a) >>> m.format 'i' >>> buf = io.BytesIO(bytearray(5*8)) >>> buf.readinto(m) buf.readinto() calls PyObject_AsWriteBuffer(), which requests a simple buffer from the memoryview, thus casting the 'i' data type to the implied type 'B'. The consumer can see that a cast has occurred because the new buffer's format field is NULL. This seems fine for the one-dimensional case. Numpy currently also allows such casts for multidimensional contiguous and non-contiguous arrays. See below for the examples; I don't want to distract from the main point of the post, which is this: I'm seeking a clear specification for the Python documentation that determines under what circumstances casts to 'B' should succeed. I'll formulate the points as statements for clarity, but in fact they are also questions: 1) An exporter of a C-contiguous array with ndim <= 1 MUST honor a PyBUF_SIMPLE request, setting format, shape and strides to NULL and itemsize to 1. As a corner case, an array with ndim = 0, format = "L" (or other) would also morph into a buffer of unsigned bytes. test_ctypes currently makes use of this. 2) An exporter of a C-contiguous buffer with ndim > 1 MUST honor a PyBUF_SIMPLE request, setting format, shape, and strides to NULL and itemsize to 1. 3) An exporter of a buffer that is not C-contiguous MUST raise BufferError in response to a PyBUF_SIMPLE request. Why am I looking for such rigid rules? The problem with memoryview is that it has to act as a re-exporter itself. For several reasons (performance of chained memoryviews, garbage collection, early release, etc.) it has been decided that the new memoryview object has a managed buffer that takes a snapshot of the original exporter's buffer (See: http://bugs.python.org/issue10181). Now, since getbuffer requests to the memoryview object cannot be redirected to the original object, strict rules are needed for memory_getbuf(). Could you agree with these rules? Point 2) isn't clear from the PEP itself. I assumed it because Numpy currently allows it, and it appears harmless. Stefan Krah Examples: ========= Cast a multidimensional contiguous array: ----------------------------------------- I think itemsize in the result should be 1. [_testbuffer.ndarray is from http://hg.python.org/features/pep-3118#memoryview] >>> from _testbuffer import * >>> from numpy import * >>> from _testbuffer import ndarray as pyarray >>> >>> exporter = ndarray(shape=[3,4], dtype="L") # Issue a PyBUF_SIMPLE request to 'exporter' and act as a re-exporter: >>> x = pyarray(exporter, getbuf=PyBUF_SIMPLE) >>> x.len 96 >>> x.shape () >>> x.strides () >>> x.format '' >>> x.itemsize # I think this should be 1, not 8. 8 Cast a multidimensional non-contiguous array: --------------------------------------------- This is clearly not right, since y.buf points to a location that the consumer cannot handle without shape and strides. >>> nd = ndarray(buffer=bytearray(96), shape=[3,4], dtype="L") [182658 refs] >>> exporter = nd[::-1, ::-2] [182661 refs] >>> exporter array([[0, 0], [0, 0], [0, 0]], dtype=uint64) [182659 refs] >>> y = pyarray(exporter, getbuf=PyBUF_SIMPLE) [182665 refs] >>> y.len 48 [182666 refs] >>> y.strides () [182666 refs] >>> y.shape () [182666 refs] >>> y.format '' [182666 refs] >>> y.itemsize 8 [182666 refs]

2 3

A question about dtype syntax
by Marquette Jean-Baptiste 01 Sep '11

01 Sep '11

Hi all, I have this piece of code: Stats = [CatBase, round(stats.mean(Data.Ra), 5), round(stats.mean(Data.Dec), 5), len(Sep), round(stats.mean(Sep),4), round(stats.stdev(Sep),4)] print Stats if First: StatsAll = np.array(np.asarray(Stats), dtype=('a11, f8, f8, i4, f8, f8')) First = False else: StatsAll = np.vstack((StatsAll, np.asarray(Stats))) print len(StatsAll) This yields the error: ['bs3000k.cat', 280.60341, -7.09118, 9480, 0.2057, 0.14] Traceback (most recent call last): File "/Users/marquett/workspace/Distort/src/StatsSep.py", line 40, in <module> StatsAll = np.array(np.asarray(Stats), dtype=('a11, f8, f8, i4, f8, f8')) ValueError: could not convert string to float: bs3000k.cat What's wrong ? Thanks for your help Cheers JB

4 9

non-uniform discrete sampling with given probabilities (w/ and w/o replacement)
by Christopher Jordan-Squire 31 Aug '11

31 Aug '11

In numpy, is there a way of generating a random integer in a specified range where the integers in that range have given probabilities? So, for example, generating a random integer between 1 and 3 with probabilities [0.1, 0.2, 0.7] for the three integers? I'd like to know how to do this without replacement, as well. If the probabilities are uniform, there are a number of ways, including just shuffling the data and taking the first however-many elements of the shuffle. But this doesn't apply with non-uniform probabilities. Similarly, one could try arbitrary-sampling-method X (such as inverse-cdf sampling) and then rejecting repeats. But that is clearly sub-optimal if the number of samples desired is near the same order of magnitude as the total population, or if the probabilities are very skewed. (E.g. a weighted sample of size 2 without replacement from [0,1,2] with probabilities [0.999,.00005, 0.00005] will take a long time if you just sample repeatedly until you have two distinct samples.) I know parts of what I want can be done in scipy.statistics using a discrete_rv or with the python standard library's random package. I would much prefer to do it only using numpy because the eventual application shouldn't have a scipy dependency and should use the same random seed as numpy.random. (For more background, what I want is to create a function like sample in R, where I can give it an array-like of doo-hickeys and another array-like of probabilities associated with each doo-hickey, and then generate a random sample of doo-hickeys with those probabilities. One step for that is generating ints, to use as indices, with the same probabilities. I'd like a version of this to be in numpy/scipy, but it doesn't really belong in scipy since it doesn't -Chris JS

3 5

Question on LinAlg Inverse Algorithm
by Mark Janikas 31 Aug '11

31 Aug '11

Hello All, Last week I posted a question involving the identification of linear dependent columns of a matrix... but now I am finding an interesting result based on the linalg.inv() function... sometime I am able to invert a matrix that has linear dependent columns and other times I get the LinAlgError()... this suggests that there is some kind of random component to the INV method. Is this normal? Thanks much ahead of time, MJ Mark Janikas Product Developer ESRI, Geoprocessing 380 New York St. Redlands, CA 92373 909-793-2853 (2563) mjanikas(a)esri.com<mailto:mjanikas@esri.com>

5 11

Numpy performance boost
by Dieter Weber 31 Aug '11

31 Aug '11

Hi, just wanted to show an example of how python3 + numpy compares with just python3 and many other languages and language implementations: http://shootout.alioth.debian.org/u64q/performance.php?test=mandelbrot#about The python3 program using numpy is #6 and you find it with the "interesting alternative" programs on the bottom because it was disqualified for doing things differently. It is 6.3x slower than the fastest program and well ahead of all other interpreted languages. Thanks to all contributors for making numpy such a great piece of software! Greetings, Dieter

3 3

Model Predictive Control package
by Davide 31 Aug '11

31 Aug '11

Dear List, Does anybody knows if there is a python package for simulating LTI dynamic systems controlled with a model predictive controller? I am writing some code which does the job, but the math is not super-easy and i would not like to reinvent the wheel and loose to much time. I will soon publish such codes somewhere, i.e Github, so anyone interested can pick it up. Cheers, Davide

1 0

NumPy-Discussion August 2011