[petsc-dev] [petsc-users] How to compute RARt with A and R as distributed (MPI) matrices ?

[Franck Houssen]

>From a user point of view, PETSc provides transpose, RARt and PtAP.

If one among most (rookie, not so experienced) users wants to compute NtMN, he could start with M and:
1. N,  so he could call RARt (with A=M and R= Nt) after having transposed N.
2. N,  so he could call PtAP (with A=M and P= N ).
3. Nt, so he could call RARt (with A=M and R= Nt).
4. Nt, so he could call PtAP (with A=M and P= N ) after having transposed Nt.

>From a user (theoretical, philosophical - name it as you like) point of view, the approaches/patterns 1, 2, 3 and 4 are identical, right ?
But, the implementation is not (this is not such a problem - there may have good reasons for that, one can not support always everything).

If the (possibly unlucky) user, had started with N and M as:
a. sequentiel  matrices, then, 1 is OK, 2 is KO, 3 is OK, 4 is KO.
b. distributed matrices, then, 1 is KO, 2 is OK, 3 is KO, 4 is OK.
There is no way, a user could understand why a pattern works with (sequential) N1 and M1, but not with (unluckyly distributed) N2 and M2... Unless there is something in the doc and/or a FAQ to guide the user that says "this kind of Mat is not supported by this kind of method because the implementation is not yet ready for that". The user get stuck in a dead end with no clue to get out.

Anyway, no big deal.

> > I believe adding a short note in MatRARt/PtAP pages and
> > MatRARt/PtAPNumeric/Symbolic would be welcome to say "this method is meant
> > to do/support this but not that”.
> 
> An operation like RARt (or PtAP) may have a specialized implementation, or be
> available though composition of operatotions, i.e. PtAP = P^t * (A * P) or
> compute transpose of P (explicitly) and then call the triple matrix product.
> It’s not possible to keep track of all the possible cases.

No, it is not possible: that's where the doc comes in.