Is it fair to say that MatPreallocateInitialize and MatPreallocateSet
need to be in the same scope for this to work? If so, that will never
work with my code. (Unless I can call repeatedly without having it
reinitialize things.)
I already had everything setup to work with my own bookkeeping (STL
vector and binary_search, in my tests STL set took too much memory).
But I need get the diagonal columns!
So be it. I'm willing to hack. May I assume PETSc distributes
columns with the following algorithm:
int diag_cols = total_cols / processors;
if (mpi_rank() < total_cols % processors)
diag_cols += 1;
Then assume the offsets follow the MPI rank? (i.e. higher rank has
later columns)
Thanks,
-Andrew
On Wed, Jul 2, 2008 at 3:34 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> On Jul 2, 2008, at 3:25 PM, Andrew Colombi wrote:
>
>> Is MatGetOwnershipRangeColumn a new function?
>
> Looks like it.
>
>> because I don't think it
>> I have it. I'm using petsc-2.3.3-p13. I searched petsc/include for
>> RangeColumn and didn't find anything.
>>
>> I could switch over to MatPreallocateInitialize, but I'm not entirely
>> sure how it works. The docs reference a chapter in the user manual
>> about performance,
>>
>> "See the chapter in the users manual on performance for more details"
>>
>> But I can't find anything in the PDF manual that mentions
>> PreallocateInitialize or PreallocateSet (I searched the text for that
>> string).
>
> The manual pages for these macros are pretty extensive; they are not
> discussed in the manual.
>>
>>
>> My first question is, how does MatPreallocateInitialize allocate dnz
>> and onz? It says we should not malloc or free them, but it takes
>> PetscInt *dnz not PetscInt **dnz? My understanding of how these
>> utilities work is that you do your normal matrix assembly, but use
>> MatPreallocateSet instead of MatSetValues. Then it just does some
>> counting and when you're done dnz and onz are correct. Am I correct?
>
> Yes. See src/dm/da/src/utils/fdda.c for many examples of its use.
>
> Barry
>
>>
>> Thanks,
>> -Andrew
>>
>> On Wed, Jul 2, 2008 at 1:52 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>>
>>> MatGetOwnershipRangeColumn() will give you the needed information. If you
>>> are coding from C/C++
>>> you might consider using the macro MatPreallocateInitialize() and
>>> friends. These utilities make it
>>> somewhat easier, in some circumstances, to count the entries in the
>>> diagonal and off-diagonal blocks.
>>>
>>> Barry
>>>
>>> On Jul 2, 2008, at 11:58 AM, Andrew Colombi wrote:
>>>
>>>> Hmm... this has not been my experience. When A is not square I'm
>>>> finding the that the diagonal portion is also not square. Instead it
>>>> assigns columns to the diagonal portion by evenly distributing the columns
>>>> among processors. For example, a 4 by 8 would look like:
>>>>
>>>> (please excuse the embedded HTML, I need it to specify a fixed-width
>>>> font)
>>>>
>>>> c0 c1 c2 c3 c4 c5 c6 c7
>>>> p0 r0 X X X X
>>>> p0 r1 X X X X
>>>> p1 r2 X X X X
>>>> p1 r3 X X X X
>>>>
>>>>
>>>> X marks the diagonal portions for each processor. I have attempted to
>>>> verify this using the program below. If you run this code it will complain
>>>> that processor 1 must allocate new memory for "(0, 3)" (which I am
>>>> interpreting to mean the local (0, 3), otherwise known as (2, 3)). Further
>>>> more you can stop the complaint by incrementing the off-diagonal
>>>> preallocation, and decrementing the diagonal preallocation. Note, that
>>>> specifying 0 is some special value that isn't 0. If you specify 0 you won't
>>>> get any malloc errors.
>>>>
>>>> Mat A;
>>>> int rows = 4;
>>>> int cols = 8;
>>>>
>>>> if (mpi_rank() == 0)
>>>> {
>>>> MatCreateMPIAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, rows,
>>>> cols,
>>>> 2, PETSC_NULL, 1, PETSC_NULL, &A);
>>>> MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR);
>>>>
>>>> static const int R = 2;
>>>> int set_rows[R] = {0, 1};
>>>> static const int C = 2;
>>>> int set_cols[C] = {0, 1};
>>>> double vals[C*R] = {1., 1., 1., 1.};
>>>> MatSetValues(A, R, set_rows, C, set_cols, vals, ADD_VALUES);
>>>> }
>>>> else
>>>> {
>>>> MatCreateMPIAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE,rows,
>>>> cols,
>>>> 2, PETSC_NULL, 1, PETSC_NULL, &A);
>>>> MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR);
>>>>
>>>> static const int R = 2;
>>>> int set_rows[R] = {2,3};
>>>> static const int C = 2;
>>>> int set_cols[C] = {2,3};
>>>> double vals[C*R] = {1., 1., 1., 1.};
>>>> MatSetValues(A, R, set_rows, C, set_cols, vals, ADD_VALUES);
>>>> }
>>>>
>>>> MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
>>>> MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
>>>>
>>>> On Wed, Jul 2, 2008 at 6:35 AM, Matthew Knepley <knepley at gmail.com>
>>>> wrote:
>>>>>
>>>>> By definition, the diagonal portion is square. Thus you retrieve the
>>>>> local rows [start, end), and then the diagonal block is
>>>>>
>>>>> [start, end) x [start, end)
>>>>>
>>>>> which you can do with
>>>>>
>>>>>
>>>>> http://www-unix.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/docs/manualpages/Mat/MatGetOwnershipRange.html
>>>>>
>>>>> Matt
>>>>>
>>>>> On Wed, Jul 2, 2008 at 2:52 AM, Andrew Colombi <acolombi at gmail.com>
>>>>> wrote:
>>>>>>
>>>>>> I'm letting PETSC_DECIDE how to distribute the rows and columns of my
>>>>>> non-square matrix, and I want to retrieve the columns it's decided to
>>>>>> make the "Diagonal" for a particular process. How can I do that? I'm
>>>>>> probably missing something obvious here, but I checked the docs.
>>>>>>
>>>>>> In case some context helps (or inspires a work-around):
>>>>>>
>>>>>> I want this information to count the number of d_nnz and o_nnz to
>>>>>> preallocate perfectly. This shouldn't be hard to do, but I need to
>>>>>> know the boundary between off-diagonal and diagonal to do it.
>>>>>>
>>>>>> Thanks,
>>>>>> -Andrew
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> 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
>>>>>
>>>>>
>>>>
>>>
>>
>
>
