aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authortradke <tradke@7842ec3a-9562-4be5-9c5b-06ba18f2b668>2000-10-12 11:45:59 +0000
committertradke <tradke@7842ec3a-9562-4be5-9c5b-06ba18f2b668>2000-10-12 11:45:59 +0000
commit87adfe63f5d5e8373345231a5331b14abfd14cab (patch)
tree91d3d5b13a778de67793c8c01432c3969701821e
parent52f2259ab2b7445288e2e7d7688cd8e50fa39fd1 (diff)
This commit was generated by cvs2svn to compensate for changes in r2, which
included commits to RCS files with non-trunk default branches. git-svn-id: http://svn.cactuscode.org/arrangements/CactusPUGHIO/IOHDF5Util/trunk@3 7842ec3a-9562-4be5-9c5b-06ba18f2b668
Diffstat
-rw-r--r--README7
-rw-r--r--interface.ccl5
-rw-r--r--param.ccl43
-rw-r--r--schedule.ccl15
-rw-r--r--src/DumpUtils.c430
-rw-r--r--src/DumpVar.c1034
-rw-r--r--src/ParseGeometry.c347
-rw-r--r--src/ParseVars.c212
-rw-r--r--src/RecoverVar.c809
-rw-r--r--src/Startup.c147
-rw-r--r--src/ioHDF5UtilGH.h217
-rw-r--r--src/make.code.defn5
-rw-r--r--src/make.configuration.defn20
13 files changed, 3291 insertions, 0 deletions
diff --git a/README b/README
new file mode 100644
index 0000000..ee94091
--- /dev/null
+++ b/README
@@ -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 (&parameterSize, 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, &gtype, &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
generated by cgit v1.2.3 (git 2.39.1) at 2024-06-08 23:02:41 +0000