numpy-discussion — Discussion list for all users of Numerical Python

Re: [Numpy-discussion] performance problem

[Gerard V. <ger...@gr...>]

On Mon, 30 Jan 2006 01:52:53 -0700
Travis Oliphant <oli...@ie...> wrote:

> Gerard Vermeulen wrote:
> 
> >coercion issue snipped
> >  
> >
> 
> In current SVN, I think improved on the current state by only calling a 
> scalar a signed integer if it is actually negative (previously only it's 
> data-type was checked and all Python integers get converted to 
> PyArray_LONG data-types which are signed integers).
> 
> This fixes the immediate problem, I think.
> 
> What are opinions on this scalar-coercion rule?
> 

Hmm, this is a consequence of your rule:

>>> from numpy import *; core.__version__
'0.9.5.2024'
>>> a = arange(3, dtype=uint32)
>>> a-3
array([4294967293, 4294967294, 4294967295], dtype=uint32)
>>> a+-3
array([-3, -2, -1], dtype=int64)
>>> (a-3) == (a+-3)
array([False, False, False], dtype=bool)

Do you think that the speed-up justifies this? I don't think so.


All my performance issues are discovered while writing demo examples
which transfer data between a Python wrapped C++-library and Numeric,
numarray, or numpy.  In that state of mind, it surprises me when an
uint32 ndarray gets coerced to an int64 ndarray.

I rather prefer numarray's approach which raises an overflow error for the
>>> a+-3
above.

Agreed that sometimes a silent coercion is a good thing, but somebody
who has started with an ubyte ndarray is likely to want an ubyte array in
the end.

I don't want to start a flame war, happy to write
a - uint32(3)
for numpy specific code.

Gerard

Re: [Numpy-discussion] Numerical Mathematics Consortium vs. scipy and numpy

[Christopher B. <Chr...@no...>]

Alan G Isaac wrote:
> On Mon, 30 Jan 2006, Christopher Barker apparently wrote:=20
>> how does one get the open source community to be=20
>> represented in this kind of thing?=20
>=20
> http://www.nmconsortium.org/FldRte/?id=3D70&page=3DDownloads

Were you pointing us to any particular document there?

Full membership (that is, with voting rights) requires a $2500.00 a y=
ear=20
membership fee. I think it'll be a bit tricky for an open source proj=
ect=20
to raise that kind of cash. The result is that you have to buy your w=
ay=20
into the group, which favors commercial entities.

Granted, they have operating expenses, and $2500 isn't all that much,=
=20
but it's probably too much for open-source groups to get representati=
on.

An associate membership is $250./year, which is far more manageable, =
but=20
then you don't get voting rights:

"""
The Associate Membership with no voting rights will have full access =
to=20
the consortium=92s output and be engaged in some of the consortium=
=92s=20
activities and is typically for the user community, both academia and=
=20
industry ($200 annual fee).
"""

The SciPy community is not a "user" group, it's a developer group. th=
ere=20
are other open source projects worth of representation, such as the G=
NU=20
Scientific Library.

Oh well.

-Chris


--=20
Christopher Barker, Ph.D.
Oceanographer
                                     =09=09
NOAA/OR&R/HAZMAT         (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

Chr...@no...

Re: [Numpy-discussion] performance problem

[Travis O. <oli...@ee...>]

Gerard Vermeulen wrote:

>0.9.5.2021
>  
>
>>>>t1 = Timer('a <<= 8', 'import numarray as NX; a = NX.ones(10**6, NX.UInt32)')
>>>>t2 = Timer('a <<= 8', 'import numpy as NX; a = NX.ones(10**6, NX.UInt32)')
>>>>t1.timeit(100)
>>>>        
>>>>
>0.21813011169433594
>  
>
>>>>t2.timeit(100)
>>>>        
>>>>
>1.1523458957672119
>  
>
While ultimately, this slow-down was related to a coercion issue, I did 
still wonder about the extra dereference in the 1-d loop when one of the 
inputs is a scalar.  So, I added a patch that checks for that case and 
defines a different loop.

It seemed to give a small performance boost on my system.  I'm wondering 
if such special-case coding is wise in general.  Are there other ways to 
get C-compilers to produce faster code on modern machines?

-Travis

Re: [Numpy-discussion] Verbose output from setup.py --help intentional?

[Pearu P. <pe...@sc...>]

On Mon, 30 Jan 2006, Matthew Brett wrote:

> Hi,
>
> Sorry if this is low priority, but I noticed that the output from
> setup.py --help is more than 400 lines. Also the output of (e.g)
> setup.py --help install and setup.py install --help differ (also more
> than 400 lines).

This issue is fixed in numpy SVN.

Pearu

Re[2]: [Numpy-discussion] Numerical Mathematics Consortium vs. scipy and numpy

[Alan G I. <ai...@am...>]

On Mon, 30 Jan 2006, Christopher Barker apparently wrote: 
> how does one get the open source community to be 
> represented in this kind of thing? 

http://www.nmconsortium.org/FldRte/?id=70&page=Downloads

hth,
Alan Isaac

Re: [Numpy-discussion] Numerical Mathematics Consortium vs. scipy and numpy

[Christopher B. <Chr...@no...>]

Zdene(k Hur=E1k wrote:

> the other day I ran across the following link ---->=20
> http://www.nmconsortium.org/.=20

Thanks for the pointer, it looks like a great project, let's hope it=
=20
gets some traction.

A few thoughts that are NumPy related:

-- If nothing else, it's a pretty good  list of stuff that "should" b=
e=20
included in Stock NumPy (and/or the core of SciPy).

-- I notice they have:

Matrix Operations, like: Elementwise Addition, etc.

I don't, however, see any N-d array stuff: Too bad. (they do have=
=20
"Multi-Dimensional Arrays Of Matrices" listed under present additiona=
l=20
topics. Why array of matrices, rather than arrays of scalars?

-- I'm personally not thrilled with the linear-algebra focus: for=
=20
instance the matlab style for the operators:

*  matrix multiply
=2E* elementwise addition.

I'd rather see it the other way around.

I wonder if we could get a Pythonista to join? In general,=20
unfortunately, it looks to be quite commercially-focused: how does on=
e=20
get the open source community to be represented in this kind of thing=
?

http://www.nmconsortium.org/client-area/registration/

-Chris



--=20
Christopher Barker, Ph.D.
Oceanographer
                                     =09=09
NOAA/OR&R/HAZMAT         (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

Chr...@no...

[Numpy-discussion] (some) generalized eigenvalues in numpy

[Sven S. <sve...@gm...>]

Hi,
in econometrics/statistics it's relatively common to solve a special
generalized eigenvalue problem, where the matrices are real, symmetric,
and positive definite. It bothered me to depend on scipy just because of
that. Also, scipy.linalg.eig is much more general and returns complex
types which is inconvenient for sorting etc. So I've written the
workaround function below, but it's probably not very efficient (at
least I hope the algebra is right).

I'd be happy to hear your comments, suggestions; first with respect to
the function itself, but also whether you agree that the implemented
functionality would be useful enough to include in official numpy (in a
more professional way than I have done).

cheers,
Sven


def geneigsympos(A, B):
    """ Solves symmetric-positive-definite generalized
    eigenvalue problem Az=lBz.

    Takes two real-valued symmetric matrices A and B (B must also be
    positive-definite) and returns the corresponding generalized (but
    also real-valued) eigenvalues and eigenvectors. Return format: as
    in scipy.linalg.eig, tuple (l, Z); l is directly taken from eigh
    output (a 1-dim array of length A.shape[0] ?), and Z is an array
    or matrix (depending on type of input A) with the corresponding
    eigenvectors in columns (hopefully).
    Eigenvalues are ordered ascending (thanks to eigh).
    """
    from numpy import asmatrix, asarray, linalg
    #fixme: input checks on the matrices
    LI = asmatrix(linalg.cholesky(B)).I
    C = LI * asmatrix(A) * LI.T
    evals, evecs = linalg.eigh(C)
    if type(A) == type(asarray(A)): output = "array"
    # A was passed as numpy-array
    else: output = "matrix"
    #but the evecs need to be transformed back:
    evecs = LI.T * asmatrix(evecs)
    if output == "array": return evals, asarray(evecs)
    else:   return asmatrix(evals), evecs

Re: [Numpy-discussion] ATLAS binaries - where are we?

[Travis O. <oli...@ie...>]

Matthew Brett wrote:

>Hi,
>
>I wonder if the forthcoming release of numpy is a good stimulus to
>revisit the issue of ATLAS binaries, as for:
>
>http://www.scipy.org/download/atlasbinaries/
>
>and
>
>http://www.kevinsheppard.com/research/matlabatlas/matlab_atlas.aspx
>
>It seems to me that it will be very important to have ATLAS binaries,
>at very least for windows, which install seamlessly with NumPy.
>
>I am very happy to contribute builds, and can probably make a little
>windows compile farm at my lab (they _love_ windows where I work).  
>Would that be useful?  How would that best be integrated with the
>distribution?
>  
>
Absolutely useful...   I suggest making a page on the new.scipy.org Wiki 
and uploading them there (to upload to the Wiki I think you create an 
attachment using attachment:<filename> Then, when the page is rendered 
you can click on the resulting link to upload the file).

We could easily place links to this site for downloading. 

The old SciPy site already had such a page of pre-compiled ATLAS builds.

Thanks for offering.

-Travis

Re: [Numpy-discussion] Problems with creating arrays from array descriptors and lists

[Travis O. <oli...@ie...>]

Francesc Altet wrote:

>El dl 30 de 01 del 2006 a les 02:23 -0700, en/na Travis Oliphant va
>escriure:
>  
>
>>>In [20]:numpy.array([1,1], dtype=[('f1', 'i8'),('f2','i4')])
>>>---------------------------------------------------------------------------
>>>exceptions.TypeError                                 Traceback (most
>>>recent call last)
>>>
>>>/home/faltet/<ipython console>
>>>
>>>TypeError: expected a readable buffer object
>>>
>>>Sorry, I cannot check with latest SVN because it segfaults.
>>> 
>>>
>>>      
>>>
>>It shouldn't.  Latest SVN works fine.
>>    
>>
>
>True.
>
>  
>
>>But, your example does not work (Records are only settable with tuples).
>>
>>numpy.array((1,1), ....)  would work (returning an 0-d array)
>>    
>>
>
>Why is this? I think that allowing lists as well would be highly useful.
>In fact, I'm sure that many users are more used to provide lists instead
>of tuples. Moreover, allowing tuples and not lists can be a bit
>misleading, IMO.
>  
>
You can provide lists.   It's just that fields must be entered in the 
array constructor as tuples to fill them in.  This was the easiest to 
implement (it fit in naturally with the code already written from 
Numeric, and does not seem to me to be a real limitation and in fact can 
make things more readable.

Thus,   you can do

numpy.array([[obj1, obj2, obj3]], dtype=<withfields>)

as long as the objects are tuples.

-Travis

[Numpy-discussion] ATLAS binaries - where are we?

[Matthew B. <mat...@gm...>]

Hi,

I wonder if the forthcoming release of numpy is a good stimulus to
revisit the issue of ATLAS binaries, as for:

http://www.scipy.org/download/atlasbinaries/

and

http://www.kevinsheppard.com/research/matlabatlas/matlab_atlas.aspx

It seems to me that it will be very important to have ATLAS binaries,
at very least for windows, which install seamlessly with NumPy.

I am very happy to contribute builds, and can probably make a little
windows compile farm at my lab (they _love_ windows where I work). =20
Would that be useful?  How would that best be integrated with the
distribution?

Best,

Matthew

[Numpy-discussion] Numerical Mathematics Consortium vs. scipy and numpy

[<hu...@co...>]

Hello all,

the other day I ran across the following link ----> 
http://www.nmconsortium.org/. Are you, Python numerics guys, aware of this?
I admit I have absolutely no idea about the actual usefullness and
relevancy of this consortium and their drafts. But it might be useful just
to be aware... Important players in scientific computation are involved.

Best regards,
Zdenek Hurak

Re: [Numpy-discussion] Segfault with dodgy type descriptors

[Travis O. <oli...@ie...>]

Ed Schofield wrote:

>Hi Travis,
>
>Here's a segfault from the crazy input department:
>
>  
>
>>>>cs = [('A', 1), ('B', 2)]
>>>>a = numpy.array(cs, dtype=(object,5))
>>>>        
>>>>
>
>Using the descriptor dtype=(object,x) for other values of x gives
>varying results, but all of them should probably raise a ValueError.
>  
>
These data-types are technically o.k. because (object, 5) means a 
data-type which is an array of 5 objects.

The problem was that your constructor was not creating 5 objects and so 
an error should have been raised.   It is now.

>On a related note, it seems that it's now possible again to use an array
>as a dtype.  You changed this behaviour in November after I suggested
>that this should be illegal.  (See
>http://www.scipy.net/pipermail/scipy-dev/2005-November/004126.html).  My
>argument then was that calling array(a, b.dtype) is clearer and safer
>than using array(a, b), which can be confused with array((a,b)).  Is it
>an oversight that arrays can again be used as dtypes?
>  
>
I think this is an oversight, but I'm not exactly sure how to fix it.   
Basically during design of the new data-types,  I allowed the 
possibility that any object with a dtype attribute could be considered a 
data-type.   I wasn't thinking about arrays (which of course have that 
attribute as well), I was thinking about record objects with nested 
records and making it relatively easy to specify record data-types using 
objects with a dtype attribute.  

But, I could be persuaded (again), that it is easy enough to request the 
data-type attribute should you every actually need it, then have it 
specially-coded in the descriptor converter.   In fact, I don't see any 
real benefit to it at all.  My thinking about record objects changed 
during development and this could be just a holdover.

-Travis

[Numpy-discussion] Problems with creating arrays from array descriptors and lists

[Francesc A. <fa...@ca...>]

Hi,

I think this worked some time ago, but I cannot tell for sure:

In [19]:numpy.__version__
Out[19]:'0.9.5.1980'

In [20]:numpy.array([1,1], dtype=3D[('f1', 'i8'),('f2','i4')])
---------------------------------------------------------------------------
exceptions.TypeError                                 Traceback (most
recent call last)

/home/faltet/<ipython console>

TypeError: expected a readable buffer object

Sorry, I cannot check with latest SVN because it segfaults.

--=20
>0,0<   Francesc Altet     http://www.carabos.com/
V   V   C=E1rabos Coop. V.   Enjoy Data
 "-"

Re: [Numpy-discussion] performance problem

[Travis O. <oli...@ie...>]

Gerard Vermeulen wrote:

>coercion issue snipped
>  
>

In current SVN, I think improved on the current state by only calling a 
scalar a signed integer if it is actually negative (previously only it's 
data-type was checked and all Python integers get converted to 
PyArray_LONG data-types which are signed integers).

This fixes the immediate problem, I think.

What are opinions on this scalar-coercion rule?

-Travis

[Numpy-discussion] Verbose output from setup.py --help intentional?

[Matthew B. <mat...@gm...>]

Hi,

Sorry if this is low priority, but I noticed that the output from
setup.py --help is more than 400 lines. Also the output of (e.g)
setup.py --help install and setup.py install --help differ (also more
than 400 lines).

Best,

Matthew

Re: [Numpy-discussion] performance problem

[Travis O. <oli...@ie...>]

Gerard Vermeulen wrote:

>Hi Travis, thanks for your incredible quick fixes.
>
>numpy is about 5 times slower than numarray (on my numarray-friendly bi-pentium):
>
>  
>
This is a coercion issue. 

NumPy is exercising the most complicated section of the ufunc code 
requiring buffered casting because you are asking it to combine a uint32 
array with a signed scalar (so the only commensurable type is a signed 
scalar of type int64.  This is then coerced back into the unsigned array).

Look at the data-type of b = a << 8.

(a<<8).dtype.name
'int64'


Try

a<<=val  where val = uint32(8) is in the header

You should see more commensurate times.


We can of course discuss the appropriateness of the coercion rules.   
Basically scalars don't cause coercion unless they are of a different 
kind but as of now signed and unsigned integers are considered to be of 
different kinds.   I think that there is a valid point to be made that 
all scalar integers should be treated the same since Python only has one 
way to enter an integer.  

Right now, this is not what's done, but it could be changed rather easily.

-Travis

Re: [Numpy-discussion] Segfaults when trying to build an array using the description in array protocol

[Francesc A. <fa...@ca...>]

El dg 29 de 01 del 2006 a les 23:59 -0700, en/na Travis Oliphant va
escriure:
> I think you can report anonymously.

I've just done it. However, I'd like to be informed about the follow-ups
of the ticket by e-mail. I've put my e-mail address in the ticket,
although I'm not sure if I'll get e-mail on this.

Cheers,

--=20
>0,0<   Francesc Altet     http://www.carabos.com/
V   V   C=E1rabos Coop. V.   Enjoy Data
 "-"

[Numpy-discussion] Segfaults when trying to build an array using the description in array protocol

[Francesc A. <fa...@ca...>]

Hi,

Apparently something wrong went in latest version of numpy SVN:

In [4]:numpy.__version__
Out[4]:'0.9.5.2021'

In [5]:numpy.empty((1,), dtype=3D'i8')
Out[5]:array([-5193749191705487200], dtype=3Dint64)

In [6]:numpy.empty((1,), dtype=3D[('p', 'i8')])
Segmentation fault

This usually worked well with older versions:

In [7]:numpy.__version__
Out[7]:'0.9.5.1980'

In [8]:numpy.empty((1,), dtype=3D[('p', 'i8')])
Out[8]:array([(-5192306666809247592L,)], dtype=3D(void,8))

BTW, which is the recommended way to report that sort of things, the
list or the Trac interface for SVN numpy? I'd like to start using the
second, but I haven't manage to create a user for reporting issues.

Thanks,

--=20
>0,0<   Francesc Altet     http://www.carabos.com/
V   V   C=E1rabos Coop. V.   Enjoy Data
 "-"

Re: [Numpy-discussion] another performance problem

[Travis O. <oli...@ie...>]

Gerard Vermeulen wrote:

>Hi Travis,
>
>max and min are really slow in numpy:
>
>  
>
>>>>from timeit import Timer
>>>>import Numeric; Numeric.__version__
>>>>        
>>>>
>'23.1'
>  
>
>>>>import numarray; numarray.__version__
>>>>        
>>>>
>'1.2.3'
>  
>
>>>>import numpy; numpy.__version__
>>>>        
>>>>
>'0.9.5.2021'
>  
>
>>>>t1 = Timer('a = max(b)', 'import Numeric as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>>>t2 = Timer('a = max(b)', 'import numarray as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>>>t3 = Timer('a = max(b)', 'import numpy as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>>>t1.timeit(10)
>>>>        
>>>>
>0.13748002052307129
>  
>
>>>>t2.timeit(10)
>>>>        
>>>>
>0.64549708366394043
>  
>
>>>>t3.timeit(10)
>>>>        
>>>>
>4.5433549880981445
>  
>
You shouldn't be using max(b)!!!

You should use b.max()  because  max is a Python function that is using 
the generic sequence interface to do the comparisons.  Very likely it is 
slower because in numpy you get an array scalar (which don't yet have 
their own math defined and so use the full ufunc machinery do to 
computations).

At some point the array scalars will have their own math defined and 
should be as fast as the Python scalars, so this particular "slowness" 
will go away.  But, you use maximum.reduce() instead anyway.

-Travis

[Numpy-discussion] another performance problem

[Gerard V. <ger...@gr...>]

Hi Travis,

max and min are really slow in numpy:

>>> from timeit import Timer
>>> import Numeric; Numeric.__version__
'23.1'
>>> import numarray; numarray.__version__
'1.2.3'
>>> import numpy; numpy.__version__
'0.9.5.2021'
>>> t1 = Timer('a = max(b)', 'import Numeric as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>> t2 = Timer('a = max(b)', 'import numarray as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>> t3 = Timer('a = max(b)', 'import numpy as N; b=N.cos(N.arange(0, 10**5, 1, N.Float))')
>>> t1.timeit(10)
0.13748002052307129
>>> t2.timeit(10)
0.64549708366394043
>>> t3.timeit(10)
4.5433549880981445
>>>

Gerard

[Numpy-discussion] Segfault with dodgy type descriptors

[Ed S. <sch...@ft...>]

Hi Travis,

Here's a segfault from the crazy input department:

>>> cs = [('A', 1), ('B', 2)]
>>> a = numpy.array(cs, dtype=(object,5))

Using the descriptor dtype=(object,x) for other values of x gives
varying results, but all of them should probably raise a ValueError.

On a related note, it seems that it's now possible again to use an array
as a dtype.  You changed this behaviour in November after I suggested
that this should be illegal.  (See
http://www.scipy.net/pipermail/scipy-dev/2005-November/004126.html).  My
argument then was that calling array(a, b.dtype) is clearer and safer
than using array(a, b), which can be confused with array((a,b)).  Is it
an oversight that arrays can again be used as dtypes?


-- Ed

[Numpy-discussion] performance problem

[Gerard V. <ger...@gr...>]

Hi Travis, thanks for your incredible quick fixes.

numpy is about 5 times slower than numarray (on my numarray-friendly bi-pentium):

>>> from timeit import Timer
>>> import numarray; print numarray.__version__
1.5.0
>>> import numpy; print numpy.__version__
0.9.5.2021
>>> t1 = Timer('a <<= 8', 'import numarray as NX; a = NX.ones(10**6, NX.UInt32)')
>>> t2 = Timer('a <<= 8', 'import numpy as NX; a = NX.ones(10**6, NX.UInt32)')
>>> t1.timeit(100)
0.21813011169433594
>>> t2.timeit(100)
1.1523458957672119
>>>

Numeric-23.1 is about as fast as numarray for inplace left shifts.

Gerard

Re: [Numpy-discussion] bug triggered by array.astype()

[Travis O. <oli...@ie...>]

Gerard Vermeulen wrote:

>>>>from numpy import *
>>>>core.__version__
>>>>        
>>>>
>'0.9.5.2019'
>  
>
>>>>v = linspace(0, 2*pi)
>>>>c = (255<<24)*cos(v) + (255<<16)*cos(v+2*pi/3) + (255<<8)*cos(v+4*pi/3) + 255
>>>>c.astype(UInt32)
>>>>        
>>>>
>Traceback (most recent call last):
>  File "<stdin>", line 1, in ?
>SystemError: Objects/longobject.c:257: bad argument to internal function
>  
>
This behaves like numarray now in numpy SVN where-in (255<<24) is 
converted to an int64 array during the multiplication step so that 
object arrays are not constructed.

-Travis

Re: [Numpy-discussion] bug triggered by array.astype()

[Travis O. <oli...@ie...>]

Gerard Vermeulen wrote:

>On Sun, 29 Jan 2006 01:50:57 -0700
>Travis Oliphant <oli...@ie...> wrote:
>
>[ .. ]
>
>  
>
>>In SVN, now a Overflow error is raised instead in this example because 
>>negative long integer objects are trying to be converted to unsigned 
>>integers (this is a Python error being raised).  We could check for this 
>>and instead do what the c-compiler would do in converting from signed to 
>>unsigned objects. 
>>
>>Is that preferrable?
>>
>>    
>>
>
>Well, the current SVN does now almost what Numeric does:
>
>  
>
>>>>from Numeric import *
>>>>v = arange(0, 10, 1, Float32)
>>>>c = (255<<24)*cos(v)+ (255<<16)*cos(v+2*pi/3) + (255<<8)*cos(v+4*pi/3) + 255
>>>>c.astype(UInt32)
>>>>        
>>>>
>Traceback (most recent call last):
>  File "<stdin>", line 1, in ?
>OverflowError: long int too large to convert to int
>  
>
>
>However, numarray does:
>
>  
>
>>>>from numarray import *
>>>>v = arange(0, 10, 1, Float32)
>>>>c = (255<<24)*cos(v)+ (255<<16)*cos(v+2*pi/3) + (255<<8)*cos(v+4*pi/3) + 255
>>>>c.astype(UInt32)
>>>>        
>>>>
>array([4269801984, 2294853120, 2504994432,   65861376, 1514949120,
>       1225005568, 4103764992, 3209541888, 3659448448,  398681088], type=UInt32)
>  
>
>
>This is more user friendly, but may cost more CPU cycles.
>  
>

The issue here is really that in Numeric (and numpy)  c is an object 
array while in numarray c is a float array.  Because Numeric and numpy 
support object arrays, long integers have been converted to objects. 

I think what should be changed is that long integers should only be 
converted to objects if they are bigger than the maximum integer-type on 
the platform.   That is what numarray allows and makes sense. 

In other words, the root of this difference is what array(255<<24) 
returns.  In Numeric and numpy it returns an object array.  In numarray 
it returns a UInt64 array.

-Travis

Re: [Numpy-discussion] bug triggered by array.astype()

[Gerard V. <ger...@gr...>]

On Sun, 29 Jan 2006 01:50:57 -0700
Travis Oliphant <oli...@ie...> wrote:

[ .. ]

> In SVN, now a Overflow error is raised instead in this example because 
> negative long integer objects are trying to be converted to unsigned 
> integers (this is a Python error being raised).  We could check for this 
> and instead do what the c-compiler would do in converting from signed to 
> unsigned objects. 
> 
> Is that preferrable?
>

Well, the current SVN does now almost what Numeric does:

>>> from Numeric import *
>>> v = arange(0, 10, 1, Float32)
>>> c = (255<<24)*cos(v)+ (255<<16)*cos(v+2*pi/3) + (255<<8)*cos(v+4*pi/3) + 255
>>> c.astype(UInt32)
Traceback (most recent call last):
  File "<stdin>", line 1, in ?
OverflowError: long int too large to convert to int
>>>

However, numarray does:

>>> from numarray import *
>>> v = arange(0, 10, 1, Float32)
>>> c = (255<<24)*cos(v)+ (255<<16)*cos(v+2*pi/3) + (255<<8)*cos(v+4*pi/3) + 255
>>> c.astype(UInt32)
array([4269801984, 2294853120, 2504994432,   65861376, 1514949120,
       1225005568, 4103764992, 3209541888, 3659448448,  398681088], type=UInt32)
>>>

This is more user friendly, but may cost more CPU cycles.

Gerard