Re: [eigen] inconsistent cwise() support

[Benoit Jacob <jacob.benoit.1@xxxxxxxxx>]

2009/11/18 Gael Guennebaud <gael.guennebaud@xxxxxxxxx>:
>
>
> On Wed, Nov 18, 2009 at 1:03 PM, Jitse Niesen <jitse@xxxxxxxxxxxxxxxxx>
> wrote:
>>
>> On Tue, 17 Nov 2009, Gael Guennebaud wrote:
>>
>>> I'm thinking that once we'll have a true array support, all this
>>> array/cwise stuff is going to be messy. Here I'll assume that a true Array
>>> class supporting the standard global math function (e.g., std::sin) is
>>> really needed.
>>
>> I agree that we need an easy way to apply functions like sin
>> coefficient-wise. However, the Array class has the disadvantage that A*B has
>> different meaning depending on whether A and B are Matrix or Array.
>>
>> To me, that seems a serious issue, but perhaps I'm alone in that opinion..
>
> Actually I think that's one of the main goals of the class Array.
>
>> Nevertheless, let me discuss how to mitigate this. When the Array class is
>> proposed, I thought that it was to support the different languages used in
>> different domains (e.g., and rather crudely, theoretical linear algebra for
>> Matrix, data manipulation for Array). A program in the first domain would
>> use exclusively Matrix, and a program in the second domain would use
>> exclusively Array. In that case, confusion is minimal, so I think that would
>> be good programming style.
>>
>> However, now it's envisaged that Matrix and Array are freely mixed, and if
>> cwise() is abolished it's even necessary to do so. So, the next best thing
>> is to make sure that all your variables are Matrix (supposing you're in the
>> theoretical linear algebra domain). As Gael said, this will lead to code
>> like the following for componentwise multiplication = Hadamard product:
>>
>>  MatrixXd A = ...
>>  MatrixXd B = ...
>>  MatrixXd prod = (A.array() * B.array()).matrix();
>
> Well, to be precise, the .matrix() here would be useless because Matrix
> would have a ctor from Array and vice versa. So:
>
> MatrixXd prod = A.array() * B.array();
>
> will work. Note that because Array::operator* would expect a template
> ArrayBase object, in spite of this Matrix from Array ctor, the following
> still won't work:
>
> MatrixXd prod = A.array() * B;
>
> That let me think such a ctor is safe.
>
>
>> Not only is that cumbersome, but I think that it's more difficult to
>> understand than A.cwise() * B.
>>
>>> I think it is also important to acknowledge that the current cwise() API
>>> is odd enough to make it hard to learn for new comers. They often have some
>>> difficulties to determine where .cwise() is needed (e.g., m1.cwise() *
>>> m2.cwise()), and how it propagates (though it does not, of course). For
>>> instance, if you look at "m1.cwise() * m2" and you know, or did not
>>> understand the prefix concept of .cwise(), then it seems that .cwise
>>> propagate to m2.
>>
>> I agree with this, cwise() is difficult, because of what you call the
>> prefix concept: in A.cwise().max(B), cwise() modifies max, while C++
>> programmers will think it modifies A. Perhaps a resolution would be to add a
>> cwisemax() function to MatrixBase, so that we would replace A.cwise().max(B)
>> to A.cwisemax(B). Now it's clear that the cwise prefix modifies max.
>>
>> With componentwise multiplication, that would mean that we replace
>> A.cwise() * B by A.cwisemul(B). It's a pity that we lose the multiplication
>> sign, but I think that we win on clarity: it's hard to misinterpret what
>> A.cwisemul(B) would mean.
>
> he he, this is exactly what we had at the very beginning of Eigen2. If I
> remember well the major reason, not to say the unique one, that motivated us
> to move to the cwise() prefix was to be able to overload operators. Indeed,
> if we list all the cwise operators we currently have together with their
> possible respective method names, one can see that some are very hard to
> translate (especially the +=, -=, and *= operators):
>
> + scalar <=> cwiseAdd(scalar)
> - scalar <=> cwiseSub(scalar)
> * mat <=> cwiseMul(mat)
> < mat <=> cwiseLessThan(mat)   or cwiseLT
> <= mat <=> cwiseLessEqual(mat)  or cwiseLE
> < scalar <=> cwiseLessThan(scalar)  or cwiseLT
> <= scalar <=> cwiseLessEqual(scalar)  or cwiseLE
> // same for >, >=, ==, !=
>
> += scalar   <=>  cwiseAddInPlace(scalar)
> -= scalar   <=>  cwiseSubInPlace(scalar)
> *= mat   <=>  cwiseMulInPlace(scalar)
>
> And we could also imagine to support all kind of bitwise operations: |, &,
> ~, ^, ! that we could name bitwiseXxx.
>
> Now, I also acknowledge this drawback is mitigated by the introduction of an
> Array class, because then we can expect that the need to apply such
> operators on matrices will be significantly reduced.
>
> Nevertheless, I think it is still important to be able to view matrix as
> array, and vice versa. Two options:
>
> 1 - We introduce easy to use .array() and .matrix() methods, and then we'll
> have both my proposal and your proposal, so again an infinite ways to do the
> same thing that is not very nice IMO.
>
> 2 - We assume such a conversion is worth it only if you use the
> view-as-array object multiple times, and so we enforce to use named wrapper
> objects, e.g.:
> MatrixXf m = ...;
> ArrayWrapper<MatrixXf> a(m);
> // take care that "a" and "m" reference the same data.
> a += 2;
>
>
> Finally, there is another issue that I forgot to mention in my previous
> email and that go in favor of your proposal. All those cwise
> operators/functions are not only useful for dense objects but also for
> sparse matrices and certainly all other special storage objects (especially
> abs(), abs2(), etc.). Since we want to preserve a consistent API as well as
> to reduce code duplication, this means that with my proposal the Array class
> should be able to potentially wrap any kind of matrices.... No need to say
> that significantly complicate the implementation of the Array class (e.g.,
> some will support only unary ops, some require to extend Array with other
> overloads, etc.). The implementation of your proposal is far much simpler
> because 1) Array would be reserved to dense storage, and 2) the cwise*
> functions won't be hidden into wrapper classes. (sorry if these last
> sentences does not make much sense to you, but believe me it's much simpler
> ;) )

Yes, it's much better. I am 100% for letting Array be plainly for dense storage.

>
> So after all, I'm going to like your proposal, and will wait for other
> opinions.
>
> As a side note, if we go for your proposal, we can also think about
> implementing in MatrixBase only a subset of the cwise ops available in
> Array.

Yes that's what I was going to reply: I like Jitse's proposal too but
we really don't have to provide .cwiseEverything(), only the most
common functions,
  matrix.cwiseProduct(matrix2)
  matrix.cwiseSum(scalar)
  matrix.cwiseAbs()
  matrix.cwiseAbs2()
Just that will be enough to cover 95% of cases and I don't even think
that we need the += and *= operators (but it's of course a possibility
to add a few more functions if you think it's needed).

Notice how i prefer to name functions according to what they are
returning, when they leave their arguments constant: Product instead
of Add.

On the other hand I don't think that we need such functions
matrixOperation() in the Array class. Doing matrix operations on an
Array is far, far less common, so I think that just having .matrix()
will be good enough.

Finally, I am not in favor of generalizing the use of global functions
everywhere. I understand that it makes sense to have global std::sin
and ei_sin for Array, to make that usable as any scalar type. And I
understand that for some functions, it is pleasant to have a global
version, such as dot(). But I don't think that it should become the
new rule that every function has a global version. (Referring e.g. to
above cwisemul(A,B)).

Benoit

>
> Gael.
>
>
>>
>> With the addition of Array, that leads to the following API. The idea is
>> that a module would use either the left or the right column, but not mix
>> expressions.
>>
>>                      Matrix A,B          Array A,B
>> matrix multiply:      A * B               A.matrixmul(B)
>> c-wise multiply:      A.cwisemul(B)       A * B
>> matrix exponential:   A.exp()             A.matrixexp()
>> c-wise exponential:   A.cwiseexp()        A.exp()
>>
>> The discussion on whether to introduce global functions instead of (or in
>> parallel of) member functions seems to be orthogonal to this, but just for
>> illustration, this is how it would work with global functions.
>>
>>                      Matrix A,B          Array A,B
>> matrix multiply:      A * B               matrixmul(A,B)
>> c-wise multiply:      cwisemul(A,B)       A * B
>> matrix exponential:   exp(A)              matrixexp(A)
>> c-wise exponential:   cwiseexp(A)         exp(A)
>>
>> Cheers,
>> Jitse
>>
>>
>
>