1c4762a1bSJed Brown static char help[] = "Tests I/O of vector and string attribute for HDF5 format\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscvec.h> 4c4762a1bSJed Brown #include <petscviewerhdf5.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown int main(int argc,char **args) 7c4762a1bSJed Brown { 8c4762a1bSJed Brown PetscErrorCode ierr; 9c4762a1bSJed Brown Vec u; 10c4762a1bSJed Brown PetscViewer viewer; 11c4762a1bSJed Brown char *attrReadVal, attrWriteVal[20]={"Hello World!!"}; 12c4762a1bSJed Brown 13c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 14c4762a1bSJed Brown 15c4762a1bSJed Brown /* PART 1: Generate vector, then write it in the given data format */ 16*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecCreate(PETSC_COMM_WORLD,&u)); 17*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscObjectSetName((PetscObject)u, "Test_Vec")); 18*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecSetSizes(u,PETSC_DECIDE,10)); 19*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecSetFromOptions(u)); 20*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecSet(u,0.)); 21c4762a1bSJed Brown 22c4762a1bSJed Brown /* write vector and attribute*/ 23*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5Open(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_WRITE,&viewer)); 24*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecView(u,viewer)); 25*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"Attribute value written: '%s'\n\n",attrWriteVal)); 26*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5WriteAttribute(viewer,"Test_Vec","Test_Attr",PETSC_STRING,attrWriteVal)); 27c4762a1bSJed Brown 28*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 29*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&u)); 30c4762a1bSJed Brown 31c4762a1bSJed Brown /* PART 2: Read in attribute */ 32*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5Open(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_READ,&viewer)); 33*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5ReadAttribute(viewer,"Test_Vec","Test_Attr",PETSC_STRING,NULL,&attrReadVal)); 34*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"Attribute value read: '%s'\n\n",attrReadVal)); 35*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscFree(attrReadVal)); 36c4762a1bSJed Brown 37*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 38c4762a1bSJed Brown ierr = PetscFinalize(); 39c4762a1bSJed Brown return ierr; 40c4762a1bSJed Brown } 41c4762a1bSJed Brown 42c4762a1bSJed Brown /*TEST 43c4762a1bSJed Brown 44c4762a1bSJed Brown build: 45c4762a1bSJed Brown requires: hdf5 46c4762a1bSJed Brown 47c4762a1bSJed Brown test: 48c4762a1bSJed Brown 49c4762a1bSJed Brown TEST*/ 50