Actual source code: ex233.c

  1: static char help[] = "Tests MatMPI{AIJ,BAIJ,SBAIJ}SetPreallocationCSR\n\n";

  3: #include <petscmat.h>

  5: int main(int argc,char **argv)
  6: {
  7:   PetscInt       ia[3]={0,2,4};
  8:   PetscInt       ja[4]={0,1,0,1};
  9:   PetscScalar    c[16]={0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15};
 10:   PetscInt       ia2[5]={0,4,8,12,16};
 11:   PetscInt       ja2[16]={0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3};
 12:   PetscScalar    c2[16]={0,1,4,5,2,3,6,7,8,9,12,13,10,11,14,15};
 13:   PetscMPIInt    size,rank;
 14:   Mat            ssbaij;
 15:   PetscBool      rect = PETSC_FALSE;

 17:   PetscInitialize(&argc,&argv,NULL,help);
 18:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 19:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 21:   if (rank) {
 22:     PetscInt i;
 23:     for (i = 0; i < 3; i++) ia[i] = 0;
 24:     for (i = 0; i < 5; i++) ia2[i] = 0;
 25:   }
 26:   PetscOptionsGetBool(NULL,NULL,"-rect",&rect,NULL);
 27:   MatCreate(PETSC_COMM_WORLD,&ssbaij);
 28:   MatSetBlockSize(ssbaij,2);
 29:   if (rect) {
 30:     MatSetType(ssbaij,MATMPIBAIJ);
 31:     MatSetSizes(ssbaij,4,6,PETSC_DECIDE,PETSC_DECIDE);
 32:   } else {
 33:     MatSetType(ssbaij,MATMPISBAIJ);
 34:     MatSetSizes(ssbaij,4,4,PETSC_DECIDE,PETSC_DECIDE);
 35:   }
 36:   MatSetFromOptions(ssbaij);
 37:   MatMPIAIJSetPreallocationCSR(ssbaij,ia2,ja2,c2);
 38:   MatMPIBAIJSetPreallocationCSR(ssbaij,2,ia,ja,c);
 39:   MatMPISBAIJSetPreallocationCSR(ssbaij,2,ia,ja,c);
 40:   MatViewFromOptions(ssbaij,NULL,"-view");
 41:   MatDestroy(&ssbaij);
 42:   PetscFinalize();
 43:   return 0;
 44: }

 46: /*TEST

 48:   test:
 49:     filter: grep -v type | sed -e "s/\.//g"
 50:     suffix: aijbaij_csr
 51:     nsize: 2
 52:     args: -mat_type {{aij baij}} -view -rect {{0 1}}

 54:   test:
 55:     filter: sed -e "s/\.//g"
 56:     suffix: sbaij_csr
 57:     nsize: 2
 58:     args: -mat_type sbaij -view -rect {{0 1}}

 60: TEST*/