Mailman 3 April 2016 - NumPy-Discussion

mixed mode arithmetic
by Neal Becker July 11, 2023

July 11, 2023

I've been browsing the numpy source. I'm wondering about mixed-mode arithmetic on arrays. I believe the way numpy handles this is that it never does mixed arithmetic, but instead converts arrays to a common type. Arguably, that might be efficient for a mix of say, double and float. Maybe not. But for a mix of complex and a scalar type (say, CDouble * Double), it's clearly suboptimal in efficiency. So, do I understand this correctly? If so, is that something we should improve?

4 6

Invalid value encoutered : how to prevent numpy.where to do this?
by Eric Emsellem Feb. 18, 2023

Feb. 18, 2023

Dear all, I have a code using lots of "numpy.where" to make some constrained calculations as in: data = arange(10) result = np.where(data == 0, 0., 1./data) # or data1 = arange(10) data2 = arange(10)+1.0 result = np.where(data1 > data2, np.sqrt(data1-data2), np.sqrt(data2-data2)) which then produces warnings like: /usr/bin/ipython:1: RuntimeWarning: invalid value encountered in sqrt or for the first example: /usr/bin/ipython:1: RuntimeWarning: divide by zero encountered in divide How do I avoid these messages to appear? I know that I could in principle use numpy.seterr. However, I do NOT want to remove these warnings for other potential divide/multiply/sqrt etc errors. Only when I am using a "where", to in fact avoid such warnings! Note that the warnings only happen once, but since I am going to release that code, I would like to avoid the user to get such messages which are irrelevant here (because I am testing, with the where, when NOT to divide by zero or take a sqrt of a negative number). thanks! Eric

5 4

Re: [Numpy-discussion] example reading binary Fortran file
by Neil Martinsen-Burrell July 22, 2021

July 22, 2021

David Froger <david.froger.info <at> gmail.com> writes: > Hy,My question is about reading Fortran binary file (oh no this question > again...) I've posted this before, but I finally got it cleaned up for the Cookbook. For this purpose I use a subclass of file that has methods for reading unformatted Fortran data. See http://www.scipy.org/Cookbook/FortranIO/FortranFile. I'd gladly see this in numpy or scipy somewhere, but I'm not sure where it belongs. > program makeArray > implicit none > integer,parameter:: nx=10,ny=20 > real(4),dimension(nx,ny):: ux,uy,p > integer :: i,j > open(11,file='uxuyp.bin',form='unformatted') > do i = 1,nx > do j = 1,ny > ux(i,j) = real(i*j) > uy(i,j) = real(i)/real(j) > p (i,j) = real(i) + real(j) > enddo > enddo > write(11) ux,uy > write(11) p > close(11) > end program makeArray When I run the above program compiled with gfortran on my Intel Mac, I can read it back with:: >>> import numpy as np >>> from fortranfile import FortranFile >>> f=FortranFile('uxuyp.bin', endian='<') >>> uxuy = f.readReals(prec='f') # 'f' for default reals >>> len(uxuy) 400 >>> ux = np.array(uxuy[:200]).reshape((20,10)).T >>> uy = np.array(uxuy[200:]).reshape((20,10)).T >>> p = f.readReals('f').reshape((20,10)).T >>> ux array([[ 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20.], [ 2., 4., 6., 8., 10., 12., 14., 16., 18., 20., 22., 24., 26., 28., 30., 32., 34., 36., 38., 40.], [ 3., 6., 9., 12., 15., 18., 21., 24., 27., 30., 33., 36., 39., 42., 45., 48., 51., 54., 57., 60.], [ 4., 8., 12., 16., 20., 24., 28., 32., 36., 40., 44., 48., 52., 56., 60., 64., 68., 72., 76., 80.], [ 5., 10., 15., 20., 25., 30., 35., 40., 45., 50., 55., 60., 65., 70., 75., 80., 85., 90., 95., 100.], [ 6., 12., 18., 24., 30., 36., 42., 48., 54., 60., 66., 72., 78., 84., 90., 96., 102., 108., 114., 120.], [ 7.,Proxy-Connection: keep-alive Cache-Control: max-age=0 14., 21., 28., 35., 42., 49., 56., 63., 70., 77., 84., 91., 98., 105., 112., 119., 126., 133., 140.], [ 8., 16., 24., 32., 40., 48., 56., 64., 72., 80., 88., 96., 104., 112., 120., 128., 136., 144., 152., 160.], [ 9., 18., 27., 36., 45., 54., 63., 72., 81., 90., 99., 108., 117., 126., 135., 144., 153., 162., 171., 180.], [ 10., 20., 30., 40., 50., 60., 70., 80., 90., 100., 110., 120., 130., 140., 150., 160., 170., 180., 190., 200.]]) >>> uy array([[ 1. , 0.5 , 0.33333334, 0.25 , 0.2 , 0.16666667, 0.14285715, 0.125 , 0.11111111, 0.1 , 0.09090909, 0.08333334, 0.07692308, 0.07142857, 0.06666667, 0.0625 , 0.05882353, 0.05555556, 0.05263158, 0.05 ], [ 2. , 1. , 0.66666669, 0.5 , 0.40000001, 0.33333334, 0.2857143 , 0.25 , 0.22222222, 0.2 , 0.18181819, 0.16666667, 0.15384616, 0.14285715, 0.13333334, 0.125 , 0.11764706, 0.11111111, 0.10526316, 0.1 ], [ 3. , 1.5 , 1. , 0.75 , 0.60000002, 0.5 , 0.42857143, 0.375 , 0.33333334, 0.30000001, 0.27272728, 0.25 , 0.23076923, 0.21428572, 0.2 , 0.1875 , 0.17647059, 0.16666667, 0.15789473, 0.15000001], [ 4. , 2. , 1.33333337, 1. , 0.80000001, 0.66666669, 0.5714286 , 0.5 , 0.44444445, 0.40000001, 0.36363637, 0.33333334, 0.30769232, 0.2857143 , 0.26666668, 0.25 , 0.23529412, 0.22222222, 0.21052632, 0.2 ], [ 5. , 2.5 , 1.66666663, 1.25 , 1. , 0.83333331, 0.71428573, 0.625 , 0.55555558, 0.5 , 0.45454547, 0.41666666, 0.38461539, 0.35714287, 0.33333334, 0.3125 , 0.29411766, 0.27777779, 0.2631579 , 0.25 ], [ 6. , 3. , 2. , 1.5 , 1.20000005, 1. , 0.85714287, 0.75 , 0.66666669, 0.60000002, 0.54545456, 0.5 , 0.46153846, 0.42857143, 0.40000001, 0.375 , 0.35294119, 0.33333334, 0.31578946, 0.30000001], [ 7. , 3.5 , 2.33333325, 1.75 , 1.39999998, 1.16666663, 1. , 0.875 , 0.77777779, 0.69999999, 0.63636363, 0.58333331, 0.53846157, 0.5 , 0.46666667, 0.4375 , 0.41176471, 0.3888889 , 0.36842105, 0.34999999], [ 8. , 4. , 2.66666675, 2. , 1.60000002, 1.33333337, 1.14285719, 1. , 0.8888889 , 0.80000001, 0.72727275, 0.66666669, 0.61538464, 0.5714286 , 0.53333336, 0.5 , 0.47058824, 0.44444445, 0.42105263, 0.40000001], [ 9. , 4.5 , 3. , 2.25 , 1.79999995, 1.5 , 1.28571427, 1.125 , 1. , 0.89999998, 0.81818181, 0.75 , 0.69230771, 0.64285713, 0.60000002, 0.5625 , 0.52941179, 0.5 , 0.47368422, 0.44999999], [ 10. , 5. , 3.33333325, 2.5 , 2. , 1.66666663, 1.42857146, 1.25 , 1.11111116, 1. , 0.90909094, 0.83333331, 0.76923078, 0.71428573, 0.66666669, 0.625 , 0.58823532, 0.55555558, 0.52631581, 0.5 ]]) >>> p array([[ 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21.], [ 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22.], [ 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23.], [ 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24.], [ 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25.], [ 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26.], [ 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27.], [ 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28.], [ 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28., 29.], [ 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28., 29., 30.]]) Note that you have to provide the shape information for ux and uy because fortran writes them together as a stream of 400 numbers. -Neil

10 10

Behavior of .reduceat()
by Jaime Fernández del Río May 23, 2016

May 23, 2016

Two of the oldest issues in the tracker (#834 <https://github.com/numpy/numpy/issues/834> and #835 <https://github.com/numpy/numpy/issues/835>) are about how .reduceat() handles its indices parameter. I have been taking a look at the source code, and it would be relatively easy to modify, the hardest part being to figure out what the exact behavior should be. Current behavior is that np.ufunc.reduceat(x, ind) returns [np.ufunc.reduce(a[ind[i]:ind[i+1]] for i in range(len(ind))] with a couple of caveats: 1. if ind[i] >= ind[i+1], then a[ind[i]] is returned, rather than a reduction over an empty slice. 2. an index of len(ind) is appended to the indices argument, to be used as the endpoint of the last slice to reduce over. 3. aside from this last case, the indices are required to be strictly inbounds, 0 <= index < len(x), or an error is raised The proposed new behavior, with some optional behaviors, would be: 1. if ind[i] >= ind[i+1], then a reduction over an empty slice, i.e. the ufunc identity, is returned. This includes raising an error if the ufunc does not have an identity, e.g. np.minimum. 2. to fully support the "reduction over slices" idea, some form of out of bounds indices should be allowed. This could mean either that: 1. only index = len(x) is allowed without raising an error, to allow computing the full reduction anywhere, not just as the last entry of the return, or 2. allow any index in -len(x) <= index <= len(x), with the usual meaning given to negative values, or 3. any index is allowed, with reduction results clipped to existing values (and the usual meaning for negative values). 3. Regarding the appending of that last index of len(ind) to indices, we could: 1. keep appending it, or 2. never append it, since you can now request it without an error being raised, or 3. only append it if the last index is smaller than len(x). My thoughts on the options: - The minimal, more conservative approach would go with 2.1 and 3.1. And of course 1, if we don't implement that none of this makes sense. - I kind of think 2.2 or even 2.3 are a nice enhancement that shouldn't break too much stuff. - 3.2 I'm not sure about, probably hurts more than it helps at this point, although in a brand new design you probably would either not append the last index or also prepend a zero, as in np.split. - And 3.3 seems too magical, probably not a good idea, only listed it for completeness. Any other thoughts or votes on what, if anything, should we implement, and what the deprecation of current behavior should look like? Jaime -- (\__/) ( O.o) ( > <) Este es Conejo. Copia a Conejo en tu firma y ayúdale en sus planes de dominación mundial.

3 2

Numpy arrays shareable among related processes (PR #7533)
by Matěj Týč May 18, 2016

May 18, 2016

Dear Numpy developers, I propose a pull request https://github.com/numpy/numpy/pull/7533 that features numpy arrays that can be shared among processes (with some effort). Why: In CPython, multiprocessing is the only way of how to exploit multi-core CPUs if your parallel code can't avoid creating Python objects. In that case, CPython's GIL makes threads unusable. However, unlike with threading, sharing data among processes is something that is non-trivial and platform-dependent. Although numpy (and certainly some other packages) implement some operations in a way that GIL is not a concern, consider another case: You have a large amount of data in a form of a numpy array and you want to pass it to a function of an arbitrary Python module that also expects numpy array (e.g. list of vertices coordinates as an input and array of the corresponding polygon as an output). Here, it is clear GIL is an issue you and since you want a numpy array on both ends, now you would have to copy your numpy array to a multiprocessing.Array (to pass the data) and then to convert it back to ndarray in the worker process. This contribution would streamline it a bit - you would create an array as you are used to, pass it to the subprocess as you would do with the multiprocessing.Array, and the process can work with a numpy array right away. How: The idea is to create a numpy array in a buffer that can be shared among processes. Python has support for this in its standard library, so the current solution creates a multiprocessing.Array and then passes it as the "buffer" to the ndarray.__new__. That would be it on Unixes, but on Windows, there has to be a a custom pickle method, otherwise the array "forgets" that its buffer is that special and the sharing doesn't work. Some of what has been said in the pull request & my answer to that: * ... I do see some value in providing a canonical right way to construct shared memory arrays in NumPy, but I'm not very happy with this solution, ... terrible code organization (with the global variables): * I understand that, however this is a pattern of Python multiprocessing and everybody who wants to use the Pool and shared data either is familiar with this approach or has to become familiar with[2, 3]. The good compromise is to have a separate module for each parallel calculation, so global variables are not a problem. * Can you explain why the ndarray subclass is needed? Subclasses can be rather annoying to get right, and also for other reasons. * The shmarray class needs the custom pickler (but only on Windows). * If there's some way to we can paper over the boilerplate such that users can use it without understanding the arcana of multiprocessing, then yes, that would be great. But otherwise I'm not sure there's anything to be gained by putting it in a library rather than referring users to the examples on StackOverflow [1] [2]. * What about telling users: "You can use numpy with multiprocessing. Remeber the multiprocessing.Value and multiprocessing.Aray classes? numpy.shm works exactly the same way, which means that it shares their limitations. Refer to an example: <link to numpy doc>." Notice that although those SO links contain all of the information, it is very difficult to get it up and running for a newcomer like me few years ago. * This needs tests and justification for custom pickling methods, which are not used in any of the current examples. ... * I am sorry, but don't fully understand that point. The custom pickling method of shmarray has to be there on Windows, but users don't have to know about it at all. As noted earlier, the global variable is the only way of using standard Python multiprocessing.Pool with shared objects. [1]: http://stackoverflow.com/questions/10721915/shared-memory-objects-in-python… [2]: http://stackoverflow.com/questions/7894791/use-numpy-array-in-shared-memory… [3]: http://stackoverflow.com/questions/1675766/how-to-combine-pool-map-with-arr…

8 22

Help with bit arrays
by Henrique Almeida April 29, 2016

April 29, 2016

Hello, what's the current status on numpy for loading bit-arrays ? I'm currently unable to correctly load black and white (1-bit) TIFF images. Code example follows: from PIL import Image import numpy from matplotlib import pyplot img = Image.open('oi-00.tiff') a = numpy.array(img) ^ does not work for 1-bit TIFF images PIL source shows that it incorrectly uses typestr == '|b1'. I tried to change this to '|t1', but I get : TypeError: data type "|t1" not understood My goal is to make the above code to work for black and white TIFF images the same way it works for grayscale images. Any help ?

3 8

Possible negative impact of future change to None comparison
by Benjamin Root April 28, 2016

April 28, 2016

Working on closing out some bug reports at work, and I ran into one about comparisons to 'None' will result in elementwise object comparison in the future. Now, I totally get the idea behind the change, and I am not here to argue that decision. However, I have come across a situation where the change might have an unexpected negative consequence. Consider the following: p = Pool(min(len(tiles), maxprocs)) res = p.map(_wrap_tilesum, izip(repeat(args.varname), tiles, repeat(request_start), repeat(args.timeLen), repeat(args.srcdir), repeat(args.tarredInputs), repeat(args.dataset))) p.close() p.join() (tiles, tile_start_dates, tile_end_dates, tile_lons, tile_lats) = zip(*res) if None in tiles: logging.critical("At least one tile was invalid. Aborting") raise Exception("Invalid data retrieved!") Essentailly, in the nominal case, "tiles" would be a list of numpy arrays. However, my error handling is such that if one of my subprocesses errors out, then it returns a None instead of a numpy array. So, all I am doing is testing to see if any of the items in the "tiles" list is None. I have zero desire to actually compare None with the elements in the arrays that happens to be in the list. Of course, I can rewrite this if statement as `any(tile is None for tile in tiles)`, but that isn't my point. Is `if None in tiles:` an unreasonable idiom? Cheers! Ben Root

1 0

Is this a known error in Ipython with numpy?
by Elliot Hallmark April 28, 2016

April 28, 2016

Hello, I haven't worked hard yet to create a minimal runnable (reproduce-able) example but I wanted to check if this sounds familiar to anyone. I have a pretty involved program that resizes arrays in place with arr.resize. When I run it with python it completes and gives the expected result. When I run it in Ipython, I get the following error: ``` ---> 43 self._buffer.resize((count,)+self._dim) 44 45 ValueError: cannot resize an array that references or is referenced by another array in this way. Use the resize function ``` It was consistently doing this on the same array, after resizing two others before hand, even after rebooting. But trying to track it down it goes away even if I undo everything I did to try and track it down. Does this sound familiar?

3 3

nan version of einsum
by Chris Barber April 28, 2016

April 28, 2016

Is there any interest in a nan-ignoring version of einsum a la nansum, nanprod, etc? Any idea how difficult it would be to implement? - Chris

1 1

(no subject)
by Saumyajit Dey April 27, 2016

April 27, 2016

Hi, This is Saumyajit Dey and I am looking forward to start contributing to NumPy. I have never contributed to any open source projects before so I would want to know some tips and guidelines to start contributing. Regards, Saumyajit Saumyajit Dey Junior Undergraduate Student: Department of Computer Science and Engineering National Institute of Technology Warangal (NITW), India Cell: +91-8885847028

4 5