diff options
Diffstat (limited to 'src/Reduction.c')
-rw-r--r-- | src/Reduction.c | 184 |
1 files changed, 184 insertions, 0 deletions
diff --git a/src/Reduction.c b/src/Reduction.c index b48be12..16da833 100644 --- a/src/Reduction.c +++ b/src/Reduction.c @@ -38,6 +38,8 @@ int pugh_MinVal_Array (cGH *GH, int proc, int nDims, int dims [], void *inField, int inType, CCTK_REAL *outVal); int pugh_MaxVal_Array (cGH *GH, int proc, int nDims, int dims [], void *inField, int inType, CCTK_REAL *outVal); +int pugh_Sum_Array (cGH *GH, int proc, int nDims, int dims [], + void *inField, int inType, CCTK_REAL *outVal); int pugh_Norm1_Array (cGH *GH, int proc, int nDims, int dims [], void *inField, int inType, CCTK_REAL *outVal); int pugh_Norm2_Array (cGH *GH, int proc, int nDims, int dims [], @@ -291,6 +293,129 @@ void pugh_ReduceArrayMaxVal (cGH *GH, int proc, int nDims, int dims [], /*@@ + @routine pugh_ReduceArraySum + @author Thomas Radke + @date 19 Aug 1999 + @desc + Registered PUGH reduction routine for computing the "sum" + 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 + @history + @endhistory + @var GH + @vdesc Pointer to CCTK grid hierarchy + @vtype 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 int * + @vio in + @endvar + @var nArrays + @vdesc number of input arrays + @vtype int + @vio in + @endvar + @var inArrays + @vdesc field of input arrays + @vtype void ** + @vio in + @endvar + @var inType + @vdesc (common) variable type of input arrays + @vtype int + @vio in + @endvar + @var nOutVals + @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 +@@*/ +void pugh_ReduceArraySum (cGH *GH, int proc, int nDims, int dims [], + int nArrays, void *inArrays [], int inType, + int nOutVals, void *outVals, int outType) +{ + int i, extent; + pGH *pughGH; + CCTK_REAL intermediate_result; + + if (nOutVals != 1) { + CCTK_WARN (1, "pugh_ReduceArraySum accepts only one output value per " + "input array"); + return; + } + + switch (outType) { + case CCTK_VARIABLE_REAL: + extent = sizeof (CCTK_REAL); break; + case CCTK_VARIABLE_INT: + extent = sizeof (CCTK_INT); break; + case CCTK_VARIABLE_CHAR: + extent = sizeof (CCTK_CHAR); break; + case CCTK_VARIABLE_COMPLEX: + extent = sizeof (CCTK_COMPLEX); break; + default: + CCTK_WARN (1, "pugh_ReduceArraySum: Unsupported return type"); + return; + } + + pughGH = (pGH *) GH->extensions [PUGH_GHExtension]; + + for (i = 0; i < nArrays; i++) { + + /* do the reduction on input array [i] */ + if (pugh_Sum_Array (GH, proc, nDims, dims, inArrays [i], inType, + &intermediate_result) == 0) { + + /* type-cast the result to the requested datatype */ + if (proc < 0 || proc == pughGH->myproc) { + void *outVal; + + outVal = (void *) (i*extent + (char *) outVals); + + if (outType == CCTK_VARIABLE_CHAR) + *(CCTK_CHAR *) outVal = (CCTK_CHAR) intermediate_result; + else if (outType == CCTK_VARIABLE_INT) + *(CCTK_INT *) outVal = (CCTK_INT) intermediate_result; + else if (outType == CCTK_VARIABLE_REAL) + *(CCTK_REAL *) outVal = intermediate_result; + else + CCTK_WARN (1, "Unsupported output type in pugh_ReduceArraySum"); + } + } + } +} + + +/*@@ @routine pugh_ReduceArrayNorm1 @author Thomas Radke @date 19 Aug 1999 @@ -1134,6 +1259,65 @@ int pugh_MaxVal_Array (cGH *GH, int proc, int nDims, int dims [], /*@@ + @routine pugh_Sum_Array + @date Sep 14 1999 + @author Thomas Radke + @desc Returns the sum of a distributed array. + @enddesc +@@*/ +int pugh_Sum_Array (cGH *GH, int proc, int nDims, int dims [], + void *inArray, int inType, CCTK_REAL *outVal) +{ + int i, nPoints; + pGH *pughGH; +#ifdef MPI + CCTK_REAL localOutVal; +#endif + + pughGH = (pGH *) GH->extensions [PUGH_GHExtension]; + + *outVal = 0.0; + + nPoints = dims [0]; + for (i = 1; i < nDims; i++) + nPoints *= dims [i]; + + for (i = 0; i < nPoints; i++) + if (inType == CCTK_VARIABLE_CHAR) + *outVal += ((CCTK_CHAR *) inArray) [i]; + else if (inType == CCTK_VARIABLE_INT) + *outVal += ((CCTK_INT *) inArray) [i]; + else if (inType == CCTK_VARIABLE_REAL) + *outVal += ((CCTK_REAL *) inArray) [i]; + else if (inType == CCTK_VARIABLE_COMPLEX) { + CCTK_WARN (1, "Don't know how to reduce complex variables !!!"); + return (1); + } else { + CCTK_WARN (1, "Unknown variable type in pugh_Sum_Array"); + return (1); + } + +#ifdef MPI + + if (pughGH->nprocs > 1) { + + /* *outVal contains now the local sum */ + localOutVal = *outVal; + if (proc < 0) + CACTUS_MPI_ERROR (MPI_Allreduce (&localOutVal, outVal, 1, PUGH_MPI_REAL, + MPI_SUM, pughGH->PUGH_COMM_WORLD)); + else + CACTUS_MPI_ERROR (MPI_Reduce (&localOutVal, outVal, 1, PUGH_MPI_REAL, + MPI_SUM, proc, pughGH->PUGH_COMM_WORLD)); + } + +#endif + + return (0); +} + + +/*@@ @routine pugh_Norm1_Array @date Aug 19 1999 @author Thomas Radke |