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
26: PetscFunctionBeginUser;
27: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
28: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
29: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
31: PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
32: PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
33: PetscCall(PetscOptionsGetInt(NULL, NULL, "-P", &P, NULL));
34: PetscCall(PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL));
35: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_width", &stencil_width, NULL));
36: PetscCall(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: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL));
45: stencil_type = (DMDAStencilType)st;
47: PetscCall(PetscOptionsHasName(NULL, NULL, "-oned", &flg2));
48: PetscCall(PetscOptionsHasName(NULL, NULL, "-twod", &flg2));
49: PetscCall(PetscOptionsHasName(NULL, NULL, "-threed", &flg3));
51: PetscCall(PetscOptionsHasName(NULL, NULL, "-binary", &isbinary));
52: #if defined(PETSC_HAVE_HDF5)
53: PetscCall(PetscOptionsHasName(NULL, NULL, "-hdf5", &ishdf5));
54: #endif
55: PetscCall(PetscOptionsHasName(NULL, NULL, "-mpiio", &mpiio));
56: if (flg2) {
57: PetscCall(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: PetscCall(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: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da));
62: }
63: PetscCall(DMSetFromOptions(da));
64: PetscCall(DMSetUp(da));
66: PetscCall(DMCreateGlobalVector(da, &global1));
67: PetscCall(PetscObjectSetName((PetscObject)global1, "Test_Vec"));
68: PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rdm));
69: PetscCall(PetscRandomSetFromOptions(rdm));
70: PetscCall(VecSetRandom(global1, rdm));
71: if (isbinary) {
72: if (mpiio) PetscCall(PetscOptionsSetValue(NULL, "-viewer_binary_mpiio", ""));
73: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer));
74: #if defined(PETSC_HAVE_HDF5)
75: } else if (ishdf5) {
76: PetscCall(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: PetscCall(VecView(global1, viewer));
80: PetscCall(PetscViewerDestroy(&viewer));
82: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Global vector written to temp file is \n"));
83: PetscCall(VecView(global1, PETSC_VIEWER_STDOUT_WORLD));
85: if (flg2) {
86: PetscCall(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: PetscCall(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: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da2));
91: }
92: PetscCall(DMSetFromOptions(da2));
93: PetscCall(DMSetUp(da2));
95: if (isbinary) {
96: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer));
97: #if defined(PETSC_HAVE_HDF5)
98: } else if (ishdf5) {
99: PetscCall(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: PetscCall(DMCreateGlobalVector(da2, &global2));
104: PetscCall(PetscObjectSetName((PetscObject)global2, "Test_Vec"));
105: PetscCall(VecLoad(global2, viewer));
106: PetscCall(PetscViewerDestroy(&viewer));
108: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Global vector read from temp file is \n"));
109: PetscCall(VecView(global2, PETSC_VIEWER_STDOUT_WORLD));
110: PetscCall(VecAXPY(global2, mone, global1));
111: PetscCall(VecNorm(global2, NORM_MAX, &norm));
112: if (norm != 0.0) {
113: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "ex23: Norm of difference %g should be zero\n", (double)norm));
114: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " Number of processors %d\n", size));
115: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", M, N, P, dof));
116: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n", stencil_width, (int)stencil_type, (int)pt));
117: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " dimension %d\n", 1 + (int)flg2 + (int)flg3));
118: }
120: PetscCall(PetscRandomDestroy(&rdm));
121: PetscCall(DMDestroy(&da));
122: PetscCall(DMDestroy(&da2));
123: PetscCall(VecDestroy(&global1));
124: PetscCall(VecDestroy(&global2));
125: PetscCall(PetscFinalize());
126: return 0;
127: }