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