Using SUBSET_NONZERO_PATTERN with MatAXPY

[Barry Smith]

    David,

     Thank you reporting the problem. It is a bug for the BAIJ format.

    Satish,

      Could you please make a patch for 2.3.3

     In the  file src/mat/impls/baij/seq/baij.c in the function  
MatAXPY_SeqBAIJ() locate

PetscErrorCode MatAXPY_SeqBAIJ(Mat Y,PetscScalar a,Mat X,MatStructure  
str)
{
   Mat_SeqBAIJ    *x  = (Mat_SeqBAIJ *)X->data,*y = (Mat_SeqBAIJ *)Y- 
 >data;
   PetscErrorCode ierr;
   PetscInt       i,bs=Y->rmap.bs,j,bs2;
   PetscBLASInt   one=1,bnz = PetscBLASIntCast(x->nz);

   PetscFunctionBegin;
   if (str == SAME_NONZERO_PATTERN) {
     PetscScalar alpha = a;
     BLASaxpy_(&bnz,&alpha,x->a,&one,y->a,&one);
   } else if (str == SUBSET_NONZERO_PATTERN) { /* nonzeros of X is a  
subset of Y's */
     if (y->xtoy && y->XtoY != X) {
       ierr = PetscFree(y->xtoy);CHKERRQ(ierr);
       ierr = MatDestroy(y->XtoY);CHKERRQ(ierr);
     }
     if (!y->xtoy) { /* get xtoy */
       ierr = MatAXPYGetxtoy_Private(x->mbs,x->i,x->j,PETSC_NULL, y- 
 >i,y->j,PETSC_NULL, &y->xtoy);CHKERRQ(ierr);
       y->XtoY = X;
++++++++++++++++++++++++add the next line
       ierr = PetscObjectReference((PetscObject)X);CHKERRQ(ierr);
++++++++++++++++++++++++++++
recompile in that directory with "make lib shared"

    Please let us know if you continue to have trouble at petsc-maint at mcs.anl.gov

    Barry


On Mar 21, 2008, at 10:50 AM, David Knezevic wrote:

> Hello,
>
> I've got a SeqBAIJ matrix A which I am reassembling at each time- 
> step in order to solve a time-dependent PDE. I initially set up the  
> sparsity pattern of A, and then assemble A at each time step by  
> doing operations like the following:
>
> ierr = MatZeroEntries(A);CHKERRQ(ierr);
> ierr = MatAXPY(A,const1,B1,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
> ierr = MatAXPY(A,const2,B2,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
>
> This works fine. However, the sparsity patterns of B1, B2 etc should  
> all be subsets of the sparsity pattern of A (which should be  
> retained by MatZeroEntries). But changing DIFFERENT_NONZERO_PATTERN  
> to SUBSET_NONZERO_PATTERN in the above code causes an error. In  
> particular, the code runs fine for the first time step, but on the  
> second time step I get the error:
>
> [0]PETSC ERROR: --------------------- Error Message  
> ------------------------------------
> [0]PETSC ERROR: Invalid argument!
> [0]PETSC ERROR: Wrong type of object: Parameter # 3!
> ...
> [0]PETSC ERROR: MatAXPY() line 34 in src/mat/utils/axpy.c
>
> I found a description of the same problem here:
> http://ganymed.iwr.uni-heidelberg.de/pipermail/dealii/2007/002066.html
> but I haven't found an answer. Assistance would be most appreciated.
>
> Best regards,
> David
>
>