xref: /petsc/src/dm/tutorials/ex9.c (revision 21e3ffae2f3b73c0bd738cf6d0a809700fc04bb0)
1 static char help[] = "Demonstrates HDF5 vector input/output\n\n";
2 
3 #include <petscsys.h>
4 #include <petscdm.h>
5 #include <petscdmda.h>
6 #include <petscviewerhdf5.h>
7 
8 int main(int argc, char **argv)
9 {
10   PetscViewer viewer;
11   DM          da;
12   Vec         global, local, global2;
13   PetscMPIInt rank;
14   PetscBool   flg;
15   PetscInt    ndof;
16 
17   /*
18     Every PETSc routine should begin with the PetscInitialize() routine.
19     argc, argv - These command line arguments are taken to extract the options
20                  supplied to PETSc and options supplied to MPI.
21     help       - When PETSc executable is invoked with the option -help,
22                  it prints the various options that can be applied at
23                  runtime.  The user can use the "help" variable place
24                  additional help messages in this printout.
25   */
26   PetscFunctionBeginUser;
27   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
28   /* Get number of DOF's from command line */
29   PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "DMDA VecView/VecLoad example", "");
30   {
31     ndof = 1;
32     PetscCall(PetscOptionsBoundedInt("-ndof", "Number of DOF's in DMDA", "", ndof, &ndof, NULL, 1));
33   }
34   PetscOptionsEnd();
35 
36   /* Create a DMDA and an associated vector */
37   PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 100, 90, PETSC_DECIDE, PETSC_DECIDE, ndof, 1, NULL, NULL, &da));
38   PetscCall(DMSetFromOptions(da));
39   PetscCall(DMSetUp(da));
40   PetscCall(DMCreateGlobalVector(da, &global));
41   PetscCall(DMCreateLocalVector(da, &local));
42   PetscCall(VecSet(global, -1.0));
43   PetscCall(DMGlobalToLocalBegin(da, global, INSERT_VALUES, local));
44   PetscCall(DMGlobalToLocalEnd(da, global, INSERT_VALUES, local));
45   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
46   PetscCall(VecScale(local, rank + 1));
47   PetscCall(DMLocalToGlobalBegin(da, local, ADD_VALUES, global));
48   PetscCall(DMLocalToGlobalEnd(da, local, ADD_VALUES, global));
49 
50   /* Create the HDF5 viewer for writing */
51   PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "hdf5output.h5", FILE_MODE_WRITE, &viewer));
52   PetscCall(PetscViewerSetFromOptions(viewer));
53 
54   /* Write the Vec without one extra dimension for BS */
55   PetscCall(PetscViewerHDF5SetBaseDimension2(viewer, PETSC_FALSE));
56   PetscCall(PetscObjectSetName((PetscObject)global, "noBsDim"));
57   PetscCall(VecView(global, viewer));
58 
59   /* Write the Vec with one extra, 1-sized, dimension for BS */
60   PetscCall(PetscViewerHDF5SetBaseDimension2(viewer, PETSC_TRUE));
61   PetscCall(PetscObjectSetName((PetscObject)global, "bsDim"));
62   PetscCall(VecView(global, viewer));
63 
64   PetscCall(PetscViewerDestroy(&viewer));
65   PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD));
66   PetscCall(VecDuplicate(global, &global2));
67 
68   /* Create the HDF5 viewer for reading */
69   PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "hdf5output.h5", FILE_MODE_READ, &viewer));
70   PetscCall(PetscViewerSetFromOptions(viewer));
71 
72   /* Load the Vec without the BS dim and compare */
73   PetscCall(PetscObjectSetName((PetscObject)global2, "noBsDim"));
74   PetscCall(VecLoad(global2, viewer));
75 
76   PetscCall(VecEqual(global, global2, &flg));
77   if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error: Vectors are not equal\n"));
78 
79   /* Load the Vec with one extra, 1-sized, BS dim and compare */
80   PetscCall(PetscObjectSetName((PetscObject)global2, "bsDim"));
81   PetscCall(VecLoad(global2, viewer));
82   PetscCall(PetscViewerDestroy(&viewer));
83 
84   PetscCall(VecEqual(global, global2, &flg));
85   if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error: Vectors are not equal\n"));
86 
87   /* clean up and exit */
88   PetscCall(VecDestroy(&local));
89   PetscCall(VecDestroy(&global));
90   PetscCall(VecDestroy(&global2));
91   PetscCall(DMDestroy(&da));
92 
93   PetscCall(PetscFinalize());
94   return 0;
95 }
96 
97 /*TEST
98 
99       build:
100          requires: hdf5
101 
102       test:
103          nsize: 4
104 
105       test:
106          nsize: 4
107          suffix: 2
108          args: -ndof 2
109          output_file: output/ex9_1.out
110 
111 TEST*/
112