Actual source code: ex135.c

  1: static const char help[] = "Test parallel assembly of SBAIJ matrices\n\n";

  3: #include <petscmat.h>

  7: PetscErrorCode Assemble(MPI_Comm comm,PetscInt n,const MatType mtype)
  8: {
  9:   Mat A;
 10:   PetscInt first,last,i;
 12:   PetscMPIInt rank,size;

 15:   MatCreate(PETSC_COMM_WORLD,&A);
 16:   MatSetSizes(A, PETSC_DECIDE,PETSC_DECIDE,n,n);
 17:   MatSetType(A,MATMPISBAIJ);
 18:   MPI_Comm_size(comm,&size);
 19:   MPI_Comm_rank(comm,&rank);
 20:   if (rank < size-1) {
 21:     MatMPISBAIJSetPreallocation(A,1,1,PETSC_NULL,1,PETSC_NULL);
 22:   } else {
 23:     MatMPISBAIJSetPreallocation(A,1,2,PETSC_NULL,0,PETSC_NULL);
 24:   }
 25:   MatGetOwnershipRange(A,&first,&last);
 26:   MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);
 27:   last--;
 28:   for (i=first; i<=last; i++){
 29:     MatSetValue(A,i,i,2.,INSERT_VALUES);
 30:     if (i != n-1) {MatSetValue(A,i,n-1,-1.,INSERT_VALUES);}
 31:   }
 32:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 33:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 34:   MatDestroy(&A);
 35:   return(0);
 36: }

 40: int main(int argc,char *argv[])
 41: {
 43:   MPI_Comm       comm;
 44:   PetscInt       n = 6;

 46:   PetscInitialize(&argc,&argv,NULL,help);
 47:   comm = PETSC_COMM_WORLD;
 48:   PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
 49:   Assemble(comm,n,MATMPISBAIJ);
 50:   PetscFinalize();
 51:   return 0;
 52: }