[petsc-users] questions on BAIJ matrix

[Barry Smith]

> On Jun 8, 2017, at 2:56 PM, Xiangdong <epscodes at gmail.com> wrote:
> 
> 
> On Thu, Jun 8, 2017 at 3:17 PM, Hong <hzhang at mcs.anl.gov> wrote:
> Xiangdong:
> MatCreateMPIBAIJWithArrays() is obviously buggy, and not been tested.
> 
> 
> 1) In the remark of the function MatCreateMPIBAIJWithArrays, it says " bs - the block size, only a block size of 1 is supported". Why must the block size be 1? Is this a typo? 
> 
> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateMPIBAIJWithArrays.html
> 
> It seems only bs=1 was implemented. I would not trust it without a test example. 
> 
> 2) In the Line 4040 of the implemention of MatCreateMPIBAIJWithArrays, would the matrix type be matmpibaij instead of matmpiSbaij?
> 
> This is an error. It should be matmpibaij.
> 
> http://www.mcs.anl.gov/petsc/petsc-current/src/mat/impls/baij/mpi/mpibaij.c.html#MatCreateMPIBAIJWithArrays
> 
> 4031: PetscErrorCode  MatCreateMPIBAIJWithArrays(MPI_Comm comm,PetscInt bs,PetscInt m,PetscInt n,PetscInt M,PetscInt N,const PetscInt i[],const PetscInt j[],const PetscScalar a[],Mat *mat)
> 4032: {
> 
> 4036:   if (i[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"i (row indices) must start with 0");
> 4037:   if (m < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"local number of rows (m) cannot be PETSC_DECIDE, or negative");
> 4038:   MatCreate(comm,mat);
> 4039:   MatSetSizes(*mat,m,n,M,N);
> 4040:   MatSetType(*mat,MATMPISBAIJ);
>  
>  It should be MATMPIBAIJ.
> 
> 3) I want to create a petsc matrix M equivalent to the sum of two block csr matrix/array (M1csr, M2csr). What is the best way to achieve it? I am thinking of created two petsc baij matrix (M1baij and M2baij) by calling MatCreateMPIBAIJWithArrays twice and then call MATAXPY to get the sum M=M1baij + M2baij. Is there a better way to do it?
> 
> This is an approach. However  MatCreateMPIBAIJWithArrays() needs to be fixed, tested and implemented with requested bs. What bs do you need?
> 
> Why does each bs need to be implemented separately? In the mean time, I modifed the implementation of MatCreateMPIBAIJWithArrays() a little bit to create a baij matrix with csr arrays.
> 
> MatCreate(comm,mat);
> MatSetSizes(*mat,m,n,M,N);
> MatSetType(*mat,MATMPIBAIJ);
> MatMPIBAIJSetPreallocationCSR(*mat,bs,i,j,a);
> MatSetOption(*mat,MAT_ROW_ORIENTED,PETSC_FALSE);
> MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);
> MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);
> 
> I just set the type to MATMPIBAIJ and delete the line MatSetOption before preallocation (otherwise I get error at runtime complaining using set options before preallocation) and it works fine. The only thing missing is that setting mat_row_oriented to be petsc_false has no effect on the final matrix, which I do not know how to fix.

   row_oriented is only used in MatSetValues[Blocked] for the array being passed in to set values. It has nothing to do with the actual storage of the matrix, for BAIJ the actual storage is always column oriented.

   Is this ok, why do you want to set row_oriented?

  Barry



> 
>  
> 
> Why not use MatCreate(), MatSetValuses() (set a block values at time) to create two MPIBAIJ matrices, then call MATAXPY. Since petsc MPIBAIJ matrix has different internal data structure than csr, 
> "The i, j, and a arrays ARE copied by MatCreateMPIBAIJWithArrays() into the internal format used by PETSc;", so this approach would give similar performance.
> 
> I will try this option as well. Thanks for your suggestions.
> 
> Xiangdong 
> 
> Hong
> 
>