[petsc-dev] Bug in MatGetSubMatrix

[hong at aspiritech.org]

Stefano :

> There's a bug in MatGetSubMatrix; the case MAT_REUSE_MATRIX is not covered
> by the interface when  isrow == iscol and isrow is a stride 1 IS
> representing the whole set of rows. A tentative patch is attached.
>

Actually, current code is correct. The  case MAT_REUSE_MATRIX is covered:
          if (n == rend-rstart) {
            /* special case grabbing all rows; NEED to do a global
reduction to make sure all processes are doing this */
            if (cll == MAT_INITIAL_MATRIX) {
              *newmat = mat;
              ierr    =
PetscObjectReference((PetscObject)mat);CHKERRQ(ierr);
            }
            PetscFunctionReturn(0);
          }

Since 'newmat' is same as 'mat' in this case, we do NOT create a newmat,
instead, assign it to the address of 'mat'  in the case MAT_INITIAL_MATRIX,
i.e.,
'*newmat' = mat are the same matrix.

Whenever 'mat' is changed, so does '*newmat' ,
MatGetSubMatrix(C,isrow,NULL,MAT_REUSE_MATRIX,&A) does not need to do
anything.

I created a test ex181.c, attached below.

Hong
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20150309/70f4bf48/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: ex181.c
Type: text/x-csrc
Size: 2399 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20150309/70f4bf48/attachment.bin>