From 2a7f167722503f0ef2a08268d7a8ed24555be1de Mon Sep 17 00:00:00 2001 From: tradke Date: Tue, 15 Apr 2003 15:30:11 +0000 Subject: Added new reduction operator "count" which returns the number of grid points of a grid variable which a reduction operator is working on. This closes PR CactusPUGH/1477. Added new reduction operators "L3Norm" and "L4Norm" to compute the L3 and L4 norm resp. Also added synonyms for some existing reduction operators: "average" = "mean", "norm1" = "L1Norm", "norm2" = "L2Norm", "norm_inf" = "LinfNorm". git-svn-id: http://svn.cactuscode.org/arrangements/CactusPUGH/PUGHReduce/trunk@40 d60812e6-3970-4df4-986e-c251b06effeb --- src/ReductionCount.c | 234 +++++++++++++++++++++++++++++++ src/ReductionNorm3.c | 375 ++++++++++++++++++++++++++++++++++++++++++++++++++ src/ReductionNorm4.c | 373 +++++++++++++++++++++++++++++++++++++++++++++++++ src/Startup.c | 28 +++- src/make.code.defn | 4 +- src/pugh_reductions.h | 18 ++- 6 files changed, 1019 insertions(+), 13 deletions(-) create mode 100644 src/ReductionCount.c create mode 100644 src/ReductionNorm3.c create mode 100644 src/ReductionNorm4.c diff --git a/src/ReductionCount.c b/src/ReductionCount.c new file mode 100644 index 0000000..ee4c627 --- /dev/null +++ b/src/ReductionCount.c @@ -0,0 +1,234 @@ + /*@@ + @file ReductionCount.c + @date Thu Apr 3 11:54:53 1997 + @author Thomas Radke + @desc + Defines the PUGH reduction operator "count" to get the number of + grid points of an arbitrary array. + @enddesc + @version $Id$ + @@*/ + +#include +#include + +#include "pugh_reductions.h" + +static const char *rcsid = "$Id$"; + +CCTK_FILEVERSION(CactusPUGH_PUGHReduce_ReductionCount_c) + +/* local function prototypes */ +static int ReductionCount (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]); + + +/*@@ + @routine PUGH_ReductionCountArrays + @author Thomas Radke + @date 19 Aug 1999 + @desc + Registered PUGH reduction routine for computing the counts + of a set of arrays. + The arrays are described by their dimensions and variable type. + For the number of output values only 1 is accepted. + Type casting of the result is provided by specifying the + requested output datatype. The intermediate reduction value + is always computed as a CCTK_REAL value internally. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var nDims + @vdesc number of dimensions of input arrays + @vtype int + @vio in + @endvar + @var dims + @vdesc dimensions of input arrays + @vtype const int * + @vio in + @endvar + @var nArrays + @vdesc number of input arrays + @vtype int + @vio in + @endvar + @var inarrays + @vdesc field of input arrays + @vtype const void *const + @vio in + @endvar + @var inType + @vdesc (common) variable type of input arrays + @vtype int + @vio in + @endvar + @var num_outvals + @vdesc number of values per output array + @vtype int + @vio in + @endvar + @var outvals + @vdesc pointer to buffer holding the output values + @vtype void * + @vio in + @endvar + @var outtype + @vdesc (common) variable type of output arrays + @vtype int + @vio in + @endvar +@@*/ +int PUGH_ReductionCountArrays (const cGH *GH, + int proc, + int num_dims, + const int dims[/* num_dims */], + int num_inarrays, + const void *const inarrays[/* num_inarrays */], + int intype, + int num_outvals, + void *outvals /* [num_outvals] */, + int outtype) +{ + return (PUGH_ReductionArrays (GH, proc, num_dims, dims, + intype, num_inarrays, inarrays, + outtype, num_outvals, outvals, + ReductionCount)); +} + + +/*@@ + @routine PUGH_ReductionCountGVs + @author Thomas Radke + @date 23 Apr 1999 + @desc + PUGH reduction routine to be registered for calculating the count. + It just goes through the list of variables and calls the appropriate + grouptype reduction routine. + If the number of requested output values equals the size of the + variable, global reduction is done element-wise. + Otherwise num_outvals must be 1, and global reduction is done on the + results of the local reductions. + @enddesc + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var num_outvals + @vdesc number of requested output elements + @vtype int + @vio in + @endvar + @var outtype + @vdesc type of return value + @vtype int + @vio in + @endvar + @var outvals + @vdesc array for the result values to be stored + @vtype void * + @vio in + @endvar + @var num_invars + @vdesc number of variables in the list + @vtype int + @vio in + @endvar + @var invars + @vdesc list of variables (given as their global indices) + @vtype const int * + @vio in + @endvar +@@*/ +int PUGH_ReductionCountGVs (const cGH *GH, + int proc, + int num_outvals, + int outtype, + void *outvals /* [num_outvals] */, + int num_invars, + const int invars[/* num_invars */]) +{ + return (PUGH_ReductionGVs (GH, proc, num_invars, invars, + outtype, num_outvals, outvals, + ReductionCount)); +} + + +/*****************************************************************************/ +/* local functions */ +/*****************************************************************************/ +/*@@ + @routine ReductionCount + @date Aug 19 1999 + @author Thomas Radke + @desc + Returns the number of grid points of a distributed array + @enddesc +@@*/ +static int ReductionCount (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]) +{ + int i; + + + /* avoid compiler warnings about unused parameters */ + (void) (GH + 0); + (void) (proc + 0); + (void) (num_dims + 0); + (void) (from + 0); + (void) (to + 0); + (void) (iterator + 0); + (void) (points_per_dim + 0); + (void) (intypes + 0); + (void) (inarrays + 0); + (void) (num_outvals + 0); + + /* assign the return value */ + if (proc < 0 || proc == CCTK_MyProc (GH)) + { + for (i = 0; i < num_inarrays; i++) + { + outvals[i] = num_points; + } + } + + return (0); +} diff --git a/src/ReductionNorm3.c b/src/ReductionNorm3.c new file mode 100644 index 0000000..2ccb67b --- /dev/null +++ b/src/ReductionNorm3.c @@ -0,0 +1,375 @@ +/*@@ + @file ReductionNorm3.c + @date Tue Apr 15 2003 + @author Thomas Radke + @desc + Defines the PUGH reduction operator to get the norm + of an arbitrary array which is defined as: + norm3 = $\sqrt[3]{\Sigma |a_i|^3 / np}$ + @enddesc + @version $Id$ + @@*/ + +#include +#include +#include + +#include "pugh_reductions.h" + +static const char *rcsid = "$Id$"; + +CCTK_FILEVERSION(CactusPUGH_PUGHReduce_ReductionNorm3_c) + +/* local function prototypes */ +static int ReductionNorm3 (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]); + + +/*@@ + @routine PUGH_ReductionNorm3Arrays + @author Thomas Radke + @date Tue Apr 15 2003 + @desc + Registered PUGH reduction routine for computing the "norm3" + of a set of arrays. + The arrays are described by their dimensions and variable type. + If the number of requested output values equals the number of + elements in an input array, global reduction is done element-wise. + Otherwise num_outvals must be 1, and global reduction is done on the + results of the local reductions. + Type casting of the result is provided by specifying the + requested output datatype. The intermediate reduction value + is always computed as a CCTK_REAL value internally. + @enddesc + @history + @endhistory + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var num_dims + @vdesc number of dimensions of input arrays + @vtype int + @vio in + @endvar + @var dims + @vdesc dimensions of input arrays + @vtype const int * + @vio in + @endvar + @var num_arrays + @vdesc number of input arrays + @vtype int + @vio in + @endvar + @var inarrays + @vdesc field of input arrays + @vtype const void *const + @vio in + @endvar + @var intype + @vdesc (common) variable type of input arrays + @vtype int + @vio in + @endvar + @var num_outvals + @vdesc number of values per output array + @vtype int + @vio in + @endvar + @var outvals + @vdesc pointer to buffer holding the output values + @vtype void * + @vio in + @endvar + @var outtype + @vdesc (common) variable type of output arrays + @vtype int + @vio in + @endvar +@@*/ +int PUGH_ReductionNorm3Arrays (const cGH *GH, + int proc, + int num_dims, + const int dims[/* num_dims */], + int num_inarrays, + const void *const inarrays[/* num_inarrays */], + int intype, + int num_outvals, + void *outvals /* [num_outvals] */, + int outtype) +{ + return (PUGH_ReductionArrays (GH, proc, num_dims, dims, + intype, num_inarrays, inarrays, + outtype, num_outvals, outvals, + ReductionNorm3)); +} + + +/*@@ + @routine PUGH_ReductionNorm3GVs + @author Thomas Radke + @date Tue Apr 15 2003 + @desc + PUGH reduction routine to be registered for calculating the norm3. + It just goes through the list of variables and calls the appropriate + grouptype reduction routine. + Global reduction is always done on the results of the local + reductions, so num_outvals must be 1. + @enddesc + @history + @endhistory + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var nOutVars + @vdesc number of requested output elements + @vtype int + @vio in + @endvar + @var outtype + @vdesc type of return value + @vtype int + @vio in + @endvar + @var outvals + @vdesc array for the result values to be stored + @vtype void * + @vio in + @endvar + @var num_invars + @vdesc number of variables in the list + @vtype int + @vio in + @endvar + @var invars + @vdesc list of variables (given as their global indices) + @vtype const int * + @vio in + @endvar +@@*/ +int PUGH_ReductionNorm3GVs (const cGH *GH, + int proc, + int num_outvals, + int outtype, + void *outvals /* [num_outvals] */, + int num_invars, + const int invars[/* num_invars */]) +{ + return (PUGH_ReductionGVs (GH, proc, num_invars, invars, + outtype, num_outvals, outvals, + ReductionNorm3)); +} + + +/*****************************************************************************/ +/* local functions */ +/*****************************************************************************/ +/*@@ + @routine ReductionNorm3 + @date Tue Apr 15 2003 + @author Thomas Radke + @desc Returns the "norm3" of a distributed array with + 'num_points' elements. Global reduction is done element-wise + (num_outvals == 1) or on the results of the local reductions. + + "norm3" is defined as $\sqrt{\Sigma |a_i|^3 / np}$. + @enddesc +@@*/ +static int ReductionNorm3 (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]) +{ + int i, total_outvals; + const pGH *pughGH; +#ifdef CCTK_MPI + CCTK_REAL *local_outvals; +#endif + + +/* macros to complete the ITERATE_ARRAY macro */ +#define INITIAL_REDUCTION_VALUE(array) 0 + + + for (i = total_outvals = 0; i < num_inarrays; i++) + { + switch (intypes[i]) + { + case CCTK_VARIABLE_CHAR: +#define CUBE_ABS(x) ((x) * (x) * (x)) +#define REDUCTION_OPERATION(norm3, scalar) norm3 += CUBE_ABS (scalar) + ITERATE_ARRAY (CCTK_BYTE, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_BYTE, outvals, num_outvals, total_outvals); + break; + + case CCTK_VARIABLE_INT: +#undef CUBE_ABS +#define CUBE_ABS(x) ((x) < 0 ? -((x) * (x) * (x)) : (x) * (x) * (x)) + ITERATE_ARRAY (CCTK_INT, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_INT2 + case CCTK_VARIABLE_INT2: + ITERATE_ARRAY (CCTK_INT2, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT2, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_INT4 + case CCTK_VARIABLE_INT4: + ITERATE_ARRAY (CCTK_INT4, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_INT8 + case CCTK_VARIABLE_INT8: + ITERATE_ARRAY (CCTK_INT8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT8, outvals, num_outvals, total_outvals); + break; +#endif + + case CCTK_VARIABLE_REAL: + ITERATE_ARRAY (CCTK_REAL, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_REAL4 + case CCTK_VARIABLE_REAL4: + ITERATE_ARRAY (CCTK_REAL4, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL8 + case CCTK_VARIABLE_REAL8: + ITERATE_ARRAY (CCTK_REAL8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL8, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL16 + case CCTK_VARIABLE_REAL16: + ITERATE_ARRAY (CCTK_REAL16, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL16, outvals, num_outvals, total_outvals); + break; +#endif + + case CCTK_VARIABLE_COMPLEX: +#undef REDUCTION_OPERATION +#define REDUCTION_OPERATION(norm3, scalar) norm3 += CUBE_ABS ((scalar).Re) +\ + CUBE_ABS ((scalar).Im) + ITERATE_ARRAY (CCTK_COMPLEX, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_REAL4 + ITERATE_ARRAY (CCTK_COMPLEX8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL8 + ITERATE_ARRAY (CCTK_COMPLEX16, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL8, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL16 + ITERATE_ARRAY (CCTK_COMPLEX32, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL16, outvals, num_outvals, total_outvals); + break; +#endif + + default: + CCTK_WARN (1, "PUGH_ReductionNorm3: Unknown variable type"); + return (-1); + } + } + + pughGH = PUGH_pGH (GH); + +#ifdef CCTK_MPI + /* do the global reduction from the processors' local reduction results */ + if (pughGH->nprocs > 1) + { + local_outvals = malloc (total_outvals * sizeof (CCTK_REAL)); + + /* outvals[] contains now the local sum */ + memcpy (local_outvals, outvals, total_outvals * sizeof (CCTK_REAL)); + if (proc < 0) + { + CACTUS_MPI_ERROR (MPI_Allreduce (local_outvals, outvals, total_outvals, + PUGH_MPI_REAL, MPI_SUM, pughGH->PUGH_COMM_WORLD)); + } + else + { + CACTUS_MPI_ERROR (MPI_Reduce (local_outvals, outvals, total_outvals, + PUGH_MPI_REAL, MPI_SUM, proc, pughGH->PUGH_COMM_WORLD)); + } + free (local_outvals); + } + +#endif + + /* finally assign the return value */ + if (proc < 0 || proc == pughGH->myproc) + { + for (i = 0; i < total_outvals; i++) + { + outvals[i] = pow (outvals[i] / num_points, 1.0 / 3.0); + } + } + + return (0); +} diff --git a/src/ReductionNorm4.c b/src/ReductionNorm4.c new file mode 100644 index 0000000..1476eb3 --- /dev/null +++ b/src/ReductionNorm4.c @@ -0,0 +1,373 @@ +/*@@ + @file ReductionNorm4.c + @date Tue Apr 15 2003 + @author Thomas Radke + @desc + Defines the PUGH reduction operator to get the norm + of an arbitrary array which is defined as: + norm4 = $\sqrt[4]{\Sigma (a_i^4) / np}$ + @enddesc + @version $Id$ + @@*/ + +#include +#include +#include + +#include "pugh_reductions.h" + +static const char *rcsid = "$Id$"; + +CCTK_FILEVERSION(CactusPUGH_PUGHReduce_ReductionNorm4_c) + +/* local function prototypes */ +static int ReductionNorm4 (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]); + + +/*@@ + @routine PUGH_ReductionNorm4Arrays + @author Thomas Radke + @date Tue Apr 15 2003 + @desc + Registered PUGH reduction routine for computing the "norm4" + of a set of arrays. + The arrays are described by their dimensions and variable type. + If the number of requested output values equals the number of + elements in an input array, global reduction is done element-wise. + Otherwise num_outvals must be 1, and global reduction is done on the + results of the local reductions. + Type casting of the result is provided by specifying the + requested output datatype. The intermediate reduction value + is always computed as a CCTK_REAL value internally. + @enddesc + @history + @endhistory + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var num_dims + @vdesc number of dimensions of input arrays + @vtype int + @vio in + @endvar + @var dims + @vdesc dimensions of input arrays + @vtype const int * + @vio in + @endvar + @var num_arrays + @vdesc number of input arrays + @vtype int + @vio in + @endvar + @var inarrays + @vdesc field of input arrays + @vtype const void *const + @vio in + @endvar + @var intype + @vdesc (common) variable type of input arrays + @vtype int + @vio in + @endvar + @var num_outvals + @vdesc number of values per output array + @vtype int + @vio in + @endvar + @var outvals + @vdesc pointer to buffer holding the output values + @vtype void * + @vio in + @endvar + @var outtype + @vdesc (common) variable type of output arrays + @vtype int + @vio in + @endvar +@@*/ +int PUGH_ReductionNorm4Arrays (const cGH *GH, + int proc, + int num_dims, + const int dims[/* num_dims */], + int num_inarrays, + const void *const inarrays[/* num_inarrays */], + int intype, + int num_outvals, + void *outvals /* [num_outvals] */, + int outtype) +{ + return (PUGH_ReductionArrays (GH, proc, num_dims, dims, + intype, num_inarrays, inarrays, + outtype, num_outvals, outvals, + ReductionNorm4)); +} + + +/*@@ + @routine PUGH_ReductionNorm4GVs + @author Thomas Radke + @date Tue Apr 15 2003 + @desc + PUGH reduction routine to be registered for calculating the norm4. + It just goes through the list of variables and calls the appropriate + grouptype reduction routine. + Global reduction is always done on the results of the local + reductions, so num_outvals must be 1. + @enddesc + @history + @endhistory + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype const cGH * + @vio in + @endvar + @var proc + @vdesc processor that should receive the result of operation + (negative value means all processors receive the result) + @vtype int + @vio in + @endvar + @var nOutVars + @vdesc number of requested output elements + @vtype int + @vio in + @endvar + @var outtype + @vdesc type of return value + @vtype int + @vio in + @endvar + @var outvals + @vdesc array for the result values to be stored + @vtype void * + @vio in + @endvar + @var num_invars + @vdesc number of variables in the list + @vtype int + @vio in + @endvar + @var invars + @vdesc list of variables (given as their global indices) + @vtype const int * + @vio in + @endvar +@@*/ +int PUGH_ReductionNorm4GVs (const cGH *GH, + int proc, + int num_outvals, + int outtype, + void *outvals /* [num_outvals] */, + int num_invars, + const int invars[/* num_invars */]) +{ + return (PUGH_ReductionGVs (GH, proc, num_invars, invars, + outtype, num_outvals, outvals, + ReductionNorm4)); +} + + +/*****************************************************************************/ +/* local functions */ +/*****************************************************************************/ +/*@@ + @routine ReductionNorm4 + @date Tue Apr 15 2003 + @author Thomas Radke + @desc Returns the "norm4" of a distributed array with + 'num_points' elements. Global reduction is done element-wise + (num_outvals == 1) or on the results of the local reductions. + + "norm4" is defined as $\sqrt{\Sigma (a_i^4) / np}$. + @enddesc +@@*/ +static int ReductionNorm4 (const cGH *GH, + int proc, + int num_dims, + const int from[/* dim */], + const int to[/* dim */], + int iterator[/* dim */], + const int points_per_dim[/* dim */], + int num_points, + int num_inarrays, + const int intypes[/* num_inarrays */], + const void *const inarrays[/* num_inarrays */], + int num_outvals, + CCTK_REAL outvals[/*num_inarrays*num_outvals*/]) +{ + int i, total_outvals; + const pGH *pughGH; +#ifdef CCTK_MPI + CCTK_REAL *local_outvals; +#endif + + +/* macros to complete the ITERATE_ARRAY macro */ +#define INITIAL_REDUCTION_VALUE(array) 0 +#define POWER4(x) ((x) * (x) * (x) * (x)) + + + for (i = total_outvals = 0; i < num_inarrays; i++) + { + switch (intypes[i]) + { + case CCTK_VARIABLE_CHAR: +#define REDUCTION_OPERATION(norm4, scalar) norm4 += POWER4 (scalar) + ITERATE_ARRAY (CCTK_BYTE, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_BYTE, outvals, num_outvals, total_outvals); + break; + + case CCTK_VARIABLE_INT: + ITERATE_ARRAY (CCTK_INT, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_INT2 + case CCTK_VARIABLE_INT2: + ITERATE_ARRAY (CCTK_INT2, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT2, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_INT4 + case CCTK_VARIABLE_INT4: + ITERATE_ARRAY (CCTK_INT4, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_INT8 + case CCTK_VARIABLE_INT8: + ITERATE_ARRAY (CCTK_INT8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_INT8, outvals, num_outvals, total_outvals); + break; +#endif + + case CCTK_VARIABLE_REAL: + ITERATE_ARRAY (CCTK_REAL, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_REAL4 + case CCTK_VARIABLE_REAL4: + ITERATE_ARRAY (CCTK_REAL4, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL8 + case CCTK_VARIABLE_REAL8: + ITERATE_ARRAY (CCTK_REAL8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL8, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL16 + case CCTK_VARIABLE_REAL16: + ITERATE_ARRAY (CCTK_REAL16, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL16, outvals, num_outvals, total_outvals); + break; +#endif + + case CCTK_VARIABLE_COMPLEX: +#undef REDUCTION_OPERATION +#define REDUCTION_OPERATION(norm4, scalar) norm4 += POWER4 ((scalar).Re) + \ + POWER4 ((scalar).Im) + ITERATE_ARRAY (CCTK_COMPLEX, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL, outvals, num_outvals, total_outvals); + break; + +#ifdef CCTK_REAL4 + ITERATE_ARRAY (CCTK_COMPLEX8, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL4, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL8 + ITERATE_ARRAY (CCTK_COMPLEX16, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL8, outvals, num_outvals, total_outvals); + break; +#endif + +#ifdef CCTK_REAL16 + ITERATE_ARRAY (CCTK_COMPLEX32, num_dims, inarrays[i], + from, to, iterator, points_per_dim, + CCTK_REAL16, outvals, num_outvals, total_outvals); + break; +#endif + + default: + CCTK_WARN (1, "PUGH_ReductionNorm4: Unknown variable type"); + return (-1); + } + } + + pughGH = PUGH_pGH (GH); + +#ifdef CCTK_MPI + /* do the global reduction from the processors' local reduction results */ + if (pughGH->nprocs > 1) + { + local_outvals = malloc (total_outvals * sizeof (CCTK_REAL)); + + /* outvals[] contains now the local sum */ + memcpy (local_outvals, outvals, total_outvals * sizeof (CCTK_REAL)); + if (proc < 0) + { + CACTUS_MPI_ERROR (MPI_Allreduce (local_outvals, outvals, total_outvals, + PUGH_MPI_REAL, MPI_SUM, pughGH->PUGH_COMM_WORLD)); + } + else + { + CACTUS_MPI_ERROR (MPI_Reduce (local_outvals, outvals, total_outvals, + PUGH_MPI_REAL, MPI_SUM, proc, pughGH->PUGH_COMM_WORLD)); + } + free (local_outvals); + } + +#endif + + /* finally assign the return value */ + if (proc < 0 || proc == pughGH->myproc) + { + for (i = 0; i < total_outvals; i++) + { + outvals[i] = pow (outvals[i] / num_points, 1.0/4.0); + } + } + + return (0); +} diff --git a/src/Startup.c b/src/Startup.c index 2422b58..e40f2b9 100644 --- a/src/Startup.c +++ b/src/Startup.c @@ -38,21 +38,37 @@ int PUGHReduce_Startup(void); int PUGHReduce_Startup (void) { /* Register the reduction operators provided by PUGH */ - CCTK_RegisterReductionOperator (PUGH_ReductionMinValGVs, "minimum"); - CCTK_RegisterReductionOperator (PUGH_ReductionMaxValGVs, "maximum"); - CCTK_RegisterReductionOperator (PUGH_ReductionSumGVs, "sum"); CCTK_RegisterReductionOperator (PUGH_ReductionAvgGVs, "average"); + CCTK_RegisterReductionOperator (PUGH_ReductionAvgGVs, "mean"); + CCTK_RegisterReductionOperator (PUGH_ReductionCountGVs, "count"); + CCTK_RegisterReductionOperator (PUGH_ReductionMaxValGVs, "maximum"); + CCTK_RegisterReductionOperator (PUGH_ReductionMinValGVs, "minimum"); CCTK_RegisterReductionOperator (PUGH_ReductionNorm1GVs, "norm1"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm1GVs, "L1Norm"); CCTK_RegisterReductionOperator (PUGH_ReductionNorm2GVs, "norm2"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm2GVs, "L2Norm"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm3GVs, "norm3"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm3GVs, "L3Norm"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm4GVs, "norm4"); + CCTK_RegisterReductionOperator (PUGH_ReductionNorm4GVs, "L4Norm"); CCTK_RegisterReductionOperator (PUGH_ReductionNormInfGVs, "norm_inf"); + CCTK_RegisterReductionOperator (PUGH_ReductionSumGVs, "sum"); - CCTK_RegisterReductionArrayOperator (PUGH_ReductionMinValArrays, "minimum"); - CCTK_RegisterReductionArrayOperator (PUGH_ReductionMaxValArrays, "maximum"); - CCTK_RegisterReductionArrayOperator (PUGH_ReductionSumArrays, "sum"); CCTK_RegisterReductionArrayOperator (PUGH_ReductionAvgArrays, "average"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionAvgArrays, "mean"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionCountArrays, "count"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionMaxValArrays, "maximum"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionMinValArrays, "minimum"); CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm1Arrays, "norm1"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm1Arrays, "L1Norm"); CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm2Arrays, "norm2"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm2Arrays, "L2Norm"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm3Arrays, "norm3"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm3Arrays, "L3Norm"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm4Arrays, "norm4"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionNorm4Arrays, "L4Norm"); CCTK_RegisterReductionArrayOperator (PUGH_ReductionNormInfArrays, "norm_inf"); + CCTK_RegisterReductionArrayOperator (PUGH_ReductionSumArrays, "sum"); return (0); } diff --git a/src/make.code.defn b/src/make.code.defn index b53a2d3..7e73b26 100644 --- a/src/make.code.defn +++ b/src/make.code.defn @@ -2,4 +2,6 @@ # $Header$ # Source files in this directory -SRCS = ReductionAvg.c ReductionMax.c ReductionMin.c ReductionNorm1.c ReductionNorm2.c ReductionNormInf.c ReductionSum.c Startup.c Reduction.c +SRCS = ReductionAvg.c ReductionCount.c ReductionMax.c ReductionMin.c \ + ReductionNorm1.c ReductionNorm2.c ReductionNorm3.c ReductionNorm4.c \ + ReductionNormInf.c ReductionSum.c Startup.c Reduction.c diff --git a/src/pugh_reductions.h b/src/pugh_reductions.h index 6bb56c1..b65cdfb 100644 --- a/src/pugh_reductions.h +++ b/src/pugh_reductions.h @@ -133,21 +133,27 @@ extern "C" { #endif -int PUGH_ReductionMinValGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); -int PUGH_ReductionMaxValGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); -int PUGH_ReductionSumGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionAvgGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionCountGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionMaxValGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionMinValGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNorm1GVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNorm2GVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionNorm3GVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionNorm4GVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNormInfGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionSumGVs (REDUCTION_OPERATOR_REGISTER_ARGLIST); -int PUGH_ReductionMinValArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); -int PUGH_ReductionMaxValArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); -int PUGH_ReductionSumArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionAvgArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionCountArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionMaxValArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionMinValArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNorm1Arrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNorm2Arrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionNorm3Arrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionNorm4Arrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); int PUGH_ReductionNormInfArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); +int PUGH_ReductionSumArrays (REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST); typedef int (*reduction_fn_t) (const cGH *GH, int proc, -- cgit v1.2.3