/*@@ @file RestoreFile.c @date Thu Jun 18 16:34:59 1998 @author Tom Goodale @desc Routines to restore variables from a given file @enddesc @history @hauthor Gabrielle Allen @hdate 19 Oct 1998 @hdesc Changed names ready for thorn_IO @endhistory @version $Id$ @@*/ #include #include #include "cctk.h" #include "cctk_Parameters.h" #include "CactusBase/IOUtil/src/ioGH.h" #include "CactusPUGH/PUGH/src/include/pugh.h" #include "ioHDF5GH.h" /* MPI tag base */ #define STARTUPBASE 1001 typedef struct { cGH *GH; int ioproc; int ioproc_every; int unchunked; } iterate_info_t; typedef struct { iterate_info_t *it_info; int element_size; int iohdf5_type; #ifdef CCTK_MPI MPI_Datatype mpi_type; #endif } recover_info_t; /* local function prototypes */ static herr_t processDataset (hid_t group, const char *datasetname, void *arg); /*@@ @routine IOHDF5i_RecoverVariables @date Fri Jun 19 09:19:48 1998 @author Tom Goodale @desc Reads in data from an open HDF5 file. Code was originally in StartupReader. @enddesc @calledby IOHDF5_Recover @history @endhistory @var GH @vdesc Pointer to CCTK grid hierarchy @vtype cGH @vio in @endvar @var file @vdesc HDF5 file info structure @vtype fileinfo_t * @vio in @endvar @@*/ int IOHDF5i_RecoverVariables (cGH *GH, fileinfo_t *file) { DECLARE_CCTK_PARAMETERS #ifdef CCTK_MPI pGH *pughGH; /* PUGH extension handle */ int proc; /* looper */ CCTK_INT var_info [3]; /* communication buffer for MPI */ #endif iterate_info_t info; /* info to pass down to the iterator routine */ info.GH = GH; info.ioproc = file->ioproc; info.unchunked = file->unchunked; info.ioproc_every = file->ioproc_every; #ifdef CCTK_MPI /* Get the handle for PUGH extensions */ pughGH = PUGH_pGH (GH); /* Now process the datasets. All IO processors read the datasets from their checkpoint file, verify their contents and communicate them to the non-IO processors. */ /* At first the code for the IO processors. This holds also for the single processor case. */ if (CCTK_MyProc (GH) == file->ioproc) { #endif /* iterate over all datasets starting from "/" in the HDF5 file */ CACTUS_IOHDF5_ERROR (H5Giterate (file->fid, "/", NULL, processDataset, &info)); #ifdef CCTK_MPI /* To signal completion to the non-IO processors an invalid variable index is broadcasted. */ var_info [0] = -1; for (proc = 1; proc < file->ioproc_every; proc++) for (proc = file->ioproc + 1; proc < file->ioproc + file->ioproc_every && proc < CCTK_nProcs (GH); proc++) CACTUS_MPI_ERROR (MPI_Send (var_info, 3, PUGH_MPI_INT, proc, STARTUPBASE, pughGH->PUGH_COMM_WORLD)); } else { /* And here the code for non-IO processors: */ MPI_Datatype mpi_type; MPI_Status ms; int index, timelevel, npoints; /* They don't know how many datasets there are, because the IO processors could skip some on the fly during their consistency checks. The IO Processor sends the index of the variable to be processed next. So, all non-IO processors execute a loop where the termination condition is when an invalid index was received. */ while (1) { /* receive the next variable index from my IO processor */ CACTUS_MPI_ERROR (MPI_Recv (var_info, 3, PUGH_MPI_INT, file->ioproc, STARTUPBASE, pughGH->PUGH_COMM_WORLD, &ms)); index = var_info [0]; timelevel = var_info [1]; npoints = var_info [2]; /* check for termination condition */ if (index < 0) break; switch (CCTK_VarTypeI (index)) { case CCTK_VARIABLE_CHAR: mpi_type = PUGH_MPI_CHAR; break; case CCTK_VARIABLE_INT: mpi_type = PUGH_MPI_INT; break; case CCTK_VARIABLE_REAL: mpi_type = PUGH_MPI_REAL; break; case CCTK_VARIABLE_COMPLEX: mpi_type = pughGH->PUGH_mpi_complex; break; default: CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Unsupported datatype %d", CCTK_VarTypeI (index)); continue; } /* receive following data from my IO processor */ CACTUS_MPI_ERROR (MPI_Recv (GH->data [index][timelevel], npoints,mpi_type, file->ioproc, STARTUPBASE,pughGH->PUGH_COMM_WORLD,&ms)); } } #endif /* MPI */ return (0); } /**************************** local routines ****************************/ /* local routine GetCommonAttributes() reads in the next dataset's attributes and verifies them: * checks if there is a variable with the name given by the name attribute * verifies that this variable still belongs to the same group * checks the group data info: - group type - variable type - ntimelevels - sizes (rank, dimensions) according to chunking mode If there is a mismatch a warning (warning level 2) is printed and value of 1 is returned to indicate that this dataset should be ignored. If successful, the global variable index, the group type and the timelevel to restore are stored in {*index, *gtype, *timelevel}, and 0 is returned. */ static int GetCommonAttributes (cGH *GH, hid_t dataset, const char *datasetname, int unchunked, int *index, int *grouptype, int *timelevel, int is_group) { pGExtras *extras; cGroup groupdata; int i, flag, *dims; hid_t datatype, dataspace; hsize_t rank_stored, *dims_stored; int grouptype_stored, ndims_stored, numtimelevels_stored; char *groupname, groupname_stored [128], fullvarname [128]; /* decompose the datasetname, ignore the iteration number */ if (sscanf (datasetname, "%[^@]@%d@%d", fullvarname, &i, timelevel) != 3) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Cannot parse datasetname '%s'", datasetname); return (1); } /* check if there is a matching variable */ *index = CCTK_VarIndex (fullvarname); if (*index < 0) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "No matching variable found for '%s'", fullvarname); return (1); } /* read and verify the group name */ READ_ATTRIBUTE (dataset, "groupname", H5T_C_S1, groupname_stored); groupname = CCTK_GroupNameFromVarI (*index); if (! CCTK_Equals (groupname_stored, groupname)) { CCTK_WARN (2, "Groupnames don't match"); return (1); } free (groupname); /* verify group type, variable type, dims, sizes and ntimelevels */ READ_ATTRIBUTE (dataset, "grouptype", H5T_NATIVE_INT, &grouptype_stored); READ_ATTRIBUTE (dataset, "ntimelevels", H5T_NATIVE_INT,&numtimelevels_stored); if (CCTK_GroupData (CCTK_GroupIndex (groupname_stored), &groupdata) != 0) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Cannot get group data for '%s'", fullvarname); return (1); } if (groupdata.grouptype != grouptype_stored) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Group types don't match for '%s'", fullvarname); return (1); } if (groupdata.numtimelevels != numtimelevels_stored) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Number of timelevels don't match for '%s'", fullvarname); return (1); } /* open the first chunk to determine datatype, dims and sizes if the dataset is a chunk group */ if (is_group) CACTUS_IOHDF5_ERROR (dataset = H5Dopen (dataset, "chunk0")); CACTUS_IOHDF5_ERROR (datatype = H5Dget_type (dataset)); /* The CCTK variable type defines do not correlate with the HDF5 defines so compare them explicitely here. */ flag = (H5Tget_class (datatype) == H5T_FLOAT && groupdata.vartype == CCTK_VARIABLE_REAL) || (H5Tget_class (datatype) == H5T_INTEGER && groupdata.vartype == CCTK_VARIABLE_INT) || (H5Tget_class (datatype) == H5T_INTEGER && groupdata.vartype == CCTK_VARIABLE_CHAR) || (H5Tget_class (datatype) == H5T_COMPOUND && groupdata.vartype == CCTK_VARIABLE_COMPLEX); CACTUS_IOHDF5_ERROR (H5Tclose (datatype)); if (! flag) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Variable types don't match for '%s'", fullvarname); return (1); } /* verify the dims and sizes */ CACTUS_IOHDF5_ERROR (dataspace = H5Dget_space (dataset)); CACTUS_IOHDF5_ERROR (ndims_stored = H5Sget_simple_extent_ndims (dataspace)); dims_stored = (hsize_t *) malloc (ndims_stored * sizeof (hsize_t)); CACTUS_IOHDF5_ERROR (rank_stored = H5Sget_simple_extent_dims (dataspace, dims_stored, NULL)); /* scalars are stored as rank 0 in HDF5 but have rank 1 in Cactus */ if (rank_stored == 0) rank_stored = 1; CACTUS_IOHDF5_ERROR (H5Sclose (dataspace)); flag = 0; if (groupdata.dim != rank_stored) flag = 1; switch (groupdata.grouptype) { case CCTK_SCALAR: break; case CCTK_ARRAY: case CCTK_GF: extras = ((pGA ***) PUGH_pGH (GH)->variables)[*index][*timelevel]->extras; dims = unchunked ? extras->nsize : extras->lnsize; for (i = 0; i < groupdata.dim; i++) if (dims [groupdata.dim - i - 1] != dims_stored [i]) flag = 1; break; } free (dims_stored); if (flag) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Variable sizes don't match for '%s'", fullvarname); return (1); } if (! CCTK_QueryGroupStorageI (GH, CCTK_GroupIndexFromVarI (*index))) { CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, "Can't read into variable '%s': no storage", fullvarname); return (1); } /* close the first chunk if the dataset is a chunk group */ if (is_group) CACTUS_IOHDF5_ERROR (H5Dclose (dataset)); *grouptype = groupdata.grouptype; return (0); } static int IOHDF5_RestoreGS (hid_t dataset, int index, int timelevel, recover_info_t *rec_info) { void *data; #ifdef CCTK_MPI int proc; CCTK_INT var_info [3]; pGH *pughGH; #endif data = CCTK_VarDataPtrI (rec_info->it_info->GH, timelevel, index); /* read the data into the local variable ... */ CACTUS_IOHDF5_ERROR (H5Dread (dataset, rec_info->iohdf5_type, H5S_ALL, H5S_ALL, H5P_DEFAULT, data)); #ifdef CCTK_MPI /* ... and communicate it for the MPI case */ /* Get the handle for PUGH extensions */ pughGH = PUGH_pGH (rec_info->it_info->GH); /* set the variable's index and the timelevel */ var_info [0] = index; var_info [1] = timelevel; var_info [2] = 1; /* send info and data to the non-IO processors */ for (proc = rec_info->it_info->ioproc + 1; proc < rec_info->it_info->ioproc + rec_info->it_info->ioproc_every && proc < CCTK_nProcs (rec_info->it_info->GH); proc++) { CACTUS_MPI_ERROR (MPI_Send (var_info, 3, PUGH_MPI_INT, proc, STARTUPBASE, pughGH->PUGH_COMM_WORLD)); CACTUS_MPI_ERROR (MPI_Send (data, 1, rec_info->mpi_type, proc, STARTUPBASE, pughGH->PUGH_COMM_WORLD)); } #endif return (0); } static int IOHDF5_RestoreGA (hid_t dataset, int index, int timelevel, recover_info_t *rec_info) { #ifdef CCTK_MPI int proc, npoints; CCTK_INT var_info [3]; pGH *pughGH; void *buffer, *data; hid_t filespace, memspace; pGExtras *extras; #endif /* single processor case is easy: just read the whole dataset */ if (CCTK_nProcs (rec_info->it_info->GH) == 1) { CACTUS_IOHDF5_ERROR (H5Dread (dataset, rec_info->iohdf5_type, H5S_ALL, H5S_ALL, H5P_DEFAULT, rec_info->it_info->GH->data [index][timelevel])); return (0); } #ifdef CCTK_MPI /* Get the handle for PUGH extensions */ pughGH = PUGH_pGH (rec_info->it_info->GH); /* Get the pGExtras pointer as a shortcut */ extras = ((pGA ***) pughGH->variables) [index][timelevel]->extras; /* allocate memory for the biggest chunk */ npoints = extras->rnpoints [rec_info->it_info->ioproc + 1]; for (proc = 2; proc < rec_info->it_info->ioproc_every; proc++) if (npoints < extras->rnpoints [rec_info->it_info->ioproc + proc]) npoints = extras->rnpoints [rec_info->it_info->ioproc + proc]; buffer = malloc (npoints * rec_info->element_size); /* set the variable's index and timelevel to restore */ var_info [0] = index; var_info [1] = timelevel; /* now loop over the group of processors associated to each IO processor */ for (proc = rec_info->it_info->ioproc; proc < rec_info->it_info->ioproc + rec_info->it_info->ioproc_every && proc < CCTK_nProcs (rec_info->it_info->GH); proc++) { /* read own data directly into variable */ if (proc == rec_info->it_info->ioproc) data = rec_info->it_info->GH->data [index][timelevel]; else data = buffer; if (! rec_info->it_info->unchunked) { hid_t chunk; char chunkname [32]; /* Chunked data is stored as separate chunk datasets within a group. So open, read and close the separate chunks here. */ sprintf (chunkname, "chunk%d", proc - rec_info->it_info->ioproc); CACTUS_IOHDF5_ERROR (chunk = H5Dopen (dataset, chunkname)); CACTUS_IOHDF5_ERROR (H5Dread (chunk, rec_info->iohdf5_type, H5S_ALL, H5S_ALL, H5P_DEFAULT, data)); CACTUS_IOHDF5_ERROR (H5Dclose (chunk)); } else { int i, dim; hssize_t *chunk_origin; hsize_t *chunk_dims; /* get the dimension of the variable */ dim = CCTK_GroupDimI (CCTK_GroupIndexFromVarI (index)); chunk_origin = (hssize_t *) malloc (dim * sizeof (hssize_t)); chunk_dims = (hsize_t *) malloc (dim * sizeof (hsize_t)); /* Unchunked data is read as one hyperslab per processor. So prepare the memspace and the filespace and read the hyperslab. */ for (i = 0; i < dim; i++) { chunk_dims [dim - 1 - i] = extras->rnsize [proc][i]; chunk_origin [dim - 1 - i] = extras->lb [proc][i]; } CACTUS_IOHDF5_ERROR (filespace = H5Dget_space (dataset)); CACTUS_IOHDF5_ERROR (memspace = H5Screate_simple (dim, chunk_dims, NULL)); CACTUS_IOHDF5_ERROR (H5Sselect_hyperslab (filespace, H5S_SELECT_SET, chunk_origin, NULL, chunk_dims, NULL)); CACTUS_IOHDF5_ERROR (H5Dread (dataset, rec_info->iohdf5_type, memspace, filespace, H5P_DEFAULT, data)); CACTUS_IOHDF5_ERROR (H5Sclose (memspace)); CACTUS_IOHDF5_ERROR (H5Sclose (filespace)); free (chunk_dims); free (chunk_origin); } /* send the index and the data to the non-IO processors */ if (proc != rec_info->it_info->ioproc) { var_info [2] = extras->rnpoints [proc]; CACTUS_MPI_ERROR (MPI_Send (var_info, 3, PUGH_MPI_INT, proc, STARTUPBASE, pughGH->PUGH_COMM_WORLD)); CACTUS_MPI_ERROR (MPI_Send (data, extras->rnpoints [proc], rec_info->mpi_type, proc, STARTUPBASE, pughGH->PUGH_COMM_WORLD)); } } free (buffer); #endif return (0); } static herr_t processDataset (hid_t group, const char *datasetname, void *arg) { DECLARE_CCTK_PARAMETERS pGH *pughGH; ioGH *ioUtilGH; ioHDF5GH *myGH; int index, gtype, timelevel; iterate_info_t *it_info = (iterate_info_t *) arg; recover_info_t rec_info; hid_t dataset; H5G_stat_t object_info; int is_group; /* skip the global attributes and GH extensions groups */ if (! strcmp (datasetname, GLOBAL_PARAMETERS_GROUP) || ! strcmp (datasetname, GHEXTENSIONS_GROUP)) return (0); /* Get the handle for PUGH, IOUtil, and IOHDF5 extensions */ pughGH = PUGH_pGH (it_info->GH); ioUtilGH = (ioGH *) it_info->GH->extensions [CCTK_GHExtensionHandle ("IO")]; myGH = (ioHDF5GH*) it_info->GH->extensions [CCTK_GHExtensionHandle("IOHDF5")]; CACTUS_IOHDF5_ERROR (H5Gget_objinfo (group, datasetname, 0, &object_info)); is_group = object_info.type == H5G_GROUP; if (is_group) CACTUS_IOHDF5_ERROR (dataset = H5Gopen (group, datasetname)); else CACTUS_IOHDF5_ERROR (dataset = H5Dopen (group, datasetname)); /* read in the dataset's attributes and verify them */ if (GetCommonAttributes (it_info->GH, dataset, datasetname, it_info->unchunked, &index, >ype, &timelevel, is_group)) { CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Ignoring dataset '%s'", datasetname); return (0); } /* if we read in initial data via the file reader interface check whether the user wants to have this variable restored */ if (ioUtilGH->do_inVars && ! ioUtilGH->do_inVars [index]) { if (verbose) { char *varname = CCTK_FullName (index); CCTK_VInfo (CCTK_THORNSTRING, "Ignoring variable '%s' for file reader " "recovery", varname); free (varname); } return (0); } if (verbose) { char *varname = CCTK_FullName (index); CCTK_VInfo (CCTK_THORNSTRING, "Restoring variable '%s' (timelevel %d)", varname, timelevel); free (varname); } rec_info.it_info = it_info; switch (CCTK_VarTypeI (index)) { case CCTK_VARIABLE_CHAR: rec_info.iohdf5_type = IOHDF5_CHAR; rec_info.element_size = sizeof (CCTK_CHAR); #ifdef CCTK_MPI rec_info.mpi_type = PUGH_MPI_CHAR; #endif break; case CCTK_VARIABLE_INT: rec_info.iohdf5_type = IOHDF5_INT; rec_info.element_size = sizeof (CCTK_INT); #ifdef CCTK_MPI rec_info.mpi_type = PUGH_MPI_INT; #endif break; case CCTK_VARIABLE_REAL: rec_info.iohdf5_type = IOHDF5_REAL; rec_info.element_size = sizeof (CCTK_REAL); #ifdef CCTK_MPI rec_info.mpi_type = PUGH_MPI_REAL; #endif break; case CCTK_VARIABLE_COMPLEX: rec_info.iohdf5_type = myGH->IOHDF5_COMPLEX; rec_info.element_size = sizeof (CCTK_COMPLEX); #ifdef CCTK_MPI rec_info.mpi_type = pughGH->PUGH_mpi_complex; #endif break; default: CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Unsupported datatype %d", CCTK_VarTypeI (index)); return (0); } /* Read in the data */ switch (gtype) { case CCTK_SCALAR: IOHDF5_RestoreGS (dataset, index, timelevel, &rec_info); break; case CCTK_GF: case CCTK_ARRAY: IOHDF5_RestoreGA (dataset, index, timelevel, &rec_info); break; default: CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Unsupported group type %d", gtype); return (1); } if (is_group) CACTUS_IOHDF5_ERROR (H5Gclose (dataset)); else CACTUS_IOHDF5_ERROR (H5Dclose (dataset)); return (0); }