Mailman 3 August 2001 - NumPy-Discussion

one-offset arrays
by Eric Nodwell 29 Aug '01

29 Aug '01

Does anyone have a simple method of using one-offset arrays? Specifically, can I define an array "a" so that a[1] refers to the first element? Without one-offset indexing, it seems to me that Python is minimally useful for numerical computations. Many, perhaps the majority, of numerical algorithms are one-indexed. Matlab for example is one-based for this reason. In fact it seems strange to me that a "high-level" language like Python should use zero-offset lists.

4 3

Re: [Numpy-discussion] one-offset arrays
by Eric Nodwell 29 Aug '01

29 Aug '01

>> Does anyone have a simple method of using one-offset arrays? >> Specifically, can I define an array "a" so that a[1] refers to the >> first element? >> > Inherit from the UserArray and redefine slicing to your hearts content. >> Without one-offset indexing, it seems to me that Python is minimally >> useful for numerical computations. Many, perhaps the majority, of >> numerical algorithms are one-indexed. Matlab for example is one-based >> for this reason. In fact it seems strange to me that a "high-level" >> language like Python should use zero-offset lists. >> > Wow, that is a pretty condescending-sounding statement, though I'm sure you > didn't mean it that way. Python is indeed being used for quite serious > numerical computations. I have been using Python for quite some time for > Numerical work and find it's zero-based indexing combined with the > leave-last-one-out slicing notation to be much more convenient. Oops, what I really need is a wet-ware (i.e. brain) macro which enforces the correct order of the pair (think,write,send)! The above unconsidered comments arose from the sequence (frustration,write,send,think,regret). ;) Obviously there is a lot of numerical work being done with python and people are very happy with it. But for me I still think it would be "minimally useful" without 1-indexed arrays. Here's why: In my view, the most important reason to prefer 1-based indexing versus 0-based indexing is compatibility. For numerical work, some of the languages which I use or have used are Matlab, Mathematica, Maple and Fortran. These are all 1-indexed. (C is by nature 0-indexed because it is so close to machine architecture, but with a little bit of not-entirely-clean pointer manipulation, you can easily make 1-indexed arrays and matrices.) Obviously Python can't be compatible with these languages in a strict sense, but like most people who do some programming work, I've built up a library of my own commonly used routines specific to my work; in general it's a trivial matter to translate numerical routines from one language to the another if translation is just a matter of substituting of one set of syntactical symbols and function names for anther. However it can be damn tricky to convert 1-indexed code to 0-indexed code or visa versa without introducing any errors- believe me! (Yes it's possible to call nearly any language from nearly any other language these days so in theory you don't have to recode, but there are lots of reasons why often recoding is the preferable route.) The other reason for choosing 1-based indexing is to keep the code as near to the standard notation as possible. This is one of the advantages of using a high level language - it approximates the way you think about things instead of the way the computer organizes them. Of course, this can go either way depending on the quantity in question: as a simple example a spatial vector (x,y,z) is conventionally labelled 1,2,3 (1-indexed), but a relativistic four-vector with time included (t,x,y,z) is conventionally labelled 0,1,2,3 (0-indexed). So ideally one would be able to choose the indexing-type case-by-case. I'm sure that computer programmers will argue vehemently that code which mixes both 0-indexed and 1-indexed arrays is hard to understand and maintain, but for experts in a particular field who are accustomed to certain ingrained notations, it is the code which breaks the conventional notation which is most error-prone. In my case, I'm dealing at the moment with crystal structures with which are associated certain conventional sets of vectors and tensors - all 1-indexed by convention. I find it a delicate enough task to always get the correct vector or tensor without having to remember that d[2] is actually d3. Defining d1,d2,d3 is not convenient because often the index itself needs to be calculated. I guess if I understood the reason for 0-indexed lists and tuples in Python I would be happier. In normal, everyday usage, sets, collections and lists are 1-indexed (the first item, the second item, the third item, and so on). Python is otherwise such an elegant and natural language. Why the ugly exception of making the user conform to the underlying mechanism of an array being an address plus an offset? All this is really neither here nor there, since this debate, at least as far as Python is concerned, was probably settled 10 years ago and I'm sure nobody wants to hear anything more about it at this point. As you point out, I can define my own array type with inheritance. I will also need my own range command and several other functions which haven't occured to me yet. I was hoping that there would be a standard module to implement this. By the way, what is leave-last-one-out slicing? Is it a[:-1] or is it a[0,...] or is it something else? Eric

4 5

invert() with complex arrays
by Rob 29 Aug '01

29 Aug '01

Is there any way to get invert() to work with complex arrays? Thanks, Rob. -- The Numeric Python EM Project www.members.home.net/europax

2 3

Re: MLab.py status: do you use it?
by Serge Rey 29 Aug '01

29 Aug '01

On Tue, Aug 28, 2001 at 12:09:50PM -0700, numpy-discussion-request(a)lists.sourceforge.net wrote: -snip- > I certainly could write up a test routine. I will start work on that. In > the meantime, before I put too much effort into it, I'd really like some > feedback on what people want MLab.py to be, if they see using it at all. > > As I said above, I'm not sure trying to give NumPy a MATLAB feel is a > useful goal, so what is left is having all those handy functions. > Perhaps we could shift away from MLab.py, and turn it into Utilities.py, > which could use MATLAB (and other user experience) as guidance as to > what utilities to include. > > I'd really like folks' feedback on this before I put too much work into > something no one uses, or will use. i use MLab.py fairly often but see it more as a collection of useful functions. personally, i agree with moving it into Utilities.py, as that is in line with how i use it. otoh, MLab.py does provide a nice entry point for folks coming from Matlab into numpy. so, perhaps any eventual move to Utilities.py could be accompanied with some documentation to remind users of those linkages. i'd offer to help with the testing, but my copy of matlab is back in the states (i'm in australia on sabbatical) so i'd have to go from memory (dangerous). let me know if i can help you. cheers, serge -- Sergio J. Rey http://typhoon.sdsu.edu/rey.html "NT is secure.... as long as you don't remove the shrink wrap." - G. Myers

1 0

Newbie Installation problems
by Culhane, Aedin 28 Aug '01

28 Aug '01

Dear Numpy Discussion List Users, I would really appreciate if someone could help me install numerical python on Win98. I have downloaded Numerical Python twice now (older version and 20.1.0) but I am unable to get python to import numeric. I am running python2.1. I need it to be able to calculate basis statistics (mean, std dev, etc). Will numerical python do this easily or should I write my own scripts? Thanks a million for your help, Aedin Culhane

4 3

Gaussian Quadrature routine Numpy-ization :)
by Rob 28 Aug '01

28 Aug '01

I finally got it to work, but the Numpy-ized version runs slower than the plain Python one. I think that I can transpose the NodeCord matrix once in the program and feed that in, rather than the scratch matrix that is generated here. Evidently transpose is not very fast even for smal matrices. Here is my test program: from Numeric import * Qpnt=array(([20,21,22],[23,24,25],[26,27,28])) NodeCord=array(([1,2,3],[4,5,6],[7,8,9])) TrngleNode=array((1,2,0)) #the original routine def ComputeGaussQuadPoint(QuadPoint,TrngleNode,Qpnt,NodeCord): SrcPointCol=zeros((3)) SrcPointCol[0] = Qpnt[QuadPoint,0]*NodeCord[TrngleNode[0],0]\ + Qpnt[QuadPoint,1]*NodeCord[TrngleNode[1],0]\ + Qpnt[QuadPoint,2]*NodeCord[TrngleNode[2],0] SrcPointCol[1] = Qpnt[QuadPoint,0]*NodeCord[TrngleNode[0],1]\ + Qpnt[QuadPoint,1]*NodeCord[TrngleNode[1],1]\ + Qpnt[QuadPoint,2]*NodeCord[TrngleNode[2],1] SrcPointCol[2] = Qpnt[QuadPoint,0]*NodeCord[TrngleNode[0],2]\ + Qpnt[QuadPoint,1]*NodeCord[TrngleNode[1],2]\ + Qpnt[QuadPoint,2]*NodeCord[TrngleNode[2],2] return SrcPointCol #the yet-to-be-faster routine def Quad(QuadPoint, TrngleNode, Qpnt,NodeCord): s = Qpnt[QuadPoint,:] c= take(NodeCord, TrngleNode) SrcPointCol= add.reduce(s * transpose(c),1) return SrcPointCol QuadPoint=1 print "The Correct:" print ComputeGaussQuadPoint(QuadPoint,TrngleNode,Qpnt,NodeCord) print "The New" print Quad(QuadPoint,TrngleNode,Qpnt,NodeCord) -- The Numeric Python EM Project www.members.home.net/europax

2 3

manual update
by Paul F. Dubois 27 Aug '01

27 Aug '01

I have updated the Numeric manual on the web site. I added some to the MLab section to reflect the axis arguments that were added some time ago but never reflected in the manual. The manual also documents a change to the behavior of the attribute "masked" in package MA. The document corresponds to CVS now rather than the last release. I hope to release Numeric-20.2 sometime soon but want more exercise of the MA and MLab changes first. I changed MLab.mean to return a floating point result when the argument is integer. I'm sorry if this breaks anything but the prior behavior has no mathematical meaning and then, worse, was being used by std to get the wrong standard deviation too. I changed std to work for axis > 0, as it was previously broken. I don't own Matlab and so maintaining a package that emulates something I don't have makes me nervous. I wish one of the other developers with more knowledge in this area would pitch in. Thanks to Chris Barker and others who sent in typos and suggestions.

3 2

Can this function by Numpy-ized?
by Rob 27 Aug '01

27 Aug '01

I finally got my FEM EM code working. I profiled it and this function uses up a big hunk of time. It performs gaussian integration over a triangle. I am trying to figure out how to slice the arrays so as to push it down into the C level. Does anyone have any ideas? Thanks, Rob. ps. it looks to be intractible to me. Maybe I need to look at writing a C extension. I've never done that before. ##*************************************************************************** ##Prototype: void ComputeGaussQuadPoint(int QuadPoint, int *a.TrngleNode, ## double *SrcPointCol) ##Description: To compute the coordinates of 7-point Gauss nodes of ## a triangular patch. ##Input value: ## int QuadPoint --- node index, it can be from 0 to 6. ## int *a.TrngleNode --- the three nodes of a tringular patch. ## double *SrcPointCol --- where to store the results ##Return value: none ##Global value used: a.NodeCord, a.Qpnt ##Global value modified: none ##Subroutines called: none ##Note: Not very sure. **************************************************************************** def ComputeGaussQuadPoint(QuadPoint,TrngleNode,a): SrcPointCol=zeros((3),Float) SrcPointCol[0] = a.Qpnt[QuadPoint,0]*a.NodeCord[TrngleNode[0],0]\ + a.Qpnt[QuadPoint,1]*a.NodeCord[TrngleNode[1],0]\ + a.Qpnt[QuadPoint,2]*a.NodeCord[TrngleNode[2],0] SrcPointCol[1] = a.Qpnt[QuadPoint,0]*a.NodeCord[TrngleNode[0],1]\ + a.Qpnt[QuadPoint,1]*a.NodeCord[TrngleNode[1],1]\ + a.Qpnt[QuadPoint,2]*a.NodeCord[TrngleNode[2],1] SrcPointCol[2] = a.Qpnt[QuadPoint,0]*a.NodeCord[TrngleNode[0],2]\ + a.Qpnt[QuadPoint,1]*a.NodeCord[TrngleNode[1],2]\ + a.Qpnt[QuadPoint,2]*a.NodeCord[TrngleNode[2],2] return SrcPointCol -- The Numeric Python EM Project www.members.home.net/europax

3 4

HTML manual display problems in Mozilla
by Edward C. Jones 25 Aug '01

25 Aug '01

Mozilla .0.9.3 is not correctly displaying the HTML version of the manual. The problem appears to be in footnotes. Netscape 4.76 had similar problems with the previous version of the manual. I use RedHat 7.1 Linux on a PC. Thanks, Ed Jones

1 0

RE: [Numpy-discussion] Subarray with with arbitrary index?
by Paul F. Dubois 24 Aug '01

24 Aug '01

I'm afraid that you have to do put using one-d indices. But you do *not* have to try to ravel the source. I.e., the first arg is just the name of the array. >>> from Numeric import * >>> x=array([[1,2,3],[4,5,6]]) >>> x array([[1, 2, 3], [4, 5, 6]]) >>> put(x,[0,4],[100,200]) >>> x array([[100, 2, 3], [ 4, 200, 6]]) >>> -----Original Message----- From: Huaiyu Zhu [mailto:huaiyu_zhu@yahoo.com] Sent: Friday, August 17, 2001 11:36 PM To: Paul F. Dubois Cc: John J. Lee Subject: RE: [Numpy-discussion] Subarray with with arbitrary index? Thanks, John and Paul. That is what I was looking for. It did not occur to me to look for verbs put and take, rather than words line sub- index, slice and so on. Maybe puting some of these words in the manual could help people doing a search? Now that this made the most costly part of my program about seven times faster, other problems become more prominent. One of such question is: How do we do it on more than one axis? Suppose a is a 2d array. Then put(a[1,:], b, c) works, but put(a[:,1], b, c) complains about the first argument not a continuous array. Doing transpose does not help. So do I guess it right that this is implemented only in the representation of a linear array? If so, there would be no hope of using put(a, ([2, 4], [1,2]), v) or even more exotic ones like using += on an arbitray subgrid? Regards, Huaiyu On Fri, 17 Aug 2001, Paul F. Dubois wrote: > John is right: > >>> a=Numeric.arange(8) > >>> b=Numeric.array([2,3,5]) > >>> c=Numeric.arange(3)+100 > >>> Numeric.put(a,b,c) > >>> print a > [ 0 1 100 101 4 102 6 7] > > Thanks for pointing out that I had left allclose out of the Numeric part of > the manual. I did it in the MA part and then forgot. I'm fixing it now. BTW: > There are changenotes at source forge that are sometimes ahead of the > manual. > > -----Original Message----- > From: numpy-discussion-admin(a)lists.sourceforge.net > [mailto:numpy-discussion-admin@lists.sourceforge.net]On Behalf Of John > J. Lee > Sent: Friday, August 17, 2001 6:36 AM > To: Huaiyu Zhu > Cc: numpy-discussion(a)lists.sourceforge.net > Subject: Re: [Numpy-discussion] Subarray with with arbitrary index? > > > On Thu, 16 Aug 2001, Huaiyu Zhu wrote: > > > Hi, > > > > Is it possible to assign to a subarray with arbitrary index? > > > > Suppose I have three arrays > > > > a = arange(8) > > b = array([2, 3, 5]) > > c = arange(3)+100 > > > > I want a function f, such that calling f(a, b, c) would change a to > > > > [0 1 100 101 4 102 6 7] > > f = Numeric.put > f(a, b, c) > > put used to be in Python, but it's been in C since some release 17.x.x, I > think. > > I have a sinking feeling that I must have missed something (no interpreter > here to check it works)... > > BTW, a week ago I noticed that I had reinvented the wheel in rewriting, in > an uglier and less efficient form, Numeric.allclose (hope I got the name > right). As far as I can see, it isn't listed in the manual. Did I miss > it? All it would need is the docstring copying over. > > > John > > > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)lists.sourceforge.net > http://lists.sourceforge.net/lists/listinfo/numpy-discussion >

3 2