[Numpy-discussion] missing from contributor list?

2016-11-02 Thread Sturla Molden
Why am I missing from the contributor hist here?

https://github.com/numpy/numpy/blob/master/numpy/_build_utils/src/apple_sgemv_fix.c


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Accelerate or OpenBLAS for numpy / scipy wheels?

2016-07-01 Thread Sturla Molden
Ralf Gommers  wrote:

> Thanks Sturla, interesting details as always. You didn't state your
> preference by the way, do you have one?

I use Accelerate because it is the easier for me to use when building
SciPy. But that is from a developer's perspective.

As you know, Accelerate breaks a common (ab)use of multiprocessing on POSIX
systems. While the bug is strictly speaking in multiprocessing (but
partially fixed in Python 3.4 and later), it is still a nasty surprise to
many users. E.g. a call to np.dot never returns, and there is no error
message indicating why. That speaks against using it in the wheels.

Accelerate, like MKL and FFTW, has nifty FFTs. If we start to use MKL and
Accelerate for numpy.fft (which I sometimes have fantacies about), that
would shift the balance the other way, in favour of Accelerate.

Speed wise Accelerate wins for things like dot product of two vectors or
multiplication of a vector and a matrix. For general matrix multiplication
the performance is about the same, except when matrices are very small and
Accelerate can benefit from the tiny GCD overhead. But then the Python
overhead probably dominates, so they are going to be about equal anyway. 

I am going to vote ± 0. I am really not sure which will be the better for
the binary wheels. They seem about equal to me right now. There are pros
and cons with either.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Accelerate or OpenBLAS for numpy / scipy wheels?

2016-07-01 Thread Sturla Molden

On 29/06/16 21:55, Nathaniel Smith wrote:


Speed is important, but it's far from the only consideration, especially
since differences between the top tier libraries are usually rather
small.


It is not even the more important consideration. I would say that 
correctness matters most. Everything else comes second.



Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Accelerate or OpenBLAS for numpy / scipy wheels?

2016-07-01 Thread Sturla Molden

On 29/06/16 21:55, Nathaniel Smith wrote:


Accelerate is closed, so when we hit bugs then there's often nothing we
can do except file a bug with apple and hope that it gets fixed within a
year or two. This isn't hypothetical -- we've hit cases where accelerate
gave wrong answers. Numpy actually carries some scary code right now to
work around one of these bugs by monkeypatching (!) accelerate using
dynamic linker trickiness. And, of course, there's the thing where
accelerate totally breaks multiprocessing.


Yes, those are the cons.



Apple has said that they don't consider this a bug.


Theoretically they are right, but from a practical perspective...


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Accelerate or OpenBLAS for numpy / scipy wheels?

2016-06-29 Thread Sturla Molden
Ralf Gommers  wrote:

> For most routines performance seems to be comparable, and both are much
> better than ATLAS. When there's a significant difference, I have the
> impression that OpenBLAS is more often the slower one (example:
>  href="https://github.com/xianyi/OpenBLAS/issues/533";>https://github.com/xianyi/OpenBLAS/issues/533).

Accelerate is in general better optimized for level-1 and level-2 BLAS than
OpenBLAS. There are two reasons for this: 

First, OpenBLAS does not use AVX for these kernels, but Accelerate does.
This is the more important difference. It seems the OpenBLAS devs are now
working on this.

Second, the thread pool in OpenBLAS is not as scalable on small tasks as
the "Grand Central Dispatch" (GCD) used by Accelerate. The GCD thread-pool
used by Accelerate is actually quite unique in having a very tiny overhead:
It takes only 16 extra opcodes (IIRC) for running a task on the global
parallel queue instead of the current thread. (Even if my memory is not
perfect and it is not exactly 16 opcodes, it is within that order of
magnitude.) GCD can do this because the global queues and threadpool is
actually built into the kernel of the OS. On the other hand, OpenBLAS and
MKL depends on thread pools managed in userspace, for which the scheduler
in the OS have no special knowledge. When you need fine-grained parallelism
and synchronization, there is nothing like GCD. Even a user-space spinlock
will have bigger overhead than a sequential queue in GCD. With a userspace
threadpool all threads are scheduled on a round robin basis, but with GCD
the scheduler has special knowledge about the tasks put on the queues, and
executes them as fast as possible. Accelerate therefore has an unique
advantage when running level-1 and 2 BLAS routines, with which OpenBLAS or
MKL probably never can properly compete. Programming with GCD can actually
often be counter-intuitive to someone used to deal with OpenMP, MPI or
pthreads. For example it is often better to enqueue a lot of small tasks
instead of splitting up the computation into large chunks of work. When
parallelising a tight loop, a chunk size of 1 can be great on GCD but is
likely to be horrible on OpenMP and anything else that has userspace
threads.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Integers to integer powers, let's make a decision

2016-06-05 Thread Sturla Molden
Charles R Harris  wrote:

>1. Integers to negative integer powers raise an error.
>2. Integers to integer powers always results in floats.

2

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Changing FFT cache to a bounded LRU cache

2016-05-31 Thread Sturla Molden
Joseph Martinot-Lagarde  wrote:
 
> The problem with FFTW is that its license is more restrictive (GPL), and
> because of this may not be suitable everywhere numpy.fft is.

A lot of us use NumPy linked with MKL or Accelerate, both of which have
some really nifty FFTs. And the license issue is hardly any worse than
linking with them for BLAS and LAPACK, which we do anyway. We could extend
numpy.fft to use MKL or Accelerate when they are available.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Changing FFT cache to a bounded LRU cache

2016-05-31 Thread Sturla Molden
Lion Krischer  wrote:

> I added a slightly more comprehensive benchmark to the PR. Please have a
> look. It tests the total time for 100 FFTs with and without cache. It is
> over 30 percent faster with cache which it totally worth it in my
> opinion as repeated FFTs of the same size are a very common use case.

All the calls to trancendental functions are stored in the cache. Without a
cache, we get excessive calls to sin(x) and cos(x) whenever FFTs of the
same size are repeated. This can indeed matter at lot.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-24 Thread Sturla Molden
Antoine Pitrou  wrote:

> When writing C code to interact with buffer-providing objects, you
> usually don't bother with memoryviews at all.  You just use a Py_buffer
> structure.

I was taking about "typed memoryviews" which is a Cython abstraction for a
Py_buffer struct. I was not taking about Python memoryviews, which is
something else. When writing Cython code that interacts with a buffer we
usually use typed memoryviews, not Py_buffer structs directly.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-17 Thread Sturla Molden
Matěj Týč  wrote:

> Does it mean
> that if you pass the numpy array to the child process using Queue, no
> significant amount of data will flow through it? 

This is what my shared memory arrayes do.

> Or I shouldn't pass it
> using Queue at all and just rely on inheritance? 

This is what David Baddeley's shared memory arrays do.

> Finally, I assume that
> passing it as an argument to the Process class is the worst option,
> because it will be pickled and unpickled.
 
My shared memory arrays only pickles the metadata, and can be used in this
way.

> Or maybe you refer to modules s.a. joblib that use this functionality
> and expose only a nice interface?

Joblib creates "share memory" by memory mapping a temporary file, which is
back by RAM on Libux (tempfs). It is backed by a physical file on disk on
Mac and Windows. In this resepect, joblib is much better on Linux than Mac
or Windows.


> And finally, cow means that returning large arrays still involves data
> moving between processes, whereas the shm approach has the workaround
> that you can preallocate the result array by the parent process, where
> the worker process can write to.

My shared memory arrays need no workaround dor this. They also allow shared
memory arrays to be returned to the parent process. No preallocation is
needed.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: numpy.random.random_seed

2016-05-17 Thread Sturla Molden
Stephan Hoyer  wrote:
> I have recently encountered several use cases for randomly generate random
> number seeds:
> 
> 1. When writing a library of stochastic functions that take a seed as an
> input argument, and some of these functions call multiple other such
> stochastic functions. Dask is one such example [1].
> 
> 2. When a library needs to produce results that are reproducible after
> calling numpy.random.seed, but that do not want to use the functions in
> numpy.random directly. This came up recently in a pandas pull request [2],
> because we want to allow using RandomState objects as an alternative to
> global state in numpy.random. A major advantage of this approach is that it
> provides an obvious alternative to reusing the private numpy.random._mtrand
> [3].


What about making numpy.random a finite state machine, and keeping a stack
of RandomState seeds? That is, something similar to what OpenGL does for
its matrices? Then we get two functions, numpy.random.push_seed and
numpy.random.pop_seed.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-17 Thread Sturla Molden
Matěj Týč  wrote:

>  - Parallel processing of HUGE data, and

This is mainly a Windows problem, as copy-on-write fork() will solve this
on any other platform. I am more in favor of asking  Microsoft to fix their
broken OS. 

Also observe that the usefulness of shared memory is very limited on
Windows, as we in practice never get the same base address in a spawned
process. This prevents sharing data structures with pointers and Python
objects. Anything more complex than an array cannot be shared.

What this means is that shared memory is seldom useful for sharing huge
data, even on Windows. It is only useful for this on Unix/Linux, where base
addresses can stay they same. But on non-Windows platforms, the COW will in
99.99% of the cases be sufficient, thus make shared memory superfluous
anyway. We don't need shared memory to scatter large data on Linux, only
fork.

As I see it. shared memory is mostly useful as a means to construct an
inter-process communication (IPC) protocol. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-13 Thread Sturla Molden
Feng Yu  wrote:

> Also, did you checkout http://zeromq.org/blog:zero-copy ?
> ZeroMQ is a dependency of Jupyter, so it is quite available.

ZeroMQ is great, but it lacks some crucial features. In particular it does
not support IPC on Windows. Ideally one should e.g. use Unix doman sockets
on Linux and named pipes on Windows. Most MPI implementations seems to
prefer shared memory over these mechanisms, though. Also I am not sure
about ZeroMQ and asynch i/o. I would e.g. like to use IOCP on Windows, GCD
on Mac, and a threadpool plus epoll on Linux. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-13 Thread Sturla Molden
Niki Spahiev  wrote:

> Apparently next Win10 will have fork as part of bash integration.

It is Interix/SUA rebranded "Subsystem for Linux". It remains to be seen
how long it will stay this time. Also a Python built for this subsystem
will not run on the Win32 subsystem, so there is no graphics. Also it will
not be installed by default, just like SUA.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-12 Thread Sturla Molden
Feng Yu  wrote:
 
> In most (half?) situations the result can be directly write back via
> preallocated shared array before works are spawned. Then there is no
> need to pass data back with named segments.

You can work around it in various ways, this being one of them.

Personally I prefer a parallel programming style with queues – either to
scatter arrays to workers and collecting arrays from workers, or to chain
workers together in a pipeline (without using coroutines). But exactly how
you program is a matter of taste. I want to make it as inexpensive as
possible to pass a NumPy array through a queue. If anyone else wants to
help improve parallel programming with NumPy using a different paradigm,
that is fine too. I just wanted to clarify why I stopped working on shared
memory arrays.

(As for the implementation, I am also experimenting with platform dependent
asynchronous I/O (IOCP, GCD or kqueue, epoll) to pass NumPy arrays though a
queue as inexpensively and scalably as possible. And no, there is no public
repo, as I like to experiment with my pet project undisturbed before I let
it out in the wild.)


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-12 Thread Sturla Molden
Niki Spahiev  wrote:

> Apparently next Win10 will have fork as part of bash integration.

That would be great. The lack of fork on Windows is very annoying.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-12 Thread Sturla Molden
Antoine Pitrou  wrote:

> Can you define "expensive"?

Slow enough to cause complaints on the Cython mailing list.

> You're assuming this is the cost of "buffer acquisition", while most
> likely it's the cost of creating the memoryview object itself.

Constructing a typed memoryview from a typed memoryview or a slice is fast.
Numerical code doing this intensively is still within 80-90% of the speed
of plain c code using pointer arithmetics.

 
> Buffer acquisition itself only calls a single C callback and uses a
> stack-allocated C structure. It shouldn't be "expensive".

I don't know the reason, only that buffer acquisition from NumPy arrays
with typed memoryviews is very expensive compared to assigning a typed
memoryview to another or slicing a typed memoryview.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Allan Haldane  wrote:

> You probably already know this, but I just wanted to note that the
> mpi4py module has worked around pickle too. They discuss how they
> efficiently transfer numpy arrays in mpi messages here:
> http://pythonhosted.org/mpi4py/usrman/overview.html#communicating-python-objects-and-array-data

Unless I am mistaken, they use the PEP 3118 buffer interface to support
NumPy as well as a number of other Python objects. However, this protocol
makes buffer aquisition an expensive operation. You can see this in Cython
if you use typed memory views. Assigning a NumPy array to a typed
memoryview (i,e, buffer acqisition) is slow. They are correct that avoiding
pickle means we save some memory. It also avoids creating and destroying
temporary Python objects, and associated reference counting. However,
because of the expensive buffer acquisition, I am not sure how much faster
their apporach will be. I prefer to use the NumPy C API, and bypass any
unneccesary overhead. The idea is to make IPC of NumPy arrays fast, and
then we cannot have an expensive buffer acquisition in there.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Feng Yu  wrote:

> 1. If we are talking about shared memory and copy-on-write
> inheritance, then we are using 'fork'. 

Not available on Windows. On Unix it only allows one-way communication,
from parent to child.

> 2. Picking of inherited shared memory array can be done minimally by
> just picking the array_interface and the pointer address. It is
> because the child process and the parent share the same address space
> layout, guarenteed by the fork call.

Again, not everyone uses Unix. 

And on Unix it is not trival to pass data back from the child process. I
solved that problem with Sys V IPC (pickling the name of the segment).

> 6. If we are to define a set of operations I would recommend take a
> look at OpenMP as a reference -- It has been out there for decades and
> used widely. An equiavlant to the 'omp parallel for' construct in
> Python will be a very good starting point and immediately useful.

If you are on Unix, you can just use a context manager. Call os.fork in
__enter__ and os.waitpid in __exit__. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Joe Kington  wrote:

> You're far better off just
> communicating between processes as opposed to using shared memory.

Yes.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Benjamin Root  wrote:

> Oftentimes, if one needs to share numpy arrays for multiprocessing, I would
> imagine that it is because the array is huge, right? 

That is a case for shared memory, but what. i was taking about is more
common than this. In order for processes to cooperate, they must
communicate. So we need a way to pass around NumPy arrays quickly.
Sometimes we want to use shared memory because of the size of the data, but
more often it is just used as a form of inexpensive IPC.

> So, the pickling
> approach would copy that array for each process, which defeats the purpose,
> right?

I am not sure what you mean. When I made shared memory arrays I used named
segments, and made sure only the name of the segments were pickled, not the
contents of the buffers.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Allan Haldane  wrote:

> That's interesting. I've also used multiprocessing with numpy and didn't
> realize that. Is this true in python3 too?

I am not sure. As you have noticed, pickle is faster by to orders of
magnitude on Python 3. But several microseconds is also a lot, particularly
if we are going to do this often during a computation.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
Elliot Hallmark  wrote:
> Strula, this sounds brilliant!  To be clear, you're talking about
> serializing the numpy array and reconstructing it in a way that's faster
> than pickle?

Yes. We know the binary format of NumPy arrays. We don't need to invoke the
machinery of pickle to serialize an array and write the bytes to some IPC
mechanism (pipe, tcp socket, unix socket, shared memory). The choise of IPC
mechanism might not even be relevant, and could even be deferred to a
library like ZeroMQ. The point is that if multiple peocesses are to
cooperate efficiently, we need a way to let them communicate NumPy arrays
quickly. That is where using multiprocessing hurts today, and shared memory
does not help here.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Numpy arrays shareable among related processes (PR #7533)

2016-05-11 Thread Sturla Molden
I did some work on this some years ago. I have more or less concluded that
it was a waste of effort. But first let me explain what the suggested
approach do not work. As it uses memory mapping to create shared memory
(i.e. shared segments are not named), they must be created ahead of
spawning processes. But if you really want this to work smoothly, you want
named shared memory (Sys V IPC or posix shm_open), so that shared arrays
can be created in the spawned processes and passed back.

Now for the reason I don't care about shared memory arrays anymore, and
what I am currently working on instead:

1. I have come across very few cases where threaded code cannot be used in
numerical computing. In fact, multithreading nearly always happens in the
code where I write pure C or Fortran anyway. Most often it happens in
library code that are already multithreaded (Intel MKL, Apple Accelerate
Framework, OpenBLAS, etc.), which means using it requires no extra effort
from my side. A multithreaded LAPACK library is not less multithreaded if I
call it from Python.

2. Getting shared memory right can be difficult because of hierarchical
memory and false sharing. You might not see it if you only have a multicore
CPU with a shared cache. But your code might not scale up on computers with
more than one physical processor. False sharing acts like the GIL, except
it happens in hardware and affects your C code invisibly without any
explicit locking you can pinpoint. This is also why MPI code tends to scale
much better than OpenMP code. If nothing is shared there will be no false
sharing.

3. Raw C level IPC is cheap – very, very cheap. Even if you use pipes or
sockets instead of shared memory it is cheap. There are very few cases
where the IPC tends to be a bottleneck. 

4. The reason IPC appears expensive with NumPy is because multiprocessing
pickles the arrays. It is pickle that is slow, not the IPC. Some would say
that the pickle overhead is an integral part of the IPC ovearhead, but i
will argue that it is not. The slowness of pickle is a separate problem
alltogether.

5. Share memory does not improve on the pickle overhead because also NumPy
arrays with shared memory must be pickled. Multiprocessing can bypass
pickling the RawArray object, but the rest of the NumPy array is pickled.
Using shared memory arrays have no speed advantage over normal NumPy arrays
when we use multiprocessing.

6. It is much easier to write concurrent code that uses queues for message
passing than anything else. That is why using a Queue object has been the
popular Pythonic approach to both multitreading and multiprocessing. I
would like this to continue.

I am therefore focusing my effort on the multiprocessing.Queue object. If
you understand the six points I listed you will see where this is going:
What we really need is a specialized queue that has knowledge about NumPy
arrays and can bypass pickle. I am therefore focusing my efforts on
creating a NumPy aware queue object.

We are not doing the users a favor by encouraging the use of shared memory
arrays. They help with nothing.


Sturla Molden



Matěj  Týč  wrote:
> Dear Numpy developers,
> I propose a pull request https://github.com/numpy/numpy/pull/7533 that
> features numpy arrays that can be shared among processes (with some
> effort).
> 
> Why:
> In CPython, multiprocessing is the only way of how to exploit
> multi-core CPUs if your parallel code can't avoid creating Python
> objects. In that case, CPython's GIL makes threads unusable. However,
> unlike with threading, sharing data among processes is something that
> is non-trivial and platform-dependent.
> 
> Although numpy (and certainly some other packages) implement some
> operations in a way that GIL is not a concern, consider another case:
> You have a large amount of data in a form of a numpy array and you
> want to pass it to a function of an arbitrary Python module that also
> expects numpy array (e.g. list of vertices coordinates as an input and
> array of the corresponding polygon as an output). Here, it is clear
> GIL is an issue you and since you want a numpy array on both ends, now
> you would have to copy your numpy array to a multiprocessing.Array (to
> pass the data) and then to convert it back to ndarray in the worker
> process.
> This contribution would streamline it a bit - you would create an
> array as you are used to, pass it to the subprocess as you would do
> with the multiprocessing.Array, and the process can work with a numpy
> array right away.
> 
> How:
> The idea is to create a numpy array in a buffer that can be shared
> among processes. Python has support for this in its standard library,
> so the current solution creates a multiprocessing.Array and then
> passes it as the "buffer" to the ndarray.__new__. That would be it on
> Unixes, but on Windows, there has to be a a custom pick

Re: [Numpy-discussion] Make as_strided result writeonly

2016-01-25 Thread Sturla Molden

On 25/01/16 18:06, Sebastian Berg wrote:


That said, I guess I could agree with you in the regard that there are
so many *other* awful ways to use as_strided, that maybe it really is
just so bad, that improving one thing doesn't actually help anyway ;).


That is roughly my position on this, yes. :)

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Make as_strided result writeonly

2016-01-25 Thread Sturla Molden

On 23/01/16 22:25, Sebastian Berg wrote:


Do you agree with this, or would it be a major inconvenience?


I think any user of as_strided should be considered a power user. This 
is an inherently dangerous function, that can easily segfault the 
process. Anyone who uses as_strided should be assumed to have taken all 
precautions.


-1

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Support of @= in numpy?

2015-12-27 Thread Sturla Molden
Charles R Harris  wrote:

> In any case, we support the `@` operator in 1.10, but not the `@=`
> operator. The `@=` operator is tricky to have true inplace matrix
> multiplication, as not only are elements on the left overwritten, but the
> dimensions need to be preserved.

As long as we use BLAS, we can never have true inplace matrix
multiplication because Fortran prohibits pointer aliasing. We therefore
need to make a temporary copy before we call BLAS.

But as for preserving dimensions: This should be allowed if the array is
square. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] performance solving system of equations in numpy and MATLAB

2015-12-17 Thread Sturla Molden

On 16/12/15 20:47, Derek Homeier wrote:


Getting around 30 s wall time here on a not so recent 4-core iMac, so that 
would seem to fit
(iirc Accelerate should actually largely be using the same machine code as MKL).


Yes, the same kernels, but not the same threadpool. Accelerate uses the 
GCD, MKL uses Intel TBB and Intel OpenMP (both of them). GCD scales 
better than TBB, even in Intel's own benchmarks. However, GCD uses a 
kernel threadpool (accesible via kqueue) which is not fork-safe, whereas 
MKL's threadpool is fork-safe (but will leak memory on fork).


Sturla




___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] performance solving system of equations in numpy and MATLAB

2015-12-17 Thread Sturla Molden

On 17/12/15 12:06, Francesc Alted wrote:


Pretty good.  I did not know that OpenBLAS was so close in performance
to MKL.


MKL, OpenBLAS and Accelerate are very close in performance, except for 
level-1 BLAS where Accelerate and MKL are better than OpenBLAS.


MKL requires the number of threads to be a multiple of four to achieve 
good performance, OpenBLAS and Accelerate do not. It e.g. matters if you 
have an online data acquisition and DSP system and want to dedicate one 
processor to take care of i/o tasks. In this case OpenBLAS and 
Accelerate are likely to perform better than MKL.



Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Fast vectorized arithmetic with ~32 significant digits under Numpy

2015-12-12 Thread Sturla Molden
"Thomas Baruchel"  wrote:

> While this is obviously the most relevant answer for many users because
> it will allow them to use Numpy arrays exactly
> as they would have used them with native types, the wrong thing is that
> from some point of view "true" vectorization
> will be lost.

What does "true vectorization" mean anyway?


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Memory mapping and NPZ files

2015-12-11 Thread Sturla Molden
Mathieu Dubois  wrote:

> The point is precisely that, you can't do memory mapping with Npz files 
> (while it works with Npy files).

The operating system can memory map any file. But as npz-files are
compressed, you will need to uncompress the contents in your memory mapping
to make sense of it. I would suggest you use PyTables instead of npz-files.
It allows on the fly compression and uncompression (via blosc) and will
probably do what you want.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Memory mapping and NPZ files

2015-12-10 Thread Sturla Molden
Mathieu Dubois  wrote:

> Does it make sense?

No. Memory mapping should just memory map, not do all sorts of crap.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] When to stop supporting Python 2.6?

2015-12-07 Thread Sturla Molden
Charles R Harris  wrote:

> As a strawman proposal, how about dropping moving to 2.7 and 3.4 minimum
> supported version next fall, say around numpy 1.12 or 1.13 depending on how
> the releases go.
> 
> I would like to here from the scipy folks first.

Personally I would be in favor of this, because 2.7 and 3.4 are the minimum
versions anyone should consider to use. However, for SciPy which heavily
depends on Python code, the real improvement will be when we can bump the
minimum Python version to 3.5 and write x @ y instead of dot(x,y). I am not
sure of bumping the minimum version to 3.4 before that is worth it or not.
But certainly dropping 2.6 might be a good thing already now, so we can
start to use bytes, bytearray, memoryview, etc.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Where is Jaime?

2015-12-07 Thread Sturla Molden
Charles R Harris  wrote:

> The cash economy is nothing to sniff at ;) It is big in NYC and other
> places with high taxes and bureaucratic meddling. Cash was one of the great
> inventions.

Yeah, there is a Sicilian New Yorker called "Gambino" who has been
advertising "protection from ISIS" in European newspapers lately. From what
I read his father was big at selling protection for cash, and now he is
taking up his father's business and selling protection from ISIS. To prove
his value, he claimed ISIS is so afraid of his organisation that Sicily is
a place they never dare visit. Presumably Gambino's business model depends
on a cash based economy, or at least it did.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] f2py, numpy.distutils and multiple Fortran source files

2015-12-04 Thread Sturla Molden

On 03/12/15 22:07, David Verelst wrote:


Can this workflow be incorporated into |setuptools|/|numpy.distutils|?
Something along the lines as:


Take a look at what SciPy does.

https://github.com/scipy/scipy/blob/81c096001974f0b5efe29ec83b54f725cc681540/scipy/fftpack/setup.py

Multiple Fortran files are compiled into a static library using 
"add_library", which is subsequently linked to the extension module.



Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] future of f2py and Fortran90+

2015-12-04 Thread Sturla Molden

On 03/12/15 22:38, Eric Firing wrote:


Right, but for each function that requires writing two wrappers, one in
Fortran and a second one in cython.


Yes, you need two wrappers for each function, one in Cython and one in 
Fortran 2003. That is what fwrap is supposed to automate, but it has 
been abandonware for quite a while.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Compilation problems npy_float64

2015-11-07 Thread Sturla Molden
Johan  wrote:

> Hello, I searched the forum, but couldn't find a post related to my 
> problem.  I am installing scipy via pip in cygwin environment

I think I introduced this error when moving a global variable from the
Cython module to a C++ module. The name collision with math.h was silent on
Linux, Mac, and Windows (MinGW and MSVC) -- or not even present --, and
thus went under the radar. But it eventually showed up on SunOS, and now
also on Cygwin. :-(

My apologies.

Anyhow, it should be gone now. Try SciPy master.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] isfortran compatibility in numpy 1.10.

2015-11-01 Thread Sturla Molden
Charles R Harris  wrote:

>1. Return `a.flags.f_contiguous`. This differs for 1-D arrays, but is
>most consistent with the name isfortran.

If the idea is to determine if an array can safely be passed to Fortran,
this is the correct one.

>2. Return `a.flags.f_contiguous and a.ndim > 1`, which would be backward
>compatible.

This one is just wrong.

A compromize might be to raise an exception in the case of a.ndim<2. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-10-05 Thread Sturla Molden

On 02/10/15 13:05, Daπid wrote:


Have you tried asking Python-dev for help with this? Hopefully they
would have some weight there.


It seems both GCC dev and Apple (for GCD and Accelerate) has taken a 
similar stance on this. There is tiny set of functions the POSIX 
standard demands should work on both sides of a fork without exec, but 
OpenMP, GCD, BLAS or LAPAPCK are not included. As long as there is no 
bug, it is hard to convince them to follow Intel and allow fork-based 
multiprocessing.


As it stands now, using Intel compilers and MKL is the only way to make 
this work, but Intel's development tools are not freeware.


Sturla



___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Let's move forward with the current governance document.

2015-10-05 Thread Sturla Molden
Nathaniel Smith  wrote:

> Are you planning to go around vetoing things 

I don't consider myself qualified.

> for ridiculous reasons and causing havoc?

That would be unpolite.

> And if not, then who is it that you're worried about?

I am not sure :)

I just envisioned a Roman patron shouting veto or a US senator
filibustering. Expulsion would be the appropriate recation, yes :-)


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Let's move forward with the current governance document.

2015-10-05 Thread Sturla Molden
Nathaniel Smith  wrote:

> Thanks Chuck! It looks like it's just wording tweaks / clarifications
> at this point, so nothing we need to discuss in detail on the list. If
> anyone wants to watch the sausage being made, then the link is above
> :-), and we'll continue the discussion in the PR unless anything
> substantive comes up.

Anyone has a veto? That reminds me of something that happened in the senate
of Rome; they only had a small number of vetoers, sometimes only one or
two, and even that caused havoc. I think it should be better clarified how
much contribution is needed before someone can be considered to have veto
rights. It would e.g. be ridiculous if I were to begin and veto stuff, as
my contributions are minute... OMG.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-10-02 Thread Sturla Molden
Sturla Molden  wrote:

> Cython actually requires that there is a shared address space, and it
> invokes something that strictly speaking has undefined behavior under the
> OpenMP standard. So thus, a prange block in Cython is expected to work
> correctly on a laptop with a multicore processor, but it is not expected to
> work correctly on a cluster.

OpenMP does not guarrantee that dereferencing a pointer in a parallel block
will dereference the same object across all processors. It only guarrantees
that the value of a shared object can be synchronized. There are many who
use OpenMP and think only in terms of threads that do this incorrectly.
Cython is actually among those.

S.M.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-10-02 Thread Sturla Molden
Sturla Molden  wrote:

> OpenMP has a flush pragma for synchronizing shared variables. This means
> that OpenMP is not restricted to shared memory hardware. A "pragma omp
> flush" can just as well invoke some IPC mechanism, even network
> communication. 

By the way, while this is the case for C and Fortran, it is certainly not
the case for Cython. In a Cython prange block, a shared variable is
accessed by dereferencing its address. This requires shared memory. Pure
OpenMP in C does not, because shared variables are not accessed through
pointers, but are rather normal variables that are synchronized with a
pragma.

Cython actually requires that there is a shared address space, and it
invokes something that strictly speaking has undefined behavior under the
OpenMP standard. So thus, a prange block in Cython is expected to work
correctly on a laptop with a multicore processor, but it is not expected to
work correctly on a cluster.

IIRC, Intel's cluster OpenMP is based on MPI, which means the compiler will
internally translate code with OpenMP pragmas into equivalent code that
calls MPI functions. A program written for OpenMP can then run on any
cluster that provides an MPI implementation.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-10-02 Thread Sturla Molden
Juha Jeronen  wrote:

> Mm. I've quite often run MPI locally (it's nice for multicore scientific 
> computing on Python), but I had no idea that OpenMP had cluster 
> implementations. Thanks for the tip.

Intel has been selling one, I think there are others too. 

OpenMP has a flush pragma for synchronizing shared variables. This means
that OpenMP is not restricted to shared memory hardware. A "pragma omp
flush" can just as well invoke some IPC mechanism, even network
communication. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 01/10/15 02:32, Juha Jeronen wrote:



Sounds good. Out of curiosity, are there any standard fork-safe
threadpools, or would this imply rolling our own?


You have to roll your own.

Basically use pthreads_atfork to register a callback that shuts down the 
threadpool before a fork and another callback that restarts it. Python's 
threading module does not expose the pthreads_atfork function, so you 
must call it from Cython.


I believe there is also a tiny atfork module in PyPI.


So maybe it would be better, at least at first, to make a pure-Cython
version with no attempt at multithreading?


I would start by making a pure Cython version that works correctly. The 
next step would be to ensure that it releases the GIL. After that you 
can worry about parallel processing, or just tell the user to use 
threads or joblib.



Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 01/10/15 02:20, Juha Jeronen wrote:


Then again, the matter is further complicated by considering codes that
run on a single machine, versus codes that run on a cluster.Threads
being local to each node in a cluster,


You can run MPI programs on a single machine and you get OpenMP 
implementations for clusters. Just pick an API and stick with it.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 30/09/15 18:20, Chris Barker wrote:


python threads with nogil?


This is often the easiest, particularly if we construct a fork-safe 
threadpool.



Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 30/09/15 15:57, Nathan Goldbaum wrote:


Basically, just try to compile a simple OpenMP test program in a
subprocess. If that succeeds, then great, we can add -fopenmp as a
compilation flag. If not, don't do that.


Unless you use GCC on Linux, it will be more complex than that. I.e. do 
you need -fopenmp, -openmp or /openmp? And which libraries do you need 
to link, if any? Link GOMP? Link some other OpenMP runtime libraries? 
Link pthreads? There is three different pthread implementations to 
choose between on Windows, depending on GCC version (MinGW, MinGW-w64,
and Cygwin need different pthread libs for OpenMP). It gets messy. The 
only clean way is to add the correct flags and libraries into the 
compiler classes in numpy.distutils, and then let distutils and bento 
query those.



STurla






___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 30/09/15 18:20, Chris Barker wrote:


We'd need a run-time check.


We need to amend the compiler classes in numpy.distutils with OpenMP 
relevant information (compiler flags and libraries).


The OpenMP support libraries must be statically linked.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 30/09/15 18:20, Nathaniel Smith wrote:


This is incorrect -- the only common implementation of BLAS that uses
*OpenMP* threads is OpenBLAS,


MKL and ACML also use OpenMP.

Sturla




___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

2015-09-30 Thread Sturla Molden

On 30/09/15 11:27, Daπid wrote:


Is there a nice way to ship both versions? After all, most
implementations of BLAS and friends do spawn OpenMP threads, so I don't
think it would be outrageous to take advantage of it in more places;


Some do, others don't.

ACML uses OpenMP.

GotoBLAS uses OpenMP.

Intel MKL uses Intel TBB and OpenMP (both of them).

OpenBLAS will by default use an internal threadpool. It can be 
configured to use OpenMP instead.


ATLAS uses its own threadpool.

Apple Accelerate Framework uses a kernel thread-pool called the Grand 
Central Dispatch (GCD). A library called libdispatch uses kqueue to 
access the GCD.




There are two principal problems with using OpenMP in NumPy:

One is that the GNU OpenMP threadpool is not fork-safe, and can break 
multiprocessing on some platforms (e.g. when using Python 2.7 on Linux). 
Anything that uses GCD has this nasty effect on Apple and FreeBSD as 
well. Note that the problem is actually in multiprocessing. It is not 
present on Windows (there is no fork system call) and it is avoidable 
even on Linux with Python 3.4 or later. Also the default builds of NumPy 
and SciPy on MacOSX uses GCD by default.


The other problem is that NumPy with its iterator objects and gufuncs is 
not particularly fit for multithreading. There will be a lot of 
contention for the iterator object. Without a major redesign, one thread 
would do useful work while most of the others would be busy-waiting on a 
spinlock and fighting for iterator access. Not to mention that the 
iterator object would be false shared between the processors, which 
would trash the performance if you have more than one CPU, even when 
compared to using a single thread. This means that for multithreading to 
be useful in NumPy, the loops will have to be redesigned so that the 
work sharing between the threads is taken care of ahead of creating the 
iterator, i.e. allowing us to use one iterator per thread. Each iterator 
would then iterate over a slice of the original array. This is ok, but 
we cannot simply do this by adding OpenMP pragmas in the C code.


Given the two problems mentioned above, it would likely be better to use 
a fork-safe threadpool instead of OpenMP.



Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Governance model request

2015-09-22 Thread Sturla Molden

On 22/09/15 14:31, Perry Greenfield wrote:


I’ve also stayed out of this until now. I’m surprised and disheartened at the 
amount of suspicion and distrust directed towards Travis.


I have no idea where this distrust comes from. Nobody has invested so 
much of their time in NumPy. Without Travis there would not even be a NumPy.




So much of this distrust appears based on what Continuum might do rather than 
what the actual record is.


Being CEO of Continuum should not disqualify him in any way.



I agree with Travis that virtually everyone has their own interests.


Even Guido and Linus have employers. Should we distrust Guido as Python 
BDFL because some day Dropbox will be evil? It is just nonsense.



Sturla










___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Governance model request

2015-09-22 Thread Sturla Molden

On 20/09/15 20:20, Travis Oliphant wrote:


1 - define a BDFL for the council.  I would nominate chuck Harris

2 - limit the council to 3 people.  I would nominate chuck, nathaniel,
and pauli.

3 - add me as a permanent member of the steering council.



I have stayed out of this governance debate until now.

Personally I would prefer if you were BDFL.



Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] [Python-ideas] Should our default random number generator be secure?

2015-09-20 Thread Sturla Molden

On 20/09/15 21:48, Sturla Molden wrote:


This is where a small subset of C++ would be handy. Making an uint128_t
class with overloaded operators is a nobrainer. :-)


Meh... The C++ version of PCG already has this.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] [Python-ideas] Should our default random number generator be secure?

2015-09-20 Thread Sturla Molden

On 19/09/15 18:06, Robert Kern wrote:


That said, we'd
probably end up doing a significant amount of rewriting so that we will
have a C implementation of software-uint128 arithmetic.


This is where a small subset of C++ would be handy. Making an uint128_t 
class with overloaded operators is a nobrainer. :-)


Not that we will ever use C++ in NumPy, but we do have a fair amount of 
C++ in SciPy.



Sturla



___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] [Python-ideas] Should our default random number generator be secure?

2015-09-18 Thread Sturla Molden

On 14/09/15 10:26, Robert Kern wrote:


I want fast, multiple independent streams on my
current hardware first, and PCG gives that to me.


DCMT is good for that as well.

It should be possible to implement a pluggable design of NumPy's mtrand. 
Basically call a function pointer instead of rk_double. That way any 
(P)RNG can be plugged in. Then we could e.g. allow a Python callable, a 
Cython class or an f2py function as callback. We could ship multiple 
PRNGs with NumPy, and we could allow users to supply their own.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] [Python-ideas] Should our default random number generator be secure?

2015-09-18 Thread Sturla Molden

On 14/09/15 10:34, Antoine Pitrou wrote:


If Numpy wanted to switch to a different generator, and if Numba wanted
to remain compatible with Numpy, one of the PCG functions would be an
excellent choice (also for CPU performance, incidentally).


Is Apache license ok in NumPy?

(Not sure, thus asking. The PCG suite is Apache licensed.)





___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] [Python-ideas] Should our default random number generator be secure?

2015-09-18 Thread Sturla Molden

On 14/09/15 10:34, Antoine Pitrou wrote:


Currently we don't provide those APIs on the GPU, since MT is much too
costly there.

If Numpy wanted to switch to a different generator, and if Numba wanted
to remain compatible with Numpy, one of the PCG functions would be an
excellent choice (also for CPU performance, incidentally).


We have discussed allowing plugable PRNGs in NumPy on multiple occations.

I don't think NumPy should change its default PRNG. The Mersenne Twister 
MT19937 is state of the art of numerical work. It has excellent 
numerical accuracy, a huge period, and is very fast. For most users of 
NumPy, the MT19937 is the first and last word that needs to be said 
about pseudorandom number generators. We have the default that most 
users of NumPy will ever want, unless severe flaws are discovered in the 
algorithm (not very likely).


But there are occations where it does not work well. For example for 
parallel simulations it is (sometimes) nice to use DCMT instead of a 
vanilla MT19937. There are a lot of other generators that NumPy should 
consider as well, including the one you mention and Marsaglia's 
generators. We need to refactor randomkit to use a plugable entropy 
source instead of its internal MT19937. This is actually very easy to 
implement. In that case a user could just provide a Cython class or 
callback function (C, Fortran or Python) for generating a random ulong32.


The discussion on python-ideas refers to cryptography. That is not 
equally relevant for NumPy. But with a plugable design we can let NumPy 
users use os.urandom as well.



Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] NPY_DOUBLE not declared

2015-08-17 Thread Sturla Molden
Why not do as it says instead?


#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION

/* Or: NPY_NO_DEPRECATED_API NPY_API_VERSION */

#include 
#include 



Sturla


On 17/08/15 13:11, Florian Lindner wrote:
> Hello,
>
> I try to converse a piece of C code to the new NumPy API to get rid of the
> deprecation warning.
>
> #warning "Using deprecated NumPy API, disable it by " "#defining
> NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION"
>
> As a first step I replaced arrayobject.h by numpy/npy_math.h:
>
> #include 
> #include 
>
> But this gives errors that NPY_DOUBLE is not declared.
>
> http://docs.scipy.org/doc/numpy/reference/c-api.dtype.html#c.NPY_DOUBLE
> gives no information where NPY_DOUBLE is declared, so I used the standard
> npy_math.h header.
>
> src/action/PythonAction.cpp:90:51: error: 'NPY_DOUBLE' was not declared in
> this scope
> PyArray_SimpleNewFromData(1, sourceDim, NPY_DOUBLE,
> sourceValues);
>
>
> Including numpy/npy_common.h does not change it either.
>
> Thanks,
> Florian
>


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-08-03 Thread Sturla Molden
On 03/08/15 20:51, Chris Barker wrote:

> well, IIUC, np.int  is the python integer type, which is
> a C long in all the implemtations of cPython that I know about -- but is
> that a guarantee?in the future as well?

It is a Python int on Python 2.

On Python 3 dtype=np.int means the dtype will be C long, because a 
Python int has no size limit. But np.int aliases Python int. And 
creating an array with dype=int therefore does not create an array of 
Python int, it creates an array of C long. To actually get dtype=int we 
have to write dtype=object, which is just crazy.


Sturla



___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-08-03 Thread Sturla Molden
On 03/08/15 18:25, Chris Barker wrote:

> 2) The vagaries of the standard C types:  int, long, etc (spelled
> np.intc, which is a int32 on my machine, anyway)
>  [NOTE: is there a C long dtype? I can't find it at the moment...]

There is, it is called np.int. This just illustrates the problem...


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-03 Thread Sturla Molden
Matthew Brett  wrote:

> Sure, but to avoid confusion, maybe move the discussion of image
> indexing order to another thread?
> 
> I think this thread is about memory layout, which is a different issue.

It is actually a bit convoluted and not completely orthogonal. Memory
layout does not matter for 2d ndexing, i.e. (x,y) vs. (row, column), if you
are careful when iterating, but it does matter for Nd indexing. There is a
reason to prefer (x,y,z,t,r) in column major order or (recording, time,
slice, row, column) in row major order. Otherwise you can get very
inefficient memory traversals. Then if you work with visualization
libraries that expects (x,y,z) and column major order, e.g. ITK, VTK and
OpenGL, this is really what you want to use. And the choise of indexing
(x,y,z) cannot be seen as independent of the memory layout. Remember, it is
not just a matter of mapping coordinates to pixels. The data sets are so
large in MRI processing that memory layout does matter.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-03 Thread Sturla Molden
SMRUTI RANJAN SAHOO  wrote:

> well its really great idea. i can help on python but i don't have any
> knowledge on fortran.

I have been thinking in these lines too. But I have always thought it would
be too much work for very little in return, and it might not interop
properly with libraries written for NumPy (though PEP3118 might have
changed that). I am not sure using Fortran in addition to Cython is a good
idea, but it might be. At least if we limit the number of dimenstions to,
say, 4 or less, ot would be easy to implement most of the code in
vectorized Fortran.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-03 Thread Sturla Molden
Juan Nunez-Iglesias  wrote:

> The short version is that you'll save yourself a lot of pain by starting to
> think of your images as (plane, row, column) instead of (x, y, z). 

There are several things to consider here.

1. The vertices in computer graphics (OpenGL) are (x,y,z).

2. OpenGL rotation matrices and projection matrice are stored in column
major order.

3. OpenGL frame buffers are indexed (x,y) in column major order with (0,0)
being lower left.

4. ITK and VTK depends on OpenGL and are thus using column major order.

5. Those who use Matlab or Fortran in addition to Python prefer column
major order.

6. BLAS and LAPACK use column major order.

7. The common notation in image prorcessing (as opposed to computer
graphics in geberal) is indexing (row, column), in row major order, with
(0,0) being upper left.

All in all, this is a strong case for prefering column major order and the
common mathematical notation (x,y,z). 

Also notice how the ususal notation in image pricessing differs from
OpenGL.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-02 Thread Sturla Molden
On 02/08/15 22:28, Bryan Van de Ven wrote:
> And to eliminate the order kwarg, use functools.partial to patch the zeros 
> function (or any others, as needed):

This will probably break code that depends on NumPy, like SciPy and 
scikit-image. But if NumPy is all that matters, sure go ahead and monkey 
patch. Otherwise keep the patched functions in another namespace.

:-)

Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-02 Thread Sturla Molden
On 02/08/15 22:14, Kang Wang wrote:
> Thank you all for replying!
>
> I did a quick test, using python 2.6.6, and the original numpy package
> on my Linux computer without any change.
> ==
> x = np.zeros((2,3),dtype=np.int32,order='F')
> print "x.strides ="
> print x.strides
>
> y = x + 1
> print "y.strides ="
> print y.strides
> ==
>
> Output:
> 
> x.strides =
> (4, 8)
> y.strides =
> (12, 4)
> 

Update NumPy. This is the behavior I talked about that has changed.

Now NumPy does this:


In [21]: x = np.zeros((2,3),dtype=np.int32,order='F')

In [22]: y = x + 1

In [24]: x.strides
Out[24]: (4, 8)

In [25]: y.strides
Out[25]: (4, 8)



Sturla





___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Change default order to Fortran order

2015-08-02 Thread Sturla Molden
On 02/08/15 15:55, Kang Wang wrote:

> Can anyone provide any insight/help?

There is no "default order". There was before, but now all operators 
control the order of their return arrays from the order of their input 
array. The only thing that makes C order "default" is the keyword 
argument to np.empty, np.ones and np.zeros. Just monkey patch those 
functions and it should be fine.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-08-02 Thread Sturla Molden
On 31/07/15 09:38, Julian Taylor wrote:

> A long is only machine word wide on posix, in windows its not.

Actually it is the opposite. A pointer is 64 bit on AMD64, but the 
native integer and pointer offset is only 32 bit. But it does not matter 
because it is int that should be machine word sized, not long, which it 
is on both platforms.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-08-01 Thread Sturla Molden
Chris Barker - NOAA Federal  wrote:

> Which is part of the problem with C -- if two types happen to be the
> same, the compiler is perfectly happy.

That int and long int be the same is not more problematic than int and
signed int be the same.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-07-31 Thread Sturla Molden
Chris Barker - NOAA Federal  wrote:

> Turns out I was passing in numpy arrays that I had typed as "np.int".
> It worked OK two years ago when I was testing only on 32 bit pythons,
> but today I got a bunch of failed tests on 64 bit OS-X -- a np.int is
> now a C long!

It has always been C long. It is the C long that varies between platforms.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-07-31 Thread Sturla Molden
Chris Barker  wrote:

> What about Fortan -- I've been out of that loop for ages -- does
> semi-modern Fortran use well defined integer types?

Modern Fortran is completely sane.

INTEGER without kind number (Fortran 77) is the fastest integer on the CPU.
On AMD64 that is 32 bit, because it is designed to use a 64 bit pointer
with a 32 bit offset. (That is also why Microsoft decided to use a 32 bit
long, because it by definition is the fastest integer of at least 32 bits.
One can actually claim that the C standard is violated with a 64 bit long
on AMD64.) Because of this we use a 32 bit interger in BLAS and LAPACK
linked to NumPy and SciPy.

The function KIND (Fortran 90) allows us to query the kind number of a
given variable, e.g. to find out the size of INTEGER and REAL.

The function SELECTED_INT_KIND (Fortran 90) returns the kind number of
smallest integer with a specified range. 

The function SELECTED_REAL_KIND (Fortran 90) returns the kind number of
smallest float with a given range and precision. THe returned kind number
can be used for REAL and COMPLEX.

KIND, SELECTED_INT_KIND and SELECTED_REAL_KIND will all return compile-time
constants, and can be used to declare other variables if the return value
is stored in a variable with the attribute PARAMETER. This allows te
programmer to get the REAL, COMPLEX or INTEGER the algorithm needs
numerically, without thinking about how big they need to be in bits.

ISO_C_BINDING is a Fortran 2003 module which contains kind numbers
corresponding to all C types, including size_t and void*, C structs, an
attribute for using pass-by-value semantics, controlling the C name to
avoid name mangling, as well as functions for converting between C and
Fortran pointers. It allows portable interop between C and Fortran (either
calling C from Fortran or calling Fortran from C). 

ISO_FORTRAN_ENV is a Fortran 2003 and 2008 module. In F2003 it contain kind
numbers for integers with specified size: INT8, INT16, INT32, and INT64. In
F2008 it also contains kind numbers for IEEE floating point types: REAL32,
REAL64, and REAL128. The kind numbers for floating point types can also be
used to declare complex numbers.

So with modern Fortran we have a completely portable and unambiguous type
system.

C11/C++11 is sane as well, but not quite as sane as that of modern Fortran.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Shared memory check on in-place modification.

2015-07-27 Thread Sturla Molden
On 27/07/15 22:10, Anton Akhmerov wrote:
> Hi everyone,
>
> I have encountered an initially rather confusing problem in a piece of
> code that attempted to symmetrize a matrix: `h += h.T`
> The problem of course appears due to `h.T` being a view of `h`, and
> some elements being overwritten during the __iadd__ call.

Here is another example

 >>> a = np.ones(10)
 >>> a[1:] += a[:-1]
 >>> a
array([ 1.,  2.,  3.,  2.,  3.,  2.,  3.,  2.,  3.,  2.])

I am not sure I totally dislike this behavior. If it could be made 
constent it could be used to vectorize recursive algorithms. In the case 
above I would prefer the output to be:

array([ 1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9.,  10.])

It does not happen because we do not enforce that the result of one 
operation is stored before the next two operands are read. The only way 
to speed up recursive equations today is to use compiled code.


Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-07-27 Thread Sturla Molden
Chris Barker  wrote:

> 32 bits on all (most) 32 bit platforms
> 64 bits on 64 bit Linux and OS-X
> 32 bits on 64 bit Windows (also if compiled by cygwin??)

sizeof(long) is 8 on 64-bit Cygwin. This is to make sure it is inconsistent
with MSVC and MinGW-w64, and make sure there will always be ABI mismatches
unless the headerfiles are modified accordingly.

OTOH, it is one only sane 64-bit compiler on Windows. You can actually take
code written for 64 bit Linux or OSX and expect that it will work
correctly.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-07-27 Thread Sturla Molden
Chris Barker  wrote:

> we get away with np.float, because every OS/compiler that gets any regular
> use has np.float == a c double, which is always 64 bit.

Not if we are passing an array of np.float to a ac routine that expects
float*, e.g. in OpenGL, BLAS or LAPACK. That will for sure give crazy
results, just hang, or segfault.

I got away with pisting a PR with a "bugfix" which supposedly should fix a
case of precision loss in a SciPy routine, because I thought np.float was
np.float32 and not np.float64 (which it is). But it did make me feel rather
stupid.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Proposal: Deprecate np.int, np.float, etc.?

2015-07-24 Thread Sturla Molden
Julian Taylor  wrote:

> I don't see the issue. They are just aliases so how is np.float worse 
> than just float?

I have burned my fingers on it.

Since np.double is a C double I assumed np.float is a C float. It is not.

np.int has the same problem by being a C long. Pure evil. Most users of
NumPy probably expect the np.foobar dtype to map to the corresponding
foobar C type. This is actually inconsistent and plain dangerous.

It would be much better if dtype=float meant Python float, dtype=np.float
meant C float, dtype=int meant Python int, and dtype=np.int meant C int.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Rationale for integer promotion rules

2015-07-17 Thread Sturla Molden
Matthew Brett  wrote:

>>>  Furthermore, adding int64
>>> and uint64 returns float64.
>> 
>> This is a grievous kluge, on the grounds that no-one is really sure
>> *what* to do in this case.
> 
> It doesn't seem unreasonable to me : casting int64 to uint64 or uint64
> to int64 could lead to disastrous truncation.  float64 can exactly
> represent integers +/- 2**53 and will have some defined relative error
> above that.

We could promote to Python int on Python 3 and long on Python 2. 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] future of f2py and Fortran90+

2015-07-14 Thread Sturla Molden
Eric Firing  wrote:

> I'm curious: has anyone been looking into what it would take to enable 
> f2py to handle modern Fortran in general?  And into prospects for 
> getting such an effort funded?

No need. Use Cython and Fortran 2003 ISO C bindings. That is the only
portable way to interop between Fortran and C (including CPython) anyway.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] floats for indexing, reshape - too strict ?

2015-07-02 Thread Sturla Molden
Antoine Pitrou  wrote:

> I don't think relaxing type checking here makes any good.

I agee. NumPy should do the same as Python in this case.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Clarification sought on Scipy Numpy version requirements.

2015-06-21 Thread Sturla Molden
Charles R Harris  wrote:

> Ralf, I cannot compile Scipy 0.13.3 on my system, it seems to fail her

> _decomp_update.pyx:60:0: 'cython_blas.pxd' not found

Do you have a clean SciPy 0.13.3 source tree?

cython_blas.pxd is introduced in 0.16, and should be in 0.13 at all.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] I can't tell if Numpy is configured properly with show_config()

2015-06-19 Thread Sturla Molden
Elliot Hallmark  wrote:

> And I can't help but wonder if there is further configuration I need
> to make numpy faster, or if this is just a difference between out
> machines

Try to build NumPy with Intel MKL or OpenBLAS instead. 

ATLAS is only efficient on the host computer on which it is built, and even
there it is not very fast (but far better than the reference BLAS). 

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Clarification sought on Scipy Numpy version requirements.

2015-06-19 Thread Sturla Molden
Charles R Harris  wrote:

> I'm looking to change some numpy deprecations into errors as well as remove
> some deprecated functions. The problem I see is that
> SciPy claims to support Numpy >= 1.5 and Numpy 1.5 is really, really, old.
> So the question is, does "support" mean compiles with earlier versions
> of Numpy ?

It means there is a Travis CI build with NumPy 1.6.2. So any change to the
SciPy source code must compile with NumPy 1.6 and any later version of
NumPy. 

There is no Travis CI build with NumPy 1.5. I don't think we know for sure
if it is really compatible with the current SciPy.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Python 3 and isinstance(np.int64(42), int)

2015-06-18 Thread Sturla Molden
Nathaniel Smith  wrote:

> In py3,
> 'int' is an arbitrary width integer bignum, like py2 'long', which is
> fundamentally different from int32 and int64 in both semantics and
> implementation.

Only when stored in an ndarray.

An array scalar object does not need to care about the exact number of bits
used for storage as long as the storage is large enough, which a python int
always is.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] How to limit cross correlation window width in Numpy?

2015-06-18 Thread Sturla Molden
Mansour Moufid  wrote:

> The cross-correlation of two arrays of lengths m and n is of length
> m + n - 1, where m is usually much larger than n.

He is thinking about the situation where m == n and m is much larger than
maxlag.

Truncating the input arrays would also throw away data.

This is about correlating two long signals, not about correlating a signal
with a much shorter template.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] How to limit cross correlation window width in Numpy?

2015-06-17 Thread Sturla Molden
On 17/06/15 04:38, Honi Sanders wrote:

> I have now implemented this functionality in numpy.correlate() and 
> numpy.convolve(). https://github.com/bringingheavendown/numpy. The files that 
> were edited are:
> numpy/core/src/multiarray/multiarraymodule.c
> numpy/core/numeric.py
> numpy/core/tests/test_numeric.py
> Please look over the code, my design decisions, and the unit tests I have 
> written. This is my first time contributing, so I am not confident about any 
> of these and welcome feedback.

I'll just repeat here what I already said on Github.

I think this stems from the need to compute cross-correlograms as used 
in statistical signal analysis, whereas numpy.correlate and 
scipy.signal.correlate are better suited for matched filtering.

I think the best solution would be to add a function called 
scipy.signal.correlogram, which would return a cross-correlation and an 
array of time lags. It could take minlag and maxlag as optional arguments.

Adding maxlag and minlag arguments to numpy.convolve makes very little 
sense, as far as I am concerned.

Sturla









___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Aternative to PyArray_SetBaseObject in NumPy 1.6?

2015-06-16 Thread Sturla Molden
Eric Moore  wrote:
> You have to do it by hand in numpy 1.6.  For example see
>  href="https://github.com/scipy/scipy/blob/master/scipy/signal/lfilter.c.src#L285-L292";>https://github.com/scipy/scipy/blob/master/scipy/signal/lfilter.c.src#L285-L292

Thank you :)


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


[Numpy-discussion] Aternative to PyArray_SetBaseObject in NumPy 1.6?

2015-06-14 Thread Sturla Molden
What would be the best alternative to PyArray_SetBaseObject in NumPy 1.6?

Purpose: Keep alive an object owning data passed to
PyArray_SimpleNewFromData.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] changing ValueError to KeyError for bad field access

2015-06-12 Thread Sturla Molden
On 12/06/15 19:46, Allan Haldane wrote:

>  >>> a = np.ones(3, dtype=[('a', 'f4'), ('b', 'f4')])
>  >>> a['c']
>  KeyError: 'c'

> Any opinions?

Sounds good to me. But it will probably break someones code.


Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Verify your sourceforge windows installer downloads

2015-05-28 Thread Sturla Molden
Pauli Virtanen  wrote:

> Is it possible to host them on github? I think there's an option to add
> release notes and (apparently) to upload binaries if you go to the
> "Releases" section --- there's one for each tag.

And then Sourceforge will put up tainted installers "for the benefit of
NumPy users". :)

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Verify your sourceforge windows installer downloads

2015-05-28 Thread Sturla Molden
David Cournapeau  wrote:
> IMO, this really begs the question on whether we still want to use
> sourceforge at all. At this point I just don't trust the service at all
> anymore.

Here is their lame excuse:

https://sourceforge.net/blog/gimp-win-project-wasnt-hijacked-just-abandoned/

It probably means this:

If NumPy installers are moved away from Sourceforge, they will set up a
mirror and load the mirrored installers with all sorts of crapware. It is
some sort of racket the mob couldn't do better.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Verify your sourceforge windows installer downloads

2015-05-28 Thread Sturla Molden
Julian Taylor  wrote:

> It has been reported that sourceforge has taken over the gimp
> unofficial windows downloader page and temporarily bundled the
> installer with unauthorized adware:
> https://plus.google.com/+gimp/posts/cxhB1PScFpe

WTF?

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Strategy for OpenBLAS

2015-05-26 Thread Sturla Molden
Matthew Brett  wrote:

> I am getting the impression that OpenBLAS is looking like the most
> likely medium term solution for open-source stack builds of numpy and
> scipy on Linux and Windows at least.

I think you right. 

OpenBLAS might even be a long-term solution. We should also consider that
GotoBLAS (and GotoBLAS2) powered some of the World's most expensive
superomputers for a decade. It is not like this is untested software. 

The remaining test errors on Windows are also due to MSVC and MinGW-w64
differences, not due to OpenBLAS itself, and those are not relevant on
Linux.

On Apple, I am not sure which is better. Accelerate is faster in some
corner cases (level-1 BLAS with AVX, operations on very small matrices),
but it has issues with multiprocessing (GCD's threadpool is not forksafe).
Apart from that OpenBLAS and Accelerate are about equivalent in
performance. I have built OpenBLAS on OSX with clang and gfortran, it works
like a charm. So it might be worth cobsidering for binary wheels on OSX as
well.

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Backwards-incompatible improvements to numpy.random.RandomState

2015-05-24 Thread Sturla Molden
On 24/05/15 10:22, Antony Lee wrote:

> Comments, and help for writing tests (in particular to make sure
> backwards compatibility is maintained) are welcome.

I have one comment, and that is what makes random numbers so special? 
This applies to the rest of NumPy too, fixing a bug can sometimes change 
the output of a function.

Personally I think we should only make guarantees about the data types, 
array shapes, and things like that, but not about the values. Those who 
need a particular version of NumPy for exact reproducibility should 
install the version of Python and NumPy they need. That is why virtual 
environments exist.

I am sure a lot will disagree with me on this. So please don't take this 
as flamebait.


Sturla



___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Backwards-incompatible improvements to numpy.random.RandomState

2015-05-24 Thread Sturla Molden
On 24/05/15 20:04, Nathaniel Smith wrote:

> I'm not sure what you're envisioning as needing a deprecation cycle? The
> neat thing about random is that we already have a way for users to say
> that they want replicability -- the use of an explicit seed --

No, this is not sufficient for random numbers. Random sampling and 
ziggurat generators are examples. If we introduce a change (e.g. a 
bugfix) that will affect the number of calls to the entropy source, just 
setting the seed will in general not be enough to ensure backwards 
compatibility. That is e.g. the case with using ziggurat samplers 
instead of the current transcendental transforms for normal, exponential 
and gamma distributions. While ziggurat is faster (and to my knowledge) 
more accurate, it will also make a different number of calls to the 
entropy source, and hence the whole sequence will be affected, even if 
you do set a random seed.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Backwards-incompatible improvements to numpy.random.RandomState

2015-05-24 Thread Sturla Molden
On 24/05/15 17:13, Anne Archibald wrote:
> Do we want a deprecation-like approach, so that eventually people who
> want replicability will specify versions, and everyone else gets bug
> fixes and improvements? This would presumably take several major
> versions, but it might avoid people getting unintentionally trapped on
> this version.
>
> Incidentally, bug fixes are complicated: if a bug fix uses more or fewer
> raw random numbers, it breaks repeatability not just for the call that
> got fixed but for all successive random number generations.


If a function has a bug, changing it will change the output of the 
function. This is not special for random numbers. If not retaining the 
old erroneous output means we break-backwards compatibility, then no 
bugs can ever be fixed, anywhere in NumPy. I think we need to clarify 
what we mean by backwards compatibility for random numbers. What 
guarantees should we make from one version to another?


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] binary wheels for numpy?

2015-05-18 Thread Sturla Molden
On 18/05/15 21:57, Chris Barker wrote:
> On Sun, May 17, 2015 at 9:23 PM, Matthew Brett  > wrote:
>
> I believe OpenBLAS does run-time selection too.
>
>
> very cool! then an excellent option if we can get it to work (make that
> you can get it to work, I'm not doing squat in this effort other than
> nudging...)


Carl Kleffner has built binary wheels for NumPy and SciPy with OpenBLAS 
configured for run-time hardware detection. I don't remember at the top 
of my head where you can download them for testing. IIRC there remaining 
test failures were not related to OpenBLAS.

Sturla


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] binary wheels for numpy?

2015-05-18 Thread Sturla Molden
On 18/05/15 06:09, Chris Barker wrote:

> IIUC, The Intel libs have the great advantage of run-time selection of
> hardware specific code -- yes? So they would both work and give high
> performance on most machines (all?).

OpenBLAS can also be built for dynamic architecture with hardware 
auto-detection. IIRC you build with DYNAMIC_ARCH=1 instead of specifying 
TARGET.

Apple Accelerate Framework does this as well.


Sturla










___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] binary wheels for numpy?

2015-05-17 Thread Sturla Molden
On 17/05/15 20:54, Ralf Gommers wrote:

> I suspect; OpenBLAS seems like the way to go (?).

I think OpenBLAS is currently the most promising candidate to replace 
ATLAS. But we need to build OpenBLAS with MinGW gcc, due to AT&T syntax 
in the assembly code. I am not sure if the old toolchain is good enough, 
or if we will need Carl Kleffner's binaries.

Sturla






___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] binary wheels for numpy?

2015-05-17 Thread Sturla Molden
Matthew Brett  wrote:

> Yes, unfortunately we can't put MKL binaries on pypi because of the
> MKL license - see

I believe we can, because we asked Intel for permission. From what I heard
the response was positive.

 But it doesn't mean we should. :-)

Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Automatic number of bins for numpy histograms

2015-04-18 Thread Sturla Molden
Jaime Fernández del Río  wrote:

> I think we have an explicit rule against C++, although I may be wrong.

Currently there is Python, C and Cython in NumPy.

SciPy also has C++ and Fortran code.


Sturla

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] OS X wheels: speed versus multiprocessing

2015-04-06 Thread Sturla Molden
On 07/04/15 02:41, Nathaniel Smith wrote:

> Sure, but in some cases accelerate reduces speed by a factor of infinity
> by hanging, and OpenBLAS may or may not give wrong answers (but
> quickly!) since apparently they don't do regression tests, so we have to
> pick our poison.

OpenBLAS is safer on Mac than Windows (no MinGW related errors on Mac) 
so we should try it and see what happens.

GotoBLAS2 used to be great so it can't be that bad :-)




___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


  1   2   3   4   5   6   7   8   9   >