Actual source code: ex207.c

  1: static char help[] = "Test MatCreateRedundantMatrix for a BAIJ matrix.\n\
  2:                       Contributed by Lawrence Mitchell, Feb. 21, 2017\n\n";

  4: #include <petscmat.h>
  5: int main(int argc, char **args)
  6: {
  7:   Mat         A, B;
  8:   Vec         diag;
  9:   PetscMPIInt size, rank;

 12:   PetscInitialize(&argc, &args, (char *)0, help);
 13:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 14:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);

 16:   MatCreate(PETSC_COMM_WORLD, &A);
 17:   MatSetSizes(A, 2, 2, PETSC_DETERMINE, PETSC_DETERMINE);
 18:   MatSetBlockSize(A, 2);
 19:   MatSetType(A, MATBAIJ);
 20:   MatSetUp(A);

 22:   MatCreateVecs(A, &diag, NULL);
 23:   VecSet(diag, 1.0);
 24:   MatDiagonalSet(A, diag, INSERT_VALUES);
 25:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 26:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 27:   MatView(A, PETSC_VIEWER_STDOUT_WORLD);

 29:   MatCreateRedundantMatrix(A, size, MPI_COMM_NULL, MAT_INITIAL_MATRIX, &B);
 30:   if (rank == 0) MatView(B, PETSC_VIEWER_STDOUT_SELF);

 32:   MatDestroy(&A);
 33:   MatDestroy(&B);
 34:   VecDestroy(&diag);
 35:   PetscFinalize();
 36:   return 0;
 37: }

 39: /*TEST

 41:    test:

 43:    test:
 44:       suffix: 2
 45:       nsize: 3

 47: TEST*/