Actual source code: ex23.c


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

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

  7: int main(int argc, char **argv)
  8: {
  9:   PetscMPIInt     size;
 10:   PetscInt        N = 6, m = PETSC_DECIDE, n = PETSC_DECIDE, p = PETSC_DECIDE, M = 8, dof = 1, stencil_width = 1, P = 5, pt = 0, st = 0;
 11:   PetscBool       flg2, flg3, native = PETSC_FALSE;
 12:   DMBoundaryType  bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE, bz = DM_BOUNDARY_NONE;
 13:   DMDAStencilType stencil_type = DMDA_STENCIL_STAR;
 14:   DM              da;
 15:   Vec             global1, global2, global3, global4;
 16:   PetscScalar     mone = -1.0;
 17:   PetscReal       norm;
 18:   PetscViewer     viewer;
 19:   PetscRandom     rdm;

 22:   PetscInitialize(&argc, &argv, (char *)0, help);
 23:   PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL);
 24:   PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL);
 25:   PetscOptionsGetInt(NULL, NULL, "-P", &P, NULL);
 26:   PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL);
 27:   PetscOptionsGetInt(NULL, NULL, "-stencil_width", &stencil_width, NULL);
 28:   PetscOptionsGetInt(NULL, NULL, "-periodic", &pt, NULL);
 29:   PetscOptionsGetBool(NULL, NULL, "-native", &native, NULL);
 30:   if (pt == 1) bx = DM_BOUNDARY_PERIODIC;
 31:   if (pt == 2) by = DM_BOUNDARY_PERIODIC;
 32:   if (pt == 3) {
 33:     bx = DM_BOUNDARY_PERIODIC;
 34:     by = DM_BOUNDARY_PERIODIC;
 35:   }
 36:   if (pt == 4) bz = DM_BOUNDARY_PERIODIC;

 38:   PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL);
 39:   stencil_type = (DMDAStencilType)st;

 41:   PetscOptionsHasName(NULL, NULL, "-one", &flg2);
 42:   PetscOptionsHasName(NULL, NULL, "-two", &flg2);
 43:   PetscOptionsHasName(NULL, NULL, "-three", &flg3);
 44:   if (flg2) {
 45:     DMDACreate2d(PETSC_COMM_WORLD, bx, by, stencil_type, M, N, m, n, dof, stencil_width, 0, 0, &da);
 46:   } else if (flg3) {
 47:     DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, stencil_type, M, N, P, m, n, p, dof, stencil_width, 0, 0, 0, &da);
 48:   } else {
 49:     DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da);
 50:   }
 51:   DMSetFromOptions(da);
 52:   DMSetUp(da);

 54:   DMCreateGlobalVector(da, &global1);
 55:   PetscRandomCreate(PETSC_COMM_WORLD, &rdm);
 56:   PetscRandomSetFromOptions(rdm);
 57:   DMCreateGlobalVector(da, &global2);
 58:   DMCreateGlobalVector(da, &global3);
 59:   DMCreateGlobalVector(da, &global4);

 61:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer);
 62:   if (native) PetscViewerPushFormat(viewer, PETSC_VIEWER_NATIVE);
 63:   VecSetRandom(global1, rdm);
 64:   VecView(global1, viewer);
 65:   VecSetRandom(global3, rdm);
 66:   VecView(global3, viewer);
 67:   if (native) PetscViewerPopFormat(viewer);
 68:   PetscViewerDestroy(&viewer);

 70:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer);
 71:   if (native) PetscViewerPushFormat(viewer, PETSC_VIEWER_NATIVE);
 72:   VecLoad(global2, viewer);
 73:   VecLoad(global4, viewer);
 74:   if (native) PetscViewerPopFormat(viewer);
 75:   PetscViewerDestroy(&viewer);

 77:   if (native) {
 78:     Vec       filenative;
 79:     PetscBool same;
 80:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer);
 81:     DMDACreateNaturalVector(da, &filenative);
 82:     /* DMDA "natural" Vec does not commandeer VecLoad.  The following load will only work when run on the same process
 83:      * layout, where as the standard VecView/VecLoad (using DMDA and not PETSC_VIEWER_NATIVE) can be read on a different
 84:      * number of processors. */
 85:     VecLoad(filenative, viewer);
 86:     VecEqual(global2, filenative, &same);
 87:     if (!same) {
 88:       PetscPrintf(PETSC_COMM_WORLD, "ex23: global vector does not match contents of file\n");
 89:       VecView(global2, 0);
 90:       VecView(filenative, 0);
 91:     }
 92:     PetscViewerDestroy(&viewer);
 93:     VecDestroy(&filenative);
 94:   }

 96:   VecAXPY(global2, mone, global1);
 97:   VecNorm(global2, NORM_MAX, &norm);
 98:   if (norm != 0.0) {
 99:     MPI_Comm_size(PETSC_COMM_WORLD, &size);
100:     PetscPrintf(PETSC_COMM_WORLD, "ex23: Norm of difference %g should be zero\n", (double)norm);
101:     PetscPrintf(PETSC_COMM_WORLD, "  Number of processors %d\n", size);
102:     PetscPrintf(PETSC_COMM_WORLD, "  M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", M, N, P, dof);
103:     PetscPrintf(PETSC_COMM_WORLD, "  stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n", stencil_width, (int)stencil_type, (int)bx);
104:     PetscPrintf(PETSC_COMM_WORLD, "  dimension %d\n", 1 + (int)flg2 + (int)flg3);
105:   }
106:   VecAXPY(global4, mone, global3);
107:   VecNorm(global4, NORM_MAX, &norm);
108:   if (norm != 0.0) {
109:     MPI_Comm_size(PETSC_COMM_WORLD, &size);
110:     PetscPrintf(PETSC_COMM_WORLD, "ex23: Norm of difference %g should be zero\n", (double)norm);
111:     PetscPrintf(PETSC_COMM_WORLD, "  Number of processors %d\n", size);
112:     PetscPrintf(PETSC_COMM_WORLD, "  M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", M, N, P, dof);
113:     PetscPrintf(PETSC_COMM_WORLD, "  stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n", stencil_width, (int)stencil_type, (int)bx);
114:     PetscPrintf(PETSC_COMM_WORLD, "  dimension %d\n", 1 + (int)flg2 + (int)flg3);
115:   }

117:   PetscRandomDestroy(&rdm);
118:   DMDestroy(&da);
119:   VecDestroy(&global1);
120:   VecDestroy(&global2);
121:   VecDestroy(&global3);
122:   VecDestroy(&global4);
123:   PetscFinalize();
124:   return 0;
125: }

127: /*TEST

129:    test:
130:       nsize: {{1  3}}
131:       args: -one -dof {{1 2 3}} -stencil_type {{0 1}}

133:    test:
134:       suffix: 3
135:       nsize: {{2 4}}
136:       args: -two -dof {{1 3}} -stencil_type {{0 1}}

138:    test:
139:       suffix: 4
140:       nsize: {{1 4}}
141:       args: -three -dof {{2 3}} -stencil_type {{0 1}}

143:    test:
144:       suffix: 2
145:       nsize: 2
146:       args: -two -native

148: TEST*/