Mailman 3 January 2015 - NumPy-Discussion

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

18 Feb '23

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 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

Generalize hstack/vstack --> stack; Block matrices like in matlab
by Stefan Otte 09 Jan '16

09 Jan '16

Hey, quite often I work with block matrices. Matlab offers the convenient notation [ a b; c d ] to stack matrices. The numpy equivalent is kinda clumsy: vstack([hstack([a,b]), hstack([c,d])]) I wrote the little function `stack` that does exactly that: stack([[a, b], [c, d]]) In my case `stack` replaced `hstack` and `vstack` almost completely. If you're interested in including it in numpy I created a pull request [1]. I'm looking forward to getting some feedback! Best, Stefan [1] https://github.com/numpy/numpy/pull/5057

12 24

PR added: frozen dimensions in gufunc signatures
by Jaime Fernández del Río 23 Jun '15

23 Jun '15

Hi, I have just sent a PR (https://github.com/numpy/numpy/pull/5015), adding the possibility of having frozen dimensions in gufunc signatures. As a proof of concept, I have added a `cross1d` gufunc to `numpy.core.umath_tests`: In [1]: import numpy as np In [2]: from numpy.core.umath_tests import cross1d In [3]: cross1d.signature Out[3]: '(3),(3)->(3)' In [4]: a = np.random.rand(1000, 3) In [5]: b = np.random.rand(1000, 3) In [6]: np.allclose(np.cross(a, b), cross1d(a, b)) Out[6]: True In [7]: %timeit np.cross(a, b) 10000 loops, best of 3: 76.1 us per loop In [8]: %timeit cross1d(a, b) 100000 loops, best of 3: 13.1 us per loop In [9]: c = np.random.rand(1000, 2) In [10]: d = np.random.rand(1000, 2) In [11]: cross1d(c, d) --------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-11-72c66212e40c> in <module>() ----> 1 cross1d(c, d) ValueError: cross1d: Operand 0 has a mismatch in its core dimension 0, with gufunc signature (3),(3)->(3) (size 2 is different from 3) The speed up over `np.cross` is nice, and while `np.cross` is not the best of examples, as it needs to handle more sizes, in many cases this will allow producing gufuncs that work without a Python wrapper redoing checks that are best left to the iterator, such as dimension sizes. It still needs tests, but before embarking on fully developing those, I wanted to make sure that there is an interest on this. I would also like to further enhance gufuncs providing computed dimensions, e.g. making it possible to e.g. define `pairwise_cross` with signature '(n, 3)->($m, 3)', where the $ indicates that m is a computed dimension, that would have to be calculated by a function passed to the gufunc constructor and stored in the gufunc object, based on the other core dimensions. In this case it would make $m be n*(n-1), so that all pairwise cross products between 3D vectors could be computed. The syntax with '$' is kind of crappy, so any suggestions on how to better express this in the signature are more than welcome, as well as any feedback on the merits (or lack of them) of implementing this. Jaime -- (\__/) ( O.o) ( > <) Este es Conejo. Copia a Conejo en tu firma y ayúdale en sus planes de dominación mundial.

4 4

Problem with _dotblas.pyd when using Matplotlib for 3d plot
by Geoffrey Mégardon 19 Mar '15

19 Mar '15

Hi, I am on Windows 8.1, 64-bit, using the usual version of Anaconda-64-bit for windows. I just come from the mailing list of matplotlib: my PC is not able to do 3D plot with matplotlib. Python crashes completely if I try to plot something in 3D and to show it: "Python.exe has stopped working..." The main problem in resolving my issue is that the python traceback and the faulthandler package do not succeed in catching the errors. Only the Windows crash report allows us to see what happened before the crash. The windows report tells that the file _dotblas.pyd from numpy is likely responsible for the crash. So I tried a simple: import numpy; numpy.test() numpy.test() does not work! I dont get any traceback from python. It just crashes as before: "Python.exe has stopped working..." And again, the Windows report blames _dotblas.pyd :) Here the new Windows report for the numpy.test(), I highlighted the lines which blame _dotblas.pyd: Version=1 EventType=APPCRASH EventTime=130608138227608275 ReportType=2 Consent=1 UploadTime=130608138229728384 ReportIdentifier=216b16d2-6f5c-11e4-bec3-48d22435da2b IntegratorReportIdentifier=216b16d1-6f5c-11e4-bec3-48d22435da2b NsAppName=python.exe Response.BucketId=398b7eee350a0fd8a7a96f705d3488f6 Response.BucketTable=4 Response.LegacyBucketId=85979911964 Response.type=4 Sig[0].Name=Application Name Sig[0].Value=python.exe Sig[1].Name=Application Version Sig[1].Value=0.0.0.0 Sig[2].Name=Application Timestamp Sig[2].Value=53b4679e *Sig[3].Name=Fault Module NameSig[3].Value=_dotblas.pyd* Sig[4].Name=Fault Module Version Sig[4].Value=0.0.0.0 Sig[5].Name=Fault Module Timestamp Sig[5].Value=545678cb Sig[6].Name=Exception Code Sig[6].Value=c000001d Sig[7].Name=Exception Offset Sig[7].Value=0000000000324022 DynamicSig[1].Name=OS Version DynamicSig[1].Value=6.3.9600.2.0.0.768.101 DynamicSig[2].Name=Locale ID DynamicSig[2].Value=2057 DynamicSig[22].Name=Additional Information 1 DynamicSig[22].Value=f32f DynamicSig[23].Name=Additional Information 2 DynamicSig[23].Value=f32feb95f950f918532aa47b4372840e DynamicSig[24].Name=Additional Information 3 DynamicSig[24].Value=9882 DynamicSig[25].Name=Additional Information 4 DynamicSig[25].Value=98823b6c7f579b24e92112ab827fe4a1 UI[2]=C:\Anaconda\python.exe UI[3]=python.exe has stopped working UI[4]=Windows can check online for a solution to the problem. UI[5]=Check online for a solution and close the program UI[6]=Check online for a solution later and close the program UI[7]=Close the program LoadedModule[0]=C:\Anaconda\python.exe LoadedModule[1]=C:\WINDOWS\SYSTEM32\ntdll.dll LoadedModule[2]=C:\WINDOWS\system32\KERNEL32.DLL LoadedModule[3]=C:\WINDOWS\system32\KERNELBASE.dll LoadedModule[4]=C:\Anaconda\python27.dll LoadedModule[5]=C:\WINDOWS\WinSxS\amd64_microsoft.vc90.crt_1fc8b3b9a1e18e3b_9.0.30729.8387_none_08e793bfa83a89b5\MSVCR90.dll LoadedModule[6]=C:\WINDOWS\system32\USER32.dll LoadedModule[7]=C:\WINDOWS\system32\ADVAPI32.dll LoadedModule[8]=C:\WINDOWS\system32\SHELL32.dll LoadedModule[9]=C:\WINDOWS\system32\GDI32.dll LoadedModule[10]=C:\WINDOWS\system32\msvcrt.dll LoadedModule[11]=C:\WINDOWS\SYSTEM32\sechost.dll LoadedModule[12]=C:\WINDOWS\system32\RPCRT4.dll LoadedModule[13]=C:\WINDOWS\SYSTEM32\combase.dll LoadedModule[14]=C:\WINDOWS\system32\SHLWAPI.dll LoadedModule[15]=C:\WINDOWS\system32\IMM32.DLL LoadedModule[16]=C:\WINDOWS\system32\MSCTF.dll LoadedModule[17]=C:\Anaconda\DLLs\_socket.pyd LoadedModule[18]=C:\WINDOWS\system32\WS2_32.dll LoadedModule[19]=C:\WINDOWS\system32\NSI.dll LoadedModule[20]=C:\Anaconda\DLLs\_ssl.pyd LoadedModule[21]=C:\Anaconda\DLLs\_ctypes.pyd LoadedModule[22]=C:\WINDOWS\system32\ole32.dll LoadedModule[23]=C:\WINDOWS\system32\OLEAUT32.dll LoadedModule[24]=C:\Anaconda\DLLs\_hashlib.pyd LoadedModule[25]=C:\WINDOWS\SYSTEM32\CRYPTSP.dll LoadedModule[26]=C:\WINDOWS\system32\rsaenh.dll LoadedModule[27]=C:\WINDOWS\SYSTEM32\bcrypt.dll LoadedModule[28]=C:\WINDOWS\SYSTEM32\CRYPTBASE.dll LoadedModule[29]=C:\WINDOWS\SYSTEM32\bcryptPrimitives.dll LoadedModule[30]=C:\Users\User\AppData\Roaming\Python-Eggs\faulthandler-2.4-py2.7-win-amd64.egg-tmp\faulthandler.pyd LoadedModule[31]=C:\Anaconda\lib\site-packages\numpy\core\multiarray.pyd LoadedModule[32]=C:\Anaconda\lib\site-packages\numpy\core\umath.pyd LoadedModule[33]=C:\Anaconda\lib\site-packages\numpy\core\_dotblas.pyd LoadedModule[34]=C:\Anaconda\lib\site-packages\numpy\core\libiomp5md.dll LoadedModule[35]=C:\Anaconda\lib\site-packages\numpy\core\scalarmath.pyd LoadedModule[36]=C:\Anaconda\lib\site-packages\numpy\lib\_compiled_base.pyd LoadedModule[37]=C:\Anaconda\lib\site-packages\numpy\linalg\lapack_lite.pyd LoadedModule[38]=C:\Anaconda\lib\site-packages\numpy\linalg\_umath_linalg.pyd LoadedModule[39]=C:\Anaconda\lib\site-packages\numpy\fft\fftpack_lite.pyd LoadedModule[40]=C:\Anaconda\lib\site-packages\numpy\random\mtrand.pyd LoadedModule[41]=C:\Anaconda\DLLs\_multiprocessing.pyd State[0].Key=Transport.DoneStage1 State[0].Value=1 FriendlyEventName=Stopped working ConsentKey=APPCRASH AppName=python.exe AppPath=C:\Anaconda\python.exe NsPartner=windows NsGroup=windows8 ApplicationIdentity=5B036AF1EC2E20F320DBF28D119DE93D Any idea, of what I can do? I try to reinstall anaconda, I also try conda remove numpy and conda install numpy. That does not seem to change anything. -- -- MEGARDON Geoffrey

5 4

3D array and the right hand rule
by Dieter Van Eessen 11 Feb '15

11 Feb '15

Hello, I'm a novice with respect to scientific computing using python. I've read that numpy.array isn't arranged according to the 'right-hand-rule' (right-hand-rule => thumb = +x; index finger = +y, bend middle finder = +z). This is also confirmed by an old message I dug up from the mailing list archives. (see message below) I guess this has consequences for certain algorithms (been a while, should actually revise some algebra textbooks). At least it does when I try to mentally visualize a 3D array when I'm not sure which dimensions to use... What are the consequences in using 3D arrays in for example transformation algorithms or other algorithms which expect certain shape? Should I always 'reshape' before using them or do most algorithms take this into account in the underlying algebra? Does the 'Fortran contiguous' shape always respect the 'right-hand-rule' (for 3D arrays)? And just to be sure: Is the information telling if an array is C-contiguous or Fortran-contiguous always stored within the array? kind regards, Dieter / Old message From: Anne Archibald <peridot.faceted <at> gmail.com> Subject: Re: dimension aligment <http://news.gmane.org/find-root.php?message_id=ce557a360805201104r3fd2e192y…> Newsgroups: gmane.comp.python.numeric.general <http://news.gmane.org/gmane.comp.python.numeric.general> Date: 2008-05-20 18:04:46 GMT (6 years, 35 weeks, 5 days, 4 hours and 34 minutes ago) 2008/5/20 Thomas Hrabe <thrabe <at> burnham.org>: > given a *3d* *array* > a = > numpy.*array*([[[1,2,3],[4,5,6]],[[7,8,9],[10,11,12]],[[13,14,15],[16,17,18]],[[19,20,21],[22,23,24]]]) > a.shape > returns (4,2,3) > > so I assume the first digit is the 3rd dimension, second is 2nd dim and > third is the first. > > how is the data aligned in memory now? > according to the strides it should be > 1,2,3,4,5,6,7,8,9,10,... > *right*? > > if I had an *array* of more dimensions, the first digit returned by shape > should always be the highest dim. You are basically *right*, but this is a surprisingly subtle issue for numpy. A numpy *array* is basically a block of memory and some description. One piece of that description is the type of data it contains (i.e., how to interpret each chunk of memory) for example int32, float64, etc. Another is the sizes of all the various dimensions. A third piece, which makes many of the things numpy does possible, is the "strides". The way numpy works is that basically it translates A[i,j,k] into a lookup of the item in the memory block at position i*strides[0]+j*strides[1]+k*strides[2] This means, if you have an *array* A and you want every second element (A[::2]), all numpy needs to do is *hand* you back a new *array* pointing to the same data block, but with strides[0] doubled. Similarly if you want to transpose a two-dimensional *array*, all it needs to do is exchange strides[0] and strides[1]; no data need be moved. This means, though, that if you are *hand*ed a numpy *array*, the elements can be arranged in memory in quite a complicated fashion. Sometimes this is no problem - you can always use the strides to find it all. But sometimes you need the data arranged in a particular way. numpy defines two particular ways: "C contiguous" and "FORTRAN contiguous". "C contiguous" *array*s are what you describe, and they're what numpy produces by default; they are arranged so that the *right*most index has the smallest stride. "FORTRAN contiguous" *array*s are arranged the other way around; the leftmost index has the smallest stride. (This is how FORTRAN *array*s are arranged in memory.) There is also a special case: the reshape() function changes the shape of the *array*. It has an "order" argument that describes not how the elements are arranged in memory but how you want to think of the elements as arranged in memory for the reshape operation. Anne /Old message

2 2

new mingw-w64 based numpy and scipy wheel (still experimental)
by Carl Kleffner 09 Feb '15

09 Feb '15

I took time to create mingw-w64 based wheels of numpy-1.9.1 and scipy-0.15.1 source distributions and put them on https://bitbucket.org/carlkl/mingw-w64-for-python/downloads as well as on binstar.org. The test matrix is python-2.7 and 3.4 for both 32bit and 64bit. Feedback is welcome. The wheels can be pip installed with: pip install -i https://pypi.binstar.org/carlkl/simple numpy pip install -i https://pypi.binstar.org/carlkl/simple scipy Some technical details: the binaries are build upon OpenBLAS as accelerated BLAS/Lapack. OpenBLAS itself is build with dynamic kernels (similar to MKL) and automatic runtime selection depending on the CPU. The minimal requested feature supplied by the CPU is SSE2. SSE1 and non-SSE CPUs are not supported with this builds. This is the default for 64bit binaries anyway. OpenBLAS is deployed as part of the numpy wheel. That said, the scipy wheels mentioned above are dependant on the installation of the OpenBLAS based numpy and won't work i.e. with an installed numpy-MKL. For the numpy 32bit builds there are 3 failures for special FP value tests, due to a bug in mingw-w64 that is still present. All scipy versions show up 7 failures with some numerical noise, that could be ignored (or corrected with relaxed asserts in the test code). PR's for numpy and scipy are in preparation. The mingw-w64 compiler used for building can be found at https://bitbucket.org/carlkl/mingw-w64-for-python/downloads.

7 31

Characteristic of a Matrix.
by Colin J. Williams 04 Feb '15

04 Feb '15

One of the essential characteristics of a matrix is that it be rectangular. This is neither spelt out or checked currently. The Doc description refers to a class: - *class *numpy.matrix[source] <http://github.com/numpy/numpy/blob/v1.9.1/numpy/matrixlib/defmatrix.py#L206> Returns a matrix from an array-like object, or from a string of data. A matrix is a specialized 2-D array that retains its 2-D nature through operations. It has certain special operators, such as * (matrix multiplication) and ** (matrix power). This illustrates a failure, which is reported later in the calculation: A2= np.matrix([[1, 2, -2], [-3, -1, 4], [4, 2 -6]]) Here 2 - 6 is treated as an expression. Wikipedia offers: In mathematics <http://en.wikipedia.org/wiki/Mathematics>, a *matrix* (plural *matrices*) is a rectangular <http://en.wikipedia.org/wiki/Rectangle> *array <http://en.wiktionary.org/wiki/array>*[1] <http://en.wikipedia.org/wiki/Matrix_%28mathematics%29#cite_note-1> of numbers <http://en.wikipedia.org/wiki/Number>, symbols <http://en.wikipedia.org/wiki/Symbol_%28formal%29>, or expressions <http://en.wikipedia.org/wiki/Expression_%28mathematics%29>, arranged in *rows <http://en.wiktionary.org/wiki/row>* and *columns <http://en.wiktionary.org/wiki/column>*.[2] <http://en.wikipedia.org/wiki/Matrix_%28mathematics%29#cite_note-2>[3] <http://en.wikipedia.org/wiki/Matrix_%28mathematics%29#cite_note-3> The individual items in a matrix are called its *elements* or *entries*. An example of a matrix with 2 rows and 3 columns is [image: \begin{bmatrix}1 & 9 & -13 \\20 & 5 & -6 \end{bmatrix}.]In the Numpy context, the symbols or expressions need to be evaluable. Colin W.

11 20

Views of a different dtype
by Jaime Fernández del Río 04 Feb '15

04 Feb '15

HI all, There has been some recent discussion going on on the limitations that numpy imposes to taking views of an array with a different dtype. As of right now, you can basically only take a view of an array if it has no Python objects and neither the old nor the new dtype are structured. Furthermore, the array has to be either C or Fortran contiguous. This seem to be way too strict, but the potential for disaster getting a loosening of the restrictions wrong is big, so it should be handled with care. Allan Haldane and myself have been looking into this separately and discussing some of the details over at github, and we both think that the only true limitation that has to be imposed is that the offsets of Python objects within the new and old dtypes remain compatible. I have expanded Allan's work from here: https://github.com/ahaldane/numpy/commit/e9ca367 to make it as flexible as I have been able. An implementation of the algorithm in Python, with a few tests, can be found here: https://gist.github.com/jaimefrio/b4dae59fa09fccd9638c#file-dtype_compat-py I would appreciate getting some eyes on it for correctness, and to make sure that it won't break with some weird dtype. I am also trying to figure out what the ground rules for stride and shape conversions when taking a view with a different dtype should be. I submitted a PR (gh-5508) a couple for days ago working on that, but I am not so sure that the logic is completely sound. Again, to get more eyes on it, I am going to reproduce my thoughts here on the hope of getting some feedback. The objective would be to always allow a view of a different dtype (given that the dtypes be compatible as described above) to be taken if: - The itemsize of the dtype doesn't change. - The itemsize changes, but the array being viewed is the result of slicing and transposing axes of a contiguous array, and it is still contiguous, defined as stride == dtype.itemsize, along its smallest-strided dimension, and the itemsize of the newtype exactly divides the size of that dimension. - Ideally taking a view should be a reversible process, i.e. if oldtype = arr.dtype, then doing arr.view(newtype).view(oldtype) should give you back a view of arr with the same original shape, strides and dtype. This last point can get tricky if the minimal stride dimension has size 1, as there could be several of those, e.g.: >>> a = np.ones((3, 4, 1), dtype=float)[:, :2, :].transpose(0, 2, 1) >>> a.flags.contiguous False >>> a.shape (3, 1, 2) >>> a.strides # the stride of the size 1 dimension could be anything, ignore it! (32, 8, 8) b = a.view(complex) # this fails right now, but should work >>> b.flags.contiguous False >>> b.shape (3, 1, 1) >>> b.strides # the stride of the size 1 dimensions could be anything, ignore them! (32, 16, 16) c = b.view(float) # which of the two size 1 dimensions should we expand? "In the face of ambiguity refuse the temptation to guess" dictates that last view should raise an error, unless we agree and document some default. Any thoughts? Then there is the endless complication one could get into with arrays created with as_strided. I'm not smart enough to figure when and when not those could work, but am willing to retake the discussion if someone wiser si interested. With all these in mind, my proposal for the new behavior is that taking a view of an array with a different dtype would require: 1. That the newtype and oldtype be compatible, as defined by the algorithm checking object offsets linked above. 2. If newtype.itemsize == oldtype.itemsize no more checks are needed, make it happen! 3. If the array is C/Fortran contiguous, check that the size in bytes of the last/first dimension is evenly divided by newtype.itemsize. If it does, go for it. 4. For non-contiguous arrays: 1. Ignoring dimensions of size 1, check that no stride is smaller than either oldtype.itemsize or newtype.itemsize. If any is found this is an as_strided product, sorry, can't do it! 2. Ignoring dimensions of size 1, find a contiguous dimension, i.e. stride == oldtype.itemsize 1. If found, check that it is the only one with that stride, that it is the minimal stride, and that the size in bytes of that dimension is evenly divided by newitem,itemsize. 2. If none is found, check if there is a size 1 dimension that is also unique (unless we agree on a default, as mentioned above) and that newtype.itemsize evenly divides oldtype.itemsize. Apologies for the long, dense content, but any thought or comments are very welcome. Jaime -- (\__/) ( O.o) ( > <) Este es Conejo. Copia a Conejo en tu firma y ayúdale en sus planes de dominación mundial.

6 13

NumPy-Discussion January 2015