Actual source code: ex191.c

  1: static char help[] = "Tests MatLoad() for dense matrix with uneven dimensions set in program\n\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **args)
  6: {
  7:   Mat          A;
  8:   PetscViewer  fd;
  9:   PetscMPIInt  rank;
 10:   PetscScalar *Av;
 11:   PetscInt     i;

 14:   PetscInitialize(&argc, &args, (char *)0, help);
 15:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);

 17:   MatCreateDense(PETSC_COMM_WORLD, 6, 6, 12, 12, NULL, &A);
 18:   MatDenseGetArray(A, &Av);
 19:   for (i = 0; i < 6 * 12; i++) Av[i] = (PetscScalar)i;
 20:   MatDenseRestoreArray(A, &Av);

 22:   /* Load matrices */
 23:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "ex191matrix", FILE_MODE_WRITE, &fd);
 24:   PetscViewerPushFormat(fd, PETSC_VIEWER_NATIVE);
 25:   MatView(A, fd);
 26:   MatDestroy(&A);
 27:   PetscViewerPopFormat(fd);
 28:   PetscViewerDestroy(&fd);

 30:   MatCreate(PETSC_COMM_WORLD, &A);
 31:   MatSetType(A, MATDENSE);
 32:   if (rank == 0) {
 33:     MatSetSizes(A, 4, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_DETERMINE);
 34:   } else {
 35:     MatSetSizes(A, 8, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_DETERMINE);
 36:   }
 37:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "ex191matrix", FILE_MODE_READ, &fd);
 38:   MatLoad(A, fd);
 39:   PetscViewerDestroy(&fd);
 40:   MatView(A, PETSC_VIEWER_STDOUT_WORLD);
 41:   PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO_DETAIL);
 42:   MatView(A, PETSC_VIEWER_STDOUT_WORLD);
 43:   PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
 44:   MatDestroy(&A);
 45:   PetscFinalize();
 46:   return 0;
 47: }

 49: /*TEST

 51:    test:
 52:       nsize: 2
 53:       filter: grep -v alloced

 55: TEST*/