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