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