Mailman 3 March 2014 - 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

C-API: multidimensional array indexing?
by Johann Bauer March 31, 2016

March 31, 2016

Dear experts, is there a C-API function for numpy which implements Python's multidimensional indexing? Say, I have a 2d-array PyArrayObject * M; and an index int i; how do I extract the i-th row or column M[i,:] respectively M[:,i]? I am looking for a function which gives again a PyArrayObject * and which is a view to M (no copied data; the result should be another PyArrayObject whose data and strides points to the correct memory portion of M). I searched the API documentation, Google and mailing lists for quite a long time but didn't find anything. Can you help me? Thanks, Johann

3 3

taking a 2D uneven surface slice
by Bryan Woods Dec. 3, 2014

Dec. 3, 2014

I'm trying to do something that at first glance I think should be simple but I can't quite figure out how to do it. The problem is as follows: I have a 3D grid Values[Nx, Ny, Nz] I want to slice Values at a 2D surface in the Z dimension specified by Z_index[Nx, Ny] and return a 2D slice[Nx, Ny]. It is not as simple as Values[:,:,Z_index]. I tried this: >>> values.shape (4, 5, 6) >>> coords.shape (4, 5) >>> slice = values[:,:,coords] >>> slice.shape (4, 5, 4, 5) >>> slice = np.take(values, coords, axis=2) >>> slice.shape (4, 5, 4, 5) >>> Obviously I could create an empty 2D slice and then fill it by using np.ndenumerate to fill it point by point by selecting values[i, j, Z_index[i, j]]. This just seems too inefficient and not very pythonic.

3 3

building 32 bit numpy on 64 bit linux
by David Jones May 19, 2014

May 19, 2014

I'm trying to compile 32-bit numpy on a 64 bit Centos 6 system, but fails with the message: "Broken toolchain: cannot link a simple C program" It get's the compile flags right, but not the linker: C compiler: gcc -pthread -fno-strict-aliasing -g -O2 -m32 -DNDEBUG -g -fwrapv -O3 -Wall -Wstrict-prototypes -fPIC compile options: '-Inumpy/core/src/private -Inumpy/core/src -Inumpy/core -Inumpy/core/src/npymath -Inumpy/core/src/multiarray -Inumpy/core/src/umath -Inumpy/core/src/npysort -Inumpy/core/include -I/opt/python/ia32/include/python2.7 -c' gcc: _configtest.c gcc -pthread _configtest.o -o _configtest _configtest.o: could not read symbols: File in wrong format collect2: ld returned 1 exit status I'm bulding it using a 32bit python build that I compiled on the same system. I tried: OPT="-m32" FOPT="-m32" python setup.py build and setarch x86_64 -B python setup.py build But with the same results. Someone worked around this by altering ccompiler.py, but I'm trying to find a cleaner solution. See:http://stackoverflow.com/questions/11265057/how-do-i-install-a-32-bit-v… Is there a standard way of doing this? Regards, David J.

3 13

Dates and times and Datetime64 (again)
by Sankarshan Mudkavi April 28, 2014

April 28, 2014

Hey all, It's been a while since the last datetime and timezones discussion thread was visited (linked below): http://thread.gmane.org/gmane.comp.python.numeric.general/53805 It looks like the best approach to follow is the UTC only approach in the linked thread with an optional flag to indicate the timezone (to avoid confusing applications where they don't expect any timezone info). Since this is slightly more useful than having just a naive datetime64 package and would be open to extension if required, it's probably the best way to start improving the datetime64 library. If we do wish to have full timezone support it would very likely lead to performance drops (as reasoned in the thread) and we would need to have a dedicated, maintained tzinfo package, at which point it would make much more sense to just incorporate the pytz library. (I also don't have the expertise to implement this, so I would be unable to help resolve the current logjam) I would like to start writing a NEP for this followed by implementation, however I'm not sure what the format etc. is, could someone direct me to a page where this information is provided? Please let me know if there are any ideas, comments etc. Cheers, Sankarshan

10 49

f2py links extensions to incorrect python installation on OSX / Anaconda
by Alex Goodman April 21, 2014

April 21, 2014

Hi all, I have used f2py in the past on a Linux machine with virtually no issues. However on my Mac, I get the following error when importing an f2py generated extension: Fatal Python error: PyThreadState_Get: no current thread Abort trap: 6 After doing some research I found out that the extension is linked to the wrong python installation: otool -L add.so add.so: ./add.so (compatibility version 0.0.0, current version 0.0.0) /System/Library/Frameworks/Python.framework/Versions/2.7/Python (compatibility version 2.7.0, current version 2.7.2) /usr/lib/libSystem.B.dylib (compatibility version 1.0.0, current version 169.3.0) This seems odd because I am using the f2py executable included in Anaconda 1.9.1. I can easily fix this problem by manually using install_name_tool -change on the extension to link the correct library location, but this is really cumbersome. Is there an alternative solution, such as an additional command-line argument when invoking f2py? For what it is worth, I am also using Version 14.0.2 of the Intel Fortran Compiler. Thanks, Alex -- Alex Goodman Graduate Research Assistant Department of Atmospheric Science Colorado State University

4 4

Default builds of OpenBLAS development branch are now fork safe
by Nathaniel Smith April 6, 2014

April 6, 2014

Hey all, Just a heads up: thanks to the tireless work of Olivier Grisel, the OpenBLAS development branch is now fork-safe when built with its default threading support. (It is still not thread-safe when built using OMP for threading and gcc, but this is not the default.) Gory details: https://github.com/xianyi/OpenBLAS/issues/294 Check it out - if it works you might want to consider lobbying your favorite distro to backport it. -n

10 61

Resolving the associativity/precedence debate for @
by Nathaniel Smith April 6, 2014

April 6, 2014

Hi all, After 88 emails we don't have a conclusion in the other thread (see [1] for background). But we have to come to some conclusion or another if we want @ to exist :-). So I'll summarize where the discussion stands and let's see if we can find some way to resolve this. The fundamental question is whether a chain like (a @ b @ c) should be evaluated left-to-right (left-associativity) or right-to-left (right-associativity). DATA SOURCE 1: This isn't a democratic vote, but it's useful to get a sense of people's intuitions. Counting messages in the other thread, opinion seems to be pretty evenly split: == "Votes" for right-associativity == Weak-right: [2] [3] [5] Tight-right: [4] [6] Same-right: [11] == "Votes" for left-associativity == Same-left: [7] [8] [14] [15] [16] Tight-left: [9] Weak-left: [12] There's also the "grouping" option (described in [10]), but that's received very little support (just [13]). DATA SOURCE 2: Several people have suggested that performance considerations mean that right-to-left evaluation is more common in practice than left-to-right evaluation. But, if we look at actual usage in Python code, that's not what we find: when people call dot() in chains, then they're about evenly split, and actually use the left-to-right, left-associative order slightly more often than the right-to-left, right-associative order: http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069578.html DATA SOURCE 3: And if we look at other languages, then we find: == "Votes" for right-associativity == <none> == "Votes" for left-associativity == Same-left: Matlab, Julia, IDL, GAUSS Tight-left: R And Mathematica uses the "grouping" approach. ARGUMENTS: The final outcome of this is that I need to write a piece of text that says what our (at least rough) consensus is, and lays out the reasons. So long as the "vote" is so evenly split, I can't really do this. But I can imagine what the different pieces of text might look like. THE CASE FOR LEFT-ASSOCIATIVITY: If I were writing this text in favor of left-associativity, I'd point out: - "Special cases aren't special enough to break the rules". Every single operator in Python besides ** is left-associative (and ** has very compelling arguments for right associativity). @ does not have similarly compelling arguments. If we were having this debate about "*", then it'd probably be much more lopsided towards left-associativity. So sure, there's something about @ that makes right-associativity *more* appealing than for most other operators. But not *much* more appealing -- left-associativity still comes out at least slightly ahead in all of the above measures. And there are a lot of benefits to avoiding special cases -- it gives fewer rules to memorize, fewer rules to remember, etc. So @ may be a special case, but it's not special enough. - Other languages with @ operators almost overwhelmingly use the "same-left" rule, and I've never heard anyone complain about this, so clearly nothing horrible will happen if we go this way. We have no comparable experience for right-associativity. - Given left-associativity, then there's good agreement about the appropriate precedence. If we choose right-associativity then it's much less clear (which will then make things harder for experts to remember, harder for non-experts to guess, etc.). Note that one of the votes for right-associativity even preferred the "same-right" rule, which is not even technically possible... This strikes me as a nice solid case. THE CASE FOR RIGHT-ASSOCIATIVITY: If I were writing this text in favor of right-associativity, I'd point out: - Because matrix multiplication has a tight conceptual association with function application/composition, many mathematically sophisticated users have an intuition that a matrix expression like R S x proceeds from right-to-left, with first S transforming x, and then R transforming the result. This isn't universally agreed, but at the least this intuition is more common than for other operations like 2 * 3 * 4 that everyone reads as going from left-to-right. - There might be some speed argument, if people often write things like "Mat @ Mat @ vec"? But no-one has found any evidence that people actually do write such things often. - There's been discussion of how right-associativity might maybe perhaps be nice for non-matmul applications? But I can't use those arguments [17] [18]. - ...... I got nothin'. I am fine with any outcome here. (I'm actually listed under *both* tight-right and same-left in the straw poll above ;-).) I'm totally happy to go back to Guido et al and argue for right-associativity. BUT if you all want me to do that then you need to give me some better arguments to use :-). One way to do this might be to go through the ((a @ b) @ c) and (a @ (b @ c)) examples I found (the scripts are available [19], and I can help modify them to spit out more details), look at the actual code, and demonstrate that the left-to-right ((a @ b) @ c) cases are mostly ones where evaluation order doesn't matter (i.e., they could have just as well been written the other way), and the right-to-left (a @ (b @ c)) ones are ones where right-to-left really is better than left-to-right. I have no idea if this is true, and it'll require some reading of the surrounding code to figure out what the matrix shapes are, but if it *is* true then at least that'd be something solid that right-associativity advocates could point to. WHAT NOW: If seeing this summary laid out caused you to change your mind one way or the other, then please reply and say so! If you think of some way to get more data that could favor one or the other option (running some kind of usability experiment on people? [20]), then please share! If you think of some other arguments in favor of left-associativity, then please share! If you think of some other arguments in favor of right-associativity, especially ones that are based on something besides your gut feeling, then PLEASE PLEASE share! Thanks, -n [1] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069444.html [2] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069446.html [3] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069450.html [4] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069452.html [5] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069455.html [6] https://mail.python.org/pipermail/python-ideas/2014-March/027124.html [7] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069512.html [8] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069513.html [9] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069467.html [10] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069530.html [11] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069537.html [12] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069540.html [13] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069571.html [14] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069514.html [15] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069531.html [16] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069567.html [17] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069527.html [18] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069584.html [19] https://gist.github.com/njsmith/9157645 [20] http://mail.scipy.org/pipermail/numpy-discussion/2014-March/069584.html -- Nathaniel J. Smith Postdoctoral researcher - Informatics - University of Edinburgh http://vorpus.org

12 19