Stefano,
There is currently a PCFieldSplitSetBlockSize() and a PCSPAISetBlockSize()
PetscErrorCode MatSetBlockSize(Mat mat,PetscInt bs)
{
PetscErrorCode ierr;
PetscFunctionBegin;
PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
PetscValidLogicalCollectiveInt(mat,bs,2);
ierr = PetscLayoutSetBlockSize(mat->rmap,bs);CHKERRQ(ierr);
ierr = PetscLayoutSetBlockSize(mat->cmap,bs);CHKERRQ(ierr);
PetscErrorCode PetscLayoutSetBlockSize(PetscLayout map,PetscInt bs)
{
PetscFunctionBegin;
if (bs < 0) PetscFunctionReturn(0);
if (map->n > 0 && map->n % bs) SETERRQ2(map->comm,PETSC_ERR_ARG_INCOMP,"Local
size %D not compatible with block size %D",map->n,bs);
if (map->range && map->bs > 0 && map->bs != bs)
SETERRQ2(map->comm,PETSC_ERR_ARG_INCOMP,"Cannot change block size %D to
%D",map->bs,bs);
if (map->mapping) {
PetscInt lbs;
PetscErrorCode ierr;
ierr = ISLocalToGlobalMappingGetBlockSize(map->mapping,&lbs);CHKERRQ(ierr);
if (lbs != bs) SETERRQ2(map->comm,PETSC_ERR_PLIB,"Blocksize of
localtoglobalmapping %D must match that of layout %D",lbs,bs);
}
map->bs = bs;
PetscFunctionReturn(0);
So it doesn't want to change the block size if map->range has been set (i.e.
PetscLayoutSetUp()) was called. But if the matrix is AIJ so long as !(map->n %
bs) one can actually set the block size at any time. So how about if I change
the Mat code to allow setting the block size late for AIJ? Maybe we can also
get rid of PCFieldSplitSetBlockSize() and a PCSPAISetBlockSize()?
One issue that comes up with this change is: should setting the block size
later automatically reset the block size of the two inner AIJ matrices? Or
should they remain one? I could have it update the inner block sizes at the
expense of a tiny bit more code. But one complication is that
MatSetUpMultiply_MPIAIJ() always sets the column block size of the off-diagaonl
matrix to 1; so this needs to be respected.
The PETSc design policy of "there is one way to do something" dictates that it
is best not to have PCSetBlockSize() unless absolutely necessary so I am trying
to see if it is absolutely necessary.
Barry
> On Oct 10, 2016, at 10:18 AM, Stefano Zampini <[email protected]>
> wrote:
>
> Mark,
>
> I was wondering if you could provide an API/command line customization for
> PCGAMG to explicitly provide the "block size" of the problem to the gamg
> solver, instead of using the block size of the matrix. The matrix block size
> could be always used as a fallback if this new information as not been
> provided.
>
> I'm interfacing a complicated finite element package to PETSc , where you
> cannot do any assumption on the blocksize of the problem during matrix setup,
> and (AFAIK) I cannot change the block size with a call to MatSetBlockSize
> after the matrix has been finalized without hacking (that I would like to
> avoid).
>
> Thanks
> --
> Stefano
