Actual source code: ex33.c


  2: static char help[] = "Tests VecView()/VecLoad() for DMDA vectors (this tests DMDAGlobalToNatural()).\n\n";

  4: #include <petscdm.h>
  5: #include <petscdmda.h>
  6: #include <petscviewerhdf5.h>

  8: int main(int argc, char **argv)
  9: {
 10:   PetscMPIInt     rank, size;
 11:   PetscInt        N = 6, M = 8, P = 5, dof = 1;
 12:   PetscInt        stencil_width = 1, pt = 0, st = 0;
 13:   PetscBool       flg2, flg3, isbinary, mpiio;
 14:   DMBoundaryType  bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE, bz = DM_BOUNDARY_NONE;
 15:   DMDAStencilType stencil_type = DMDA_STENCIL_STAR;
 16:   DM              da, da2;
 17:   Vec             global1, global2;
 18:   PetscScalar     mone = -1.0;
 19:   PetscReal       norm;
 20:   PetscViewer     viewer;
 21:   PetscRandom     rdm;
 22: #if defined(PETSC_HAVE_HDF5)
 23:   PetscBool ishdf5;
 24: #endif

 27:   PetscInitialize(&argc, &argv, (char *)0, help);
 28:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 29:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 31:   PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL);
 32:   PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL);
 33:   PetscOptionsGetInt(NULL, NULL, "-P", &P, NULL);
 34:   PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL);
 35:   PetscOptionsGetInt(NULL, NULL, "-stencil_width", &stencil_width, NULL);
 36:   PetscOptionsGetInt(NULL, NULL, "-periodic", &pt, NULL);
 37:   if (pt == 1) bx = DM_BOUNDARY_PERIODIC;
 38:   if (pt == 2) by = DM_BOUNDARY_PERIODIC;
 39:   if (pt == 4) {
 40:     bx = DM_BOUNDARY_PERIODIC;
 41:     by = DM_BOUNDARY_PERIODIC;
 42:   }

 44:   PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL);
 45:   stencil_type = (DMDAStencilType)st;

 47:   PetscOptionsHasName(NULL, NULL, "-oned", &flg2);
 48:   PetscOptionsHasName(NULL, NULL, "-twod", &flg2);
 49:   PetscOptionsHasName(NULL, NULL, "-threed", &flg3);

 51:   PetscOptionsHasName(NULL, NULL, "-binary", &isbinary);
 52: #if defined(PETSC_HAVE_HDF5)
 53:   PetscOptionsHasName(NULL, NULL, "-hdf5", &ishdf5);
 54: #endif
 55:   PetscOptionsHasName(NULL, NULL, "-mpiio", &mpiio);
 56:   if (flg2) {
 57:     DMDACreate2d(PETSC_COMM_WORLD, bx, by, stencil_type, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, &da);
 58:   } else if (flg3) {
 59:     DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, stencil_type, M, N, P, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, 0, &da);
 60:   } else {
 61:     DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da);
 62:   }
 63:   DMSetFromOptions(da);
 64:   DMSetUp(da);

 66:   DMCreateGlobalVector(da, &global1);
 67:   PetscObjectSetName((PetscObject)global1, "Test_Vec");
 68:   PetscRandomCreate(PETSC_COMM_WORLD, &rdm);
 69:   PetscRandomSetFromOptions(rdm);
 70:   VecSetRandom(global1, rdm);
 71:   if (isbinary) {
 72:     if (mpiio) PetscOptionsSetValue(NULL, "-viewer_binary_mpiio", "");
 73:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer);
 74: #if defined(PETSC_HAVE_HDF5)
 75:   } else if (ishdf5) {
 76:     PetscViewerHDF5Open(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer);
 77: #endif
 78:   } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Invalid Viewer : Run with -binary or -hdf5 option");
 79:   VecView(global1, viewer);
 80:   PetscViewerDestroy(&viewer);

 82:   PetscPrintf(PETSC_COMM_WORLD, "Global vector written to temp file is \n");
 83:   VecView(global1, PETSC_VIEWER_STDOUT_WORLD);

 85:   if (flg2) {
 86:     DMDACreate2d(PETSC_COMM_WORLD, bx, by, stencil_type, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, &da2);
 87:   } else if (flg3) {
 88:     DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, stencil_type, M, N, P, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, 0, &da2);
 89:   } else {
 90:     DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da2);
 91:   }
 92:   DMSetFromOptions(da2);
 93:   DMSetUp(da2);

 95:   if (isbinary) {
 96:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer);
 97: #if defined(PETSC_HAVE_HDF5)
 98:   } else if (ishdf5) {
 99:     PetscViewerHDF5Open(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer);
100: #endif
101:   } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Invalid Viewer : Run with -binary or -hdf5 option");

103:   DMCreateGlobalVector(da2, &global2);
104:   PetscObjectSetName((PetscObject)global2, "Test_Vec");
105:   VecLoad(global2, viewer);
106:   PetscViewerDestroy(&viewer);

108:   PetscPrintf(PETSC_COMM_WORLD, "Global vector read from temp file is \n");
109:   VecView(global2, PETSC_VIEWER_STDOUT_WORLD);
110:   VecAXPY(global2, mone, global1);
111:   VecNorm(global2, NORM_MAX, &norm);
112:   if (norm != 0.0) {
113:     PetscPrintf(PETSC_COMM_WORLD, "ex23: Norm of difference %g should be zero\n", (double)norm);
114:     PetscPrintf(PETSC_COMM_WORLD, "  Number of processors %d\n", size);
115:     PetscPrintf(PETSC_COMM_WORLD, "  M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", M, N, P, dof);
116:     PetscPrintf(PETSC_COMM_WORLD, "  stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n", stencil_width, (int)stencil_type, (int)pt);
117:     PetscPrintf(PETSC_COMM_WORLD, "  dimension %d\n", 1 + (int)flg2 + (int)flg3);
118:   }

120:   PetscRandomDestroy(&rdm);
121:   DMDestroy(&da);
122:   DMDestroy(&da2);
123:   VecDestroy(&global1);
124:   VecDestroy(&global2);
125:   PetscFinalize();
126:   return 0;
127: }