[Numpy-discussion] @ operator

[Nathaniel Smith]

On Thu, Sep 11, 2014 at 12:10 PM, Charles R Harris
<charlesr.harris at gmail.com> wrote:
>
> On Wed, Sep 10, 2014 at 10:08 PM, Nathaniel Smith <njs at pobox.com> wrote:
>>
>> My vote is:
>>
>> __matmul__/__rmatmul__ do the standard dispatch stuff that all __op__
>> methods do (so I guess check __array_priority__ or whatever it is we
>> always do). I'd also be okay with ignoring __array_priority__ on the
>> grounds that __numpy_ufunc__ is better, and there's no existing code
>> relying on __array_priority__ support in __matmul__.
>>
>> Having decided that we are actually going to run, they dispatch
>> unconditionally to np.newdot(a, b) (or newdot(a, b, out=a) for the
>> in-place version), similarly to how e.g. __add__ dispatches to np.add.
>>
>> newdot acts like a standard gufunc with all the standard niceties,
>> including __numpy_ufunc__ dispatch.
>>
>> ("newdot" here is intended as a placeholder name, maybe it should be
>> np.linalg.matmul or something else to be bikeshed later. I also vote
>> that eventually 'dot' become an alias for this function, but whether
>> to do that is an orthogonal discussion for later.)
>>
> If we went the ufunc route, I think we would want three of them, matxvec,
> vecxmat, and matxmat, because the best inner loops would be different in the
> three cases,

Couldn't we write a single inner loop like:

void ufunc_loop(blah blah) {
    if (arg1_shape[0] == 1 && arg2_shape[1] == 1) {
        call DOT
    } else if (arg2_shape[0] == 1) {
        call GEMV
    } else if (...) {
        ...
    } else {
        call GEMM
    }
}
?

I can't see any reason that this would be measureably slower than
having multiple ufunc loops: the checks are extremely cheap, will
usually only be done once (because most @ calls will only enter the
inner loop once), and if they are done repeatedly in the same call
then they'll come out the same way very time and thus be predictable
branches.

> but they couldn't be straight ufuncs themselves, as we don't
> need the other options, `reduce`, etc.,

Not sure what you mean -- ATM gufuncs don't support reduce anyway, but
if someone felt like implementing it then it would be cool to get
dot.reduce for free -- it is a useful/meaningful operation. Is there
some reason that supporting more ufunc features is bad?

> but they can't be exactly like the
> linalg machinery, because we do want subclasses to be able to override.

Subclasses can override ufuncs using __numpy_ufunc__.

> Hmm...
>
> The ufunc machinery has some funky aspects. For instance, there are
> hardwired checks for `__radd__` and other such operators in
> PyUFunc_GenericFunction  that allows subclasses to overide the ufunc. Those
> options should really be part of the PyUFuncObject.

Are you mentioning this b/c it's an annoying thing that talking about
ufuncs reminded you of, or is there a specific impact on __matmul__
that you see?

-n

-- 
Nathaniel J. Smith
Postdoctoral researcher - Informatics - University of Edinburgh
http://vorpus.org