"Diagonal" Columns of MPIAIJ Matrix?

[Andrew Colombi]

Is MatGetOwnershipRangeColumn a new function? 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).

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?

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
>> >
>> >
>>
>