Complete binary operators for dense matrix/vector types

This was discussed a bit on the Developer Meeting 2018. The topic is still controversial, so will need a bit of benchmarking.

TODO:

• Figure out which operators exactly are missing.
• Figure out how people would be using those operators.
  (Several users have actually implemented these operators in their own libs, those could be used to get an idea.)
• Benchmark. Figure out which kind of uses do give different performance.
  (I started with that a little here joe/dense-types-benchmark based on google benchmark lib, though there isn't much dense types specific yet)
• Decide whether to go ahead with implementation.
• Implement.

@oliver.sander You volunteered to provide the overloads in a branch. I thought about that, it is probably not a big help after all, so don't bother if you haven't started yet. For the benchmark, the overloads can just be implemented in the benchmark itself.

Hi @joe , thanks for taking over. The operators I think would be helpful are

Scalar * Vector
Vector * Scalar
Vector + Vector
Vector - Vector

Scalar * Matrix
Matrix * Scalar
Matrix + Matrix
Matrix - Matrix

Vector * Matrix
Matrix * Vector
Matrix * Matrix

I also have a long list of operations in my usercode. One problem is the fact that a scalar could be K, FieldVector<K,1>, or FieldMatrix<K,1,1>. And it happens that you have to multiply any of these scalars with a matrix a vector or any other of these scalars

I propose that FieldVector<K,1> and FieldMatrix<K,1,1> should not be treated as scalars here. I'd like to get away from that anyway.

I am going to ignore FieldVector<K,1> and FieldMatrix<K,1,1> -- yes, they should at minimum have all the operations their larger counterparts have, but for performance purposes they are just scalars.

(Edit: I do not mean to imply they should be usable as scalars -- that is an orthogonal issue)

More operators that I'd expect:

Vector / Scalar
Matrix / Scalar

The implicit conversion should take care of FieldVector<K, 1> and FieldMatrix<K, 1, 1> as long as the Matrix/Vector combinations are properly constraint in the dimensions.

Actually, implicit conversion does not help here, due to conflicting deduced types, see https://godbolt.org/z/40I_yb (or maybe I misunderstood your point)

I would also love to ignore the FieldVector<K,1> and FieldMatrix<K,1,1>, but they are returned by some dune functions, like dune-istl BlockVector::operator[] if the block-type is FieldVector<K,1>.

but they are returned by some dune functions, like dune-istl BlockVector::operator[] if the block-type is FieldVector<K,1>.

I would still ignore those cases. My mid-range goal is to get to a point where you can write BlockVector<double> anyway.

Actually, implicit conversion does not help here, due to conflicting deduced types, see https://godbolt.org/z/40I_yb (or maybe I misunderstood your point)

more like this https://godbolt.org/z/LV9RLq. As far as I understood the notes this additional constraint would be necessary for an eventual ScalarType anyway.

@janick.gerstenberger I have implemented exactly your approach. The problem with this is that you get ambiguous calls at some points and you have to implement several more specializations to solve those, or you constrain your function even more. Maybe this is what you meant.

@oliver.sander I would support this and in the meantime everyone who needs these wrapped scalar types has to implement these operations herself.

Although it is not a binary operation, I would like to have negation:

-vector
-matrix

And maybe (this could possible be discussed), operations with a diagonal matrix:

matrix * diag
matrix + diag
matrix - diag

diag * matrix
diag + matrix
diag - matrix

I'm happy with most of the proposed methods, but we should IMO avoid

vector * vector
vector * matrix

because they lead to a non-associative product because

matrix*(vector*vector) != (matrix*vector)*vector

To also support these operations we should probably introduce a way to explicitly interpret a vector as row-vector. This would have the nice side-effect that you could also consistently implement the dyadic product of two vectors.

Rudimentary implementation in !705 (merged) .

The implementation in !705 (merged) is decent now -- it now only modifies FieldMatrix and FieldVector. Any objections to merging this?

As @simon.praetorius I am also in favor of not treating FieldMatrix and FieldVectors as scalars. Mainly because it requires A LOT of code to maintain it while its benefit is not that big. Yes, it makes some operations handy, but it also some type conversions obscure. Moreover, lately I have been trying to compile all the tests with Visual Studio and this is a rather big bottleneck. VS complains that there are ambiguities in the conversion of types and I do not have too much to argue against it. A typical case that fails VS is the following

using K = std::complex<double>;
Dune::FieldMatrix<K,1,1> A;
std::complex<K> result = (A *= K(2));

In which order this shoud be computed?

• Convert first to scalar type, that is std::complex<double>, and then do the *= operation using the built-in C++ operation.
• Use the *= operation of DenseMatrix implicitly converting 2 into dobule, and then converting Dune::FieldMatrix<K,1,1> into a std::complex<K>.

I cannot really tell whichis the correct way, but I can certainly tell that there certain reduncaies that should be avoided.

mentioned in merge request !953 (closed)

mentioned in merge request !1289