Actual source code: ex35.c

  1: static char help[] = "MatLoad test for loading matrices that are created by DMCreateMatrix() and\n\
  2:                       stored in binary via MatView_MPI_DA.MatView_MPI_DA stores the matrix\n\
  3:                       in natural ordering. Hence MatLoad() has to read the matrix first in\n\
  4:                       natural ordering and then permute it back to the application ordering.This\n\
  5:                       example is used for testing the subroutine MatLoad_MPI_DA\n\n";

  7: #include <petscdm.h>
  8: #include <petscdmda.h>

 10: int main(int argc, char **argv)
 11: {
 12:   PetscInt    X = 10, Y = 8, Z = 8;
 13:   DM          da;
 14:   PetscViewer viewer;
 15:   Mat         A;

 17:   PetscFunctionBeginUser;
 18:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
 19:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp.dat", FILE_MODE_WRITE, &viewer));

 21:   /* Read options */
 22:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-X", &X, NULL));
 23:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-Y", &Y, NULL));
 24:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-Z", &Z, NULL));

 26:   /* Create distributed array and get vectors */
 27:   PetscCall(DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, X, Y, Z, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL, NULL, &da));
 28:   PetscCall(DMSetFromOptions(da));
 29:   PetscCall(DMSetUp(da));
 30:   PetscCall(DMSetMatType(da, MATMPIAIJ));
 31:   PetscCall(DMCreateMatrix(da, &A));
 32:   PetscCall(MatShift(A, X));
 33:   PetscCall(MatView(A, viewer));
 34:   PetscCall(MatDestroy(&A));
 35:   PetscCall(PetscViewerDestroy(&viewer));

 37:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp.dat", FILE_MODE_READ, &viewer));
 38:   PetscCall(DMCreateMatrix(da, &A));
 39:   PetscCall(MatLoad(A, viewer));

 41:   /* Free memory */
 42:   PetscCall(MatDestroy(&A));
 43:   PetscCall(PetscViewerDestroy(&viewer));
 44:   PetscCall(DMDestroy(&da));
 45:   PetscCall(PetscFinalize());
 46:   return 0;
 47: }