Mailman 3 August 2012 - NumPy-Discussion

Re: [Numpy-discussion] example reading binary Fortran file
by Neil Martinsen-Burrell 22 Jul '21

22 Jul '21

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 31 Mar '16

31 Mar '16

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

SVD errors
by mtrumpis＠berkeley.edu 04 Jun '13

04 Jun '13

Hello list.. I've run into two SVD errors over the last few days. Both errors are identical in numpy/scipy. I've submitted a ticket for the 1st problem (numpy ticket #990). Summary is: some builds of the lapack_lite module linking against system LAPACK (not the bundled dlapack_lite.o, etc) give a "LinAlgError: SVD did not converge" exception on my matrix. This error does occur using Mac's Accelerate framework LAPACK, and a coworker's Ubuntu LAPACK version. It does not seem to happen using ATLAS LAPACK (nor using Octave/Matlab on said Ubuntu) Just today I've come across a negative singular value cropping up in an SVD of a different matrix. This error does occur on my ATLAS LAPACK based numpy, as well as on the Ubuntu setup. And once again, it does not happen in Octave/Matlab. I'm using numpy 1.3.0.dev6336 -- don't know what the Ubuntu box is running. Here are some npy files for the two different cases: https://cirl.berkeley.edu/twiki/pub/User/MikeTrumpis/noconverge_operator.npy https://cirl.berkeley.edu/twiki/pub/User/MikeTrumpis/negsval_operator.npy Mike

6 7

example reading binary Fortran file
by David Froger 13 Feb '13

13 Feb '13

Hy, My question is about reading Fortran binary file (oh no this question again...) Until now, I was using the unpack module like that : def lread(f,fourBeginning,fourEnd,*tuple): from struct import unpack """Reading a Fortran binary file in litte-endian""" if fourBeginning: f.seek(4,1) for array in tuple: for elt in xrange(array.size): transpose(array).flat[elt] = unpack(array.dtype.char,f.read(array.itemsize))[0] if fourEnd: f.seek(4,1) After googling, I read that fopen and npfille was deprecated, and we should use numpy.fromfile and ndarray.tofile, but despite of the documentaion, the cookbook, the mailling list and google, I don't succed in making a simple example. Considering the simple Fortran code below what is the Python script to read the four arrrays? What about if my pc is little endian and the file big endian? I think it will be a good idea to put the Fortran writting-arrays code and the Python reading-array script in the cookbook and maybe a page to help people comming from Fortran to start with Python ? Best, David Froger 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

8 22

Issue tracking
by Charles R Harris 23 Oct '12

23 Oct '12

Hi All, The issue tracking discussion seems to have died. Since github issues looks to be a viable alternative at this point, I propose to turn it on for the numpy repository and start directing people in that direction. Thoughts? Chuck

9 40

how is y += x computed when y.strides = (0, 8) and x.strides=(16, 8) ?
by Sebastian Walter 17 Oct '12

17 Oct '12

Hi, I'm using numpy 1.6.1 on Ubuntu 12.04.1 LTS. A code that used to work with an older version of numpy now fails with an error. Were there any changes in the way inplace operations like +=, *=, etc. work on arrays with non-standard strides? For the script: ------- start of code ------- import numpy x = numpy.arange(6).reshape((3,2)) y = numpy.arange(2) print 'x=\n', x print 'y=\n', y u,v = numpy.broadcast_arrays(x, y) print 'u=\n', u print 'v=\n', v print 'v.strides=\n', v.strides v += u print 'v=\n', v # expectation: v = [[6,12], [6,12], [6,12]] print 'u=\n', u print 'y=\n', y # expectation: y = [6,12] ------- end of code ------- I get the output -------- start of output --------- x= [[0 1] [2 3] [4 5]] y= [0 1] u= [[0 1] [2 3] [4 5]] v= [[0 1] [0 1] [0 1]] v.strides= (0, 8) v= [[4 6] [4 6] [4 6]] u= [[0 1] [2 3] [4 5]] y= [4 6] -------- end of output -------- I would have expected that v += u performs an element-by-element += v[0,0] += u[0,0] # increments y[0] v[0,1] += u[0,1] # increments y[1] v[1,0] += u[1,0] # increments y[0] v[1,1] += u[1,1] # increments y[1] v[2,0] += u[2,0] # increments y[0] v[2,1] += u[2,1] # increments y[1] yielding the result y = [6,12] but instead one obtains y = [4, 6] which could be the result of v[2,0] += u[2,0] # increments y[0] v[2,1] += u[2,1] # increments y[1] Is this the intended behavior? regards, Sebastian

5 5

Fwd: [numpy] ENH: Initial implementation of a 'neighbor' calculation (#303)
by Nathaniel Smith 15 Oct '12

15 Oct '12

[Manual PR notification] ---------- Forwarded message ---------- From: timcera Date: Sat, Jun 9, 2012 at 10:13 PM Subject: [numpy] ENH: Initial implementation of a 'neighbor' calculation (#303) To: njsmith <njs(a)pobox.com> Each element is assigned the result of a function based on it's neighbors. Neighbors are selected based on a weight array. It uses the new pad routines to pad arrays if neighboring values are required that would be off the edge of the input array. Will be great to have the masked array settled because right now you can only sort of exclude from the neighborhood using a zero in the weight array. Zero or np.IGNORE don't affect np.sum, but functions like np.mean and np.std would give different answers. Because of this my early implementations of neighbor included an optional mask array along with the weight array, but I decided would be best to wait for the new masked arrays. This in some ways could be considered a generalization of a convolution, and comparison with existing numpy/scipy convolution results are included in the tests. The advantage to neighbor is that any function that accepts a 1-d array, and returns a single result, can be used instead of convolution only using summation. The convolution functions require the weight array to be flipped to get the same answer as neighbor. You can merge this Pull Request by running: git pull https://github.com/timcera/numpy neighbor Or you can view, comment on it, or merge it online at: https://github.com/numpy/numpy/pull/303 -- Commit Summary -- * ENH: Initial implementation of a 'neighbor' calculation where the each -- File Changes -- M numpy/lib/__init__.py (2) A numpy/lib/neighbor.py (305) A numpy/lib/tests/test_neighbor.py (278) -- Patch Links -- https://github.com/numpy/numpy/pull/303.patch https://github.com/numpy/numpy/pull/303.diff --- Reply to this email directly or view it on GitHub: https://github.com/numpy/numpy/pull/303

7 14

dtype '|S0' not understood
by David Warde-Farley 06 Oct '12

06 Oct '12

Howdy, It seems it's possible using e.g. In [25]: dtype([('foo', str)])Out[25]: dtype([('foo', '|S0')]) to get yourself a zero-length string. However dtype('|S0') results in a TypeError: data type not understood. I understand the stupidity of creating a 0-length string field but it's conceivable that it's accidental. For example, it could lead to a situation where you've created that field, are missing all the data you had meant to put in it, serialize with np.save, and upon np.load aren't able to get _any_ of your data back because the dtype descriptor is considered bogus (can you guess why I thought of this scenario?). It seems that either dtype(str) should do something more sensible than zero-length string, or it should be possible to create it with dtype('| S0'). Which should it be? David

3 3

Long-standing issue with using numpy in embedded CPython
by Yang Zhang 25 Sep '12

25 Sep '12

It turns out that there's a long-standing problem in numpy that prevents it from being used in embedded CPython environments: http://stackoverflow.com/questions/7592565/when-embedding-cpython-in-java-w… http://mail.scipy.org/pipermail/numpy-discussion/2009-July/044046.html Is there any fix or workaround for this? Thanks. -- Yang Zhang http://yz.mit.edu/

9 15

ANN: NumPy 1.7.0b1 release
by Ondřej Čertík 09 Sep '12

09 Sep '12

Hi, I'm pleased to announce the availability of the first beta release of NumPy 1.7.0b1. Sources and binary installers can be found at https://sourceforge.net/projects/numpy/files/NumPy/1.7.0b1/ Please test this release and report any issues on the numpy-discussion mailing list. The following problems are known and we'll work on fixing them before the final release: http://projects.scipy.org/numpy/ticket/2187 http://projects.scipy.org/numpy/ticket/2185 http://projects.scipy.org/numpy/ticket/2066 http://projects.scipy.org/numpy/ticket/1588 http://projects.scipy.org/numpy/ticket/2076 http://projects.scipy.org/numpy/ticket/2101 http://projects.scipy.org/numpy/ticket/2108 http://projects.scipy.org/numpy/ticket/2150 http://projects.scipy.org/numpy/ticket/2189 I would like to thank Ralf for a lot of help with creating binaries and other help for this release. Cheers, Ondrej ========================= NumPy 1.7.0 Release Notes ========================= This release includes several new features as well as numerous bug fixes and refactorings. It supports Python 2.4 - 2.7 and 3.1 - 3.2. Highlights ========== * ``where=`` parameter to ufuncs (allows the use of boolean arrays to choose where a computation should be done) * ``vectorize`` improvements (added 'excluded' and 'cache' keyword, general cleanup and bug fixes) * ``numpy.random.choice`` (random sample generating function) Compatibility notes =================== In a future version of numpy, the functions np.diag, np.diagonal, and the diagonal method of ndarrays will return a view onto the original array, instead of producing a copy as they do now. This makes a difference if you write to the array returned by any of these functions. To facilitate this transition, numpy 1.7 produces a FutureWarning if it detects that you may be attempting to write to such an array. See the documentation for np.diagonal for details. Similar to np.diagonal above, in a future version of numpy, indexing a record array by a list of field names will return a view onto the original array, instead of producing a copy as they do now. As with np.diagonal, numpy 1.7 produces a FutureWarning if it detects that you may be attemping to write to such an array. See the documentation for array indexing for details. The default casting rule for UFunc out= parameters has been changed from 'unsafe' to 'same_kind'. Most usages which violate the 'same_kind' rule are likely bugs, so this change may expose previously undetected errors in projects that depend on NumPy. Full-array boolean indexing has been optimized to use a different, optimized code path. This code path should produce the same results, but any feedback about changes to your code would be appreciated. Attempting to write to a read-only array (one with ``arr.flags.writeable`` set to ``False``) used to raise either a RuntimeError, ValueError, or TypeError inconsistently, depending on which code path was taken. It now consistently raises a ValueError. The <ufunc>.reduce functions evaluate some reductions in a different order than in previous versions of NumPy, generally providing higher performance. Because of the nature of floating-point arithmetic, this may subtly change some results, just as linking NumPy to a different BLAS implementations such as MKL can. If upgrading from 1.5, then generally in 1.6 and 1.7 there have been substantial code added and some code paths altered, particularly in the areas of type resolution and buffered iteration over universal functions. This might have an impact on your code particularly if you relied on accidental behavior in the past. New features ============ Reduction UFuncs Generalize axis= Parameter ------------------------------------------- Any ufunc.reduce function call, as well as other reductions like sum, prod, any, all, max and min support the ability to choose a subset of the axes to reduce over. Previously, one could say axis=None to mean all the axes or axis=# to pick a single axis. Now, one can also say axis=(#,#) to pick a list of axes for reduction. Reduction UFuncs New keepdims= Parameter ---------------------------------------- There is a new keepdims= parameter, which if set to True, doesn't throw away the reduction axes but instead sets them to have size one. When this option is set, the reduction result will broadcast correctly to the original operand which was reduced. Datetime support ---------------- .. note:: The datetime API is *experimental* in 1.7.0, and may undergo changes in future versions of NumPy. There have been a lot of fixes and enhancements to datetime64 compared to NumPy 1.6: * the parser is quite strict about only accepting ISO 8601 dates, with a few convenience extensions * converts between units correctly * datetime arithmetic works correctly * business day functionality (allows the datetime to be used in contexts where only certain days of the week are valid) The notes in `doc/source/reference/arrays.datetime.rst <https://github.com/numpy/numpy/blob/maintenance/1.7.x/doc/source/reference/…>`_ (also available in the online docs at `arrays.datetime.html <http://docs.scipy.org/doc/numpy/reference/arrays.datetime.html>`_) should be consulted for more details. Custom formatter for printing arrays ------------------------------------ See the new ``formatter`` parameter of the ``numpy.set_printoptions`` function. New function numpy.random.choice --------------------------------- A generic sampling function has been added which will generate samples from a given array-like. The samples can be with or without replacement, and with uniform or given non-uniform probabilities. New function isclose -------------------- Returns a boolean array where two arrays are element-wise equal within a tolerance. Both relative and absolute tolerance can be specified. The function is NA aware. Preliminary multi-dimensional support in the polynomial package --------------------------------------------------------------- Axis keywords have been added to the integration and differentiation functions and a tensor keyword was added to the evaluation functions. These additions allow multi-dimensional coefficient arrays to be used in those functions. New functions for evaluating 2-D and 3-D coefficient arrays on grids or sets of points were added together with 2-D and 3-D pseudo-Vandermonde matrices that can be used for fitting. Ability to pad rank-n arrays ---------------------------- A pad module containing functions for padding n-dimensional arrays has been added. The various private padding functions are exposed as options to a public 'pad' function. Example:: pad(a, 5, mode='mean') Current modes are ``constant``, ``edge``, ``linear_ramp``, ``maximum``, ``mean``, ``median``, ``minimum``, ``reflect``, ``symmetric``, ``wrap``, and ``<function>``. New argument to searchsorted ---------------------------- The function searchsorted now accepts a 'sorter' argument that is a permuation array that sorts the array to search. C API ----- New function ``PyArray_RequireWriteable`` provides a consistent interface for checking array writeability -- any C code which works with arrays whose WRITEABLE flag is not known to be True a priori, should make sure to call this function before writing. NumPy C Style Guide added (``doc/C_STYLE_GUIDE.rst.txt``). Changes ======= General ------- The function np.concatenate tries to match the layout of its input arrays. Previously, the layout did not follow any particular reason, and depended in an undesirable way on the particular axis chosen for concatenation. A bug was also fixed which silently allowed out of bounds axis arguments. The ufuncs logical_or, logical_and, and logical_not now follow Python's behavior with object arrays, instead of trying to call methods on the objects. For example the expression (3 and 'test') produces the string 'test', and now np.logical_and(np.array(3, 'O'), np.array('test', 'O')) produces 'test' as well. Deprecations ============ General ------- Specifying a custom string formatter with a `_format` array attribute is deprecated. The new `formatter` keyword in ``numpy.set_printoptions`` or ``numpy.array2string`` can be used instead. The deprecated imports in the polynomial package have been removed. C-API ----- Direct access to the fields of PyArrayObject* has been deprecated. Direct access has been recommended against for many releases. Expect similar deprecations for PyArray_Descr* and other core objects in the future as preparation for NumPy 2.0. The macros in old_defines.h are deprecated and will be removed in the next major release (>= 2.0). The sed script tools/replace_old_macros.sed can be used to replace these macros with the newer versions. You can test your code against the deprecated C API by #defining NPY_NO_DEPRECATED_API to the target version number, for example NPY_1_7_API_VERSION, before including any NumPy headers.

12 21

NumPy-Discussion August 2012