diff options
-rw-r--r-- | README | 7 | ||||
-rw-r--r-- | interface.ccl | 5 | ||||
-rw-r--r-- | param.ccl | 43 | ||||
-rw-r--r-- | schedule.ccl | 15 | ||||
-rw-r--r-- | src/DumpUtils.c | 430 | ||||
-rw-r--r-- | src/DumpVar.c | 1034 | ||||
-rw-r--r-- | src/ParseGeometry.c | 347 | ||||
-rw-r--r-- | src/ParseVars.c | 212 | ||||
-rw-r--r-- | src/RecoverVar.c | 809 | ||||
-rw-r--r-- | src/Startup.c | 147 | ||||
-rw-r--r-- | src/ioHDF5UtilGH.h | 217 | ||||
-rw-r--r-- | src/make.code.defn | 5 | ||||
-rw-r--r-- | src/make.configuration.defn | 20 |
13 files changed, 3291 insertions, 0 deletions
@@ -0,0 +1,7 @@ +Cactus Code Thorn IOHDF5Util +Authors : ... +CVS info : $Header$ +-------------------------------------------------------------------------- + +Purpose of the thorn: + diff --git a/interface.ccl b/interface.ccl new file mode 100644 index 0000000..dc99c55 --- /dev/null +++ b/interface.ccl @@ -0,0 +1,5 @@ +# Interface definition for thorn IOHDF5Util +# $Header$ + +implements: IOHDF5Util +inherits: IO, Hyperslab diff --git a/param.ccl b/param.ccl new file mode 100644 index 0000000..3f6b267 --- /dev/null +++ b/param.ccl @@ -0,0 +1,43 @@ +# Parameter definitions for thorn IOHDF5Util +# $Header$ + +###################### +# Hyperslab parameters +###################### +STRING origin "Default origin" +{ + .* :: "Comma separated list of positive integer values" +} "0,0,0" +STRING downsampling "Default downsampling" +{ + .* :: "Comma separated list of positive integer values" +} "1,1,1" +STRING length "Default length of the hyperslab to stream" +{ + .* :: "Comma separated list of integer values" +} "-1,-1,-1" +STRING direction "Default direction of hyperslab to stream" +{ + .* :: "Comma separated list of positive integer values" +} "0,1,2" + +CCTK_INT slabdim "default dimension of slab" +{ + 1:3 :: "dimension of slab (1,2,3)" +} 3 + + +############################################################################# +### import IOUtil parameters +############################################################################# +shares: IO + +################ +# various things +################ +USES BOOLEAN verbose "" +{ +} +USES BOOLEAN out3D_parameters "" +{ +} diff --git a/schedule.ccl b/schedule.ccl new file mode 100644 index 0000000..eafde10 --- /dev/null +++ b/schedule.ccl @@ -0,0 +1,15 @@ +# Schedule definitions for thorn IOHDF5Util +# $Header$ + +######################################################################## +### register IOHDF5Util routines +######################################################################## +schedule IOHDF5Util_Startup at STARTUP after IOUtil_Startup +{ + LANG:C +} "IOHDF5Util startup routine" + +schedule IOHDF5Util_Terminate at TERMINATE before Driver_Terminate +{ + LANG:C +} "IOHDF5Util termination routine" diff --git a/src/DumpUtils.c b/src/DumpUtils.c new file mode 100644 index 0000000..73d7fe3 --- /dev/null +++ b/src/DumpUtils.c @@ -0,0 +1,430 @@ + /*@@ + @file DumpUtils.c + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + Utility routines for dumping data into HDF5 files. + @enddesc + @version $Id$ + @@*/ + + +#include <stdlib.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "CactusBase/IOUtil/src/ioGH.h" +#include "CactusBase/IOUtil/src/ioutil_CheckpointRecovery.h" +#include "ioHDF5UtilGH.h" + +/* the rcs ID and its dummy function to use it */ +static char *rcsid = "$Id$"; +CCTK_FILEVERSION(BetaThorns_IOHDF5Util_DumpUtils_c) + + +/*@@ + @routine IOHDF5Util_DumpGH + @date Tue Oct 10 2000 + @author Thomas Radke + @desc + Dump all variables of the given GH along with + the global attributes and parameters + to the (already opened) HDF5 checkpoint file. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var timers + @vdesc array of timers for timing the checkpointing + @vtype int [] + @vio in + @endvar + @var file + @vdesc HDF5 file handle + @vtype hid_t + @vio in + @endvar + + @returntype int + @returndesc + always 0 (for success) + @endreturndesc +@@*/ +int IOHDF5Util_DumpGH (cGH *GH, + int timers[], + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + int index; + int maxdim; + int timelevel, current_timelevel; + int *old_downsample; + int old_out_single; + ioGH *ioUtilGH; + const char *thorn; + const char *impl; + ioHDF5Geo_t slab; + + + /* Get the GH extension for IOUtil */ + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + + /* disable downsampling after saving original downsampling params */ + maxdim = CCTK_MaxDim (); + old_downsample = (int *) malloc (maxdim * sizeof (int)); + for (index = 0; index < maxdim; index++) + { + old_downsample[index] = ioUtilGH->downsample[index]; + ioUtilGH->downsample[index] = 1; + slab.direction[index] = index; + slab.downsample[index] = 1; + } + + /* disable output in single precision */ + old_out_single = ioUtilGH->out_single; + ioUtilGH->out_single = 0; + + /* sync all active groups */ + for (index = CCTK_NumGroups () - 1; index >= 0; index--) + { + if (CCTK_IsImplementationActive ( + CCTK_ImpFromVarI (CCTK_FirstVarIndexI (index)))) + { + CCTK_SyncGroupI (GH, index); + } + } + + /* start CP_PARAMETERS_TIMER timer */ + if (timers) + { + CCTK_TimerResetI (timers[CP_PARAMETERS_TIMER]); + CCTK_TimerStartI (timers[CP_PARAMETERS_TIMER]); + } + + /* Great; Now start dumping away! */ + if (file >= 0) + { + /* all GH extension variables and parameter stuff which is not + specific to any dataset goes into dedicated groups */ + if (out3D_parameters) + { + IOHDF5Util_DumpParameters (GH, file); + } + + IOHDF5Util_DumpGHExtensions (GH, file); + } + + if (verbose) + { + CCTK_INFO ("Dumping variables ..."); + } + + /* stop CP_PARAMETERS_TIMER timer and start CP_VARIABLES_TIMER */ + if (timers) + { + CCTK_TimerStopI (timers[CP_PARAMETERS_TIMER]); + CCTK_TimerResetI (timers[CP_VARIABLES_TIMER]); + CCTK_TimerStartI (timers[CP_VARIABLES_TIMER]); + } + + /* ... now the variables */ + for (index = CCTK_NumVars () - 1; index >= 0; index--) + { + /* find out the thorn implementing variable with index */ + impl = CCTK_ImpFromVarI (index); + thorn = CCTK_ImplementationThorn (impl); + if (! thorn) + { + thorn = impl; + } + + /* let only variables pass which belong to an active thorn and + have storage assigned */ + if (! CCTK_IsThornActive (thorn) || + ! CCTK_QueryGroupStorageI (GH, CCTK_GroupIndexFromVarI (index))) + { + continue; + } + + if (verbose && file >= 0) + { + CCTK_VInfo (CCTK_THORNSTRING, " %s", CCTK_VarName (index)); + } + + /* get the current timelevel */ + current_timelevel = CCTK_NumTimeLevelsFromVarI (index) - 1; + if (current_timelevel > 0) + { + current_timelevel--; + } + + slab.sdim = slab.vdim = CCTK_GroupDimFromVarI (index); + + /* now dump all timelevels up to the current */ + for (timelevel = 0; timelevel <= current_timelevel; timelevel++) + { + IOHDF5Util_DumpVar (GH, index, timelevel, &slab, file, 0); + } + + } /* end of loop over all variables */ + + /* stop CP_VARIABLES_TIMER timer */ + if (timers) + { + CCTK_TimerStopI (timers[CP_VARIABLES_TIMER]); + } + + /* restore output precision flag */ + ioUtilGH->out_single = old_out_single; + + /* restore original downsampling params */ + memcpy (ioUtilGH->downsample, old_downsample, maxdim * sizeof (int)); + free (old_downsample); + + return (0); +} + + +/*@@ + @routine IOHDF5Util_DumpCommonAttributes + @date May 21 1999 + @author Thomas Radke + @desc + Add "Common" attributes to the given dataset, these are: + <ul> + <li> the variable's groupname + <li> the grouptype + <li> number of timelevels + <li> simulation time + <li> origin + <li> bounding box + <li> gridspacings + <li> global grid size + </ul> + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var index + @vdesc CCTK index of the variable to dump + @vtype int + @vio in + @endvar + @var timelevel + @vdesc timelevel of the variable to dump + @vtype int + @vio in + @endvar + @var global_shape + @vdesc global shape of the grid + @vtype CCTK_INT [] + @vio in + @endvar + @var slab + @vdesc pointer to hyperslab structure describing the dataset + @vtype ioHDF5Geo_t * + @vio in + @endvar + @var dataset + @vdesc the dataset handle where the attributes should be attached to + @vtype hid_t + @vio in + @endvar +@@*/ +void IOHDF5Util_DumpCommonAttributes (cGH *GH, + int index, + int timelevel, + CCTK_INT global_shape[], + ioHDF5Geo_t *slab, + hid_t dataset) +{ + DECLARE_CCTK_PARAMETERS + int dim; + char *groupname; + CCTK_INT attr_int; + CCTK_REAL *attr_real; + ioHDF5UtilGH *myGH; + + + /* Get the handle for IOHDF5Util extensions */ + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + /* attributes describing the variable */ + dim = CCTK_GroupDimI (CCTK_GroupIndexFromVarI (index)); + groupname = CCTK_GroupNameFromVarI (index); + WRITE_ATTRIBUTE ("groupname", groupname, dataset, + myGH->scalar_dataspace, 0, myGH->IOHDF5_STRING); + free (groupname); + attr_int = CCTK_GroupTypeFromVarI (index); + WRITE_ATTRIBUTE ("grouptype", &attr_int, dataset, + myGH->scalar_dataspace, 0, IOHDF5_INT); + attr_int = CCTK_NumTimeLevelsFromVarI (index); + WRITE_ATTRIBUTE ("ntimelevels", &attr_int, dataset, + myGH->scalar_dataspace, 0, IOHDF5_INT); + WRITE_ATTRIBUTE ("global_size", global_shape, dataset, + myGH->array_dataspace, dim, IOHDF5_INT); + WRITE_ATTRIBUTE ("time", &GH->cctk_time, dataset, + myGH->scalar_dataspace, 0, IOHDF5_REAL); + + /* attributes describing the underlying grid */ + dim = CCTK_MaxDim (); + attr_real = (CCTK_REAL *) malloc (2 * dim * sizeof (CCTK_REAL)); + CCTK_CoordRange (GH, &attr_real[0], &attr_real[dim + 0], -1, "x", "cart3d"); + CCTK_CoordRange (GH, &attr_real[1], &attr_real[dim + 1], -1, "y", "cart3d"); + CCTK_CoordRange (GH, &attr_real[2], &attr_real[dim + 2], -1, "z", "cart3d"); + + WRITE_ATTRIBUTE ("origin", attr_real, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("min_ext", attr_real, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("max_ext", attr_real + dim, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("delta", GH->cctk_delta_space, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + free (attr_real); + + /* attributes describing the hyperslab */ + if (slab) + { + attr_real = (CCTK_REAL *) malloc (4 * slab->sdim * sizeof (CCTK_REAL)); + for (dim = 0; dim < slab->sdim; dim++) + { + attr_real[dim + 0*slab->sdim] = + slab->origin[slab->direction[dim]] * + GH->cctk_delta_space[slab->direction[dim]]; + attr_real[dim + 1*slab->sdim] = + slab->origin[dim] * GH->cctk_delta_space[dim]; + attr_real[dim + 2*slab->sdim] = + (slab->origin[dim] + slab->actlen[dim]-1) * + GH->cctk_delta_space[dim] * slab->downsample[dim]; + attr_real[dim + 3*slab->sdim] = + GH->cctk_delta_space[slab->direction[dim]] * + slab->downsample[slab->direction[dim]]; + } + WRITE_ATTRIBUTE ("origin_slab", attr_real + 0*dim, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("min_ext_slab", attr_real + 1*dim, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("max_ext_slab", attr_real + 2*dim, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + WRITE_ATTRIBUTE ("delta_slab", attr_real + 3*dim, dataset, + myGH->array_dataspace, dim, IOHDF5_REAL); + free (attr_real); + } +} + + + /*@@ + @routine IOHDF5Util_DumpParameters + @date Thu Jan 20 2000 + @author Thomas Radke + @desc + Collects the parameters of all active implementations + into a single string and writes it as an attribute + attached to the CACTUS_PARAMETERS_GROUP group in the HDF5 file. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var file + @vdesc the HDF5 file handle + @vtype hid_t + @vio in + @endvar +@@*/ +void IOHDF5Util_DumpParameters (cGH *GH, + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + ioHDF5UtilGH *myGH; + char *parameters; + hid_t group; + + + if (verbose) + { + CCTK_INFO ("Dumping GH extensions ..."); + } + + /* Get the parameter string buffer */ + parameters = IOUtil_GetAllParameters (GH); + + if (parameters) + { + /* Get the handle for IOHDF5Util extensions */ + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + IOHDF5_ERROR (group = H5Gcreate (file, CACTUS_PARAMETERS_GROUP, 0)); + WRITE_ATTRIBUTE (ALL_PARAMETERS, parameters, group, + myGH->scalar_dataspace, 0, myGH->IOHDF5_STRING); + IOHDF5_ERROR (H5Gclose (group)); + + free (parameters); + } +} + + + /*@@ + @routine IOHDF5Util_DumpGHExtensions + @date Thu Jan 20 2000 + @author Thomas Radke + @desc + Attaches the current values of important flesh and IO variables + as attributes to the GLOBAL_ATTRIBUTES_GROUP group + in the HDF5 file. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var file + @vdesc the HDF5 file handle + @vtype hid_t + @vio in + @endvar +@@*/ +void IOHDF5Util_DumpGHExtensions (cGH *GH, + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + int value; + hid_t group; + ioGH *ioUtilGH; + ioHDF5UtilGH *myGH; + + + if (verbose) + { + CCTK_INFO ("Dumping GH extensions ..."); + } + + /* Get the handles for IOUtil and IOHDF5 extensions */ + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + IOHDF5_ERROR (group = H5Gcreate (file, GLOBAL_ATTRIBUTES_GROUP, 0)); + + value = CCTK_MainLoopIndex (); + WRITE_ATTRIBUTE ("main_loop_index", &value, group, + myGH->scalar_dataspace, 0, H5T_NATIVE_INT); + value = CCTK_nProcs (GH); + WRITE_ATTRIBUTE ("nprocs", &value, group, + myGH->scalar_dataspace, 0, H5T_NATIVE_INT); + WRITE_ATTRIBUTE ("ioproc_every", &ioUtilGH->ioproc_every, group, + myGH->scalar_dataspace, 0, H5T_NATIVE_INT); + WRITE_ATTRIBUTE ("unchunked", &ioUtilGH->unchunked, group, + myGH->scalar_dataspace, 0, H5T_NATIVE_INT); + WRITE_ATTRIBUTE ("cctk_time", &GH->cctk_time, group, + myGH->scalar_dataspace, 0, IOHDF5_REAL); + WRITE_ATTRIBUTE ("cctk_iteration", &GH->cctk_iteration, group, + myGH->scalar_dataspace, 0, H5T_NATIVE_INT); + + IOHDF5_ERROR (H5Gclose (group)); +} diff --git a/src/DumpVar.c b/src/DumpVar.c new file mode 100644 index 0000000..c08c281 --- /dev/null +++ b/src/DumpVar.c @@ -0,0 +1,1034 @@ +/*@@ + @file DumpVar.c + @date Fri May 21 1999 + @author Thomas Radke + @desc + Routines for writing variables into HDF5 data or checkpoint files. + These routines are used by other HDF5 IO methods. + @enddesc + @@*/ + + +#include <stdlib.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "CactusPUGH/PUGH/src/include/pugh.h" +#include "CactusPUGH/PUGHSlab/src/PUGHSlab.h" +#include "CactusBase/IOUtil/src/ioGH.h" +#include "ioHDF5UtilGH.h" + + +/* #define DEBUG_ME 1 */ + +#define IOTAGBASE 20000 /* This may break on more than 2000 processors */ + + +/* info structure describing how to dump a given CCTK variable type */ +typedef struct +{ + int iohdf5_type; /* the HDF5 type to use */ + int element_size; /* size of a single data element */ +#ifdef CCTK_MPI + MPI_Datatype mpi_type; /* the data type for MPI communication */ +#endif +} dumpInfo; + + +/* local function prototypes */ +static int IOHDF5Util_DumpGS (cGH *GH, + int index, + int timelevel, + int iohdf5_type, + int check_exisiting_objects, + hid_t file); +static int IOHDF5Util_DumpGA (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + dumpInfo *info, + int check_exisiting_objects, + hid_t file); +static int IOHDF5Util_getDumpData (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void **outme, + int *free_outme, + CCTK_INT *geom, + int element_size); +static void IOHDF5Util_procDump (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void *outme, + CCTK_INT *geom, + int proc, + int iohdf5_type, + int check_exisiting_objects, + hid_t file); +#ifdef CCTK_MPI +#ifdef HAVE_PARALLEL +static void IOHDF5Util_collectiveDump (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void *outme, + CCTK_INT *geom, + int hdf5io_type, + hid_t file); +#endif /* HAVE_PARALLEL */ +#endif /* MPI */ + + + +/*@@ + @routine IOHDF5Util_DumpVar + @date 16 Apr 1999 + @author Thomas Radke + @desc + Generic dump routine, just calls the type-specific dump routine. + @enddesc + + @calls IOHDF5Util_DumpGS + IOHDF5Util_DumpGA + + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var index + @vdesc index of the variable to be dumped + @vtype int + @vio in + @endvar + @var timelevel + @vdesc the timelevel to dump + @vtype int + @vio in + @endvar + @var slab + @vdesc pointer to the hyperslab structure + @vtype ioHDF5Geo_t * + @vio in + @endvar + @var file + @vdesc the HDF5 file handle + @vtype hid_t + @vio in + @endvar + @var check_exisiting_objects + @vdesc flag indicating if we should check for already existing objects + before these are created anew + This might happen for existing files reopened after recovery. + @vtype int + @vio in + @endvar + + @returntype int + @returndesc + -1 if variable has unknown type + or return code of either + @seeroutine IOHDF5Util_DumpGS or + @seeroutine IOHDF5Util_DumpGA + @endreturndesc +@@*/ +int IOHDF5Util_DumpVar (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + hid_t file, + int check_exisiting_objects) +{ + int gtype; + pGH *pughGH; + ioHDF5UtilGH *myGH; + ioGH *ioUtilGH; + dumpInfo info; + int retval; + + + /* Get the handle for PUGH, IOUtil, and IOHDF5Util extensions */ + pughGH = PUGH_pGH (GH); + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + switch (CCTK_VarTypeI (index)) + { + case CCTK_VARIABLE_CHAR: + info.iohdf5_type = IOHDF5_CHAR; + info.element_size = sizeof (CCTK_CHAR); +#ifdef CCTK_MPI + info.mpi_type = PUGH_MPI_CHAR; +#endif + break; + + case CCTK_VARIABLE_INT: + info.iohdf5_type = IOHDF5_INT; + info.element_size = sizeof (CCTK_INT); +#ifdef CCTK_MPI + info.mpi_type = PUGH_MPI_INT; +#endif + break; + + case CCTK_VARIABLE_REAL: + info.iohdf5_type = ioUtilGH->out_single ? IOHDF5_REAL4 : IOHDF5_REAL; + info.element_size = ioUtilGH->out_single ? + sizeof (CCTK_REAL4) : sizeof (CCTK_REAL); +#ifdef CCTK_MPI + info.mpi_type = ioUtilGH->out_single ? PUGH_MPI_REAL4 : PUGH_MPI_REAL; +#endif + break; + + case CCTK_VARIABLE_COMPLEX: + info.iohdf5_type = myGH->IOHDF5_COMPLEX; + info.element_size = sizeof (CCTK_COMPLEX); +#ifdef CCTK_MPI + info.mpi_type = pughGH->PUGH_mpi_complex; +#endif + break; + + default: + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Unsupported variable datatype '%s'", + CCTK_VarTypeName (CCTK_VarTypeI (index))); + return (-1); + } + + /* now branch to the appropriate dump routine */ + gtype = CCTK_GroupTypeFromVarI (index); + if (gtype == CCTK_SCALAR) + { + retval = IOHDF5Util_DumpGS (GH, index, timelevel, info.iohdf5_type, + check_exisiting_objects, file); + } + else if (gtype == CCTK_ARRAY || gtype == CCTK_GF) + { + retval = IOHDF5Util_DumpGA (GH, index, timelevel, slab, &info, + check_exisiting_objects, file); + } + else + { + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Invalid group type %d", gtype); + retval = -1; + } + + return (retval); +} + + +/*@@ + @routine IOHDF5Util_DumpGS + @date May 21 1999 + @author Thomas Radke + @desc + Dumps a CCTK_SCALAR type variable into a HDF5 file. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH + @vio in + @endvar + @var index + @vdesc index of the variable to be dumped + @vtype int + @vio in + @endvar + @var timelevel + @vdesc the timelevel to dump + @vtype int + @vio in + @endvar + @var iohdf5_type + @vdesc the HDF5 datatype + @vtype int + @vio in + @endvar + @var check_exisiting_objects + @vdesc flag indicating if we should check for already existing objects + before these are created anew + This might happen for existing files reopened after recovery. + @vtype int + @vio in + @endvar + @var file + @vdesc the HDF5 file to dump to + @vtype hid_t + @vio in + @endvar + + @returntype int + @returndesc + always 0 + @endreturndesc +@@*/ +static int IOHDF5Util_DumpGS (cGH *GH, + int index, + int timelevel, + int iohdf5_type, + int check_exisiting_objects, + hid_t file) +{ + ioGH *ioUtilGH; + ioHDF5UtilGH *myGH; + hid_t dataset; + char *fullname; + char *datasetname; + CCTK_INT global_shape[1] = {0}; + + + /* Get the handles for IOHDF5Util and IOUtil extensions */ + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + if (CCTK_MyProc (GH) != ioUtilGH->ioproc) + { + return (0); + } + + fullname = CCTK_FullName (index); + datasetname = (char *) malloc (strlen (fullname) + 80); + sprintf (datasetname, "%s timelevel %d at iteration %d", + fullname, timelevel, GH->cctk_iteration); + + /* if restart from recovery delete an already existing dataset + so that it can be created as anew */ + if (check_exisiting_objects) + { + IOHDF5_ERROR (H5Eset_auto (NULL, NULL)); + H5Gunlink (file, datasetname); + IOHDF5_ERROR (H5Eset_auto (myGH->print_error_fn, myGH->print_error_fn_arg)); + } + + IOHDF5_ERROR (dataset = H5Dcreate (file, datasetname, iohdf5_type, + myGH->scalar_dataspace, H5P_DEFAULT)); + IOHDF5_ERROR (H5Dwrite (dataset, iohdf5_type, H5S_ALL, H5S_ALL, H5P_DEFAULT, + CCTK_VarDataPtrI (GH, timelevel, index))); + + IOHDF5Util_DumpCommonAttributes (GH, index, timelevel, global_shape, + NULL, dataset); + + IOHDF5_ERROR (H5Dclose (dataset)); + + free (datasetname); + free (fullname); + + return (0); +} + + +/*@@ + @routine IOHDF5Util_DumpGA + @date May 21 1999 + @author Paul Walker + @desc + Dumps a CCTK_GF or CCTK_ARRAY type variable into a HDF5 file. + @enddesc + + @calls IOHDF5Util_getDumpData + IOHDF5Util_procDump + IOHDF5Util_collectiveDump + + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH + @vio in + @endvar + @var index + @vdesc index of the variable to be dumped + @vtype int + @vio in + @endvar + @var timelevel + @vdesc the timelevel to dump + @vtype int + @vio in + @endvar + @var slab + @vdesc pointer to the hyperslab structure + @vtype ioHDF5Geo_t * + @vio in + @endvar + @var info + @vdesc info structure describing + - the HDF5 datatype to store + - the size of an element of variable + - the MPI datatype to communicate + @vtype dumpInfo + @vio in + @endvar + @var check_exisiting_objects + @vdesc flag indicating if we should check for already existing objects + before these are created anew + This might happen for existing files reopened after recovery. + @vtype int + @vio in + @endvar + @var file + @vdesc the HDF5 file to dump to + @vtype hid_t + @vio in + @endvar + + @returntype int + @returndesc + 0 for success + or returncode of @seeroutine IOHDF5Util_getDumpData + @endreturndesc +@@*/ +static int IOHDF5Util_DumpGA (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + dumpInfo *info, + int check_exisiting_objects, + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + ioGH *ioUtilGH; + pGH *pughGH; + int myproc; + int nprocs; + CCTK_INT *geom; /* Lower bounds, sizes and global shape of data */ + void *outme; /* The data pointer to dump ... */ + int free_outme; /* and whether it needs freeing */ + int retval; + void *tmpd; + int incoming; + int outgoing; +#ifdef CCTK_MPI + int i, j; + MPI_Status ms; +#endif + + + /* Get the handles for PUGH and IOUtil extensions */ + pughGH = PUGH_pGH (GH); + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + + /* allocate the geometry buffer */ + geom = (CCTK_INT *) malloc (3 * slab->sdim * sizeof (CCTK_INT)); + + /* Get the pointer to the data we want to output (includes downsampling) */ + retval = IOHDF5Util_getDumpData (GH, index, timelevel, slab, &outme, + &free_outme, geom, info->element_size); + if (retval < 0) + { + return (retval); + } + + myproc = CCTK_MyProc (GH); + nprocs = CCTK_nProcs (GH); + +#ifdef CCTK_MPI +#ifdef HAVE_PARALLEL + if (ioUtilGH->unchunked) + { + IOHDF5Util_collectiveDump (GH, index, timelevel, slab, outme, geom, + info->iohdf5_type, file); + if (free_outme) + { + free (outme); + } + free (geom); + return (0); + } +#endif +#endif + + /* Dump data held on IO processor */ + if (myproc == ioUtilGH->ioproc) + { + IOHDF5Util_procDump (GH, index, timelevel, slab, outme, geom, + myproc, info->iohdf5_type, check_exisiting_objects, + file); + } +#ifdef CCTK_MPI + if (myproc == ioUtilGH->ioproc) + { + /* Dump data from all other processors */ + for (i = myproc+1; i < myproc+ioUtilGH->ioproc_every && i < nprocs; i++) + { + /* receive geometry */ + CACTUS_MPI_ERROR (MPI_Recv (geom, 3*slab->sdim, PUGH_MPI_INT, i, + 2*i + IOTAGBASE + 1, pughGH->PUGH_COMM_WORLD, + &ms)); + + incoming = geom[slab->sdim]; + for (j = 1; j < slab->sdim; j++) + { + incoming *= geom[slab->sdim + j]; + } + + /* receive data */ + if (incoming > 0) + { + tmpd = malloc (incoming * info->element_size); + CACTUS_MPI_ERROR (MPI_Recv (tmpd, incoming, info->mpi_type, i, + 2*i + IOTAGBASE, pughGH->PUGH_COMM_WORLD, + &ms)); + + IOHDF5Util_procDump (GH, index, timelevel, slab, tmpd, geom, i, + info->iohdf5_type, check_exisiting_objects, file); + free (tmpd); + } + + } /* End loop over processors */ + + } + else + { + + /* Send the geometry and data from the non-IO processors + * to the IO processors + */ + outgoing = geom[slab->sdim]; + for (i = 1; i < slab->sdim; i++) + { + outgoing *= geom[slab->sdim + i]; + } + + /* send geometry */ + CACTUS_MPI_ERROR (MPI_Send (geom, 3*slab->sdim, PUGH_MPI_INT, + ioUtilGH->ioproc, 2*myproc + IOTAGBASE + 1, + pughGH->PUGH_COMM_WORLD)); + if (outgoing > 0) + { + /* send data if outgoing>0 */ + CACTUS_MPI_ERROR (MPI_Send (outme, outgoing, info->mpi_type, + ioUtilGH->ioproc, 2*myproc + IOTAGBASE, + pughGH->PUGH_COMM_WORLD)); +#ifdef DEBUG_ME + printf ("Processor %d sent %d data points\n", myproc, outgoing); +#endif + } + } + + /* wait for every processor to catch up */ + CCTK_Barrier (GH); +#endif + + /* free allocated resources */ + if (free_outme) + { + free (outme); + } + free (geom); + + return (retval); +} + + +/**************************************************************************/ +/* local functions */ +/**************************************************************************/ +static int IOHDF5Util_getDumpData (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void **outme, + int *free_outme, + CCTK_INT *geom, + int element_size) +{ + ioGH *ioUtilGH; + pGExtras *extras=NULL; + CCTK_REAL4 *single_ptr; + int myproc, do_downsample; + int idim, sdim, vdim; /*iteration, slab dim, variable dim */ + int *hsizes,*hsizes_global,*hsizes_offset; /* geometry information */ + int sdir[3]; /* slab direction FIXME */ + int locdowns[3]; /* Holds the downsampling in the order + specified in slab->direction[] */ + int slabstart[3]; /* global start index for slab to extract */ + int ip; + void *data; + + + myproc = CCTK_MyProc (GH); + + /* get GH extension for IOUtil */ + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + + /* get the pointer to PUGH specific structures */ + extras = ((pGA***) PUGH_pGH(GH)->variables)[index][timelevel]->extras; + data = CCTK_VarDataPtrI (GH, timelevel, index); + + /* Shortcuts */ + vdim = slab->vdim; + sdim = slab->sdim; + + if (vdim>3) + { + CCTK_WARN(1,"Cannot treat GFs with dim>3, check with cactus support\n"); + return(-1); + } + + for (idim=0;idim<sdim;idim++) + locdowns[idim]=slab->downsample[slab->direction[idim]]; + + /* Set downsample flag */ + do_downsample = 0; + for (idim = 0; idim < sdim; idim++) + do_downsample |= locdowns[idim] > 1; + + /* Simple case if no downsampling and sdim==vdim */ + if (! do_downsample && sdim==vdim) + { + if (ioUtilGH->out_single) { + single_ptr = (CCTK_REAL4 *) malloc (extras->npoints*sizeof (CCTK_REAL4)); + for (ip = 0; ip < extras->npoints; ip++) + single_ptr[ip] = (CCTK_REAL4) ((CCTK_REAL *) data)[ip]; + + *outme = single_ptr; + *free_outme = 1; + } + else + { + *outme = data; + *free_outme = 0; + } + + for (idim = 0; idim < sdim; idim++) { + geom[idim + 0*sdim] = extras->lb[myproc][idim]; + geom[idim + 1*sdim] = extras->lnsize[idim]; + geom[idim + 2*sdim] = extras->nsize[idim]; + } + } + else /* Call Hyperslab else */ + { + + /* allocate array to hold size information of slab */ + hsizes = (int*)malloc(sdim*sizeof(int)); + hsizes_global= (int*)malloc(sdim*sizeof(int)); + hsizes_offset= (int*)malloc(sdim*sizeof(int)); + + for (idim=0;idim<sdim;idim++) + { + hsizes[idim]=0; + hsizes_global[idim]=0; + hsizes_offset[idim]=0; + } + + /* TEMPORARY FIX DUE TO INCONSISTENT DIR SPECIFICATION */ + /* direction vector of 1d line in 3d volume */ + if (sdim==1) { + sdir[0] = (slab->direction[0]==0) ? 1:0; + sdir[1] = (slab->direction[0]==1) ? 1:0; + sdir[2] = (slab->direction[0]==2) ? 1:0; + } + /* norm vector of 2d surface in 3d volume */ + else if (sdim==2) + { + sdir[0] = ((slab->direction[0]==1)&&(slab->direction[1]==2)) ? 1:0; + sdir[1] = ((slab->direction[0]==0)&&(slab->direction[1]==2)) ? 1:0; + sdir[2] = ((slab->direction[0]==0)&&(slab->direction[1]==1)) ? 1:0; + } + /* spanning directions for 3d */ + else if (sdim==3) + { + sdir[0] = slab->direction[0]; + sdir[1] = slab->direction[1]; + sdir[2] = slab->direction[2]; + } + + /* Get the start indeces for the slab from origin. */ + /* Slab_start from the geo structure is not a true start index, it is currently + used as the intersection point of three surfaces , that's why two entries + have to be set to zero in the case of a surface (one entry for line) FIXME */ + for (idim=0;idim<vdim;idim++) + slabstart[idim] = slab->origin[idim]; + for (idim=0;idim<sdim;idim++) + slabstart[slab->direction[idim]]=0; + + + if (Hyperslab_GetLocalHyperslab (GH, index, timelevel, sdim, slab->origin, + sdir, slab->length, locdowns, outme, + hsizes, hsizes_global, hsizes_offset)< 0) { + char *fullname = CCTK_FullName (index); + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Failed to extract hyperslab for variable '%s'", fullname); + free (fullname); + *free_outme = 0; + return(-1); + + } + + *free_outme = 1; + + for (idim = 0; idim < sdim; idim++) { + geom[idim + 0*sdim] = hsizes_offset[idim]; + geom[idim + 1*sdim] = hsizes[idim]; + geom[idim + 2*sdim] = hsizes_global[idim]; + slab->actlen[idim]= hsizes_global[idim]; + } + } +#ifdef DEBUG_ME + printf ("***** %dD \n",sdim); + printf ("Lower bound: %d", (int) geom[0]); + for (idim = 1; idim < sdim; idim++) + printf (" %d", (int) geom[idim]); + printf ("\n"); + printf ("Chunk size : %d", (int) geom[1*sdim + 0]); + for (idim = 1; idim < sdim; idim++) + printf (" %d", (int) geom[1*sdim + idim]); + printf ("\n"); + printf ("Global size: %d", (int) geom[2*sdim + 0]); + for (idim = 1; idim < sdim; idim++) + printf (" %d", (int) geom[2*sdim + idim]); + printf ("\n"); + printf ("Downsampling: "); + for (idim = 0; idim < sdim; idim++) + printf (" %d", locdowns[idim]); + printf ("\n"); + printf ("Direction: "); + for (idim = 0; idim < sdim; idim++) + printf (" %d",slab->direction[idim]); + printf ("\n"); + printf ("Slab Start: %d", (int) slabstart[0]); + for (idim = 1; idim < vdim; idim++) + printf (" %d", (int) slabstart[idim]); + printf ("\n"); + printf ("Slab intersect cntr: %d", (int) slab->origin[0]); + for (idim = 1; idim < vdim; idim++) + printf (" %d", (int) slab->origin[idim]); + printf ("\n"); + printf ("SlabDelta: "); + for (idim = 0; idim < sdim; idim++) + printf (" %f",GH->cctk_delta_space[slab->direction[idim]]*locdowns[idim]); + printf ("\n"); +#endif + + return(1); +} + + +/*@@ + @routine IOHDF5Util_procDump + @author Thomas Radke + @date May 21 1999 + @desc + All IO processors dump the data of their group into the file, + either as one chunk per processor or unchunked + (depending on parameter "unchunked"). + @enddesc + + @calls IOHDF5Util_DumpCommonAttributes + + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var index + @vdesc global index of the variable to be dumped + @vtype int + @vio in + @endvar + @var timelevel + @vdesc the timelevel to store + @vtype int + @vio in + @endvar + @var outme + @vdesc pointer to the chunk to dump + @vtype void * + @vio in + @endvar + @var geom + @vdesc bounds and sizes of the output + @vtype CCTK_INT[3*dim] + @vio in + @vcomment geom[0*dim..1*dim-1]: lower bounds + geom[1*dim..2*dim-1]: (local) data sizes + geom[2*dim..3*dim-1]: global space + @endvar + @var proc + @vdesc the processor which's chunk is to be dumped + @vtype int + @vio in + @endvar + @var hdf5io_type + @vdesc the HDF5 datatype identifier + @vtype int + @vio in + @endvar + @var file + @vdesc the HDF5 file handle + @vtype hid_t + @vio in + @endvar +@@*/ +static void IOHDF5Util_procDump (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void *outme, + CCTK_INT *geom, + int proc, + int iohdf5_type, + int check_exisiting_objects, + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + int i, sdim; + int myproc; + ioGH *ioUtilGH; + ioHDF5UtilGH *myGH; + hid_t group, dataset, memspace, filespace; + char *fullname, *datasetname, *chunkname; + hssize_t *chunk_origin; + hsize_t *chunk_dims, *file_dims; + int locpoints; + + + ioUtilGH = (ioGH *) CCTK_GHExtension (GH, "IO"); + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + sdim = slab->sdim; + myproc = CCTK_MyProc (GH); + + /* Check if there are points to output */ + locpoints = geom[sdim]; + for (i = 1; i < sdim; i++) + { + locpoints *= geom[sdim + i]; + } + if (locpoints <= 0) + { + return; + } + + chunk_origin = (hssize_t *) malloc (sdim * sizeof (hssize_t)); + chunk_dims = (hsize_t *) malloc (2*sdim * sizeof (hsize_t)); + file_dims = chunk_dims + sdim; + + /* HDF5 needs it in reverse order */ + for (i = 0; i < sdim; i++) + { + chunk_origin[i] = geom[1*sdim - 1 - i]; + chunk_dims [i] = geom[2*sdim - 1 - i]; + file_dims [i] = geom[3*sdim - 1 - i]; + } + + /* build the unique dataset name from the variable's full name, + the timelevel and the current iteration number */ + fullname = CCTK_FullName (index); + datasetname = (char *) malloc (strlen (fullname) + 80); + sprintf (datasetname, "%s timelevel %d at iteration %d", + fullname, timelevel, GH->cctk_iteration); + + /* create the memspace according to chunk dims */ + IOHDF5_ERROR (memspace = H5Screate_simple (sdim, chunk_dims, NULL)); + + /* if restart from recovery delete an already existing dataset + so that it can be created as anew */ + if (proc == myproc && check_exisiting_objects) + { + IOHDF5_ERROR (H5Eset_auto (NULL, NULL)); + H5Gunlink (file, datasetname); + IOHDF5_ERROR (H5Eset_auto (myGH->print_error_fn, myGH->print_error_fn_arg)); + } + + if (ioUtilGH->unchunked) + { + /* create the (global) filespace and set the hyperslab for the chunk */ + IOHDF5_ERROR (filespace = H5Screate_simple (sdim, file_dims, NULL)); + IOHDF5_ERROR (H5Sselect_hyperslab (filespace, H5S_SELECT_SET, + chunk_origin, NULL, chunk_dims, NULL)); + + /* the IO processor creates the dataset and adds the common attributes + when writing its own data, otherwise the dataset is reopened */ + if (proc == myproc) + { + IOHDF5_ERROR (dataset = H5Dcreate (file, datasetname, + iohdf5_type, filespace, H5P_DEFAULT)); + IOHDF5Util_DumpCommonAttributes (GH, index, timelevel, &geom[2*sdim], + slab, dataset); + } + else + { + IOHDF5_ERROR (dataset = H5Dopen (file, datasetname)); + } + + /* write the data */ + IOHDF5_ERROR (H5Dwrite (dataset, iohdf5_type, memspace, filespace, + H5P_DEFAULT, outme)); + + /* and close the file dataspace */ + IOHDF5_ERROR (H5Sclose (filespace)); + } + else + { + /* the IO processor creates the chunk group and adds common attributes */ + if (proc == myproc) + { + IOHDF5_ERROR (group = H5Gcreate (file, datasetname, 0)); + IOHDF5Util_DumpCommonAttributes (GH, index, timelevel, &geom[2*sdim], + slab, group); + IOHDF5_ERROR (H5Gclose (group)); + } + + /* now the chunk dataset for each processor is created within the group */ + chunkname = (char *) malloc (strlen (datasetname) + 20); + sprintf (chunkname, "%s/chunk%d", datasetname, proc - myproc); + /* create the chunk dataset and dump the chunk data */ + IOHDF5_ERROR (dataset = H5Dcreate (file, chunkname, + iohdf5_type, memspace, H5P_DEFAULT)); + IOHDF5_ERROR (H5Dwrite (dataset, iohdf5_type, H5S_ALL, H5S_ALL, + H5P_DEFAULT, outme)); + + /* add the "origin" attribute for the chunk */ + WRITE_ATTRIBUTE ("chunk_origin", geom, dataset, myGH->array_dataspace, sdim, + IOHDF5_INT); + + free (chunkname); + } + + /* close the dataset and the memspace */ + IOHDF5_ERROR (H5Dclose (dataset)); + IOHDF5_ERROR (H5Sclose (memspace)); + + /* free allocated resources */ + free (datasetname); + free (fullname); + free (chunk_origin); + free (chunk_dims); +} + + +#ifdef CCTK_MPI +#ifdef HAVE_PARALLEL +/*@@ + @routine IOHDF5Util_collectiveDump + @author Thomas Radke + @date May 21 1999 + @desc + All processors dump their data into an unchunked file. + @enddesc + + @calls IOHDF5Util_DumpCommonAttributes + + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var index + @vdesc global index of the variable to be dumped + @vtype int + @vio in + @endvar + @var timelevel + @vdesc the timelevel to store + @vtype int + @vio in + @endvar + @var outme + @vdesc pointer to the chunk to dump + @vtype void * + @vio in + @endvar + @var geom + @vdesc bounds and sizes of the output + @vtype CCTK_INT[3*dim] + @vio in + @vcomment geom[0*dim..1*dim-1]: lower bounds + geom[1*dim..2*dim-1]: (local) data sizes + geom[2*dim..3*dim-1]: global space + @endvar + @var hdf5io_type + @vdesc the HDF5 datatype identifier + @vtype int + @vio in + @endvar + @var file + @vdesc the HDF5 file handle + @vtype hid_t + @vio in + @endvar +@@*/ +static void IOHDF5Util_collectiveDump (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + void *outme, + CCTK_INT *geom, + int hdf5io_type, + hid_t file) +{ + DECLARE_CCTK_PARAMETERS + int i, dim; + ioHDF5GH *myGH; + hid_t dataset, memspace, filespace; + char *name, datasetname[128]; + hssize_t *chunk_origin; + hsize_t *chunk_dims, *file_dims; + + + myGH = (ioHDF5GH *) CCTK_GHExtension (GH, "IOHDF5Util"); + + /* get the dimension of the variable */ + sdim = slab->sdim; + + chunk_origin = (hssize_t *) malloc (sdim * sizeof (hssize_t)); + chunk_dims = (hsize_t *) malloc (2*sdim * sizeof (hsize_t)); + file_dims = chunk_dims + sdim; + + /* HDF5 needs it in reverse order */ + for (i = 0; i < sdim; i++) + { + chunk_origin[i] = geom[1*sdim - 1 - i]; + chunk_dims [i] = geom[2*sdim - 1 - i]; + file_dims [i] = geom[3*sdim - 1 - i]; + } + + /* build the unique dataset name */ + name = CCTK_FullName (index); + sprintf (datasetname, "%s@%d@%d", name, GH->cctk_iteration, timelevel); + free (name); + + /* create the memspace according to chunk dims */ + IOHDF5_ERROR (memspace = H5Screate_simple (sdim, chunk_dims, NULL)); + + /* create the (global) filespace and set the hyperslab for the chunk */ + IOHDF5_ERROR (filespace = H5Screate_simple (sdim, file_dims, NULL)); + IOHDF5_ERROR (H5Sselect_hyperslab (filespace, H5S_SELECT_SET, + chunk_origin, NULL, chunk_dims, NULL)); + + /* the IO processor creates the dataset and adds the common attributes + when writing its own data, otherwise the dataset is reopened */ + /* if restart from recovery delete an already existing group + so that it can be created as anew */ + if (check_exisiting_objects) + { + IOHDF5_ERROR (H5Eset_auto (NULL, NULL)); + H5Gunlink (file, datasetname); + IOHDF5_ERROR (H5Eset_auto (myGH->print_error_fn, myGH->print_error_fn_arg)); + } + IOHDF5_ERROR (dataset = H5Dcreate (file, datasetname, + hdf5io_type, filespace, H5P_DEFAULT)); + if (CCTK_MyProc (GH) == 0) + { + IOHDF5Util_DumpCommonAttributes (GH, index, timelevel, &geom[2*sdim], + slab, dataset); + } + + /* write the data */ + IOHDF5_ERROR (H5Dwrite (dataset, hdf5io_type, memspace, + filespace, H5P_DEFAULT, outme)); + + /* and close the file dataspace */ + IOHDF5_ERROR (H5Sclose (filespace)); + + /* close the dataset and the memspace */ + IOHDF5_ERROR (H5Dclose (dataset)); + IOHDF5_ERROR (H5Sclose (memspace)); + + free (chunk_origin); + free (chunk_dims); +} +#endif /* HAVE_PARALLEL */ +#endif /* MPI */ diff --git a/src/ParseGeometry.c b/src/ParseGeometry.c new file mode 100644 index 0000000..83c0af0 --- /dev/null +++ b/src/ParseGeometry.c @@ -0,0 +1,347 @@ + +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "cctk.h" +#include "cctk_Arguments.h" +#include "cctk_Parameters.h" +#include "cctk_GNU.h" +#include "util_String.h" + +#include "ioHDF5UtilGH.h" + +/* Do not apply the USE marcro, parser will get confused + by curly braces in quoted strings */ + +/* CCTK_NO_AUTOUSE_MACRO */ + +/*$#define HEAVYDEBUG$*/ + + +int GeometryParserH5stream(const char *before, char **outname, ioHDF5Geo_t *geo) +{ + DECLARE_CCTK_PARAMETERS + + regmatch_t pmatch[3],gmatch[6]; + int matched, ierr=0,retval=0, verb=0, deb=0; + + const char *argument; + char *varname=NULL, *geo_s=NULL; + const char *token=NULL; + char *dim_s=NULL, *ori_s=NULL, *dir_s=NULL; + char *len_s=NULL, *down_s=NULL; + + int dim = 0, idim, iidim; + int index; + + char info[8000]; + + /* Debugging switches */ +#if 0 + if (CCTK_Equals(h5verbose,"debug")) { + deb =1; + verb =1; + } + if (CCTK_Equals(h5verbose,"yes")) + verb =1; +#endif + + sprintf(info,"\n\nHDF5stream GeometryParser \nargument: >%s<\n",before); + + if((matched = CCTK_RegexMatch(before, + "([A-Za-z][A-Za-z0-9_]*::[A-Za-z][A-Za-z0-9_]*)\\[?(.*)?\\]?", 3, pmatch)) != 0) { + if (matched<0) CCTK_WARN(1,"Error matching in GeometryParser"); +#ifdef HEAVYDEBUG + printf("matched %d rm_so/rm_eo: %d %d; %d %d\n", + matched, + (int)pmatch[1].rm_so,(int)pmatch[1].rm_eo, + (int)pmatch[2].rm_so,(int)pmatch[2].rm_eo); +#endif + + if(pmatch[1].rm_so != -1 && + (pmatch[1].rm_eo-pmatch[1].rm_so > 0)) + { + varname = (char*) malloc((int)(pmatch[1].rm_eo-pmatch[1].rm_so+1) + *sizeof(char)); + strncpy(varname,before+pmatch[1].rm_so, + (int)(pmatch[1].rm_eo-pmatch[1].rm_so)); + varname[(int)(pmatch[1].rm_eo-pmatch[1].rm_so)]='\0'; + + *outname = varname; + sprintf(info,"%sOUTNAME : >%s< \n",info,*outname); + + if ((index = CCTK_VarIndex(varname))>=0) { + geo->vdim = CCTK_GroupDimI(CCTK_GroupIndexFromVarI(index)); + } + else + { + sprintf(info,"%sOUTNAME : no appropriate gridfunction found:>%s< \n", + info,*outname); + geo->vdim = -1; + } + } + else + { + CCTK_WARN(1,"No variable name found in IOStreamedHDF5::out_vars"); + retval = -1; + return(retval); + } + + if(pmatch[2].rm_so != -1 && + (pmatch[2].rm_eo-pmatch[2].rm_so > 0)) + { + geo_s = (char*) malloc((int)(pmatch[2].rm_eo-pmatch[2].rm_so+1) + *sizeof(char)); + strncpy(geo_s,before+pmatch[2].rm_so, + (int)(pmatch[2].rm_eo-pmatch[2].rm_so)); + geo_s[(int)(pmatch[2].rm_eo-pmatch[2].rm_so)]='\0'; + sprintf(info,"%sGEOMETRY: >%s< \n",info,geo_s); + } else { + sprintf(info,"%sGEOMETRY: No geometry specified -> use default\n",info); + } + } + else { + retval=-1; + return(retval); + } + + /* Pattern match the hyperslab string geo_s*/ + if (geo_s) { + if((matched = CCTK_RegexMatch(geo_s, + "\\{(.*)\\},\\{(.*)\\},\\{(.*)\\},\\{(.*)\\},\\{(.*)\\}", + 6, gmatch)) != 0) + { + + /* SLAB DIMENSION */ + if(gmatch[1].rm_so != -1 && + (gmatch[1].rm_eo-gmatch[1].rm_so > 0)) + { + dim_s = (char*) malloc((int)(gmatch[1].rm_eo-gmatch[1].rm_so+1) + *sizeof(char)); + strncpy(dim_s,geo_s+gmatch[1].rm_so, + (int)(gmatch[1].rm_eo-gmatch[1].rm_so)); + dim_s[(int)(gmatch[1].rm_eo-gmatch[1].rm_so)]='\0'; + + if (deb) + sprintf(info,"%sDIMENSION: >%s< \n",info,dim_s); + + dim = atoi(dim_s); + geo->sdim = dim; + } + else { + if (deb) + sprintf(info,"%s","DIMENSION: No dimension\n"); + retval = -1; + } + free(dim_s); + + /* DIRECTION */ + ierr = -1; + if(gmatch[2].rm_so != -1 && + (gmatch[2].rm_eo-gmatch[2].rm_so > 0)) + { + dir_s = (char*) malloc((int)(gmatch[2].rm_eo-gmatch[2].rm_so+1) + *sizeof(char)); + strncpy(dir_s,geo_s+gmatch[2].rm_so, + (int)(gmatch[2].rm_eo-gmatch[2].rm_so)); + dir_s[(int)(gmatch[2].rm_eo-gmatch[2].rm_so)]='\0'; + + if (deb) + sprintf(info,"%sDIRECTION: >%s< \n",info,dir_s); + + idim = 0; + argument = dir_s; + while((token = Util_StrSep(&argument,","))) { + geo->direction[idim++]=atoi(token); + } + geo->direction[idim] = atoi(argument); + if (idim==dim-1) ierr = 0; + } + if (ierr<0) { + sprintf(info,"%sDIRECTION: dimension not consistent: >%s< (Use default)\n", + info,dir_s); + } + free(dir_s); + + /* ORIGIN */ + ierr = -1; + if(gmatch[3].rm_so != -1 && + (gmatch[3].rm_eo-gmatch[3].rm_so > 0)) + { + ori_s = (char*) malloc((int)(gmatch[3].rm_eo-gmatch[3].rm_so+1) + *sizeof(char)); + strncpy(ori_s,geo_s+gmatch[3].rm_so, + (int)(gmatch[3].rm_eo-gmatch[3].rm_so)); + ori_s[(int)(gmatch[3].rm_eo-gmatch[3].rm_so)]='\0'; + + if (deb) + sprintf(info,"%sDIRECTION: >%s< \n",info,ori_s); + + idim=0; + argument = ori_s; + while((token = Util_StrSep(&argument,","))) { + geo->origin[idim++]=atoi(token); + } + geo->origin[idim] = atoi(argument); + if (idim==(geo->vdim-1)) ierr = 0; + else if (idim==0) { + for (iidim=1;iidim<geo->vdim;iidim++) + geo->origin[iidim] = geo->origin[0]; + ierr = 0; + } + } + if (ierr<0) { + sprintf(info,"%sORIGIN: dimension not consistent: >%s< \n", + info,ori_s); + ierr = -1; + } + free(ori_s); + + /* LENGTH */ + ierr = -1; + if(gmatch[4].rm_so != -1 && + (gmatch[4].rm_eo-gmatch[4].rm_so > 0)) + { + len_s = (char*) malloc((int)(gmatch[4].rm_eo-gmatch[4].rm_so+1) + *sizeof(char)); + strncpy(len_s,geo_s+gmatch[4].rm_so, + (int)(gmatch[4].rm_eo-gmatch[4].rm_so)); + len_s[(int)(gmatch[4].rm_eo-gmatch[4].rm_so)]='\0'; + + if (deb) + sprintf(info,"%sLENGTH: >%s< \n",info,len_s); + + idim = 0; + argument = len_s; + while((token = Util_StrSep(&argument,","))) { + geo->length[idim++]=atoi(token); + } + geo->length[idim] = atoi(argument); + if (idim==dim-1) ierr = 0; + else if (idim==0) { + for (iidim=1;iidim<dim;iidim++) + geo->length[iidim] = geo->length[0]; + ierr = 0; + } + } + if (ierr<0) { + sprintf(info,"%sLENGTH: dimension not consistent: >%s< (Use Default)\n", + info,len_s); + } + free(len_s); + + + /* DOWNSAMPLING */ + ierr = -1; + if(gmatch[5].rm_so != -1 && + (gmatch[5].rm_eo-gmatch[5].rm_so > 0)) { + + down_s = (char*) malloc((int)(gmatch[5].rm_eo-gmatch[5].rm_so+1) + *sizeof(char)); + strncpy(down_s,geo_s+gmatch[5].rm_so, + (int)(gmatch[5].rm_eo-gmatch[5].rm_so)); + down_s[(int)(gmatch[5].rm_eo-gmatch[5].rm_so)]='\0'; + + if (deb) + sprintf(info,"%sDOWNSAMPLING: >%s< \n",info,down_s); + + + idim=0; + argument = down_s; + while((token = Util_StrSep(&argument,","))) { + geo->downsample[idim++]=atoi(token); + } + geo->downsample[idim] = atoi(argument); + if (idim==dim-1) ierr =0; + else if (idim==0) { + for (iidim=1;iidim<dim;iidim++) + geo->downsample[iidim] = geo->downsample[0]; + ierr = 0; + } + } + if (ierr<0) { + sprintf(info,"%sDOWNSAMPLING: dimension not consistent: >%s< \n", + info,down_s); + ierr =-1; + } + free(down_s); + } + if (geo_s) free(geo_s); + } + + if (verb) { + sprintf(info, "%sGeometry Data: \n",info); + sprintf(info, "%s Argument/Slab dimension: %d / %d \n", + info,geo->vdim,geo->sdim); + sprintf(info, "%s Origin: ",info); + for (idim=0;idim<geo->vdim;idim++) + sprintf(info,"%s %d ",info,geo->origin[idim]); + sprintf(info,"%s\n Downs : ",info); + for (idim=0;idim<geo->sdim;idim++) + sprintf(info,"%s %d ",info,geo->downsample[idim]); + sprintf(info,"%s\n Length: ",info); + for (idim=0;idim<geo->sdim;idim++) + sprintf(info,"%s %d ",info,geo->length[idim]); + sprintf(info,"%s\n Dirs : ",info); + for (idim=0;idim<geo->sdim;idim++) + sprintf(info,"%s %d ",info,geo->direction[idim]); + sprintf(info,"%s\n\n",info); + + printf("%s",info); + } + + USE_CCTK_PARAMETERS + + return(retval); +} + +void IOHDF5Util_DefaultGeo (ioHDF5Geo_t *geo) +{ + DECLARE_CCTK_PARAMETERS + int idim; + const char *argument, *token; + + geo->vdim = -1; + geo->sdim = slabdim; + + /* Origin */ + idim=0; + argument = origin; + while((token = Util_StrSep(&argument,","))) { + geo->origin[idim++]=atoi(token); + } + geo->downsample[idim] = atoi(argument); + + /* Direction */ + idim=0; + argument = direction; + while((token = Util_StrSep(&argument,","))) { + geo->direction[idim++]=atoi(token); + } + geo->downsample[idim] = atoi(argument); + + /* Downsample */ + idim=0; + argument = origin; + while((token = Util_StrSep(&argument,","))) { + geo->downsample[idim++]=atoi(token); + } + geo->downsample[idim] = atoi(argument); + + /* Length */ + idim=0; + argument = origin; + while((token = Util_StrSep(&argument,","))) { + geo->length[idim++]=atoi(token); + } + geo->length[idim] = atoi(argument); + + for (idim=0;idim<IOHDF5_MAXDIM;idim++) { + geo->direction[idim] = idim; + geo->origin[idim] = 0; + geo->length[idim] =-1; + geo->downsample[idim] = 1; + } + USE_CCTK_PARAMETERS +} diff --git a/src/ParseVars.c b/src/ParseVars.c new file mode 100644 index 0000000..d3ca305 --- /dev/null +++ b/src/ParseVars.c @@ -0,0 +1,212 @@ + /*@@ + @file GHExtension.c + @date Tue 9th Feb 1999 + @author Gabrielle Allen + @desc + IOUtil GH extension stuff. + @enddesc + @history + @endhistory + @@*/ + +/*#define DEBUG_IO*/ + +#include <stdlib.h> +#include <string.h> +#include <stdio.h> + +#include "cctk.h" +#include "util_String.h" +#include "cctk_Parameters.h" + +#include "ioHDF5UtilGH.h" + + +static char *rcsid = "$Header$"; +CCTK_FILEVERSION(BetaThorns_IOStreamedHDF5_ParseVars_c) + +/* =============================================================== + utility routines used by other IO thorns + ===============================================================*/ + + /*@@ + @routine IOUtil_ParseVarsForOutput + @date Sat March 6 1999 + @author Gabrielle Allen + @desc + Sets each flag in the do_output[] do_output to true + if var[i] could be found in the list of variable names. + var_list my also contain group names of variables as well as the + special keyword "all" which indicates that output is requested + on all variables. + @enddesc + @history + @endhistory + @var var_list + @vdesc list of variables and/or group names + @vtype const char * + @vio in + @endvar + @var do_output + @vdesc do_output of flags indicating output was requested for var[i] + @vtype char [] + @vio out + @endvar +@@*/ + +int IOHDF5Util_ParseVarsForOutput (const char *var_list, + char do_output[], + ioHDF5Geo_t geo_output[]) +{ + char *outname=NULL; + char *before=NULL; + char *after=NULL; + char *splitstring; + int first,groupnum,i,last,index,varnum; + ioHDF5Geo_t geo_tmp; + int ierr=0; +int GeometryParserH5stream(const char *before, char **outname, ioHDF5Geo_t *geo); + + + /* First initialise every variable to no output */ + for (i=0; i<CCTK_NumVars(); i++) + { + do_output[i] = 0; + } + + /* Initialize geometry structure with default */ + + IOHDF5Util_DefaultGeo (&geo_tmp); + + splitstring = (char *) var_list; + + /* split string at blanks */ + while(Util_SplitString(&before,&after,splitstring," ")==0) { + +#ifdef DEBUG_IO + printf(" String is -%s-\n",splitstring); + printf(" Split is -%s- and -%s-\n",before,after); +#endif + + if (CCTK_Equals(before," ")) continue; + + if (strlen(before) > 0) + { + /* Extract geometry information */ + ierr=GeometryParserH5stream(before, &outname, &geo_tmp); + if ((ierr<0)||(outname==NULL)) { + CCTK_WARN(1,"GeometryParserH5stream failed."); + return(-1); + } + + /* Look for any special tokens */ + if (CCTK_Equals(outname,"all")) + { + int ilab; + for (ilab=0;ilab<CCTK_NumVars();ilab++) { + geo_output[ilab]= geo_tmp; + do_output[ilab] = 1; + } + } + else + { + /* See if this name is implementation::variable */ + varnum = CCTK_VarIndex(outname); + if ( varnum < 0 ) + { + /* See if this name is implementation::group */ + groupnum = CCTK_GroupIndex(outname); + if (groupnum < 0) + { + char *msg; + msg = (char *)malloc((100+strlen(outname))*sizeof(char)); + sprintf(msg,"Ignoring <%s> in IO string (invalid token)",outname); + CCTK_WARN(2,msg); + free(msg); + } + else + { + /* We have a group so now need all the variables in the group */ + first = CCTK_FirstVarIndexI(groupnum); + last = first+CCTK_NumVarsInGroupI(groupnum)-1; + for (index=first;index<=last;index++) { + geo_output[index]=geo_tmp; + do_output[index] =1; + } + } + } + else + { + geo_output[varnum]= geo_tmp; + do_output[varnum] = 1; + } + } + } + if(before) + { + free(before); + before = NULL; + } + if(outname) + { + free(outname); + outname = NULL; + } + if(splitstring != var_list) + { + free(splitstring); + } + + splitstring = after; + + } + + if (strlen(splitstring)>0) { + + /* Extract geometry information */ + ierr=GeometryParserH5stream(splitstring, &outname, &geo_tmp); + if (ierr<0) printf("GeometryParser failed: >%s<\n",before); + + /* Special cases */ + if (CCTK_Equals(outname,"all")) { + int ilab; + for (ilab=0;ilab<CCTK_NumVars();ilab++) { + geo_output[ilab]=geo_tmp; + do_output [ilab]=1; + } + } else { + + varnum = CCTK_VarIndex(outname); + if (varnum < 0) { + groupnum = CCTK_GroupIndex(outname); + if (groupnum < 0) { + char *msg; + msg = (char *)malloc((100+strlen(outname))*sizeof(char)); + sprintf(msg,"Ignoring %s in IO string (invalid token)",outname); + CCTK_WARN(2,msg); + free(msg); + } else { + /* We have a group so now need all the variables in the group */ + first = CCTK_FirstVarIndexI(groupnum); + last = first+CCTK_NumVarsInGroupI(groupnum)-1; + for (index=first;index<=last;index++) { + geo_output[index]=geo_tmp; + do_output[index] =1; + } + } + } + else { + geo_output[varnum]=geo_tmp; + do_output[varnum] =1; + } + } + } + + if (before) free(before); + if (after) free(after); + if(splitstring != var_list) + { + free(splitstring); + } + return(0); +} diff --git a/src/RecoverVar.c b/src/RecoverVar.c new file mode 100644 index 0000000..fe1e015 --- /dev/null +++ b/src/RecoverVar.c @@ -0,0 +1,809 @@ + /*@@ + @file RecoverVar.c + @date Thu Jun 18 16:34:59 1998 + @author Tom Goodale + @desc + Routines to recover variables from a given HDF5 data or + checkpoint file. + These routines are used by other HDF5 IO methods. + @enddesc + @version $Id$ + @@*/ + + +#include <stdlib.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "CactusBase/IOUtil/src/ioGH.h" +#include "CactusBase/IOUtil/src/ioutil_CheckpointRecovery.h" +#include "CactusPUGH/PUGH/src/include/pugh.h" +#include "ioHDF5UtilGH.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 IOHDF5Util_RecoverVariables + @date Fri Jun 19 09:19:48 1998 + @author Tom Goodale + @desc + Reads in data from an open HDF5 file. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype cGH * + @vio in + @endvar + @var fileinfo + @vdesc pointer to info structure describing the HDF5 file + @vtype fileinfo_t * + @vio in + @endvar +@@*/ +int IOHDF5Util_RecoverVariables (cGH *GH, + fileinfo_t *fileinfo) +{ + 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 = fileinfo->ioproc; + info.unchunked = fileinfo->unchunked; + info.ioproc_every = fileinfo->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) == fileinfo->ioproc) + { +#endif /* CCTK_MPI */ + + /* iterate over all datasets starting from "/" in the HDF5 file */ + IOHDF5_ERROR (H5Giterate (fileinfo->file, "/", 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 < fileinfo->ioproc_every; proc++) + for (proc = fileinfo->ioproc + 1; + proc < fileinfo->ioproc + fileinfo->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, fileinfo->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, fileinfo->ioproc, STARTUPBASE, + pughGH->PUGH_COMM_WORLD, &ms)); + } + } +#endif /* CCTK_MPI */ + + return (0); +} + + +/* NOTE: Although we could read the GH extensions from multiple recovery files + in parallel, this is done only on by processor 0 here. + Broadcasting the GH extensions is found faster than + sending it in a loop from each IO processor to all the non IOs + (don't have subcommunicators yet) */ +int IOHDF5Util_RecoverGHextensions (cGH *GH, + fileinfo_t *fileinfo) +{ + hid_t group; + CCTK_REAL4 real4Buffer; + CCTK_INT4 int4Buffer[3]; + + + if (CCTK_MyProc (GH) == 0) + { + /* all the important global attributes and GH extensions + are stored in the GLOBAL_ATTRIBUTES_GROUP group */ + group = H5Gopen (fileinfo->file, GLOBAL_ATTRIBUTES_GROUP); + int4Buffer[0] = group >= 0; + + if (int4Buffer[0]) + { + + READ_ATTRIBUTE (group, "cctk_iteration", IOHDF5_INT4, &int4Buffer[1]); + READ_ATTRIBUTE (group, "main_loop_index", IOHDF5_INT4, &int4Buffer[2]); + READ_ATTRIBUTE (group, "cctk_time", IOHDF5_REAL4, &real4Buffer); + + IOHDF5_ERROR (H5Gclose (group)); + + } + else + { + CCTK_WARN (1, "Can't find global attributes group. " + "Is this really a Cactus HDF5 datafile ?"); + } + } + +#ifdef CCTK_MPI + /* Broadcast the GH extensions to all processors */ + /* NOTE: We have to use MPI_COMM_WORLD here + because PUGH_COMM_WORLD is not yet set up at parameter recovery time. + We also assume that PUGH_MPI_INT4 is a compile-time defined datatype. */ + CACTUS_MPI_ERROR (MPI_Bcast (int4Buffer, 3, PUGH_MPI_INT4, 0,MPI_COMM_WORLD)); + if (int4Buffer[0]) + { + CACTUS_MPI_ERROR (MPI_Bcast (&real4Buffer, 1, PUGH_MPI_REAL4, 0, + MPI_COMM_WORLD)); + } +#endif + + if (int4Buffer[0]) + { + GH->cctk_time = real4Buffer; + GH->cctk_iteration = int4Buffer[1]; + CCTK_SetMainLoopIndex ((int) int4Buffer[2]); + } + + /* Return 0 for success otherwise negative */ + return (int4Buffer[0] ? 0 : -1); +} + + +/* NOTE: Although we could read the parameters from multiple recovery files + in parallel, this is done only on by processor 0 here. + Broadcasting the complete parameter string is found faster than + sending it in a loop from each IO processor to all the non IOs + (don't have subcommunicators yet) */ +int IOHDF5Util_RecoverParameters (fileinfo_t *fileinfo) +{ + DECLARE_CCTK_PARAMETERS + hid_t group, attr, atype; + char *parameters; + CCTK_INT4 parameterSize; + cGH *GH = NULL; /* There's no cGH set up yet so we pass + a NULL pointer to CCTK_MyProc() */ + + /* To make the compiler happy */ + group = attr = 0; + parameters = NULL; + + if (CCTK_MyProc (GH) == 0) + { + if (verbose) + { + CCTK_VInfo (CCTK_THORNSTRING, "Recovering parameters from checkpoint " + "file '%s'", fileinfo->filename); + } + + /* the parameters are stored in the CACTUS_PARAMETERS_GROUP group + in the attribute ALL_PARAMETERS */ + group = H5Gopen (fileinfo->file, CACTUS_PARAMETERS_GROUP); + if (group > 0) + { + attr = H5Aopen_name (group, ALL_PARAMETERS); + } + + if (group > 0 && attr > 0) + { + IOHDF5_ERROR (atype = H5Aget_type (attr)); + parameterSize = H5Tget_size (atype); + parameters = (char *) malloc (parameterSize + 1); + IOHDF5_ERROR (H5Aread (attr, atype, parameters)); + parameters[parameterSize] = 0; + IOHDF5_ERROR (H5Tclose (atype)); + } + else + { + CCTK_WARN (1, "Can't find global parameters. " + "Is this really a Cactus HDF5 checkpoint file ?"); + } + + if (attr > 0) + { + IOHDF5_ERROR (H5Aclose (attr)); + } + if (group > 0) + { + IOHDF5_ERROR (H5Gclose (group)); + } + } + +#ifdef CCTK_MPI + /* Broadcast the parameter buffer size to all processors */ + /* NOTE: We have to use MPI_COMM_WORLD here + because PUGH_COMM_WORLD is not yet set up at parameter recovery time. + We also assume that PUGH_MPI_INT4 is a compile-time defined datatype. */ + CACTUS_MPI_ERROR (MPI_Bcast (¶meterSize, 1, PUGH_MPI_INT4, 0, + MPI_COMM_WORLD)); +#endif + + if (parameterSize > 0) + { +#ifdef CCTK_MPI + if (CCTK_MyProc (GH) != 0) + { + parameters = (char *) malloc (parameterSize + 1); + } + + CACTUS_MPI_ERROR (MPI_Bcast (parameters, parameterSize + 1, PUGH_MPI_CHAR, + 0, MPI_COMM_WORLD)); +#endif + + IOUtil_SetAllParameters (parameters); + + free (parameters); + } + + /* Return positive value for success otherwise negative */ + return (parameterSize > 0 ? 1 : -1); +} + + +/************************************************************************/ +/* 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, "%[^ ] timelevel %d at iteration %d", + fullvarname, timelevel, &i) != 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) + { + IOHDF5_ERROR (dataset = H5Dopen (dataset, "chunk0")); + } + 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); + + 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 */ + IOHDF5_ERROR (dataspace = H5Dget_space (dataset)); + IOHDF5_ERROR (ndims_stored = H5Sget_simple_extent_ndims (dataspace)); + dims_stored = (hsize_t *) malloc (ndims_stored * sizeof (hsize_t)); + 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; + } + 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) + { + IOHDF5_ERROR (H5Dclose (dataset)); + } + + *grouptype = groupdata.grouptype; + + return (0); +} + + +static int IOHDF5Util_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 ... */ + 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 /* CCTK_MPI */ + + return (0); +} + + +static int IOHDF5Util_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) + { + 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); + IOHDF5_ERROR (chunk = H5Dopen (dataset, chunkname)); + IOHDF5_ERROR (H5Dread (chunk, rec_info->iohdf5_type, H5S_ALL, + H5S_ALL, H5P_DEFAULT, data)); + 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]; + } + + IOHDF5_ERROR (filespace = H5Dget_space (dataset)); + IOHDF5_ERROR (memspace = H5Screate_simple (dim, chunk_dims, NULL)); + IOHDF5_ERROR (H5Sselect_hyperslab (filespace, H5S_SELECT_SET, + chunk_origin, NULL, chunk_dims, NULL)); + + IOHDF5_ERROR (H5Dread (dataset, rec_info->iohdf5_type, memspace, + filespace, H5P_DEFAULT, data)); + + IOHDF5_ERROR (H5Sclose (memspace)); + 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 /* CCTK_MPI */ + + return (0); +} + + +static herr_t processDataset (hid_t group, + const char *datasetname, + void *arg) +{ + DECLARE_CCTK_PARAMETERS + pGH *pughGH; + ioGH *ioUtilGH; + ioHDF5UtilGH *h5UtilGH; + 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, CACTUS_PARAMETERS_GROUP) || + ! strcmp (datasetname, GLOBAL_ATTRIBUTES_GROUP)) + { + return (0); + } + + /* Get the handle for PUGH, IOUtil, and IOHDF5Util extensions */ + pughGH = PUGH_pGH (it_info->GH); + ioUtilGH = (ioGH *) CCTK_GHExtension (it_info->GH, "IO"); + h5UtilGH = (ioHDF5UtilGH *) CCTK_GHExtension (it_info->GH, "IOHDF5Util"); + + IOHDF5_ERROR (H5Gget_objinfo (group, datasetname, 0, &object_info)); + is_group = object_info.type == H5G_GROUP; + if (is_group) + { + IOHDF5_ERROR (dataset = H5Gopen (group, datasetname)); + } + else + { + 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 = h5UtilGH->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 variable datatype '%s'", + CCTK_VarTypeName (CCTK_VarTypeI (index))); + return (0); + } + + /* Read in the data */ + switch (gtype) + { + case CCTK_SCALAR: + IOHDF5Util_RestoreGS (dataset, index, timelevel, &rec_info); + break; + case CCTK_GF: + case CCTK_ARRAY: + IOHDF5Util_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) + { + IOHDF5_ERROR (H5Gclose (dataset)); + } + else + { + IOHDF5_ERROR (H5Dclose (dataset)); + } + + return (0); +} diff --git a/src/Startup.c b/src/Startup.c new file mode 100644 index 0000000..133eb9e --- /dev/null +++ b/src/Startup.c @@ -0,0 +1,147 @@ + /*@@ + @file Startup.c + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + Startup and termination routines for IOHDF5Util. + @enddesc + @version $Id$ +@@*/ + + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "ioHDF5UtilGH.h" + + +/* local function prototypes */ +static void *IOHDF5Util_SetupGH (tFleshConfig *config, + int convergence_level, + cGH *GH); + + + /*@@ + @routine IOHDF5Util_Startup + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + The startup registration routine for IOHDF5Util. + Registers the GH extensions needed for IOHDF5Util + and the routine to set it up. + @enddesc + + @calls CCTK_RegisterGHExtensionSetupGH +@@*/ +void IOHDF5Util_Startup (void) +{ + /* check dynamic thorn dependencies */ + if (CCTK_GHExtensionHandle ("IO") < 0) + { + CCTK_WARN (1, "Thorn IOUtil was not activated. " + "No HDF5 IO methods will be registered."); + } + else if (CCTK_GHExtensionHandle ("PUGH") < 0) + { + CCTK_WARN (1, "Thorn PUGH was not activated. " + "No HDF5 IO methods will be registered."); + } + else + { + CCTK_RegisterGHExtensionSetupGH (CCTK_RegisterGHExtension ("IOHDF5Util"), + IOHDF5Util_SetupGH); + } +} + + + /*@@ + @routine IOHDF5Util_Terminate + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + The termination registration routine for IOHDF5Util. + It frees any resources kept on the GH extensions. + @enddesc + @var GH + @vdesc Pointer to CCTK GH + @vtype cGH * + @vio in + @endvar +@@*/ +void IOHDF5Util_Terminate (cGH *GH) +{ + ioHDF5UtilGH *myGH; + + + myGH = (ioHDF5UtilGH *) CCTK_GHExtension (GH, "IOHDF5Util"); + if (myGH) + { + /* close the dataspaces used to write scalar and array attributes */ + if (myGH->scalar_dataspace >= 0) + { + IOHDF5_ERROR (H5Sclose (myGH->scalar_dataspace)); + } + if (myGH->array_dataspace >= 0) + { + IOHDF5_ERROR (H5Sclose (myGH->array_dataspace)); + } + + /* close the predefined complex and string datatypes */ + if (myGH->IOHDF5_COMPLEX >= 0) + { + IOHDF5_ERROR (H5Tclose (myGH->IOHDF5_COMPLEX)); + } + if (myGH->IOHDF5_STRING >= 0) + { + IOHDF5_ERROR (H5Tclose (myGH->IOHDF5_STRING)); + } + } +} + + +/****************************************************************************/ +/* local routines */ +/****************************************************************************/ + /*@@ + @routine IOHDF5Util_SetupGH + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + Allocate a GH extension structure for IOHDF5Util and set it up. + @enddesc + + @returntype void * + @returndesc + pointer to the allocated GH extension structure + @endreturndesc +@@*/ +static void *IOHDF5Util_SetupGH (tFleshConfig *config, + int convergence_level, + cGH *GH) +{ + DECLARE_CCTK_PARAMETERS + ioHDF5UtilGH *myGH; + + + myGH = (ioHDF5UtilGH *) malloc (sizeof (ioHDF5UtilGH)); + + /* save the original error printing routine and its argument */ + IOHDF5_ERROR (H5Eget_auto (&myGH->print_error_fn, &myGH->print_error_fn_arg)); + + /* predefine dataspaces for writing scalar and array attributes */ + /* the dimension of the array dataspace is set when used */ + IOHDF5_ERROR (myGH->scalar_dataspace = H5Screate (H5S_SCALAR)); + IOHDF5_ERROR (myGH->array_dataspace = H5Screate (H5S_SIMPLE)); + + /* predefine a IOHDF5_COMPLEX datatype */ + IOHDF5_ERROR (myGH->IOHDF5_COMPLEX = + H5Tcreate (H5T_COMPOUND, sizeof (CCTK_COMPLEX))); + IOHDF5_ERROR (H5Tinsert (myGH->IOHDF5_COMPLEX, "real", + offsetof (CCTK_COMPLEX, Re), IOHDF5_REAL)); + IOHDF5_ERROR (H5Tinsert (myGH->IOHDF5_COMPLEX, "imag", + offsetof (CCTK_COMPLEX, Im), IOHDF5_REAL)); + + /* predefine a C string datatype */ + IOHDF5_ERROR (myGH->IOHDF5_STRING = H5Tcopy (H5T_C_S1)); + + return (myGH); +} diff --git a/src/ioHDF5UtilGH.h b/src/ioHDF5UtilGH.h new file mode 100644 index 0000000..de2de80 --- /dev/null +++ b/src/ioHDF5UtilGH.h @@ -0,0 +1,217 @@ + /*@@ + @header ioHDF5UtilGH.h + @date Fri Oct 6 2000 + @author Thomas Radke + @desc + The GH extensions structure for IOHDF5Util. + @version $Id$ + @@*/ + +#ifndef _IOUTILHDF5_IOUTILHDF5GH_H_ +#define _IOUTILHDF5_IOUTILHDF5GH_H_ + +/* include the HDF5 header file */ +#include <hdf5.h> + + +/*****************************************************************************/ +/* some useful macros */ +/*****************************************************************************/ +/* Check error flags from HDF5 calls */ +#define IOHDF5_ERROR(fn_call) \ + do \ + { \ + int error_code = fn_call; \ + \ + \ + if (error_code < 0) \ + { \ + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, \ + "HDF5 call '%s' returned error code %d\n", \ + #fn_call, error_code); \ + } \ + } while (0) + +/* macro for writing an attribute */ +#define WRITE_ATTRIBUTE(name, value, dataset, dataspace, dim, datatype) \ + do \ + { \ + int len; \ + hid_t attr; \ + void *val = value; \ + hsize_t arrayDim = dim; \ + \ + \ + if (H5Tget_class (datatype) == H5T_STRING) \ + { \ + len = strlen ((char *) val); \ + if (len == 0) /* HDF5 doesn't like zero-len strings */ \ + { \ + len++; \ + } \ + IOHDF5_ERROR (H5Tset_size (datatype, len)); \ + } \ + if (arrayDim > 0) \ + { \ + IOHDF5_ERROR (H5Sset_extent_simple (dataspace, 1, \ + &arrayDim, NULL)); \ + } \ + IOHDF5_ERROR (attr = H5Acreate (dataset, name, datatype, \ + dataspace, H5P_DEFAULT)); \ + IOHDF5_ERROR (H5Awrite (attr, datatype, val)); \ + IOHDF5_ERROR (H5Aclose (attr)); \ + } while (0); + +/* macro for reading an attribute */ +#define READ_ATTRIBUTE(dataset, attrname, requested_type, buffer) \ + do \ + { \ + hid_t attr, attrtype; \ + hsize_t asize = 0; \ + \ + \ + if ((attr = H5Aopen_name (dataset, attrname)) < 0) \ + { \ + CCTK_WARN (1, "Can't find " attrname " attribute"); \ + } \ + if (requested_type == H5T_C_S1) \ + { \ + IOHDF5_ERROR (attrtype = H5Aget_type (attr)); \ + IOHDF5_ERROR (asize = H5Tget_size (attrtype)); \ + if (asize + 1 >= sizeof (buffer)) \ + { \ + CCTK_WARN (1, "Can't read " attrname " attribute (too long)");\ + } \ + } \ + else \ + { \ + attrtype = requested_type; \ + } \ + if (H5Aread (attr, attrtype, buffer) < 0) \ + { \ + CCTK_WARN (1, "Can't read " attrname " attribute"); \ + } \ + if (requested_type == H5T_C_S1) \ + { \ + ((char *) buffer)[asize] = 0; \ + IOHDF5_ERROR (H5Tclose (attrtype)); \ + } \ + IOHDF5_ERROR (H5Aclose (attr)); \ + } while (0) + + +/* Mapping off CCTK datatypes to HDF5 datatypes */ +#define IOHDF5_REAL4 H5T_NATIVE_FLOAT + +#ifdef CCTK_REAL_PRECISION_16 +#define IOHDF5_REAL H5T_NATIVE_LDOUBLE +#elif CCTK_REAL_PRECISION_8 +#define IOHDF5_REAL H5T_NATIVE_DOUBLE +#elif CCTK_REAL_PRECISION_4 +#define IOHDF5_REAL H5T_NATIVE_FLOAT +#endif + +#define IOHDF5_INT (sizeof (CCTK_INT) == sizeof (int) ? \ + H5T_NATIVE_INT : H5T_NATIVE_SHORT) +#define IOHDF5_INT4 (sizeof (int) == 4 ? \ + H5T_NATIVE_INT : H5T_NATIVE_SHORT) +#define IOHDF5_CHAR H5T_NATIVE_CHAR + + +/* Constants to index timers within the timers array */ +#define CP_PARAMETERS_TIMER 0 +#define CP_VARIABLES_TIMER 1 +#define CP_TOTAL_TIMER 2 +#define RECOVERY_TIMER 3 +#define IOHDF5_NUM_TIMERS 4 + +/* names of the groups that hold global attributes and parameters */ +#define GLOBAL_ATTRIBUTES_GROUP "Global Attributes" +#define CACTUS_PARAMETERS_GROUP "Cactus Parameters" +#define ALL_PARAMETERS "All Parameters" + + +/* Geometry information structure for output variable */ +/* FIXME: allocate arrays dynamically */ +#define IOHDF5_MAXDIM 5 +typedef struct +{ + int vdim; + int sdim; + int direction[IOHDF5_MAXDIM]; + int origin[IOHDF5_MAXDIM]; + int length[IOHDF5_MAXDIM]; + int downsample[IOHDF5_MAXDIM]; + int actlen[IOHDF5_MAXDIM]; /* actual index slab length (by PUGHSlab)*/ +} ioHDF5Geo_t; + + +/* Structure describing a given recovery file */ +typedef struct +{ + int is_HDF5_file; /* flag indicating valid file info */ + hid_t file; /* HDF5 file handle */ + char *filename; /* complete file name for recovery */ + int nprocs; /* number of total processors */ + int ioproc; /* the associated IO processor */ + int ioproc_every; /* how many IO processors there are */ + int unchunked; /* whether data was written chunked or unchunked */ +} fileinfo_t; + + +/* IOHDF5Util GH extension structure */ +typedef struct +{ + /* while HDF5 error output is disabled + we keep the original error printing routine and its argument in here */ + H5E_auto_t print_error_fn; + void *print_error_fn_arg; + + /* predefined dataspaces for writing scalar and array attributes */ + hid_t scalar_dataspace, array_dataspace; + + /* predefined datatype for writing CCTK_COMPLEX types */ + hid_t IOHDF5_COMPLEX; + + /* predefined datatype for writing C string string attributes */ + hid_t IOHDF5_STRING; + +} ioHDF5UtilGH; + + +#ifdef __cplusplus +extern "C" +{ +#endif + +/* exported functions */ +int IOHDF5Util_ParseVarsForOutput (const char *var_list, + char do_output[], + ioHDF5Geo_t geo_output[]); +void IOHDF5Util_DefaultGeo (ioHDF5Geo_t *slab); +void IOHDF5Util_DumpParameters (cGH *GH, hid_t group); +void IOHDF5Util_DumpGHExtensions (cGH *GH, hid_t group); +int IOHDF5Util_DumpGH (cGH *GH, int timers[], hid_t file); +void IOHDF5Util_DumpCommonAttributes (cGH *GH, + int index, + int timelevel, + CCTK_INT global_shape[], + ioHDF5Geo_t *slab, + hid_t dataset); +int IOHDF5Util_DumpVar (cGH *GH, + int index, + int timelevel, + ioHDF5Geo_t *slab, + hid_t file, + int check_exisiting_objects); +int IOHDF5Util_RecoverParameters (fileinfo_t *filenfo); +int IOHDF5Util_RecoverGHextensions (cGH *GH, + fileinfo_t *filenfo); +int IOHDF5Util_RecoverVariables (cGH *GH, + fileinfo_t *filenfo); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif /* _IOUTILHDF5_IOUTILHDF5GH_H_ */ diff --git a/src/make.code.defn b/src/make.code.defn new file mode 100644 index 0000000..a3caec7 --- /dev/null +++ b/src/make.code.defn @@ -0,0 +1,5 @@ +# Main make.code.defn file for thorn IOHDF5Util +# $Header$ + +# Source files in this directory +SRCS = Startup.c DumpUtils.c DumpVar.c RecoverVar.c ParseVars.c ParseGeometry.c diff --git a/src/make.configuration.defn b/src/make.configuration.defn new file mode 100644 index 0000000..135342a --- /dev/null +++ b/src/make.configuration.defn @@ -0,0 +1,20 @@ +# make.configuration.defn for IOHDF5Util + +# make sure that IOHDF5Util was configured with HDF5 and PUGH + +ifeq ($(strip $(HDF5_LIBS)), ) +$(NAME): MissingHDF5 +.pseudo: MissingHDF5 +MissingHDF5: + @echo "IOHDF5Util: requires HDF5" + @echo "IOHDF5Util: Please configure with HDF5 or remove IOHDF5Util from Thornlist !" + exit 2 +endif + +ifeq ($(findstring CactusPUGH/PUGH,$(THORNS)),) +.pseudo: MissingPUGHinIOHDF5Util +MissingPUGHinIOHDF5Util: + @echo "IOHDF5Util: requires PUGH" + @echo "IOHDF5Util: Please add CactusPUGH/PUGH or remove IOHDF5Util from Thornlist !" + exit 2 +endif |