[petsc-dev] Dealing with off-processor aij-matrices not having preallocation information
Matthew Knepley
knepley at gmail.com
Thu Mar 21 09:21:51 CDT 2013
On Thu, Mar 21, 2013 at 10:08 AM, Karl Rupp <rupp at mcs.anl.gov> wrote:
> Hi,
>
>
> > 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.
>>
>>
>> If nnz is nonnegative, yes.
>>
>
> right, it's actually nz >= 0.
>
>
> 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 thought perhaps you were going to leave it to the _SeqAIJ parts to
>> store the nonew state. If the _MPIAIJ part resolves PETSC_DECIDE
>> independently from the _SeqAIJ part, then they have to *always* coincide
>> exactly, otherwise we'll have very confusing bugs. You can likely
>> forward PETSC_DECIDE down to the _SeqAIJ parts, then query a->nonew and
>> b->nonew to set the mpiaij->nonew flag. (I'm not sure it's consistent
>> now, but we want it to be consistent.)
>>
>
> There is no mpiaij->nonew flag in Mat_MPIAIJ, this is all handled in the
> two diagonal/off-processor matrices of type Mat_SeqAIJ. I don't want to
> change this ;-)
> I just want to resolve the way way a->nonew and b->nonew is set, so that
> it does not require a hack in MatMPIAIJSetPreallocation_**MPIAIJ().
>
> A closer look showed three options:
> a) just forward PETSC_DECIDE to MatSeqAIJSetPreallocation_**SeqAIJ(),
> don't modify the nonew-flags from _MPIAIJ(). This almost fixes the problem,
> the only drawback is that the number of nonzeros on the off-processor
> blocks is then also 5 instead of 2 when using PETSC_DECIDE.
>
I vote a). No new function, gets rid of a stupid optimization that does no
one any good, and is the simplest.
Matt
> b) Allow to pass default values to MatSeqAIJSetPreallocation_**SeqAIJ().
> The cause of the issue is that MatSeqAIJSetPreallocation_**SeqAIJ() does
> two things: It interprets PETSC_DECIDE and it does the preallocation. If we
> split these two tasks into two functions, i.e.
>
> - MatSeqAIJSetPreallocation_**SeqAIJ() for interpreting PETSC_DECIDE and
> dealing with nonew, then calling
> - MatSeqAIJSetPreallocation_**SeqAIJAlloc() for the actual preallocation
>
> we can directly call MatSeqAIJSetPreallocation_**SeqAIJAlloc() from the
> _MPIAIJ() part and cleanly set MAT_NEW_NONZERO_ALLOCATION_ERR for the on-
> and off-processor matrices.
>
> c) Improve the existing hack by querying the state of a->nonew and
> b->nonew before calling MatSeqAIJSetPreallocation_**SeqAIJ() and after
> that reset a user-defined state if necessary.
>
> I'll create a pull request for b), as it is consistent with the current
> behavior and appears to be the cleanest solution to this.
>
> Best regards,
> Karli
>
>
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20130321/8cce3bbf/attachment.html>
More information about the petsc-dev
mailing list