From f1d419cc71594cd0e4bb2ae426af55000712bb28 Mon Sep 17 00:00:00 2001 From: tradke Date: Sat, 16 Mar 2002 22:04:30 +0000 Subject: Allow individual reductions on variables by parsing the options string in the 'IOBasic::outScalar_vars' parameter. This closes PR CactusBase/895. git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/IOBasic/trunk@116 b589c3ab-70e8-4b4d-a09f-cba2dd200880 --- src/OutputScalar.c | 249 +++++++++++++++++++++-------------------------------- src/Startup.c | 227 ++++++++++++++++++++++++++++++++++++++++++------ src/WriteScalar.c | 209 +++++++++++++++++++++++--------------------- src/iobasicGH.h | 26 +++--- 4 files changed, 423 insertions(+), 288 deletions(-) (limited to 'src') diff --git a/src/OutputScalar.c b/src/OutputScalar.c index 1281cff..9169821 100644 --- a/src/OutputScalar.c +++ b/src/OutputScalar.c @@ -27,7 +27,6 @@ CCTK_FILEVERSION(CactusBase_IOBasic_OutputScalar_c) ******************** Internal Routines ************************ ********************************************************************/ static void CheckSteerableParameters (iobasicGH *myGH); -static void SetOutputFlag (int vindex, const char *optstring, void *arg); /*@@ @@ -38,8 +37,7 @@ static void SetOutputFlag (int vindex, const char *optstring, void *arg); Loops over all variables and outputs them if necessary @enddesc @calls IOBasic_TimeForScalarOutput - IOBasic_WriteScalarGS - IOBasic_WriteScalarGA + IOBasic_WriteScalar @var GH @vdesc Pointer to CCTK GH @@ -55,7 +53,7 @@ static void SetOutputFlag (int vindex, const char *optstring, void *arg); @@*/ int IOBasic_ScalarOutputGH (const cGH *GH) { - int i, vindex, retval; + int vindex, retval; const char *name; iobasicGH *myGH; @@ -91,15 +89,7 @@ int IOBasic_ScalarOutputGH (const cGH *GH) #endif /* Make the IO call */ - if (CCTK_GroupTypeFromVarI (vindex) == CCTK_SCALAR) - { - i = IOBasic_WriteScalarGS (GH, vindex, name); - } - else - { - i = IOBasic_WriteScalarGA (GH, vindex, name); - } - if (i == 0) + if (IOBasic_WriteScalar (GH, vindex, name) == 0) { /* Register GF as having 0D output this iteration */ myGH->outScalar_last [vindex] = GH->cctk_iteration; @@ -111,66 +101,6 @@ int IOBasic_ScalarOutputGH (const cGH *GH) } -/*@@ - @routine IOBasic_ScalarOutputVarAs - @date Sat March 6 1999 - @author Gabrielle Allen - @desc - Unconditional output of a variable using IOBasic's SCALAR - output method - @enddesc - @calls IOBasic_WriteScalarGS - IOBasic_WriteScalarGA - - @var GH - @vdesc Pointer to CCTK GH - @vtype const cGH * - @vio in - @endvar - @var fullname - @vdesc complete name of variable to output - @vtype const char * - @vio in - @endvar - @var alias - @vdesc alias name of variable to output (used to generate output filename) - @vtype const char * - @vio in - @endvar - - @returntype int - @returndesc - return code of @seeroutine IOBasic_WriteScalarGS or - @seeroutine IOBasic_WriteScalarGA - @endreturndesc -@@*/ -int IOBasic_ScalarOutputVarAs (const cGH *GH, const char *fullname, const char *alias) -{ - int vindex, retval; - - - vindex = CCTK_VarIndex (fullname); - -#ifdef IOBASIC_DEBUG - printf("\nIn IOBasic_ScalarOutputVarAs\n-------------------\n"); - printf(" Fullname = -%s-\n",fullname); - printf(" Alias = -%s-\n",alias); - printf(" Index = %d\n",vindex); -#endif - - if (CCTK_GroupTypeFromVarI (vindex) == CCTK_SCALAR) - { - retval = IOBasic_WriteScalarGS (GH, vindex, alias); - } - else - { - retval = IOBasic_WriteScalarGA (GH, vindex, alias); - } - - return (retval); -} - - /*@@ @routine IOBasic_TimeForScalarOutput @date Sat March 6 1999 @@ -198,45 +128,37 @@ int IOBasic_ScalarOutputVarAs (const cGH *GH, const char *fullname, const char * @@*/ int IOBasic_TimeForScalarOutput (const cGH *GH, int vindex) { - int return_type; - iobasicGH *myGH; + int retval; char *fullname; + iobasicGH *myGH; - /* Default is do not do output */ - return_type = 0; - /* Get the GH extension for IOBasic */ myGH = (iobasicGH *) CCTK_GHExtension (GH, "IOBasic"); CheckSteerableParameters (myGH); /* check if any output was requested */ - if (myGH->outScalar_every <= 0) + if (myGH->outScalar_every <= 0 || GH->cctk_iteration % myGH->outScalar_every + || myGH->scalar_reductions[vindex].num_reductions == 0) { - return (0); + retval = 0; } - - /* Check if this variable should be output */ - if (myGH->do_outScalar [vindex] && - (GH->cctk_iteration % myGH->outScalar_every) == 0) + else { /* Check if variable wasn't already output this iteration */ - if (myGH->outScalar_last [vindex] == GH->cctk_iteration) + retval = myGH->outScalar_last [vindex] != GH->cctk_iteration; + if (! retval) { fullname = CCTK_FullName (vindex); CCTK_VWarn (5, __LINE__, __FILE__, CCTK_THORNSTRING, - "Already done IOBasic scalar output for '%s' in " - "current iteration (probably via triggers)", fullname); + "Already done scalar output for '%s' in current iteration " + "(probably via triggers)", fullname); free (fullname); } - else - { - return_type = 1; - } } - return (return_type); + return (retval); } @@ -262,8 +184,7 @@ int IOBasic_TimeForScalarOutput (const cGH *GH, int vindex) @returntype int @returndesc - return code of @seeroutine IOBasic_WriteScalarGS or - @seeroutine IOBasic_WriteScalarGA + return code of @seeroutine IOBasic_WriteScalar @endreturndesc @@*/ int IOBasic_TriggerScalarOutput (const cGH *GH, int vindex) @@ -285,14 +206,7 @@ int IOBasic_TriggerScalarOutput (const cGH *GH, int vindex) #endif /* Do the Scalar output */ - if (CCTK_GroupTypeFromVarI (vindex) == CCTK_SCALAR) - { - retval = IOBasic_WriteScalarGS (GH, vindex, name); - } - else - { - retval = IOBasic_WriteScalarGA (GH, vindex, name); - } + retval = IOBasic_WriteScalar (GH, vindex, name); if (retval == 0) { @@ -304,62 +218,120 @@ int IOBasic_TriggerScalarOutput (const cGH *GH, int vindex) } -/**************************** local functions ******************************/ -/* check if steerable parameters have changed */ +/******************************************************************** + ******************** Internal Routines ************************ + ********************************************************************/ + /*@@ + @routine CheckSteerableParameters + @date Tue 31 Jul 2001 + @author Thomas Radke + @desc + Re-evaluates 'IOBasic::outScalar_every' and/or 'IO::out_every' + resp. to set myGH->outScalar_every to the frequency of scalar + output. Re-evaluates 'IOBasic::outScalar_vars' and + 'IOBasic::outScalar_reductions'. + @enddesc + @calls CCTK_ParameterQueryTimesSet + + @var myGH + @vdesc Pointer to IOBasic's GH extension + @vtype iobasicGH * + @vio in + @endvar +@@*/ static void CheckSteerableParameters (iobasicGH *myGH) { - int i, num_vars, out_old; - int times_set; - char *fullname, *msg, *oldmsg; - static int outScalar_vars_lastset = -1; + int vindex, out_old, times_set, update_reductions_list, num_vars; + iobasic_reduction_t *reduction, *next; + static int outScalar_vars_lastset = -1; + static int outScalar_reductions_lastset = -1; + iobasic_parseinfo_t info; + char *msg, *oldmsg, *fullname; DECLARE_CCTK_PARAMETERS - /* How often to output */ + /* how often to output */ out_old = myGH->outScalar_every; myGH->outScalar_every = out_every > 0 ? out_every : -1; if (outScalar_every > 0) { myGH->outScalar_every = outScalar_every; } - - if (myGH->outScalar_every != out_old) + if (myGH->outScalar_every != out_old && ! CCTK_Equals (newverbose, "none")) { - if (CCTK_Equals (newverbose, "standard") || - CCTK_Equals (newverbose, "full")) + if (myGH->outScalar_every > 0) + { + CCTK_VInfo (CCTK_THORNSTRING, "Scalar: Periodic output every %d " + "iterations", myGH->outScalar_every); + } + else { - CCTK_VInfo (CCTK_THORNSTRING, "Scalar: Output every %d iterations", - myGH->outScalar_every); + CCTK_INFO ("Scalar: Periodic output turned off"); } } - /* re-parse the 'outScalar_vars' parameter if it was changed */ + /* return if there's nothing to do */ + if (myGH->outScalar_every <= 0) + { + return; + } + + /* check if the 'outScalar_reductions' parameter if it was changed */ + times_set = CCTK_ParameterQueryTimesSet ("outScalar_reductions", + CCTK_THORNSTRING); + update_reductions_list = times_set != outScalar_reductions_lastset; + outScalar_reductions_lastset = times_set; + + /* check if the 'outScalar_vars' parameter if it was changed */ times_set = CCTK_ParameterQueryTimesSet ("outScalar_vars", CCTK_THORNSTRING); - if (times_set != outScalar_vars_lastset) + update_reductions_list |= times_set != outScalar_vars_lastset; + outScalar_vars_lastset = times_set; + + if (update_reductions_list) { + /* free old reduction lists ... */ num_vars = CCTK_NumVars (); - memset (myGH->do_outScalar, 0, num_vars); - CCTK_TraverseString (outScalar_vars, SetOutputFlag, myGH->do_outScalar, - CCTK_GROUP_OR_VAR); + for (vindex = 0; vindex < num_vars; vindex++) + { + if (myGH->scalar_reductions[vindex].num_reductions > 0) + { + myGH->scalar_reductions[vindex].num_reductions = 0; + reduction = myGH->scalar_reductions[vindex].reductions; + while (reduction) + { + next = reduction->next; + free (reduction->name); + free (reduction); + reduction = next; + } + } + } - if (myGH->outScalar_every && - (CCTK_Equals (newverbose, "standard") || - CCTK_Equals (newverbose, "full"))) + /* ... and create new ones */ + info.reduction_list = myGH->scalar_reductions; + info.reductions_string = outScalar_reductions; + if (CCTK_TraverseString (outScalar_vars, IOBasic_AssignReductionList, &info, + CCTK_GROUP_OR_VAR) < 0) + { + CCTK_WARN (1, "Failed to parse 'IOBasic::outScalar_vars' parameter"); + } + else if (! CCTK_Equals (newverbose, "none")) { msg = NULL; - for (i = 0; i < num_vars; i++) + for (vindex = 0; vindex < num_vars; vindex++) { - if (myGH->do_outScalar[i]) + if (myGH->scalar_reductions[vindex].num_reductions) { - fullname = CCTK_FullName (i); + fullname = CCTK_FullName (vindex); if (! msg) { - Util_asprintf (&msg, "Scalar: Output requested for %s", fullname); + Util_asprintf (&msg, "Scalar: Periodic output requested for '%s'", + fullname); } else { oldmsg = msg; - Util_asprintf (&msg, "%s %s", oldmsg, fullname); + Util_asprintf (&msg, "%s, '%s'", oldmsg, fullname); free (oldmsg); } free (fullname); @@ -371,26 +343,5 @@ static void CheckSteerableParameters (iobasicGH *myGH) free (msg); } } - - /* Save the last setting of 'outScalar_vars' parameter */ - outScalar_vars_lastset = times_set; - } - -} - - -/* callback for CCTK_TraverseString() to set the output flag - for the given variable */ -static void SetOutputFlag (int vindex, const char *optstring, void *arg) -{ - char *flags = (char *) arg; - - - flags[vindex] = 1; - - if (optstring) - { - CCTK_VWarn (5, __LINE__, __FILE__, CCTK_THORNSTRING, - "Optional string '%s' in variable name ignored", optstring); } } diff --git a/src/Startup.c b/src/Startup.c index 30c6bc4..3821d86 100644 --- a/src/Startup.c +++ b/src/Startup.c @@ -11,6 +11,7 @@ #include #include #include +#include #include "cctk.h" #include "cctk_Parameters.h" @@ -73,7 +74,6 @@ void IOBasic_Startup (void) @calls CCTK_RegisterIOMethod CCTK_RegisterIOMethodOutputGH - CCTK_RegisterIOMethodOutputVarAs CCTK_RegisterIOMethodTimeToOutput CCTK_RegisterIOMethodTriggerOutput @@ -119,7 +119,6 @@ static void *IOBasic_SetupGH (tFleshConfig *config, /* Register the IOBasic routines as output methods */ i = CCTK_RegisterIOMethod ("Scalar"); CCTK_RegisterIOMethodOutputGH (i, IOBasic_ScalarOutputGH); - CCTK_RegisterIOMethodOutputVarAs (i, IOBasic_ScalarOutputVarAs); CCTK_RegisterIOMethodTimeToOutput (i, IOBasic_TimeForScalarOutput); CCTK_RegisterIOMethodTriggerOutput (i, IOBasic_TriggerScalarOutput); @@ -128,8 +127,7 @@ static void *IOBasic_SetupGH (tFleshConfig *config, CCTK_RegisterIOMethodTimeToOutput (i, IOBasic_TimeForInfoOutput); CCTK_RegisterIOMethodTriggerOutput (i, IOBasic_TriggerInfoOutput); - if (CCTK_Equals (newverbose, "standard") || - CCTK_Equals( newverbose, "full")) + if (! CCTK_Equals (newverbose, "none")) { CCTK_INFO ("I/O Method 'Scalar' registered"); CCTK_INFO ("Scalar: Output of scalar quantities (grid scalars, " @@ -142,51 +140,61 @@ static void *IOBasic_SetupGH (tFleshConfig *config, i = CCTK_NumVars (); newGH->info_reductions = (iobasic_reductionlist_t *) - calloc (i, sizeof (iobasic_reductionlist_t)); - newGH->do_outScalar = (char *) malloc (i * sizeof (char)); - newGH->outInfo_last = (int *) malloc (i * sizeof (int)); - newGH->outScalar_last = (int *) malloc (i * sizeof (int)); + calloc (2 * i, sizeof (iobasic_reductionlist_t)); + newGH->scalar_reductions = newGH->info_reductions + i; + newGH->outInfo_last = (int *) malloc (2 * i * sizeof (int)); + newGH->outScalar_last = newGH->outInfo_last + i; - memset (newGH->outInfo_last, -1, i * sizeof (int)); - memset (newGH->outScalar_last, -1, i * sizeof (int)); + memset (newGH->outInfo_last, -1, 2 * i * sizeof (int)); newGH->filenameListScalar = NULL; - /* Check whether "IOBasic::outdirScalar" was set. + /* Check whether "IOBasic::outdir" was set. If so take this dir otherwise default to "IO::outdir" */ - if (CCTK_ParameterQueryTimesSet ("outdirScalar", CCTK_THORNSTRING) <= 0) + if (CCTK_ParameterQueryTimesSet ("outdir", CCTK_THORNSTRING) > 0 || + CCTK_ParameterQueryTimesSet ("outdirScalar", CCTK_THORNSTRING) > 0) { - outdirScalar = outdir; + if (CCTK_ParameterQueryTimesSet ("outdir", CCTK_THORNSTRING) <= 0) + { + CCTK_WARN (1, "Parameter 'IOBasic::outdirScalar is depricated in " + "BETA12, please use 'IOBasic::outdir' instead"); + outdir = outdirScalar; + } + } + else + { + outdir = *(const char **) + CCTK_ParameterGet ("outdir", CCTK_ImplementationThorn ("IO"), &i); } + /* skip the directory pathname if output goes into current directory */ - if (strcmp (outdirScalar, ".")) + if (strcmp (outdir, ".")) { - i = strlen (outdirScalar); - newGH->outdirScalar = (char *) malloc (i + 2); - strcpy (newGH->outdirScalar, outdirScalar); - newGH->outdirScalar[i] = '/'; - newGH->outdirScalar[i+1] = 0; + i = strlen (outdir); + newGH->outdir = (char *) malloc (i + 2); + strcpy (newGH->outdir, outdir); + newGH->outdir[i] = '/'; + newGH->outdir[i+1] = 0; } else { - newGH->outdirScalar = ""; + newGH->outdir = ""; } /* create the output dir */ - if (*newGH->outdirScalar && CCTK_MyProc (GH) == 0) + if (*newGH->outdir && CCTK_MyProc (GH) == 0) { - i = IOUtil_CreateDirectory (GH, newGH->outdirScalar, 0, 0); + i = IOUtil_CreateDirectory (GH, newGH->outdir, 0, 0); if (i < 0) { CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "IOBasic_SetupGH: Couldn't create Scalar output directory " - "'%s'", newGH->outdirScalar); + "'%s'", newGH->outdir); } else if (i >= 0 && CCTK_Equals (newverbose, "full")) { - CCTK_VInfo (CCTK_THORNSTRING, - "Scalar: Output to directory '%s'", - newGH->outdirScalar); + CCTK_VInfo (CCTK_THORNSTRING, "Scalar: Output to directory '%s'", + newGH->outdir); } } @@ -194,3 +202,170 @@ static void *IOBasic_SetupGH (tFleshConfig *config, return (newGH); } + + + /*@@ + @routine IOBasic_AssignReductionList + @date Tue 31 Jul 2001 + @author Thomas Radke + @desc + Callback routine called by CCTK_TraverseString() to set the + reduction list for a given variable. + For CCTK_GF and CCTK_ARRAY variables, it builds a chained list + of reduction operators according to the settings of 'optstring' + or the 'IOBasic::out{Info|Scalar}_reductions' parameter. + @enddesc + @calls CCTK_GroupTypeFromVarI + CCTK_ReductionHandle + + @var vindex + @vdesc index of the variable to set info output + @vtype int + @vio in + @endvar + @var optstring + @vdesc option string for this variable + @vtype const char * + @vio in + @endvar + @var arg + @vdesc user-supplied argument to callback routine + @vtype void * + @vio in + @endvar +@@*/ +void IOBasic_AssignReductionList (int vindex, const char *optstring, void *arg) +{ + const char *string_start, *string_end; + char *reduction_op, *reduction_op_list; + int reduction_handle; + iobasic_reductionlist_t *list; + iobasic_reduction_t **new_reduction; + const iobasic_parseinfo_t *info; + + + info = (const iobasic_parseinfo_t *) arg; + list = info->reduction_list + vindex; + + if (CCTK_GroupTypeFromVarI (vindex) == CCTK_SCALAR) + { + if (optstring) + { + CCTK_VWarn (5, __LINE__, __FILE__, CCTK_THORNSTRING, + "option list '%s' for variable '%s' ignored", + optstring, CCTK_VarName (vindex)); + } + + list->reductions = (iobasic_reduction_t *) + malloc (sizeof (iobasic_reduction_t)); + if (strncmp (CCTK_VarTypeName (CCTK_VarTypeI (vindex)), + "CCTK_VARIABLE_COMPLEX", 21)) + { + list->num_reductions = 1; + list->reductions->name = strdup ("scalar value"); + list->reductions->next = NULL; + } + else + { + list->num_reductions = 2; + list->reductions->name = strdup ("real part"); + list->reductions->next = (iobasic_reduction_t *) + malloc (sizeof (iobasic_reduction_t)); + list->reductions->next->name = strdup ("imag part"); + list->reductions->next->next = NULL; + } + + return; + } + + /* initialize to empty list */ + list->num_reductions = 0; + list->reductions = NULL; + + if (optstring) + { + if (strncmp (optstring, "reductions=<", 12) == 0 && + optstring[strlen (optstring) - 1] == '>') + { + reduction_op_list = strdup (optstring + 12); + reduction_op_list[strlen (reduction_op_list) - 1] = 0; + } + else + { + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOBasic_AssignReductionList: invalid syntax for option list " + "'%s'", optstring); + return; + } + } + else + { + reduction_op_list = strdup (info->reductions_string); + } + + /* now loop over all reduction operators */ + string_start = reduction_op_list; + reduction_op = (char *) malloc (strlen (string_start) + 1); + while (string_start && *string_start) + { + /* skip leading spaces */ + while (isspace ((int) *string_start)) + { + string_start++; + } + if (! *string_start) + { + break; + } + + /* advance to end of the operator string */ + string_end = string_start + 1; + while (*string_end && ! isspace ((int) *string_end)) + { + string_end++; + } + + /* copy the operator string */ + strncpy (reduction_op, string_start, string_end - string_start); + reduction_op[string_end - string_start] = 0; + string_start = string_end; + + /* get the reduction handle from the reduction operator */ + reduction_handle = CCTK_ReductionHandle (reduction_op); + if (reduction_handle < 0) + { + CCTK_VWarn (5, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOBasic_AssignReductionList: Invalid reduction operator " + "'%s'", reduction_op); + continue; + } + + /* add new reduction to end of list */ + new_reduction = &list->reductions; + while (*new_reduction) + { + if (strcmp ((*new_reduction)->name, reduction_op) == 0) + { + new_reduction = NULL; + break; + } + new_reduction = &((*new_reduction)->next); + } + if (new_reduction == NULL) + { + CCTK_VWarn (3, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOBasic_AssignReductionList: Duplicate reduction operator " + "'%s' will be ignored", reduction_op); + continue; + } + + *new_reduction = (iobasic_reduction_t *) malloc (sizeof (**new_reduction)); + (*new_reduction)->handle = reduction_handle; + (*new_reduction)->name = strdup (reduction_op); + (*new_reduction)->next = NULL; + list->num_reductions++; + } + + free (reduction_op_list); + free (reduction_op); +} diff --git a/src/WriteScalar.c b/src/WriteScalar.c index 413463c..12083de 100644 --- a/src/WriteScalar.c +++ b/src/WriteScalar.c @@ -24,6 +24,8 @@ static const char *rcsid = "$Header$"; CCTK_FILEVERSION(CactusBase_IOBasic_WriteScalar_c) +static int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias); +static int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias); static FILE *OpenScalarFile (const cGH *GH, int vindex, const char *filename, @@ -31,6 +33,75 @@ static FILE *OpenScalarFile (const cGH *GH, const char *description, const char *aliasname); + /*@@ + @routine IOBasic_WriteScalar + @date Mon Jun 19 2000 + @author Thomas Radke + @desc + Write a scalar value of a CCTK variable into an ASCII file + suitable for postprocessing by either gluplot or xgraph. + This routine just calls the appropriate output routine + according the the variable type. + @enddesc + + @calls CCTK_QueryGroupStorageI + CCTK_GroupIndexFromVarI + CCTK_GroupTypeFromVarI + IOBasic_WriteScalarGS + IOBasic_WriteScalarGA + + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var vindex + @vdesc CCTK index of the variable to output + @vtype int + @vio in + @endvar + @var alias + @vdesc alias name to use for building the output filename + @vtype const char * + @vio in + @endvar + + @returntype int + @returndesc + return code of either @seeroutine IOBasic_WriteScalarGS or + @seeroutine IOBasic_WriteScalarGA, or
+ -1 if variable has no storage assigned + @endreturndesc +@@*/ +int IOBasic_WriteScalar (const cGH *GH, int vindex, const char *alias) +{ + int retval; + char *fullname; + + + /* check if variable has storage assigned */ + if (! CCTK_QueryGroupStorageI (GH, CCTK_GroupIndexFromVarI (vindex))) + { + fullname = CCTK_FullName (vindex); + CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOBasic_WriteScalar: No scalar output for '%s' (no storage)", + fullname); + free (fullname); + retval = -1; + } + else if (CCTK_GroupTypeFromVarI (vindex) == CCTK_SCALAR) + { + retval = IOBasic_WriteScalarGS (GH, vindex, alias); + } + else + { + retval = IOBasic_WriteScalarGA (GH, vindex, alias); + } + + return (retval); +} + + /*@@ @routine IOBasic_WriteScalarGA @date Mon Jun 19 2000 @@ -41,8 +112,7 @@ static FILE *OpenScalarFile (const cGH *GH, suitable for postprocessing by either gluplot or xgraph. @enddesc - @calls CCTK_QueryGroupStorageI - CCTK_Reduce + @calls CCTK_Reduce CCTK_ReductionHandle IOUtil_RestartFromRecovery IOUtil_AdvertiseFile @@ -65,24 +135,17 @@ static FILE *OpenScalarFile (const cGH *GH, @returntype int @returndesc - 0 for success, or
- -1 if variable has no storage assigned + 0 for success @endreturndesc @@*/ -int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) +static int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) { - int ierr; - int reduction_handle; iobasicGH *myGH; FILE *file; char *filename; - char *reduction_op; - char *string_start; - char *string_end; char format_str[15]; const char *file_extension; - char *fullname; - CCTK_REAL reduction_value; + iobasic_reduction_t *reduction; union { char *non_const_ptr; @@ -102,17 +165,6 @@ int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) reductions.const_ptr = outScalar_reductions; GH_fake_const.const_ptr = GH; - /* first, check if variable has storage assigned */ - if (! CCTK_QueryGroupStorageI (GH, CCTK_GroupIndexFromVarI (vindex))) - { - fullname = CCTK_FullName (vindex); - CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOBasic_WriteScalarGA: No scalar output for '%s' (no storage)", - fullname); - free (fullname); - return (-1); - } - /* set output format */ sprintf (format_str, "%%%s\t%%%s\n", out_format, out_format); @@ -121,100 +173,71 @@ int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) /* get the GH extension handle for IOBasic */ myGH = (iobasicGH *) CCTK_GHExtension (GH, "IOBasic"); - - /* allocate string for the reduction operator */ - reduction_op = (char *) malloc (strlen (reductions.const_ptr) + 1); + reduction = myGH->scalar_reductions[vindex].reductions; /* now loop over all reduction operators */ - string_start = reductions.non_const_ptr; - while (string_start && *string_start) + while (reduction) { - /* skip leading spaces */ - while (isspace ((int) *string_start)) - { - string_start++; - } - if (! *string_start) - { - break; - } - - /* advance to end of the operator string */ - string_end = string_start + 1; - while (*string_end && ! isspace ((int) *string_end)) - { - string_end++; - } - - /* copy the operator string */ - strncpy (reduction_op, string_start, string_end - string_start); - reduction_op[string_end - string_start] = 0; - string_start = string_end; - - /* get the reduction handle from the reduction operator */ - reduction_handle = CCTK_ReductionHandle (reduction_op); - if (reduction_handle < 0) + /* do the reduction (all processors) */ + reduction->is_valid = CCTK_Reduce (GH_fake_const.non_const_ptr, 0, + reduction->handle, 1, + CCTK_VARIABLE_REAL, + &reduction->value, 1, vindex) == 0; + if (! reduction->is_valid) { CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOBasic_WriteScalarGA: Invalid reduction operator '%s'", - reduction_op); - continue; + "IOBasic_WriteInfo: Internal error in reduction '%s'", + reduction->name); } - - /* do the reduction (all processors) */ - ierr = CCTK_Reduce (GH_fake_const.non_const_ptr, 0, reduction_handle, 1, - CCTK_VARIABLE_REAL, &reduction_value, 1, vindex); - - /* dump the reduction value to file by processor 0 */ - if (ierr == 0 && CCTK_MyProc (GH) == 0) + else if (CCTK_MyProc (GH) == 0) { - + /* dump the reduction value to file by processor 0 */ if (new_filename_scheme) { /* allocate filename string buffer and build the filename */ - filename = (char *) malloc (strlen (myGH->outdirScalar) + - strlen (alias) + strlen (reduction_op) + + filename = (char *) malloc (strlen (myGH->outdir) + + strlen (alias) + strlen (reduction->name) + strlen (file_extension) + 2); - sprintf (filename, "%s%s_%s%s", myGH->outdirScalar, alias, - reduction_op, file_extension); + sprintf (filename, "%s%s_%s%s", myGH->outdir, alias, + reduction->name, file_extension); } else { /* FIXME: this check can go if we switch to the new filename scheme */ - if (strcmp (reduction_op, "minimum") == 0) + if (strcmp (reduction->name, "minimum") == 0) { file_extension = "min"; } - else if (strcmp (reduction_op, "maximum") == 0) + else if (strcmp (reduction->name, "maximum") == 0) { file_extension = "max"; } - else if (strcmp (reduction_op, "norm1") == 0) + else if (strcmp (reduction->name, "norm1") == 0) { file_extension = "nm1"; } - else if (strcmp (reduction_op, "norm2") == 0) + else if (strcmp (reduction->name, "norm2") == 0) { file_extension = "nm2"; } else { - file_extension = reduction_op; + file_extension = reduction->name; } /* allocate filename string buffer and build the filename */ - filename = (char *) malloc (strlen (myGH->outdirScalar) + + filename = (char *) malloc (strlen (myGH->outdir) + strlen (alias) + strlen (file_extension)+5); - sprintf (filename, "%s%s_%s.tl", myGH->outdirScalar, alias, + sprintf (filename, "%s%s_%s.tl", myGH->outdir, alias, file_extension); } - file = OpenScalarFile (GH, vindex, filename, reduction_op, + file = OpenScalarFile (GH, vindex, filename, reduction->name, "Reduction on Grid Arrays", alias); if (file) { /* write the data and close the file */ - fprintf (file, format_str, GH->cctk_time, reduction_value); + fprintf (file, format_str, GH->cctk_time, reduction->value); fclose (file); } else @@ -225,10 +248,9 @@ int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) } free (filename); } - } - /* free allocated resources */ - free (reduction_op); + reduction = reduction->next; + } return (0); } @@ -243,8 +265,7 @@ int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) suitable for postprocessing by either gluplot or xgraph. @enddesc - @calls CCTK_QueryGroupStorageI - CCTK_Reduce + @calls CCTK_Reduce CCTK_ReductionHandle IOUtil_RestartFromRecovery IOUtil_AdvertiseFile @@ -267,19 +288,18 @@ int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias) @returntype int @returndesc - 0 for success, or
- -1 if variable has no storage assigned + 0 for success @endreturndesc @@*/ -int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias) +static int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias) { - DECLARE_CCTK_PARAMETERS FILE *file; void *data; iobasicGH *myGH; - char *fullname, *filename; + char *filename; const char *file_extension; char format_str_real[15], format_str_int[15]; + DECLARE_CCTK_PARAMETERS /* output is done by processor 0 only */ @@ -288,17 +308,6 @@ int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias) return (0); } - /* first, check if variable has storage assigned */ - if (! CCTK_QueryGroupStorageI (GH, CCTK_GroupIndexFromVarI (vindex))) - { - fullname = CCTK_FullName (vindex); - CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOBasic_WriteScalarGS: No scalar output for '%s' (no storage)", - fullname); - free (fullname); - return (-1); - } - /* set the output format string for the desired notation */ sprintf (format_str_real, "%%%s\t%%%s\n", out_format, out_format); sprintf (format_str_int, "%%%s\t%%d\n", out_format); @@ -317,9 +326,9 @@ int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias) } /* build the output filename */ - filename = (char *) malloc (strlen (myGH->outdirScalar) + strlen (alias) + + filename = (char *) malloc (strlen (myGH->outdir) + strlen (alias) + strlen (file_extension) + 1); - sprintf (filename, "%s%s%s", myGH->outdirScalar, alias, file_extension); + sprintf (filename, "%s%s%s", myGH->outdir, alias, file_extension); /* create/reopen the file */ file = OpenScalarFile (GH, vindex, filename, "tl", "Scalar value", alias); diff --git a/src/iobasicGH.h b/src/iobasicGH.h index c75b70e..e957a00 100644 --- a/src/iobasicGH.h +++ b/src/iobasicGH.h @@ -26,25 +26,26 @@ typedef struct IOBASIC_REDUCTIONLIST iobasic_reduction_t *reductions; } iobasic_reductionlist_t; +typedef struct IOBASIC_PARSEINFO +{ + const char *reductions_string; + iobasic_reductionlist_t *reduction_list; +} iobasic_parseinfo_t; + typedef struct IOBASIC_GH { /* how often to output */ - int outScalar_every; - int outInfo_every; + int outScalar_every, outInfo_every; char info_reductions_changed; - /* flags indicating output for variable i */ - char *do_outScalar; - /* directory in which to place scalar output */ - char *outdirScalar; + char *outdir; /* The last iteration output */ - int *outInfo_last; - int *outScalar_last; + int *outInfo_last, *outScalar_last; - /* The reduction lists for info output for all variables */ - iobasic_reductionlist_t *info_reductions; + /* The reduction lists for info and scalar output for all variables */ + iobasic_reductionlist_t *info_reductions, *scalar_reductions; /* database for names of output files that were already created */ pNamedData *filenameListScalar; @@ -59,9 +60,8 @@ int IOBasic_TimeForInfoOutput (const cGH *GH, int vindex); int IOBasic_ScalarOutputGH (const cGH *GH); int IOBasic_TriggerScalarOutput (const cGH *GH, int vindex); int IOBasic_TimeForScalarOutput (const cGH *GH, int vindex); -int IOBasic_ScalarOutputVarAs (const cGH *GH, const char *vname, const char *alias); /* other function prototypes */ +void IOBasic_AssignReductionList (int vindex, const char *optstring, void *arg); int IOBasic_WriteInfo (const cGH *GH, int vindex); -int IOBasic_WriteScalarGS (const cGH *GH, int vindex, const char *alias); -int IOBasic_WriteScalarGA (const cGH *GH, int vindex, const char *alias); +int IOBasic_WriteScalar (const cGH *GH, int vindex, const char *alias); -- cgit v1.2.3