Mailman 3 August 2001 - NumPy-Discussion

simple array help?
by Pete Shinners 17 Sep '01

17 Sep '01

i have a simple array with floating values >>> print a array([ 0. , 0.2, 0.4, 0.6, 0.8]) what i need to do is change each index into another dimension with 3 indices. hmm, kind of confusing, an example will explain better >>> b = array(zip(a,a,a)) >>> print b array([[ 0. , 0. , 0. ], [ 0.2, 0.2, 0.2], [ 0.4, 0.4, 0.4], [ 0.6, 0.6, 0.6], [ 0.8, 0.8, 0.8]]) ok, this does the job, but once i start using large arrays, the back-n-forth conversion between arrays and python lists is costing me quite a bit. is there a way i can reshape the array this way without the call to python's 'zip'? (no, 'map' doesn't count either) i've tried much fooling with NewAxis, resize, and friends. but either i haven't stumbled upon the correct combination or i'm not trying the correct tools. i'd like to keep this as quick as possible, so hopefully it can be done without anything too elaborate. thanks (doh, sorry for the previous distutils related post... hit the wrong mailing list)

6 19

matrix invert() problems
by Rob 31 Aug '01

31 Aug '01

I found a slight discrepancy in my Python EM code output when I use the invert() function. It may be that the C routine that I am comparing it to actually has the bug. For now I don't know. It would be cool to have a large two dimensional array that has a known inverse. I think NIST keeps an array collection, but I don't know if they have the inverses. Rob. -- The Numeric Python EM Project www.members.home.net/europax

1 1

Re:RE: [Numpy-discussion] Cholesky factorization
by Jeff Whitaker 31 Aug '01

31 Aug '01

I've uploaded a new version of my LinearAlgebra patch that adds a Cholesky decomposition routine (cholesky_decomposition), as per Frederik's request. -Jeff On Fri, 31 Aug 2001, Fredrik Stenberg wrote: > You read my mind...... > > ;)) > > /fredriks > > > ------- Ursprungligt meddelande ------- > Fr�n: jsw(a)cdc.noaa.gov > Datum: Fri, 31 Aug 2001 08:32:01 -0700 > > This message was sent from Geocrawler.com by "jeff whitaker" > <jsw(a)cdc.noaa.gov> > > > Fredrik: I could easily add this to the > LinearAlgebra patches I've just submitted. Are > you thinking of an interface to the lapack routine > dpotrf (cholesky factorization of a real symmteric > pos definite matrix), or something else? > > -Jeff > > --------------------------------------- > Hi.. > > Has anyone written a Cholesky factorization for > general > matrices? > > Any package out there? > > The only one i know of is the python > wrappers around lapack (pylapack). > > I could write it myself, but i prefer to be lazy ;)) > > /fredriks > > > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)lists.sourceforge.net > http://lists.sourceforge.net/lists/listinfo/numpy-discussion > > > > Geocrawler.com - The Knowledge Archive > > > > -- Jeffrey S. Whitaker Phone : (303)497-6313 Meteorologist FAX : (303)497-6449 NOAA/OAR/CDC R/CDC1 Email : jsw(a)cdc.noaa.gov 325 Broadway Web : www.cdc.noaa.gov/~jsw Boulder, CO, USA 80303-3328 Office : Skaggs Research Cntr 1D-124

1 0

RE: [Numpy-discussion] one-offset arrays
by Rodgers, Kevin 31 Aug '01

31 Aug '01

Indices in Fortran can be defined at the time an array is declared, so Fortran is not necessarily "one-indexed". Kevin Rodgers Northrop Grumman Ryan Aeronautical -----Original Message----- From: Eric Nodwell [mailto:nodwell@physics.ubc.ca] Sent: Tuesday, August 28, 2001 10:01 PM To: Numpy-discussion(a)lists.sourceforge.net Subject: Re: [Numpy-discussion] one-offset arrays >> 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 _______________________________________________ Numpy-discussion mailing list Numpy-discussion(a)lists.sourceforge.net http://lists.sourceforge.net/lists/listinfo/numpy-discussion

1 0

Cholesky factorization
by Fredrik Stenberg 31 Aug '01

31 Aug '01

Hi.. Has anyone written a Cholesky factorization for general matrices? Any package out there? The only one i know of is the python wrappers around lapack (pylapack). I could write it myself, but i prefer to be lazy ;)) /fredriks

2 1

MLab.py status: do you use it?
by Chris Barker 31 Aug '01

31 Aug '01

Hi all, Paul recently requested help with testing and debugging MLab.py, and I offered to help out. He then sent me this note: "Paul F. Dubois" wrote: > MLab.py started out as a package that was supposed to give numpy a > Matlab look and feel. So useful questions are: are there really useful > things that should be added? I can't answer these questions myself. Interestingly enough, despite the fact that I came to NumPy having used MATLAB heavily for five years (and one dissertation), I have made very little use of MLab.py. I'm wondering if anyone is, indeed, using it. As for trying to give numpy a MATLAB look and feel, I question the usefulness of that. NumPy looks and feels a little different than MATLAB, and, frankly, I like NumPy's feel better for the most part (now that rich comparisons have been added, anyway). The one thing that MATLAB has that NumPy doesn't, that I really miss, is list indexing. Having to use put and take is so much less elegant! I'm sure this isn't the least bit trivial to add, however. (is it even possible, given Python's idea of slicing?) Anyway, what I do think MATLAB does provide, and MLab.py should, is a whole bunch of utility functions for various common manipulations. > Do we really have what we have correct? > Most of all, these is no test routine for it. If you could make one > following the model of the test.py in directory Test, it would be great. > The idea is to have something that does not require the runner of the > test to know anything. It just runs silently unless something is wrong. 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. -Chris -- Christopher Barker, Ph.D. ChrisHBarker(a)home.net --- --- --- http://members.home.net/barkerlohmann ---@@ -----@@ -----@@ ------@@@ ------@@@ ------@@@ Oil Spill Modeling ------ @ ------ @ ------ @ Water Resources Engineering ------- --------- -------- Coastal and Fluvial Hydrodynamics -------------------------------------- ------------------------------------------------------------------------

4 4

Subarray with with arbitrary index?
by Huaiyu Zhu 31 Aug '01

31 Aug '01

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] This would be useful in certain applications that would otherwise require sparse matrices. Huaiyu Zhu

3 3

invert() function works just fine
by Rob 31 Aug '01

31 Aug '01

My problems with the invert() function has been resolved. I found that somewhere in my program it was trying to invert [0]. So I just added a conditional statement ahead of it. Evidently the C program did the same sort of thing. So now I've shaved another 20% of execution time off of my program! Rob. -- The Numeric Python EM Project www.members.home.net/europax

1 0

LinearAlgebra patch
by Jeff Whitaker 31 Aug '01

31 Aug '01

Has anyone tried my patch to LinearAlgebra.py and lapack_litemodule.c? I've found major speedups (by a factor of 3-5) in all operations involving an svd or a eigenanalysis of a symmetric matrix. Paul's initial comment was that since it was not self-contained (it required linking external blas and lapack libs) he could not consider putting it in. I have rectified that now by including new versions of f2c_lite.c, blas_lite.c and d(z)lapack_lite.c that support the newer, faster lapack routines. I would appreciate any reports on how it works (or doesn't) on other platforms - I've tried it in on MacOS X and Solaris. Grab it from the 'Patches' section of the sourceforge project page. -Jeff -- Jeffrey S. Whitaker Phone : (303)497-6313 Meteorologist FAX : (303)497-6449 NOAA/OAR/CDC R/CDC1 Email : jsw(a)cdc.noaa.gov 325 Broadway Web : www.cdc.noaa.gov/~jsw Boulder, CO, USA 80303-3328 Office : Skaggs Research Cntr 1D-124

3 3

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

30 Aug '01

Thanks to everyone who responded to my comments about one-offset arrays in Python. I understand much better now why zero-offset arrays are a good choice for Python. For the reasons already discussed (mostly ease of translating one-offset algorithms), one-offset arrays would nevertheless be useful in certain situations, so I went ahead and wrote a class for a one-offset array. It ended up being a somewhat bigger job than I expected, and since it is rather too much material for a mailing-list submission I created a web site for it. Anyone who is interested in this topic please have a look: http://www.physics.ubc.ca/~mbelab/python/arrayone/index.html I'd love some feedback along the following lines: 1. There are 4 outstanding problems/bugs which I could identify but was unable to fix. Most of these seem to be due to limitations of the UserArray class, but they could also be due to limitations of the programmer. None of these problems make the class unusable. One of them is the issue which Travis identified about functions taking UserArrays but returning standard arrays. It would be nice if the resulting discussion led to something be done... Any suggestions for fixing the other 3 issues would be most welcome. 2. I would like this module to be as accessible as possible to others, particularly those who are just getting started with Python and may therefore be especially encumbered with one-offset baggage. How can I do this? Submit it to a Python resource such as the Vaults of Parnassus? Get people to link to my web page? Does someone with a more prominant website want to volunteer to host it? Is there any possibility of having it included eventually with Numerical Python or Scientific Python? Probably some people will consider it too trivial to include, but I found it less trivial than I expected once all the cases were covered - and there are still the 4 outstanding issues. (Yes the code sucks! Suggestions for improvement are most welcome!) Why make people reinvent the wheel? People coming from MatLab for example might be inclined to give up and go back to MatLab. Eric P.S. In case anyone is interested in the outstanding problems but for some reason doesn't want to or can't visit the web site, here are the summaries: Problem 1 --------- Description for loops don't work with variables of type arrayone Example X=arrayone((1,2,3)) for item in X: print item This just generates an error. Reason In Python, a for loop works by starting at x(0) and incrementing until an out-of-bounds error occurs. arrayone types have no 0 element. Work-around: Cast to type array in for loops. For example X=arrayone((1,2,3)) for item in array(X): print item Possible solutions: Is it possible to "overload" "for" so that it behaves differently for arrayone type? Problem 2 --------- Description Slicing an arrayone from 0 doesn't generate an error but it should. Example: X=arrayone((1,2,3)) X[0:3] This returns ArrayOne.arrayone([1, 2, 3]) instead of an error. Reason X[:3] results in a call to __getslice__ with low=0. This cannot be distinguished from X[0:3]. Therefore in order to deal correctly with the X[:] case, we have to assume that low=0 means an unspecified low. Work-around If you don't trust your code, you have to implement specific checking for this condition before slicing. Possible solution If it was possible to get access to the unparsed input (literally '0:3' for example), this could be fixed. Problem 3 --------- Description Negative subscripts return a slice offset by -1 from the expected result. Example X=arrayone((1,2,3,4)) X[-3:-2] This returns ArrayOne.arrayone([1, 2]) instead of ArrayOne.arrayone([2, 3]). Reason X[-3:-2] in the above example results in a call to __getslice__ with low=1, high=2, which cannot be distinguished from X[1:2]. Work-around Don't use negative index slicing with arrayone types. Possible solution: If had access to unparsed input, could fix this. Problem 4 --------- Description Functions which take arrayone arguments return type array. Work-around Have to cast back into arrayone, for example: Y=arrayone(somefunc(X)) Possible solution See http://www.geocrawler.com/lists/3/SourceForge/1329/0/6505926

2 2

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 26 Aug '01

26 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

For the next round on the manual
by Eric Hagemann 23 Aug '01

23 Aug '01

Type-o in the manual Check the definition for the 'convolve' function. Manual tells of mode=0 as default when code and experiment point to mode=2 being default No biggie -- but it might be nice to correct if anybody is editing the manual any time soon Cheers Eric

3 2

patch to speed up LinearAlgebra module
by Jeff Whitaker 23 Aug '01

23 Aug '01

I've uploaded a patch to LinearAlgebra.py and lapack_litemodule.c that speeds up some computations (svd, linear least-squares, and hermitian eigenanalysis) by 3-5 times when linking with external (preferably optimized) blas and lapack libs. The patch replaces calls to dgesvd, dgelss and dsyev with the newer, faster lapack routines dgesdd, dgelsd and dsyevd. Here are my timing results before and after (G4 350mhz powermac, optimized ATLAS blas lib and lapack 3.0 from netlib.org): Before: svd of 1000X1000 matrix takes 609.19 generalized inverse of 1000X1000 matrix takes 744.36 linear least-squares solution of 1000X1000 matrix takes 451.68 eigenvectors of 1000X1000 symmetric matrix takes 210.08 After: svd of 1000X1000 matrix takes 142.55 generalized inverse of 1000X1000 matrix takes 237.08 linear least-squares solution of 1000X1000 matrix takes 81.91 eigenvectors of 1000X1000 symmetric matrix takes 56.23 Note that since these newer lapack routines are not in lapack_lite, you must link external blas and lapack3 libs. -Jeff -- Jeffrey S. Whitaker Phone : (303)497-6313 Meteorologist FAX : (303)497-6449 NOAA/OAR/CDC R/CDC1 Email : jsw(a)cdc.noaa.gov 325 Broadway Web : www.cdc.noaa.gov/~jsw Boulder, CO, USA 80303-3328 Office : Skaggs Research Cntr 1D-124

1 0

Matlab_NumPy_Array?
by Karshi Hasanov 22 Aug '01

22 Aug '01

Hi all, I don't understand how to conver arrays from Matlab to Numpy. Here what I did: 1) Converted Matlab A(i,j,k) to 'A.dat' binary file 2) file_open=open('~/A.dat', 'rb') 3) read_file=file_open.read() 4) A_string=fromstring(read_file, 'f') 5) A=reshape(A_string, (i,j,k)) where i,j,k are size of A[i,j,k] I am not getting the right. What is the best/right way of convering arrays? Thanks

2 1

MA patch submitter, please upload file
by Paul F. Dubois 18 Aug '01

18 Aug '01

Would the person who submitted the patch concerning MA binary ops, please upload the file to the patch, or email it to me? You have identified a real problem whose answer is not obvious so I'd like to see your idea of how to fix it. If you log in to SourceForge before submitting bugs and patches then I can know who you are to thank you or make followup inquiries like this one.

1 0

SWIG + Numeric 20.1.0 = core dump: the fix
by Edward C. Jones 12 Aug '01

12 Aug '01

If Numeric and SWIG are used together in creating a Python extension and there are C/C++ files as part of the extension then Numeric 20.1.0 will probably break your code. This happens because changes were made to arrayobject.h between Numeric 20.0.0 and 20.1.0. The changes that need to be made in your code are explained by Konrad Hinsen: (http://www.geocrawler.com/archives/3/1329/2001/5/0/5718231/) FROM: Konrad Hinsen DATE: 05/09/2001 09:12:25 SUBJECT: [Numpy-discussion] Modified header files Recently we had a discussion about how to use NumPy arrays from extension modules with multiple source files, on various platforms. The modified header files that are attached to this message provide a (hopefully!) general solution. In fact, I propose to make them part of the official distribution, unless there are objections. If used like before, these header files give exactly the same result as the ones in NumPy 20.0.0. However, they permit to define the name of the C API pointer array and make it globally visible. Under the condition that the chosen name is unique, this should not create problems under any platform, no matter if static or dynamic linking is used. To use NumPy features in multiple source file extension modules, you have to write #define PY_ARRAY_UNIQUE_SYMBOL PyArrayXXX #include "Numeric/arrayobject.h" in the main source file (the one that contains the init function) and #define PY_ARRAY_UNIQUE_SYMBOL PyArrayXXX #define NO_IMPORT_ARRAY #include "Numeric/arrayobject.h" in all the others. The symbol you choose instead of PyArrayXXX should contain both the name of the imported module (array) and the name of the importing module (whatever your module is called) in order to be unique with a reasonably high probability. The same applies to the Ufunc module, just replace "array" by "ufunc" in the example. I have also applied the "static" correction to the Ufunc header file, there is no reason not to do it. Konrad. -- ----------------------------------------------------------------------- Konrad Hinsen | E-Mail: hinsen(a)NO-SPAM.cnrs-orleans.fr Centre de Biophysique Moleculaire (CNRS) | Tel.: +33-2.38.25.56.24 Rue Charles Sadron | Fax: +33-2.38.63.15.17 45071 Orleans Cedex 2 | Deutsch/Esperanto/English/ France | Nederlands/Francais ----------------------------------------------------------------------- I do this by modifying one header file that is included in all the ".c" and ".i" files. In the following example the header file is "spam.h". ------------------------------------------------------------ doit.py: #! /usr/bin/python import Numeric, spam spam.funct() print 'spam.funct completed' ------------------------------------------------------------ spam.h: #include "Python.h" #define PY_ARRAY_UNIQUE_SYMBOL Py_Array_API_spam #ifndef SWIG_FILE_WITH_INIT #define NO_IMPORT_ARRAY #endif #include "/usr/include/python2.1/Numeric/arrayobject.h" void funct(void); ------------------------------------------------------------ spam.i: %module spam %{ #define SWIG_FILE_WITH_INIT #include "spam.h" %} %init %{ import_array(); printf("called import_array\n"); %} void funct(void); ------------------------------------------------------------ spammodule.c #include "spam.h" void funct(void) { PyArrayObject *pao; int dims[2]; dims[0] = 100; dims[1] = 200; printf("Calling PyArray_FromDims\n"); pao = (PyArrayObject*) PyArray_FromDims(2, dims, PyArray_UBYTE); printf("Completed PyArray_FromDims\n"); Py_DECREF(pao); } ------------------------------------------------------------ compile script: swig -python spam.i gcc -c -Wall spammodule.c spam_wrap.c -I/usr/include/python2.1 ld -shared spammodule.o spam_wrap.o -o spam.so ------------------------------------------------------------ I use RedHat 7.1 Linux on a PC, Python 2.1, Numeric 20.1.0, SWIG 1.1 build 883, and gcc 2.96. Anyone is welcome to use the above example in documentation.

1 0

SWIG / Numeric core dump
by Edward C. Jones 10 Aug '01

10 Aug '01

The following Python / SWIG / Numeric program dumps core. I use RedHat Linux 7.1 on a PC. Python 2.1, Numeric 20.1.0, SWIG 1.1 build 883, gcc 2.96. What is the problem? The output looks like: > doit.py called import_array Calling PyArray_FromDims Segmentation fault (core dumped) > The program files are: ---------------------------------------------- doit.py: #! /usr/bin/python import Numeric, spam spam.funct() print 'spam.funct completed' ---------------------------------------------- spam.h: #include "Python.h" #include "/usr/include/python2.1/Numeric/arrayobject.h" void funct(void); ---------------------------------------------- spam.i %module spam %{ #include "spam.h" %} %init %{ import_array(); printf("called import_array\n"); %} void funct(void); ---------------------------------------------- spammodule.c #include "spam.h" void funct(void) { PyArrayObject *pao; int dims[2]; dims[0] = 100; dims[1] = 200; printf("Calling PyArray_FromDims\n"); pao = (PyArrayObject*) PyArray_FromDims(2, dims, PyArray_UBYTE); printf("Completed PyArray_FromDims\n"); Py_DECREF(pao); } ---------------------------------------------- swigit: swig -python spam.i gcc -c -Wall spammodule.c spam_wrap.c -I/usr/include/python2.1 ld -shared spammodule.o spam_wrap.o -o spam.so

2 1

NumPy installation instructions for Windows
by Just van Rossum 10 Aug '01

10 Aug '01

I develop a package that depends on NumPy, and I would like to add a line or two to my documentation about how to install NumPy under Windows -- just the basic set of modules, not the extensions. I understand users should grab the .exe archive. The FAQ only mentions "The .exe files are Windows installers of the sort that you are used to". Does this mean everything is explained and done by the installer itself, or is there some additional action the user should perform? Thanks! (I have no Windows box here so I can't see for myself...) Just

2 1

swig and C/C++
by Tariq Rashid 08 Aug '01

08 Aug '01

in fact ... even simpler ... can anyone point me to an example of passing Numeric arrays to and from C ... the swig examples are fine for integers / floats ... but as i'm a new to this stuff I need a bit if guidance with Numeric arrays ... any ideas / pointers would be very gratefully receieved... i'm working on wavelet stuff which i want to give back to the Numeric community.... thanks tariq (tariq_rashid at lineone.net) ---------------------------- Today's Topics: 1. NumTut (=?ks_c_5601-1987?B?wMyx1LrA?=) 2. swig and Numeric (Paul F. Dubois)

1 0

fast constructor for arrays from byte data (unpickling?)
by Hung Jung Lu 08 Aug '01

08 Aug '01

Hi, Arrays have a method called tostring() which generates the binary data. Is there an inverse function to that? That is, generating an array from the binary data string? For large matrices, pickling/unpickling is a bit too much overhead (data stored as ASCII character strings instead of binary data strings.) I know, I am talking about a factor 4 here. But there is a big difference between 1 minute loading time and 4 minute loading time. I would imagine this is a very common problem/request for people working with large matrices. And I am sure hacking the C code to provide another fast constructor for arrays from binary strings won't be too hard. The questions is: has anyone already tried it? Is it already there? (For the kludge masters: one kludge is of course to store the binary data on disk, then use cStringIO to build the pickled file and then unpickle from the cStringIO. Speed is probably OK since the pickled file lives on RAM. But that's a kludge. :) ) regards, Hung Jung __________________________________________________ Do You Yahoo!? Make international calls for as low as $.04/minute with Yahoo! Messenger http://phonecard.yahoo.com/

4 11

mmap (was: fast constructor for arrays from byte data)
by Phil Austin 08 Aug '01

08 Aug '01

Chris Barker writes: > I don't understand mmap at all. From the name, it sounds like the entire > contents of the file is mapped into memory, so the memory would get used > as soon as you set it up. If anyone knows, I'd like to hear... This might be more than you need to know, but there's a copy of the article: "Why aren't you using mmap() yet?" by Kevin Sheehan at: http://thrush.eos.ubc.ca/users/phil/mmap/mmap.pdf (note reversed page order) and a Solaris/SGI email debate in mhonarc format: http://thrush.eos.ubc.ca/users/phil/mmap/threads.html bottom line: you really need a vm implementation that supports madvise() to get the maximum benefit from mmap. (Linux 2.4 might have madvise, I haven't tried it out yet -- 2.2 didn't). The way I used to use mmap was through a numpy array class that was initiated from a pointer to a mmap'd file. This made reading/writing to the file as simple as indexing into the numpy array. I haven't been back to reimplement that class on linux since I left solaris, though. You can also get some of the benefits of mmap via netcdf, which will use mmap if it exists. Regards, Phil

2 1

swig and Numeric
by Paul F. Dubois 08 Aug '01

08 Aug '01

Does someone have an example of using SWIG to connect Numeric arrays to C code that expects double* arguments?

2 1

permutation and array subsets
by Serge Rey 08 Aug '01

08 Aug '01

i'm trying to use the results of permutation (from RandomArray.py) to reorder the rows of an array: >>> r=arange(20) >>> id=permutation(20) >>> y=r[id] Traceback (most recent call last): File "<stdin>", line 1, in ? IndexError: invalid index >>> i can use a for loop to permutate the ordering of r, but that has a speed hit associated with it. can anyone offer a suggestion for a way to do this more efficiently? thanks, serge -- Sergio J. Rey http://typhoon.sdsu.edu/rey.html

2 1

swig and Numeric
by Paul F. Dubois 07 Aug '01

07 Aug '01

Does someone have an example of using SWIG to connect Numeric arrays to C= code that expects double* arguments?

1 0

swig and Numeric
by Paul F. Dubois 07 Aug '01

07 Aug '01

Does someone have an example of using SWIG to connect Numeric arrays to C= code that expects double* arguments?

1 0

NumTut
by 이규봉 07 Aug '01

07 Aug '01

Dear who may be concerned, I want to install the package NumTut. But I can't find the directory of source distribution. Please let me know how to get the directory odf souce distribution and the subdirectory Demo. Thanks Gyoy-Bong Lee from Korea

1 0

new to list
by Rob 03 Aug '01

03 Aug '01

Hello, I'm new to the list so I thought I'd introduce myself. I use Numpy at work for some Python programs that interface with Agilent and/or Anritsu microwave equipment, using the GPIB bus. However, my main interest lately is porting some EM codes into Numeric Python. I've had varied success. See my web site www.members.home.net/europax. Right now I'm trying to Numpy-ize this FEM-MOM code. I've profiled it using the Python profile class. Using MOM for a dipole takes 5 minutes of cpu :) :) Interestingly, some of the worst routines are simple ones that really have little to do with matrix math. Anyway, I'm having fun. Its wonderful to be able to focus on algorithms rather than wading through arcane programming syntax !! Rob.

1 0

A savespace bug, perhaps?
by Nathaniel Gray 02 Aug '01

02 Aug '01

Hi everybody, I'm using Python 2.1.1 and Numpy 20.1.0 This seems to be a problem: Python 2.1.1 (#2, Jul 31 2001, 14:10:42) [GCC 2.96 20000731 (Linux-Mandrake 8.0 2.96-0.48mdk)] on linux2 Type "copyright", "credits" or "license" for more information. >>> from Numeric import * >>> from MLab import zeros >>> joe = zeros(3,'l',1) >>> fred = array([1,2,3],'l') >>> joe += fred Traceback (most recent call last): File "<stdin>", line 1, in ? TypeError: return array has incorrect type This fails whether or not fred has the spacesaver flag set. This, however, works fine: >>> joe = zeros(3, 'l') >>> fred = array([1,2,3], 'l') >>> joe += fred >>> joe array([1, 2, 3]) Am I doing something wrong or is this a bug? Please CC to me since I'm not subscribed to the list. Thanks a lot! -Nathan -- Nathaniel Gray California Institute of Technology Computation and Neural Systems --

2 1