[petsc-dev] How to check that MatMatMult is available

[Pierre Jolivet]

> On 20 Sep 2019, at 7:36 AM, Jed Brown <jed at jedbrown.org <mailto:jed at jedbrown.org>> wrote:
> 
> Pierre Jolivet via petsc-dev <petsc-dev at mcs.anl.gov <mailto:petsc-dev at mcs.anl.gov>> writes:
> 
>> Hello,
>> Given a Mat A, I’d like to know if there is an implementation available for doing C=A*B
>> I was previously using MatHasOperation(A, MATOP_MATMAT_MULT, &hasMatMatMult) but the result is not correct in at least two cases:
> 
> Do you want MATOP_MAT_MULT and MATOP_TRANSPOSE_MAT_MULT?

Ah, OK, MATMAT => MatMatMat, so you are right.
I’ll make sure MatHasOperation_Transpose and MatHasOperation_Nest return the correct value then.

>> 1) A is a MATTRANSPOSE and the underlying Mat B=A^T has a MatTransposeMatMult implementation (there is currently no MATOP_MATTRANSPOSEMAT_MULT)
>> 2) A is a MATNEST. This could be fixed in MatHasOperation_Nest, by checking MATOP_MATMAT_MULT of all matrices in the MATNEST, but this would be incorrect again if there is a single MATTRANPOSE in there
>> What is then the proper way to check that I can indeed call MatMatMult(A,…)?
>> Do I need to copy/paste all this https://www.mcs.anl.gov/petsc/petsc-current/src/mat/interface/matrix.c.html#line9801 <https://www.mcs.anl.gov/petsc/petsc-current/src/mat/interface/matrix.c.html#line9801> in user code?
> 
> Unfortunately, I don't think there is any common interface for the
> string handling, though it would make sense to add one because code of
> this sort is copied many times:
> 
>   /* dispatch based on the type of A and B from their PetscObject's PetscFunctionLists. */
>   char multname[256];
>   ierr = PetscStrncpy(multname,"MatMatMult_",sizeof(multname));CHKERRQ(ierr);
>   ierr = PetscStrlcat(multname,((PetscObject)A)->type_name,sizeof(multname));CHKERRQ(ierr);
>   ierr = PetscStrlcat(multname,"_",sizeof(multname));CHKERRQ(ierr);
>   ierr = PetscStrlcat(multname,((PetscObject)B)->type_name,sizeof(multname));CHKERRQ(ierr);
>   ierr = PetscStrlcat(multname,"_C",sizeof(multname));CHKERRQ(ierr); /* e.g., multname = "MatMatMult_seqdense_seqaij_C" */
>   ierr = PetscObjectQueryFunction((PetscObject)B,multname,&mult);CHKERRQ(ierr);
>   if (!mult) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_INCOMP,"MatMatMult requires A, %s, to be compatible with B, %s",((PetscObject)A)->type_name,((PetscObject)B)->type_name);
> 
>> Thanks,
>> Pierre
>> 
>> PS: in my case, C and B are always of type MATDENSE. Should we handle
>> this in MatMatMult and never error out for such a simple case. 
> 
> I would say yes.

Perfect, that’s all I need!

Thanks,
Pierre

>> Indeed, one can just loop on the columns of B and C by doing multiple
>> MatMult. This is what I’m currently doing in user code when
>> hasMatMatMult == PETSC_FALSE.

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20190920/d9d867e8/attachment.html>