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 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 33 /* Get number of DOF's from command line */ 34 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"DMDA VecView/VecLoad example","");CHKERRQ(ierr); 35 { 36 ndof = 1; 37 PetscOptionsBoundedInt("-ndof","Number of DOF's in DMDA","",ndof,&ndof,NULL,1); 38 } 39 ierr = PetscOptionsEnd();CHKERRQ(ierr); 40 41 /* Create a DMDA and an associated vector */ 42 ierr = DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,100,90,PETSC_DECIDE,PETSC_DECIDE,ndof,1,NULL,NULL,&da);CHKERRQ(ierr); 43 ierr = DMSetFromOptions(da);CHKERRQ(ierr); 44 ierr = DMSetUp(da);CHKERRQ(ierr); 45 ierr = DMCreateGlobalVector(da,&global);CHKERRQ(ierr); 46 ierr = DMCreateLocalVector(da,&local);CHKERRQ(ierr); 47 ierr = VecSet(global,-1.0);CHKERRQ(ierr); 48 ierr = DMGlobalToLocalBegin(da,global,INSERT_VALUES,local);CHKERRQ(ierr); 49 ierr = DMGlobalToLocalEnd(da,global,INSERT_VALUES,local);CHKERRQ(ierr); 50 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 51 ierr = VecScale(local,rank+1);CHKERRQ(ierr); 52 ierr = DMLocalToGlobalBegin(da,local,ADD_VALUES,global);CHKERRQ(ierr); 53 ierr = DMLocalToGlobalEnd(da,local,ADD_VALUES,global);CHKERRQ(ierr); 54 55 /* Create the HDF5 viewer for writing */ 56 ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD,"hdf5output.h5",FILE_MODE_WRITE,&viewer);CHKERRQ(ierr); 57 ierr = PetscViewerSetFromOptions(viewer);CHKERRQ(ierr); 58 59 /* Write the Vec without one extra dimension for BS */ 60 ierr = PetscViewerHDF5SetBaseDimension2(viewer, PETSC_FALSE); 61 ierr = PetscObjectSetName((PetscObject) global, "noBsDim");CHKERRQ(ierr); 62 ierr = VecView(global,viewer);CHKERRQ(ierr); 63 64 /* Write the Vec with one extra, 1-sized, dimension for BS */ 65 ierr = PetscViewerHDF5SetBaseDimension2(viewer, PETSC_TRUE); 66 ierr = PetscObjectSetName((PetscObject) global, "bsDim");CHKERRQ(ierr); 67 ierr = VecView(global,viewer);CHKERRQ(ierr); 68 69 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 70 ierr = MPI_Barrier(PETSC_COMM_WORLD);CHKERRQ(ierr); 71 ierr = VecDuplicate(global,&global2);CHKERRQ(ierr); 72 73 /* Create the HDF5 viewer for reading */ 74 ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD,"hdf5output.h5",FILE_MODE_READ,&viewer);CHKERRQ(ierr); 75 ierr = PetscViewerSetFromOptions(viewer);CHKERRQ(ierr); 76 77 /* Load the Vec without the BS dim and compare */ 78 ierr = PetscObjectSetName((PetscObject) global2, "noBsDim");CHKERRQ(ierr); 79 ierr = VecLoad(global2,viewer);CHKERRQ(ierr); 80 81 ierr = VecEqual(global,global2,&flg);CHKERRQ(ierr); 82 if (!flg) { 83 ierr = PetscPrintf(PETSC_COMM_WORLD,"Error: Vectors are not equal\n");CHKERRQ(ierr); 84 } 85 86 /* Load the Vec with one extra, 1-sized, BS dim and compare */ 87 ierr = PetscObjectSetName((PetscObject) global2, "bsDim");CHKERRQ(ierr); 88 ierr = VecLoad(global2,viewer);CHKERRQ(ierr); 89 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 90 91 ierr = VecEqual(global,global2,&flg);CHKERRQ(ierr); 92 if (!flg) { 93 ierr = PetscPrintf(PETSC_COMM_WORLD,"Error: Vectors are not equal\n");CHKERRQ(ierr); 94 } 95 96 /* clean up and exit */ 97 ierr = VecDestroy(&local);CHKERRQ(ierr); 98 ierr = VecDestroy(&global);CHKERRQ(ierr); 99 ierr = VecDestroy(&global2);CHKERRQ(ierr); 100 ierr = DMDestroy(&da);CHKERRQ(ierr); 101 102 ierr = PetscFinalize(); 103 return ierr; 104 } 105 106 107 /*TEST 108 109 build: 110 requires: hdf5 111 112 test: 113 nsize: 4 114 115 test: 116 nsize: 4 117 suffix: 2 118 args: -ndof 2 119 output_file: output/ex9_1.out 120 121 TEST*/ 122