[petsc-users] Matrix preallocation

[Samuel Estes]

Great. I think I've got it now. Thanks so much!

On Sun, Feb 6, 2022 at 1:48 PM Matthew Knepley <knepley at gmail.com> wrote:

> On Sun, Feb 6, 2022 at 2:39 PM Samuel Estes <samuelestes91 at gmail.com>
> wrote:
>
>> First of all, thank you so much for the detailed answers!
>> That clears up most of my confusion. Just to clarify let me make sure I
>> understand everything:
>> 1. So it seems that your advice is just to call the
>> MatXAIJSetPreallocation() routine (rather than make separate calls to other
>> preallocation routines such as MatSeqAIJPreallcation() and
>> MatMPIAIJPreallcation())? This will preallocate for Seq and MPI matrices
>> (among other types) and the decision will be made at runtime?
>>
>
> Yes.
>
>
>> 2. So the old way was just to call MatSeqAIJPreallcation() and
>> MatMPIAIJPreallcation()? I had thought that if you called the sequential
>> version for a parallel matrix or vice versa that the program would crash.
>> It seems that instead nothing happens, so by calling both you were covered
>> in either case. Is this correct? Essentially, one would always execute and
>> the other would do nothing?
>>
>
> Yes, PETSc works like Objective-C, in that if you call a routine for a
> subclass, and the object is of a different subclass, it is just ignored.
>
>
>> 3. If the dnnzu and onnzu arrays for the upper triangular parts seem
>> unnecessary for me I can just call the MatXAIJSetPreallocation() routine
>> with those values set to 'NULL'? I'm really just looking for a blend of the
>> MatSeqAIJPreallcation() and MatMPIAIJPreallcation() routines so it seems
>> that I only need dnnz and onnz arrays.
>>
>
> Yes.
>
>
>> 4. If we ignore the dnnzu and onnzu upper triangular parameters, then the
>> MatXAIJSetPreallocation() routine has the same parameters as the
>> MatMPIAIJSetPreallocation() routine without the option to give constant
>> integer values for diagonal and off-diagonal parts of the matrix so it
>> seems clear to me how this routine would work for a parallel matrix when
>> the upper triangular parameters are NULL. If the matrix is sequential then
>> does MPIXAIJSetPreallocation() just ignore the onnz parameter and just use
>> the dnnz parameter as the sole argument for determining preallocation?
>>
>
> Yes.
>
>   Thanks,
>
>     Matt
>
>
>> Thanks again for all the help. I think I understand well enough to use
>> this routine now. The questions above are mostly just to clarify and
>> double-check that I understood your previous responses correctly.
>>
>> Sam
>>
>> On Sat, Feb 5, 2022 at 10:12 AM Mark Adams <mfadams at lbl.gov> wrote:
>>
>>>
>>>
>>> On Sat, Feb 5, 2022 at 10:36 AM Matthew Knepley <knepley at gmail.com>
>>> wrote:
>>>
>>>> On Fri, Feb 4, 2022 at 11:47 PM Samuel Estes <samuelestes91 at gmail.com>
>>>> wrote:
>>>>
>>>>> Hi,
>>>>>
>>>>> I have a very basic question about matrix preallocation. I am trying
>>>>> to use the MatCreate(), MatSetFromOptions(), MatXXXXSetPreallocation()
>>>>> paradigm. I thought that I should use the MatXAIJSetPreallocation() routine
>>>>> since the code may be run with a SeqAIJ or MPIAIJ matrix but I do not
>>>>> understand all of the inputs required for the
>>>>> MatXAIJSetPreallocation routine. In particular, the dnnzu and
>>>>> onnzu variables don't quite make sense to me. Can these be NULL? I was
>>>>> basically just hoping for a routine that would preallocate for either a
>>>>> sequential or parallel matrix depending on what was given at runtime. This
>>>>> routine seems to be what I want but I don't understand it very well and the
>>>>> documentation hasn't helped me to figure it out.
>>>>>
>>>>
>>>> The example for this is here
>>>>
>>>>
>>>> https://petsc.org/main/docs/manualpages/Mat/MatMPIAIJSetPreallocation.html#MatMPIAIJSetPreallocation
>>>>
>>>> Maybe we should copy it to the XAIJ page as well. Does this help
>>>> explain the arguments?
>>>>
>>>>
>>>>> A related followup question: Is it good practice to use this function
>>>>> or should I just use the other routines like MatSeqAIJSetPreallocation()
>>>>> and MatMPIAIJSetPreallocation()?
>>>>>
>>>>> And finally my last question: if I were to use the
>>>>> MatSeqAIJSetPreallocation()/MatMPIAIJSetPreallocation() routines for
>>>>> preallocating memory, is it common to just call MatGetType() then call the
>>>>> appropriate routine depending on whether or not the matrix is parallel or
>>>>> not? I ask because when I have tested these routines out, it seems
>>>>> that MatSeqAIJSetPreallocation() works even for parallel matrices which is
>>>>> a bit confusing. I'm assuming that it just sets the diagonal part of the
>>>>> matrix?
>>>>>
>>>>
>>>> No, it definitely will not preallocate in parallel, so something else
>>>> is happening.
>>>>
>>>
>>> With only MatSeqAIJSetPreallocation in parallel it would not do any
>>> preallocation, in which case it would fall back to dynamic allocation. I'm
>>> guessing that is what is happening and it will be very slow.
>>>
>>>
>>>>
>>>>   Thanks,
>>>>
>>>>      Matt
>>>>
>>>>
>>>>> I hope that my questions were clear. Let me know if they need
>>>>> clarification and thanks in advance for the help!
>>>>>
>>>>> Sam
>>>>>
>>>>
>>>>
>>>> --
>>>> 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
>>>>
>>>> https://www.cse.buffalo.edu/~knepley/
>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>
>>>
>
> --
> 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
>
> https://www.cse.buffalo.edu/~knepley/
> <http://www.cse.buffalo.edu/~knepley/>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20220206/25e7a9bf/attachment.html>