[petsc-users] MatMPIBAIJSetPreallocationCSR i and j indices

[Barry Smith]

   I forgot that MatSetValuesBlocked() takes row oriented input even within the blocks. So if we simply call MatSetOption(mat,MAT_ROW_ORIENTED,PETSC_FALSE) on the matrix after it is created then the storage should work as I originally planned.   

  Barry

I found this discussion values(bs,nblockcols,bs,nblockrows)  very confusing. MatMPIBAIJSetPreallocationCSR() doesn't take a square block of values so I didn't see what you are talking about. You must be talking about each internal call to MatSetValuesBlocked_MPIBAIJ(B,1,&row,ncols,icols,svals,INSERT_VALUES) where it puts in one block row of values with each call.




On Sep 28, 2013, at 9:12 PM, Jed Brown <jedbrown at mcs.anl.gov> wrote:

> Barry Smith <bsmith at mcs.anl.gov> writes:
> 
>>   I don't understand the (mild) rant below, the expected format of
>>   the values is suppose to be EXACTLY what one puts in the double
>>   *array for block CSR format. That is the values are (formally) a 1
>>   dimensional array) containing each block of values (in column
>>   oriented storage, that is the first column of the values in the
>>   first block, followed by the next column in the first block, etc
>>   etc) in the same order as the j indices.  Note that this means that
>>   in Fortran the values array could be thought of as double
>>   A(bs,bs,nnz) 
> 
> No, this is not how it works.  The ordering is
> 
>  values[nblockrows][bs][nblockcols][bs]
> 
> or, in Fortran
> 
>  values(bs,nblockcols,bs,nblockrows)
> 
> where nblockrows and nblockcols are the number of rows and columns being
> inserted in this all to MatSetValuesBlocked.
> 
> this layout means that the entries are packed in exactly the same way as
> if the block indices were blown up to scalar indices and then plain
> MatSetValues was called.
> 
>>   Note that this is exactly how Matteo has his stuff organized so "it
>>   should just work". If it does not then please submit a bug report
>>   with a very small matrix that reproduces the problem.
>> 
>>   This function is exactly what it says it is, it takes in exactly a matrix in block CSR format.
>> 
>>> wondering if a layout like values[nblockrows][nblockcols][bs][bs]
>>> would not be much more convenient for users.
>> 
>>  I don't understand this. It is a dense matrix as big as the entire sparse matrix.
> 
> Lisandro means the number of block rows and cols being inserted in the
> call, not the dimension of the matrix.
> 
>> 
>>   Barry
>> 
>> 
>> On Sep 27, 2013, at 4:02 PM, Lisandro Dalcin <dalcinl at gmail.com> wrote:
>> 
>>> On 27 September 2013 21:23, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>>> Just use MatCreateMPIBAIJWithArrays() saves thinking and work.
>>> 
>>> BTW, Have you ever thought about how hard to use is this API if you
>>> want to provide the matrix values? How are users supposed to fill the
>>> "values" array for a matrix with bs>1 ? This API is certainly not
>>> Fortran-friendly (and I think it not even C-friendly). Mateo had the
>>> CSR indices as well as the block values, but I recommended to use the
>>> CSR indices to perform proper matrix preallocation, then loop over
>>> cells and call MatSetValuesBlocked(). Otherwise, he would need to
>>> create an new intermediate array and fill it the right way for
>>> MatCreateMPIBAIJWithArrays() to work as expected.
>>> 
>>> The root of the problem is how MatSetValuesBlocked() expects the array
>>> of values to be layout in memory. While I understand that the current
>>> convention make it compatible with MatSetValues(), I'm always
>>> wondering if a layout like values[nblockrows][nblockcols][bs][bs]
>>> would not be much more convenient for users.
>>> 
>>> 
>>> 
>>> -- 
>>> Lisandro Dalcin
>>> ---------------
>>> CIMEC (UNL/CONICET)
>>> Predio CONICET-Santa Fe
>>> Colectora RN 168 Km 472, Paraje El Pozo
>>> 3000 Santa Fe, Argentina
>>> Tel: +54-342-4511594 (ext 1016)
>>> Tel/Fax: +54-342-4511169