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