[petsc-dev] Dealing with off-processor aij-matrices not having preallocation information

[Karl Rupp]

Hi Jed,

 >     As the comment noted, this is quite a hack and actually breaks the
>     current implementation for GPU matrices. I think that the problem
>     can be completely resolved by not checking for PETSC_DEFAULT and
>     PETSC_DECIDE in the MPI preallocation above, but just forward d_nz
>     and o_nz to the sequential preallocation routine. This way, the
>     logic remains in MatSeqAIJSetPreallocation(), where
>     MAT_NEW_NONZERO_ALLOCATION_ERR is only set if the user provided an
>     explicit preallocation scheme.
>
>
> Hmm, the 'nonew' member is accessed in the MATMPIAIJ so we have to
> either change the code to avoid ever using it from the top level or we
> need to propagate the setting back to the parent. See also this discussion:
>
> http://lists.mcs.anl.gov/pipermail/petsc-dev/2012-January/007121.html
>
> https://bitbucket.org/petsc/petsc/commits/7827cd58

Yes, I want to preserve the behavior Barry described in 2), i.e. only 
error if the user specified an allocation pattern explicitly. I'm only 
questioning the way it is currently implemented:

If the number of nonzeros is specified as PETSC_DECIDE to
MatSeqAIJSetPreallocation_SeqAIJ(),
it does not set MAT_NEW_NONZERO_ALLOCATION_ERR. If the number of 
nonzeros is larger than zero, it does. So far, so good.

MatMPIAIJSetPreallocation_MPIAIJ(), however, first checks for the number 
of nonzeros. If it is PETSC_DECIDE, it sets defaults and then calls 
MatSeqAIJSetPreallocation_SeqAIJ() with the defaults. Because 
MatSeqAIJSetPreallocation_SeqAIJ() no longer sees PETSC_DECIDE, it sets 
MAT_NEW_NONZERO_ALLOCATION_ERR, which is then reset right after in 
MatMPIAIJSetPreallocation_MPIAIJ(). The CUSP-implementation does not 
reset the flag and thus leads to errors, even though the implementation 
of the CUSP-Preallocation routine looks perfectly fine. I rather want to 
fix the hack in MatMPIAIJSetPreallocation_MPIAIJ() than to propagate it 
to the GPU implementations as well.

I can create a pull request for the proposed changes if this is too 
complicated to be discussed via email.

Best regards,
Karli