1 #include <petsc-private/viewerimpl.h> /*I "petscsys.h" I*/ 2 #include <petscviewerhdf5.h> /*I "petscviewerhdf5.h" I*/ 3 4 typedef struct GroupList { 5 const char *name; 6 struct GroupList *next; 7 } GroupList; 8 9 typedef struct { 10 char *filename; 11 PetscFileMode btype; 12 hid_t file_id; 13 PetscInt timestep; 14 GroupList *groups; 15 } PetscViewer_HDF5; 16 17 #undef __FUNCT__ 18 #define __FUNCT__ "PetscViewerFileClose_HDF5" 19 static PetscErrorCode PetscViewerFileClose_HDF5(PetscViewer viewer) 20 { 21 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*)viewer->data; 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 ierr = PetscFree(hdf5->filename);CHKERRQ(ierr); 26 if (hdf5->file_id) H5Fclose(hdf5->file_id); 27 PetscFunctionReturn(0); 28 } 29 30 #undef __FUNCT__ 31 #define __FUNCT__ "PetscViewerDestroy_HDF5" 32 PetscErrorCode PetscViewerDestroy_HDF5(PetscViewer viewer) 33 { 34 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 35 PetscErrorCode ierr; 36 37 PetscFunctionBegin; 38 ierr = PetscViewerFileClose_HDF5(viewer);CHKERRQ(ierr); 39 if (hdf5->groups) 40 while (hdf5->groups) { 41 GroupList *tmp = hdf5->groups->next; 42 43 ierr = PetscFree(hdf5->groups->name);CHKERRQ(ierr); 44 ierr = PetscFree(hdf5->groups);CHKERRQ(ierr); 45 hdf5->groups = tmp; 46 } 47 ierr = PetscFree(hdf5);CHKERRQ(ierr); 48 PetscFunctionReturn(0); 49 } 50 51 #undef __FUNCT__ 52 #define __FUNCT__ "PetscViewerFileSetMode_HDF5" 53 PetscErrorCode PetscViewerFileSetMode_HDF5(PetscViewer viewer, PetscFileMode type) 54 { 55 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 56 57 PetscFunctionBegin; 58 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 1); 59 hdf5->btype = type; 60 PetscFunctionReturn(0); 61 } 62 63 #undef __FUNCT__ 64 #define __FUNCT__ "PetscViewerFileSetName_HDF5" 65 PetscErrorCode PetscViewerFileSetName_HDF5(PetscViewer viewer, const char name[]) 66 { 67 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 68 #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO) 69 MPI_Info info = MPI_INFO_NULL; 70 #endif 71 hid_t plist_id; 72 herr_t herr; 73 PetscErrorCode ierr; 74 75 PetscFunctionBegin; 76 ierr = PetscStrallocpy(name, &hdf5->filename);CHKERRQ(ierr); 77 /* Set up file access property list with parallel I/O access */ 78 plist_id = H5Pcreate(H5P_FILE_ACCESS); 79 #if defined(PETSC_HAVE_H5PSET_FAPL_MPIO) 80 herr = H5Pset_fapl_mpio(plist_id, PetscObjectComm((PetscObject)viewer), info);CHKERRQ(herr); 81 #endif 82 /* Create or open the file collectively */ 83 switch (hdf5->btype) { 84 case FILE_MODE_READ: 85 hdf5->file_id = H5Fopen(name, H5F_ACC_RDONLY, plist_id); 86 break; 87 case FILE_MODE_APPEND: 88 hdf5->file_id = H5Fopen(name, H5F_ACC_RDWR, plist_id); 89 break; 90 case FILE_MODE_WRITE: 91 hdf5->file_id = H5Fcreate(name, H5F_ACC_TRUNC, H5P_DEFAULT, plist_id); 92 break; 93 default: 94 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()"); 95 } 96 if (hdf5->file_id < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB, "H5Fcreate failed for %s", name); 97 H5Pclose(plist_id); 98 PetscFunctionReturn(0); 99 } 100 101 #undef __FUNCT__ 102 #define __FUNCT__ "PetscViewerCreate_HDF5" 103 PETSC_EXTERN PetscErrorCode PetscViewerCreate_HDF5(PetscViewer v) 104 { 105 PetscViewer_HDF5 *hdf5; 106 PetscErrorCode ierr; 107 108 PetscFunctionBegin; 109 ierr = PetscNewLog(v, PetscViewer_HDF5, &hdf5);CHKERRQ(ierr); 110 111 v->data = (void*) hdf5; 112 v->ops->destroy = PetscViewerDestroy_HDF5; 113 v->ops->flush = 0; 114 hdf5->btype = (PetscFileMode) -1; 115 hdf5->filename = 0; 116 hdf5->timestep = -1; 117 hdf5->groups = NULL; 118 119 ierr = PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetName_C","PetscViewerFileSetName_HDF5",PetscViewerFileSetName_HDF5);CHKERRQ(ierr); 120 ierr = PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetMode_C","PetscViewerFileSetMode_HDF5",PetscViewerFileSetMode_HDF5);CHKERRQ(ierr); 121 PetscFunctionReturn(0); 122 } 123 124 #undef __FUNCT__ 125 #define __FUNCT__ "PetscViewerHDF5Open" 126 /*@C 127 PetscViewerHDF5Open - Opens a file for HDF5 input/output. 128 129 Collective on MPI_Comm 130 131 Input Parameters: 132 + comm - MPI communicator 133 . name - name of file 134 - type - type of file 135 $ FILE_MODE_WRITE - create new file for binary output 136 $ FILE_MODE_READ - open existing file for binary input 137 $ FILE_MODE_APPEND - open existing file for binary output 138 139 Output Parameter: 140 . hdf5v - PetscViewer for HDF5 input/output to use with the specified file 141 142 Level: beginner 143 144 Note: 145 This PetscViewer should be destroyed with PetscViewerDestroy(). 146 147 Concepts: HDF5 files 148 Concepts: PetscViewerHDF5^creating 149 150 .seealso: PetscViewerASCIIOpen(), PetscViewerSetFormat(), PetscViewerDestroy(), 151 VecView(), MatView(), VecLoad(), MatLoad(), 152 PetscFileMode, PetscViewer 153 @*/ 154 PetscErrorCode PetscViewerHDF5Open(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *hdf5v) 155 { 156 PetscErrorCode ierr; 157 158 PetscFunctionBegin; 159 ierr = PetscViewerCreate(comm, hdf5v);CHKERRQ(ierr); 160 ierr = PetscViewerSetType(*hdf5v, PETSCVIEWERHDF5);CHKERRQ(ierr); 161 ierr = PetscViewerFileSetMode(*hdf5v, type);CHKERRQ(ierr); 162 ierr = PetscViewerFileSetName(*hdf5v, name);CHKERRQ(ierr); 163 PetscFunctionReturn(0); 164 } 165 166 #undef __FUNCT__ 167 #define __FUNCT__ "PetscViewerHDF5GetFileId" 168 /*@C 169 PetscViewerHDF5GetFileId - Retrieve the file id, this file ID then can be used in direct HDF5 calls 170 171 Not collective 172 173 Input Parameter: 174 . viewer - the PetscViewer 175 176 Output Parameter: 177 . file_id - The file id 178 179 Level: intermediate 180 181 .seealso: PetscViewerHDF5Open() 182 @*/ 183 PetscErrorCode PetscViewerHDF5GetFileId(PetscViewer viewer, hid_t *file_id) 184 { 185 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 186 187 PetscFunctionBegin; 188 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 189 if (file_id) *file_id = hdf5->file_id; 190 PetscFunctionReturn(0); 191 } 192 193 #undef __FUNCT__ 194 #define __FUNCT__ "PetscViewerHDF5PushGroup" 195 /*@C 196 PetscViewerHDF5PushGroup - Set the current HDF5 group for output 197 198 Not collective 199 200 Input Parameters: 201 + viewer - the PetscViewer 202 - name - The group name 203 204 Level: intermediate 205 206 .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PopGroup(),PetscViewerHDF5GetGroup() 207 @*/ 208 PetscErrorCode PetscViewerHDF5PushGroup(PetscViewer viewer, const char *name) 209 { 210 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 211 GroupList *groupNode; 212 PetscErrorCode ierr; 213 214 PetscFunctionBegin; 215 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 216 PetscValidCharPointer(name,2); 217 ierr = PetscMalloc(sizeof(GroupList), &groupNode);CHKERRQ(ierr); 218 ierr = PetscStrallocpy(name, (char**) &groupNode->name);CHKERRQ(ierr); 219 220 groupNode->next = hdf5->groups; 221 hdf5->groups = groupNode; 222 PetscFunctionReturn(0); 223 } 224 225 #undef __FUNCT__ 226 #define __FUNCT__ "PetscViewerHDF5PopGroup" 227 /*@C 228 PetscViewerHDF5PopGroup - Return the current HDF5 group for output to the previous value 229 230 Not collective 231 232 Input Parameter: 233 . viewer - the PetscViewer 234 235 Level: intermediate 236 237 .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5GetGroup() 238 @*/ 239 PetscErrorCode PetscViewerHDF5PopGroup(PetscViewer viewer) 240 { 241 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 242 GroupList *groupNode; 243 PetscErrorCode ierr; 244 245 PetscFunctionBegin; 246 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 247 if (!hdf5->groups) SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_WRONGSTATE, "HDF5 group stack is empty, cannot pop"); 248 groupNode = hdf5->groups; 249 hdf5->groups = hdf5->groups->next; 250 ierr = PetscFree(groupNode->name);CHKERRQ(ierr); 251 ierr = PetscFree(groupNode);CHKERRQ(ierr); 252 PetscFunctionReturn(0); 253 } 254 255 #undef __FUNCT__ 256 #define __FUNCT__ "PetscViewerHDF5GetGroup" 257 /*@C 258 PetscViewerHDF5GetGroup - Get the current HDF5 group for output. If none has been assigned, returns NULL. 259 260 Not collective 261 262 Input Parameter: 263 . viewer - the PetscViewer 264 265 Output Parameter: 266 . name - The group name 267 268 Level: intermediate 269 270 .seealso: PetscViewerHDF5Open(),PetscViewerHDF5PushGroup(),PetscViewerHDF5PopGroup() 271 @*/ 272 PetscErrorCode PetscViewerHDF5GetGroup(PetscViewer viewer, const char **name) 273 { 274 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5 *) viewer->data; 275 276 PetscFunctionBegin; 277 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 278 PetscValidCharPointer(name,2); 279 if (hdf5->groups) *name = hdf5->groups->name; 280 else *name = NULL; 281 PetscFunctionReturn(0); 282 } 283 284 #undef __FUNCT__ 285 #define __FUNCT__ "PetscViewerHDF5IncrementTimestep" 286 /*@C 287 PetscViewerHDF5IncrementTimestep - Increments the current timestep for the HDF5 output. Fields are stacked in time. 288 289 Not collective 290 291 Input Parameter: 292 . viewer - the PetscViewer 293 294 Level: intermediate 295 296 .seealso: PetscViewerHDF5Open(), PetscViewerHDF5SetTimestep(), PetscViewerHDF5GetTimestep() 297 @*/ 298 PetscErrorCode PetscViewerHDF5IncrementTimestep(PetscViewer viewer) 299 { 300 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 301 302 PetscFunctionBegin; 303 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 304 ++hdf5->timestep; 305 PetscFunctionReturn(0); 306 } 307 308 #undef __FUNCT__ 309 #define __FUNCT__ "PetscViewerHDF5SetTimestep" 310 /*@C 311 PetscViewerHDF5SetTimestep - Set the current timestep for the HDF5 output. Fields are stacked in time. A timestep 312 of -1 disables blocking with timesteps. 313 314 Not collective 315 316 Input Parameters: 317 + viewer - the PetscViewer 318 - timestep - The timestep number 319 320 Level: intermediate 321 322 .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5GetTimestep() 323 @*/ 324 PetscErrorCode PetscViewerHDF5SetTimestep(PetscViewer viewer, PetscInt timestep) 325 { 326 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 327 328 PetscFunctionBegin; 329 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 330 hdf5->timestep = timestep; 331 PetscFunctionReturn(0); 332 } 333 334 #undef __FUNCT__ 335 #define __FUNCT__ "PetscViewerHDF5GetTimestep" 336 /*@C 337 PetscViewerHDF5GetTimestep - Get the current timestep for the HDF5 output. Fields are stacked in time. 338 339 Not collective 340 341 Input Parameter: 342 . viewer - the PetscViewer 343 344 Output Parameter: 345 . timestep - The timestep number 346 347 Level: intermediate 348 349 .seealso: PetscViewerHDF5Open(), PetscViewerHDF5IncrementTimestep(), PetscViewerHDF5SetTimestep() 350 @*/ 351 PetscErrorCode PetscViewerHDF5GetTimestep(PetscViewer viewer, PetscInt *timestep) 352 { 353 PetscViewer_HDF5 *hdf5 = (PetscViewer_HDF5*) viewer->data; 354 355 PetscFunctionBegin; 356 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,1); 357 PetscValidPointer(timestep,2); 358 *timestep = hdf5->timestep; 359 PetscFunctionReturn(0); 360 } 361 362 #if defined(oldhdf4stuff) 363 #undef __FUNCT__ 364 #define __FUNCT__ "PetscViewerHDF5WriteSDS" 365 PetscErrorCode PetscViewerHDF5WriteSDS(PetscViewer viewer, float *xf, int d, int *dims,int bs) 366 { 367 int i; 368 PetscViewer_HDF5 *vhdf5 = (PetscViewer_HDF5*)viewer->data; 369 int32 sds_id,zero32[3],dims32[3]; 370 371 PetscFunctionBegin; 372 for (i = 0; i < d; i++) { 373 zero32[i] = 0; 374 dims32[i] = dims[i]; 375 } 376 sds_id = SDcreate(vhdf5->sd_id, "Vec", DFNT_FLOAT32, d, dims32); 377 if (sds_id < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"SDcreate failed"); 378 SDwritedata(sds_id, zero32, 0, dims32, xf); 379 SDendaccess(sds_id); 380 PetscFunctionReturn(0); 381 } 382 383 #endif 384