[petsc-users] Matrix preallocation

[Mark Adams]

Woops, I misspoke. MatXAIJSetPreallocation does not take the
integer estimate. You have to give it the nnz arrays.
You can use NULL for the upper triangular part. The two args.
Mark

On Sat, Feb 5, 2022 at 6:20 AM Mark Adams <mfadams at lbl.gov> 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?
>>
>
> Yes. The integer is ignored if you provide the array. Just make sure that
> the integer is an upper bound on the number of non-zeros in any row. If
> not, MatSetValues can be very slow. Otherwise, the performance hit is not
> bad.
>
>
>> 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.
>>
>
> That is what the "X" is for. It is basically syntactic sugar to call the
> right version of preallocate (seq, mpi, blocked versions, hypre, etc.)
>
>
>> 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.
>>
>
> Sorry,
>
>
>>
>> A related followup question: Is it good practice to use this function or
>> should I just use the other routines like MatSeqAIJSetPreallocation() and
>> MatMPIAIJSetPreallocation()?
>>
>
> Yes this is the way to go.
> TL;DR
> Before we had a sparse matrix type explosion with GPUs, you could call
> both the MPI and Seq versions and be fine. You could call the blocked
> versions I suppose to be safe.
> But now there are several GPU/device enabled matrix types and "X" is
> convenient.
>
>
>>
>> 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?
>>
>
> You could do that but no need. As I said, the old way was to call both.
>
>
>> 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?
>>
>
> There is probably a default. It should work w/o any of these calls, but
> performance would suffer if the default (10 maybe) is too small for you.
>
> Mark
>
>
>> I hope that my questions were clear. Let me know if they need
>> clarification and thanks in advance for the help!
>>
>> Sam
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20220205/6a21a6c7/attachment-0001.html>