Actual source code: ex33.c
1: static char help[] = "Tests VecView()/VecLoad() for DMDA vectors (this tests DMDAGlobalToNatural()).\n\n";
3: #include <petscdm.h>
4: #include <petscdmda.h>
5: #include <petscviewerhdf5.h>
7: int main(int argc, char **argv)
8: {
9: PetscMPIInt rank, size;
10: PetscInt N = 6, M = 8, P = 5, dof = 1;
11: PetscInt stencil_width = 1, pt = 0, st = 0;
12: PetscBool flg2, flg3, isbinary, mpiio;
13: DMBoundaryType bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE, bz = DM_BOUNDARY_NONE;
14: DMDAStencilType stencil_type = DMDA_STENCIL_STAR;
15: DM da, da2;
16: Vec global1, global2;
17: PetscScalar mone = -1.0;
18: PetscReal norm;
19: PetscViewer viewer;
20: PetscRandom rdm;
21: #if defined(PETSC_HAVE_HDF5)
22: PetscBool ishdf5;
23: #endif
25: PetscFunctionBeginUser;
26: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
27: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
28: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
30: PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
31: PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
32: PetscCall(PetscOptionsGetInt(NULL, NULL, "-P", &P, NULL));
33: PetscCall(PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL));
34: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_width", &stencil_width, NULL));
35: PetscCall(PetscOptionsGetInt(NULL, NULL, "-periodic", &pt, NULL));
36: if (pt == 1) bx = DM_BOUNDARY_PERIODIC;
37: if (pt == 2) by = DM_BOUNDARY_PERIODIC;
38: if (pt == 4) {
39: bx = DM_BOUNDARY_PERIODIC;
40: by = DM_BOUNDARY_PERIODIC;
41: }
43: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL));
44: stencil_type = (DMDAStencilType)st;
46: PetscCall(PetscOptionsHasName(NULL, NULL, "-oned", &flg2));
47: PetscCall(PetscOptionsHasName(NULL, NULL, "-twod", &flg2));
48: PetscCall(PetscOptionsHasName(NULL, NULL, "-threed", &flg3));
50: PetscCall(PetscOptionsHasName(NULL, NULL, "-binary", &isbinary));
51: #if defined(PETSC_HAVE_HDF5)
52: PetscCall(PetscOptionsHasName(NULL, NULL, "-hdf5", &ishdf5));
53: #endif
54: PetscCall(PetscOptionsHasName(NULL, NULL, "-mpiio", &mpiio));
55: if (flg2) {
56: PetscCall(DMDACreate2d(PETSC_COMM_WORLD, bx, by, stencil_type, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, &da));
57: } else if (flg3) {
58: 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));
59: } else {
60: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da));
61: }
62: PetscCall(DMSetFromOptions(da));
63: PetscCall(DMSetUp(da));
65: PetscCall(DMCreateGlobalVector(da, &global1));
66: PetscCall(PetscObjectSetName((PetscObject)global1, "Test_Vec"));
67: PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rdm));
68: PetscCall(PetscRandomSetFromOptions(rdm));
69: PetscCall(VecSetRandom(global1, rdm));
70: if (isbinary) {
71: if (mpiio) PetscCall(PetscOptionsSetValue(NULL, "-viewer_binary_mpiio", ""));
72: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer));
73: #if defined(PETSC_HAVE_HDF5)
74: } else if (ishdf5) {
75: PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "temp", FILE_MODE_WRITE, &viewer));
76: #endif
77: } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Invalid Viewer : Run with -binary or -hdf5 option");
78: PetscCall(VecView(global1, viewer));
79: PetscCall(PetscViewerDestroy(&viewer));
81: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Global vector written to temp file is \n"));
82: PetscCall(VecView(global1, PETSC_VIEWER_STDOUT_WORLD));
84: if (flg2) {
85: PetscCall(DMDACreate2d(PETSC_COMM_WORLD, bx, by, stencil_type, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, stencil_width, 0, 0, &da2));
86: } else if (flg3) {
87: 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));
88: } else {
89: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, bx, M, dof, stencil_width, NULL, &da2));
90: }
91: PetscCall(DMSetFromOptions(da2));
92: PetscCall(DMSetUp(da2));
94: if (isbinary) {
95: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer));
96: #if defined(PETSC_HAVE_HDF5)
97: } else if (ishdf5) {
98: PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "temp", FILE_MODE_READ, &viewer));
99: #endif
100: } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Invalid Viewer : Run with -binary or -hdf5 option");
102: PetscCall(DMCreateGlobalVector(da2, &global2));
103: PetscCall(PetscObjectSetName((PetscObject)global2, "Test_Vec"));
104: PetscCall(VecLoad(global2, viewer));
105: PetscCall(PetscViewerDestroy(&viewer));
107: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Global vector read from temp file is \n"));
108: PetscCall(VecView(global2, PETSC_VIEWER_STDOUT_WORLD));
109: PetscCall(VecAXPY(global2, mone, global1));
110: PetscCall(VecNorm(global2, NORM_MAX, &norm));
111: if (norm != 0.0) {
112: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "ex23: Norm of difference %g should be zero\n", (double)norm));
113: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " Number of processors %d\n", size));
114: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", M, N, P, dof));
115: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n", stencil_width, (int)stencil_type, (int)pt));
116: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " dimension %d\n", 1 + (int)flg2 + (int)flg3));
117: }
119: PetscCall(PetscRandomDestroy(&rdm));
120: PetscCall(DMDestroy(&da));
121: PetscCall(DMDestroy(&da2));
122: PetscCall(VecDestroy(&global1));
123: PetscCall(VecDestroy(&global2));
124: PetscCall(PetscFinalize());
125: return 0;
126: }