Re: [Numpy-discussion] @ operator

2014-09-09 Thread Robert Kern
On Tue, Sep 9, 2014 at 8:52 PM, Charles R Harris
 wrote:
> Hi All,
>
> I'm in the midst of implementing the '@' operator (PEP 465), and there are
> some behaviors that are unspecified by the PEP.
>
> Should the operator accept array_like for one of the arguments?

I would be mildly disappointed if it didn't.

> Does it need to handle __numpy_ufunc__, or will __array_priority__ serve?

Not sure (TBH, I don't remember what __numpy_ufunc__ does off-hand and
don't feel bothered enough to look it up).

> Do we want PyArray_Matmul in the numpy API?

Probably.

> Should a matmul function be supplied by the multiarray module?

Yes, please. It's rules are a little different than dot()'s, so we
should have a function that does it.

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


[Numpy-discussion] @ operator

2014-09-09 Thread Charles R Harris
Hi All,

I'm in the midst of implementing the '@' operator (PEP 465), and there are
some behaviors that are unspecified by the PEP.


   1. Should the operator accept array_like for one of the arguments?
   2. Does it need to handle __numpy_ufunc__, or will __array_priority__
   serve?
   3. Do we want PyArray_Matmul in the numpy API?
   4. Should a matmul function be supplied by the multiarray module?

If 3 and 4 are wanted, should they use the __numpy_ufunc__ machinery, or
will __array_priority__ serve?

Note that the type number operators, __add__ and such, currently use
__numpy_ufunc__ in combination with __array_priority__, this in addition to
the fact that they are by default using ufuncs that do the same. I'd rather
that the __*__ operators simply rely on __array_priority__.


Thoughts?

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


Re: [Numpy-discussion] SFMT (faster mersenne twister)

2014-09-09 Thread Sturla Molden
On 09/09/14 20:08, Nathaniel Smith wrote:

> There's also another reason why generator decisions should be part of
> the RandomState object itself: we may want to change the distribution
> methods themselves over time (e.g., people have been complaining for a
> while that we use a suboptimal method for generating gaussian deviates),
> but changing these creates similar backcompat problems.

Which one should we rather use? Ziggurat?

Sturla

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


Re: [Numpy-discussion] SFMT (faster mersenne twister)

2014-09-09 Thread Nathaniel Smith
On 8 Sep 2014 14:43, "Robert Kern"  wrote:
>
> On Mon, Sep 8, 2014 at 6:05 PM, Pierre-Andre Noel
>  wrote:
> >  > I think we could add new generators to NumPy though,
> >  > perhaps with a keyword to control the algorithm (defaulting to the
> > current
> >  > Mersenne Twister).
> >
> > Why not do something like the C++11 ? In , a "generator"
> > is the engine producing randomness, and a "distribution" decides what is
> > the type of outputs that you want. Here is the example on
> > http://www.cplusplus.com/reference/random/ .
> >
> >  std::default_random_engine generator;
> >  std::uniform_int_distribution distribution(1,6);
> >  int dice_roll = distribution(generator);  // generates number in
> > the range 1..6
> >
> > For convenience, you can bind the generator with the distribution (still
> > from the web page above).
> >
> >  auto dice = std::bind(distribution, generator);
> >  int wisdom = dice()+dice()+dice();
> >
> > Here is how I propose to adapt this scheme to numpy. First, there would
> > be a global generator defaulting to the current implementation of
> > Mersene Twister. Calls to numpy's "RandomState", "seed", "get_state" and
> > "set_state" would affect this global generator.
> >
> > All numpy functions associated to random number generation (i.e.,
> > everything listed on
> > http://docs.scipy.org/doc/numpy/reference/routines.random.html except
> > for "RandomState", "seed", "get_state" and "set_state") would accept the
> > kwarg "generator", which defaults to the global generator (by default
> > the current Mersene Twister).
> >
> > Now there could be other generator objects: the new Mersene Twister,
> > some lightweight-but-worse generator, or some cryptographically-safe
> > random generator. Each such generator would have "RandomState", "seed",
> > "get_state" and "set_state" methods (except perhaps the
> > criptographically-safe ones). When calling a numpy function with
> > generator=my_generator, that function uses this generator instead the
> > global one. Moreover, there would be be a function, say
> > select_default_random_engine(generator), which changes the global
> > generator to a user-specified one.
>
> I think the Python standard library's example is more instructive. We
> have new classes for each new core uniform generator. They will share
> a common superclass to share the implementation of the non-uniform
> distributions. numpy.random.RandomState will continue to be the
> current Mersenne Twister implementation, and so will the underlying
> global RandomState for all of the convenience functions in
> numpy.random. If you want the SFMT variant, you instantiate
> numpy.random.SFMT() and call its methods directly.

There's also another reason why generator decisions should be part of the
RandomState object itself: we may want to change the distribution methods
themselves over time (e.g., people have been complaining for a while that
we use a suboptimal method for generating gaussian deviates), but changing
these creates similar backcompat problems. So we need a way to say "please
give me samples using the old gaussian implementation" or whatever.

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


[Numpy-discussion] numpy ignores OPT/FOPT under Python3

2014-09-09 Thread Thomas Unterthiner
Hi!

I want to use the OPT/FOPT environment viariables to set compiler flags 
when compiling numpy. However it seems that they get ignored under 
python3. Using Ubuntu 14.04 and numpy 1.9.0, I did the following:

 >export OPT="-march=native"
 >export FOPT = "-march=native"
 > python setup.py build   # "python" executes python2.7
[...snip...]
C compiler: x86_64-linux-gnu-gcc -pthread -fno-strict-aliasing 
-march=native -fPIC
^C

obviously under python 2.7, the additional march flag works as expected. 
However, using python3:

 >export OPT="-march=native"
 >export FOPT = "-march=native"
 >python3 setup.py build
[ snip ...]
C compiler: x86_64-linux-gnu-gcc -pthread -DNDEBUG -g -fwrapv -O2 -Wall 
-Wstrict-prototypes -g -fstack-protector --param=ssp-buffer-size=4 
-Wformat -Werror=format-security -D_FORTIFY_SOURCE=2 -fPIC
^C


obviously, the flags aren't used in python 3.  Did I overlook something 
here? Do $OPT/$FOPT only work in Python 2.7 by design, is this a bug or 
did I miss something?

Cheers

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


Re: [Numpy-discussion] Generalize hstack/vstack --> stack; Block matrices like in matlab

2014-09-09 Thread josef.pktd
On Tue, Sep 9, 2014 at 8:30 AM,  wrote:

>
>
>
> On Tue, Sep 9, 2014 at 5:42 AM, Stefan Otte  wrote:
>
>> Hey,
>>
>> @Josef, I wasn't aware of `bmat` and `np.asarray(np.bmat())` does
>> basically what I want and what I'm already using.
>>
>
> I never needed any tetris or anything similar except for the matched block
> version.
>
> Just to point out two more related functions
>
> scipy.sparse also has `bmat` for sparse block matrices
> scipy.linalg and scipy.sparse have `block_diag` to complement bmat.
>
>
> What I sometimes wish for is a sparse pseudo kronecker product as
> convenience to bmat, where the first (or the second) matrix contains (0,1)
> flags where the 1's specify where to put the blocks.
> (I'm not sure what I really mean, similar to block_diag but with a
> different filling pattern.)
>

(or in analogy to kronecker, np.nonzero provides the filling pattern, but
the submatrix is multiplied by the value.)

(application: unbalanced repeated measures or panel data)

Josef
()


>
> Josef
>
>
>
>>
>> Regarding the Tetris problem: that never happened to me, but stack, as
>> Josef pointed out, can handle that already :)
>>
>> I like the idea of removing the redundant square brackets:
>> stack([[a, b], [c, d]]) --> stack([a, b], [c, d])
>> However, if the brackets are there there is no difference between
>> creating a `np.array` and stacking arrays with `np.stack`.
>>
>> If we want to get fancy and turn this PR into something bigger
>> (working our way up to a NP-complete problem ;)) then how about this.
>> I sometimes have arrays that look like:
>>   AB0
>>   0 C
>> Where 0 is a scalar but is supposed to fill the rest of the array.
>> Having something like 0 in there might lead to ambiguities though. What
>> does
>>   ABC
>>   0D0
>> mean? One could limit the "filler" to appear only on the left or the
>> right:
>>   AB0
>>   0CD
>> But even then the shape is not completely determined. So we could
>> require to have one row that only consists of arrays and determines
>> the shape. Alternatively we could have a keyword parameter `shape`:
>>   stack([A, B, 0], [0, C, D], shape=(8, 8))
>>
>> Colin, with `bmat` you can do what you're asking for. Directly taken
>> from the example:
>> >>> np.bmat('A,B; C,D')
>> matrix([[1, 1, 2, 2],
>> [1, 1, 2, 2],
>> [3, 4, 7, 8],
>> [5, 6, 9, 0]])
>>
>>
>> General question: If `bmat` already offers something like `stack`
>> should we even bother implementing `stack`? More code leads to more
>> bugs and maintenance work.
>>
>>
>> Best,
>>  Stefan
>>
>>
>>
>> On Tue, Sep 9, 2014 at 12:14 AM, cjw  wrote:
>> >
>> > On 08-Sep-14 4:40 PM, Joseph Martinot-Lagarde wrote:
>> >> Le 08/09/2014 15:29, Stefan Otte a écrit :
>> >>> Hey,
>> >>>
>> >>> quite often I work with block matrices. Matlab offers the convenient
>> notation
>> >>>
>> >>>   [ a b; c d ]
>> > This would appear to be a desirable way to go.
>> >
>> > Numpy has something similar for strings.  The above is neater.
>> >
>> > Colin W.
>> >>> to stack matrices. The numpy equivalent is kinda clumsy:
>> >>>
>> >>> vstack([hstack([a,b]), hstack([c,d])])
>> >>>
>> >>> I wrote the little function `stack` that does exactly that:
>> >>>
>> >>>   stack([[a, b], [c, d]])
>> >>>
>> >>> In my case `stack` replaced `hstack` and `vstack` almost completely.
>> >>>
>> >>> If you're interested in including it in numpy I created a pull request
>> >>> [1]. I'm looking forward to getting some feedback!
>> >>>
>> >>>
>> >>> Best,
>> >>>Stefan
>> >>>
>> >>>
>> >>>
>> >>> [1] https://github.com/numpy/numpy/pull/5057
>> >>>
>> >> The outside brackets are redundant, stack([[a, b], [c, d]]) should be
>> >> stack([a, b], [c, d])
>> >>
>> >> ___
>> >> NumPy-Discussion mailing list
>> >> NumPy-Discussion@scipy.org
>> >> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>> >
>> > ___
>> > NumPy-Discussion mailing list
>> > NumPy-Discussion@scipy.org
>> > http://mail.scipy.org/mailman/listinfo/numpy-discussion
>> ___
>> NumPy-Discussion mailing list
>> NumPy-Discussion@scipy.org
>> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>>
>
>
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Generalize hstack/vstack --> stack; Block matrices like in matlab

2014-09-09 Thread josef.pktd
On Tue, Sep 9, 2014 at 5:42 AM, Stefan Otte  wrote:

> Hey,
>
> @Josef, I wasn't aware of `bmat` and `np.asarray(np.bmat())` does
> basically what I want and what I'm already using.
>

I never needed any tetris or anything similar except for the matched block
version.

Just to point out two more related functions

scipy.sparse also has `bmat` for sparse block matrices
scipy.linalg and scipy.sparse have `block_diag` to complement bmat.


What I sometimes wish for is a sparse pseudo kronecker product as
convenience to bmat, where the first (or the second) matrix contains (0,1)
flags where the 1's specify where to put the blocks.
(I'm not sure what I really mean, similar to block_diag but with a
different filling pattern.)

Josef



>
> Regarding the Tetris problem: that never happened to me, but stack, as
> Josef pointed out, can handle that already :)
>
> I like the idea of removing the redundant square brackets:
> stack([[a, b], [c, d]]) --> stack([a, b], [c, d])
> However, if the brackets are there there is no difference between
> creating a `np.array` and stacking arrays with `np.stack`.
>
> If we want to get fancy and turn this PR into something bigger
> (working our way up to a NP-complete problem ;)) then how about this.
> I sometimes have arrays that look like:
>   AB0
>   0 C
> Where 0 is a scalar but is supposed to fill the rest of the array.
> Having something like 0 in there might lead to ambiguities though. What
> does
>   ABC
>   0D0
> mean? One could limit the "filler" to appear only on the left or the right:
>   AB0
>   0CD
> But even then the shape is not completely determined. So we could
> require to have one row that only consists of arrays and determines
> the shape. Alternatively we could have a keyword parameter `shape`:
>   stack([A, B, 0], [0, C, D], shape=(8, 8))
>
> Colin, with `bmat` you can do what you're asking for. Directly taken
> from the example:
> >>> np.bmat('A,B; C,D')
> matrix([[1, 1, 2, 2],
> [1, 1, 2, 2],
> [3, 4, 7, 8],
> [5, 6, 9, 0]])
>
>
> General question: If `bmat` already offers something like `stack`
> should we even bother implementing `stack`? More code leads to more
> bugs and maintenance work.
>
>
> Best,
>  Stefan
>
>
>
> On Tue, Sep 9, 2014 at 12:14 AM, cjw  wrote:
> >
> > On 08-Sep-14 4:40 PM, Joseph Martinot-Lagarde wrote:
> >> Le 08/09/2014 15:29, Stefan Otte a écrit :
> >>> Hey,
> >>>
> >>> quite often I work with block matrices. Matlab offers the convenient
> notation
> >>>
> >>>   [ a b; c d ]
> > This would appear to be a desirable way to go.
> >
> > Numpy has something similar for strings.  The above is neater.
> >
> > Colin W.
> >>> to stack matrices. The numpy equivalent is kinda clumsy:
> >>>
> >>> vstack([hstack([a,b]), hstack([c,d])])
> >>>
> >>> I wrote the little function `stack` that does exactly that:
> >>>
> >>>   stack([[a, b], [c, d]])
> >>>
> >>> In my case `stack` replaced `hstack` and `vstack` almost completely.
> >>>
> >>> If you're interested in including it in numpy I created a pull request
> >>> [1]. I'm looking forward to getting some feedback!
> >>>
> >>>
> >>> Best,
> >>>Stefan
> >>>
> >>>
> >>>
> >>> [1] https://github.com/numpy/numpy/pull/5057
> >>>
> >> The outside brackets are redundant, stack([[a, b], [c, d]]) should be
> >> stack([a, b], [c, d])
> >>
> >> ___
> >> NumPy-Discussion mailing list
> >> NumPy-Discussion@scipy.org
> >> http://mail.scipy.org/mailman/listinfo/numpy-discussion
> >
> > ___
> > NumPy-Discussion mailing list
> > NumPy-Discussion@scipy.org
> > http://mail.scipy.org/mailman/listinfo/numpy-discussion
> ___
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Generalize hstack/vstack --> stack; Block matrices like in matlab

2014-09-09 Thread Stefan Otte
Hey,

@Josef, I wasn't aware of `bmat` and `np.asarray(np.bmat())` does
basically what I want and what I'm already using.

Regarding the Tetris problem: that never happened to me, but stack, as
Josef pointed out, can handle that already :)

I like the idea of removing the redundant square brackets:
stack([[a, b], [c, d]]) --> stack([a, b], [c, d])
However, if the brackets are there there is no difference between
creating a `np.array` and stacking arrays with `np.stack`.

If we want to get fancy and turn this PR into something bigger
(working our way up to a NP-complete problem ;)) then how about this.
I sometimes have arrays that look like:
  AB0
  0 C
Where 0 is a scalar but is supposed to fill the rest of the array.
Having something like 0 in there might lead to ambiguities though. What does
  ABC
  0D0
mean? One could limit the "filler" to appear only on the left or the right:
  AB0
  0CD
But even then the shape is not completely determined. So we could
require to have one row that only consists of arrays and determines
the shape. Alternatively we could have a keyword parameter `shape`:
  stack([A, B, 0], [0, C, D], shape=(8, 8))

Colin, with `bmat` you can do what you're asking for. Directly taken
from the example:
>>> np.bmat('A,B; C,D')
matrix([[1, 1, 2, 2],
[1, 1, 2, 2],
[3, 4, 7, 8],
[5, 6, 9, 0]])


General question: If `bmat` already offers something like `stack`
should we even bother implementing `stack`? More code leads to more
bugs and maintenance work.


Best,
 Stefan



On Tue, Sep 9, 2014 at 12:14 AM, cjw  wrote:
>
> On 08-Sep-14 4:40 PM, Joseph Martinot-Lagarde wrote:
>> Le 08/09/2014 15:29, Stefan Otte a écrit :
>>> Hey,
>>>
>>> quite often I work with block matrices. Matlab offers the convenient 
>>> notation
>>>
>>>   [ a b; c d ]
> This would appear to be a desirable way to go.
>
> Numpy has something similar for strings.  The above is neater.
>
> Colin W.
>>> to stack matrices. The numpy equivalent is kinda clumsy:
>>>
>>> vstack([hstack([a,b]), hstack([c,d])])
>>>
>>> I wrote the little function `stack` that does exactly that:
>>>
>>>   stack([[a, b], [c, d]])
>>>
>>> In my case `stack` replaced `hstack` and `vstack` almost completely.
>>>
>>> If you're interested in including it in numpy I created a pull request
>>> [1]. I'm looking forward to getting some feedback!
>>>
>>>
>>> Best,
>>>Stefan
>>>
>>>
>>>
>>> [1] https://github.com/numpy/numpy/pull/5057
>>>
>> The outside brackets are redundant, stack([[a, b], [c, d]]) should be
>> stack([a, b], [c, d])
>>
>> ___
>> NumPy-Discussion mailing list
>> NumPy-Discussion@scipy.org
>> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
> ___
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion