NumPy-Discussion November 2011

numpy-discussion@python.org

64 participants
66 discussions

by Matthew Brett

Hi, Have I missed a fast way of doing nice float to integer conversion? By nice I mean, rounding to the nearest integer, converting NaN to 0, inf, -inf to the max and min of the integer range? The astype method and cast functions don't do what I need here: In [40]: np.array([1.6, np.nan, np.inf, -np.inf]).astype(np.int16) Out[40]: array([1, 0, 0, 0], dtype=int16) In [41]: np.cast[np.int16](np.array([1.6, np.nan, np.inf, -np.inf])) Out[41]: array([1, 0, 0, 0], dtype=int16) Have I missed something obvious? See y'all, Matthew

10 years, 2 months

Float128 integer comparison

by Matthew Brett

Hi, Continuing the exploration of float128 - can anyone explain this behavior? >>> np.float64(9223372036854775808.0) == 9223372036854775808L True >>> np.float128(9223372036854775808.0) == 9223372036854775808L False >>> int(np.float128(9223372036854775808.0)) == 9223372036854775808L True >>> np.round(np.float128(9223372036854775808.0)) == np.float128(9223372036854775808.0) True Thanks for any pointers, Best, Matthew

10 years, 2 months

float64 / int comparison different from float / int comparison

by Matthew Brett

Hi, I just ran into this confusing difference between np.float and np.float64: In [8]: np.float(2**63) == 2**63 Out[8]: True In [9]: np.float(2**63) > 2**63-1 Out[9]: True In [10]: np.float64(2**63) == 2**63 Out[10]: True In [11]: np.float64(2**63) > 2**63-1 Out[11]: False In [16]: np.float64(2**63-1) == np.float(2**63-1) Out[16]: True I believe values above 2*52 are all represented as integers in float64. http://matthew-brett.github.com/pydagogue/floating_point.html Is this this int64 issue that came up earlier in float128 comparison? Why the difference between np.float and np.float64? Thanks for any insight, Matthew

10 years, 2 months

einsum performance with many operands

by Eugenio Piasini

Hi, I was playing around with einsum (which is really cool, by the way!), and I noticed something strange (or unexpected at least) performance-wise. Here is an example: Let x and y be two NxM arrays, and W be a NxN array. I want to compute f = x^{ik} W_i^j y_{jk} (I hope the notation is clear enough) What surprises me is that it's much faster to compute the two products separately - as in ((xW)y) or (x(Wy)), rather than directly passing the three-operands expression to einsum. I did a quick benchmark, with the following script #======================================================== import numpy as np def f1(x,W,y): return np.einsum('ik,ij,jk', x,W,y) def f2(x,W,y): return np.einsum('jk,jk', np.einsum('ik,ij->jk', x, W), y) n = 30 m = 150 x=np.random.random(size=(n, m)) y=np.random.random(size=(n, m)) W=np.random.random(size=(n, n)) setup = """\ from __main__ import f1,f2,x,y,W """ if __name__ == '__main__': import timeit min_f1_time = min(timeit.repeat(stmt='f1(x,W,y)', setup=setup, repeat=10, number=10000)) min_f2_time = min(timeit.repeat(stmt='f2(x,W,y)', setup=setup, repeat=10, number=10000)) print('f1: {time}'.format(time=min_f1_time)) print('f2: {time}'.format(time=min_f2_time)) #============================================================ and I get (on a particular trial on my intel 64 bit machine running linux and numpy 1.6.1) the following: f1: 2.86719584465 f2: 0.891730070114 Of course, I know nothing of einsum's internals and there's probably a good reason for this. But still, it's quite counterintuitive for me; I just thought I'd mention it because a quick search didn't bring up anything on this topic, and AFAIK einsum's docs/examples don't cover the case with more than two operands. Anyway: I hope I've been clear enough, and that I didn't bring up some already-debated topic. Thanks in advance for any clarification, Eugenio

10 years, 2 months

Reason behind C_ABI_VERSION so high

by Sandro Tosi

Hello, in Debian we're trying to define a way to handle numpy transitions more smoothly (you can read the proposal, if interested, at http://bugs.debian.org/643873). In order to do that, we'd like to use the C_API_VERSION and C_ABI_VERSION values; while for C_API_VERSION we can see it's a quite small value, with a clear history at ./numpy/core/code_generators/cversions.txt , we don't have quite clear why C_ABI_VERSION is such a high value and how it would be incremented. Could you please shine some light on it? can we, f.e., just take 6 for API and 9 for ABI and be sure we're seeing them incremented to 7 and 10 (respectively) when needed? C_ABI_VERSION is incremented in a different way? Thanks, -- Sandro Tosi (aka morph, morpheus, matrixhasu) My website: http://matrixhasu.altervista.org/ Me at Debian: http://wiki.debian.org/SandroTosi

10 years, 2 months

NumPy nogil API

by mark florisson

Hello, First, I'd like to report a bug. It seems ndarray does not implement tp_traverse or tp_clear, so if you have a reference cycle in an ndarray with dtype object none of those objects will ever be collected. Secondly, please bear with me, I'm not a NumPy expert, but would it be possible to have NumPy export an API that can be called without the GIL? In the upcoming Cython release 0.16 (soon to be released) we will have what is called typed memoryviews [1], which allow you to obtain a typed view on any PEP 3118 buffer, and it allows you to do a lot more things without the GIL. e.g. it allows you to pass these views around, transpose them, slice them (in the same way as in NumPy but slightly more restricted, it doesn't support masks and such), index them etc without the GIL. However, there is a lot of good functionality in NumPy that is simply not accessible without going through a Python and GIL layer, only to have NumPy (possibly) release the GIL. So instead we could have an API that takes a C-level ndarray view descriptor (shape/strides/ndim (possibly suboffsets), base type description) that would do the actual work (and perhaps doesn't even need to know anything about Python) and won't need the GIL (this wouldn't work for the dtype object, of course). The Py_buffer struct comes to mind, but format strings are hard to deal with. It would be the caller's responsibility to do any necessary synchronization. This API would simply be wrapped by a Python API that gets this "view" from the PyArrayObject, and which may or may not decide to release the GIL, and call the nogil API. This wrapping API could even be written easily in Cython. As for exceptions and error handling, there are many ways we can think of doing this without requiring the GIL. One of the reasons that I think this is important is that when you're using cython.parallel [2] you don't want to hold the gil, but you do want your NumPy goodies. Cython re-implements a very small subset of NumPy to get you the core functionality, but to get back to NumPy world you have to acquire the GIL, convert your memoryview slice to an ndarray (using the buffer interface through numpy.asarray()) and then have NumPy operate on that. It's a pain to write and it's terrible for performance. Even if you forget the GIL part, there's still the (expensive and explicit) conversion. In general I think there might be many advantages to such functionality other than for Cython. There shouldn't really be a reason to tie NumPy only to the CPython platform. Anyway, what do you guys think, does this make any sense? Mark [1]: https://sage.math.washington.edu:8091/hudson/job/cython-docs/doclinks/1/src… [2]: http://docs.cython.org/src/userguide/parallelism.html

10 years, 2 months

2022

2021

2020

2019

2018

2017

2016

2015

2014

2013

2012

2011

2010

2009

2008

2007

2006

2005

2004

2003

2002

2001

2000

NumPy-Discussion November 2011