diff options
Diffstat (limited to 'src/Operator.c')
| -rw-r--r-- | src/Operator.c | 357 |
1 files changed, 91 insertions, 266 deletions
diff --git a/src/Operator.c b/src/Operator.c index 5d17ab8..08f7027 100644 --- a/src/Operator.c +++ b/src/Operator.c @@ -2,9 +2,9 @@ @file Operator.c @date Tue Apr 15 18:22:45 1997 @author Paul Walker - @desc + @desc Definition of interpolation operators for regular uniform grids. - @enddesc + @enddesc @history @date Sun Jul 04 1999 @@ -18,7 +18,7 @@ @hdesc Move all local-interpolation code from LocalInterp to here @endhistory - @version $Header$ + @version $Id$ @@*/ #include <stdlib.h> @@ -27,297 +27,122 @@ #include "cctk.h" #include "cctk_Parameters.h" -#include "CactusPUGH/PUGH/src/include/pugh.h" +#include "util_ErrorCodes.h" +#include "util_Table.h" #include "Interpolate.h" /* the rcs ID and its dummy function to use it */ static const char *rcsid = "$Header$"; CCTK_FILEVERSION(CactusBase_LocalInterp_Operator_c) -/* uncomment this to get some debugging output */ -/* #define PUGHINTERP_DEBUG 1 */ - - -/* prototypes of routines defined in this source file */ -static int LocalInterp_CheckArguments (cGH *GH, - int num_dims, - int num_points, - int num_in_arrays, - int num_out_arrays, - const int interp_coord_array_types[]); - -/******************************************************************************/ - -/*@@ - @routine LocalInterp_InterpLocal - @date Sun Jul 04 1999 - @author Thomas Radke - @desc - The interpolation operator registered with the CCTK - under the name "uniform cartesian". - - Interpolates a list of non-distributed, processor-local - input arrays to a list of output arrays (one-to-one) - at a given number of interpolation points (indicated by - their coordinates). The points are located on a coordinate - system which is assumed to be a uniform cartesian. - @enddesc - - @var GH - @vdesc Pointer to CCTK grid hierarchy - @vtype cGH * - @vio in - @endvar - @var order - @vdesc interpolation order - @vtype int - @vio in - @endvar - @var num_points - @vdesc number of points to be interpolated on this processor - @vtype int - @vio in - @endvar - @var num_dims - @vdesc dimensionality of the underlying grid - @vtype int - @vio in - @endvar - @var num_in_arrays - @vdesc number of input arrays to interpolate from - @vtype int - @vio in - @endvar - @var num_out_arrays - @vdesc number of output arrays to interpolate to - @vtype int - @vio in - @endvar - @var coord_dims - @vdesc dimensions of the underlying grid - @vtype int [num_dims] - @vio in - @endvar - @var coord_arrays - @vdesc list of grid coordinate arrays - @vtype void *const [num_dims] - @vio in - @endvar - @var coord_array_types - @vdesc CCTK data type of coordinate arrays - @vtype int [num_dims] - @vio int - @endvar - @var interp_coord_arrays - @vdesc coordinates of points to interpolate at - @vtype void *const [num_dims] - @vio in - @endvar - @var interp_coord_array_types - @vdesc CCTK data type of coordinate arrays to interpolate at - @vtype int [num_dims] - @vio out - @endvar - @var in_arrays - @vdesc list of input arrays to interpolate on - @vtype void *const [num_in_arrays] - @vio in - @endvar - @var in_array_types - @vdesc CCTK data types of input arrays - @vtype int [num_in_arrays] - @vio in - @endvar - @var out_arrays - @vdesc list of output arrays to interpolate to - @vtype void *const [num_out_arrays] - @vio out - @endvar - @var out_array_types - @vdesc CCTK data types of output arrays - @vtype int [num_out_arrays] - @vio in - @endvar - - @returntype int - @returndesc - 0 - successful interpolation - -1 - in case of any errors - @endreturndesc - @@*/ -int LocalInterp_InterpLocal (cGH *GH, - int order, - int num_points, - int num_dims, - int num_in_arrays, - int num_out_arrays, - const int coord_dims[], - const void *const coord_arrays[], - const int coord_array_types[], - const void *const interp_coord_arrays[], - const int interp_coord_array_types[], - const void *const in_arrays[], - const int in_array_types[], - void *const out_arrays[], - const int out_array_types[]) +/******************************************************************** + ******************** External Routines ************************ + ********************************************************************/ +int LocalInterp_InterpLocalUniform (int num_dims, + int table, + /***** coordinate system *****/ + const CCTK_REAL coord_origin[], + const CCTK_REAL coord_delta[], + /***** interpolation points *****/ + int num_interp_points, + int interp_coords_type_code, + const void *const interp_coords[], + /***** input arrays *****/ + int num_input_arrays, + const CCTK_INT input_array_dims[], + const CCTK_INT input_array_type_codes[], + const void *const input_arrays[], + /***** output arrays *****/ + int num_output_arrays, + const CCTK_INT output_array_type_codes[], + void *const output_arrays[]) { - int dim, point, retval; - CCTK_REAL *coords, *origin, *delta; - const CCTK_REAL *const *data; + int iterator, retval; + char key[128]; + CCTK_INT order, type, nelems; - /* check arguments */ - retval = LocalInterp_CheckArguments (GH, num_dims, num_points, - num_in_arrays, num_out_arrays, - interp_coord_array_types); - if (retval <= 0) - { - return (retval); - } - for (dim = 0; dim < num_dims; dim++) - { - if (coord_array_types[dim] != CCTK_VARIABLE_REAL) - { - CCTK_WARN (1, "Coordinates should be of type CCTK_REAL"); - return (-1); - } - if (interp_coord_array_types[dim] != CCTK_VARIABLE_REAL) - { - CCTK_WARN (1, "Interpolation coordinates should be of type CCTK_REAL"); - return (-1); - } - } - /* get the grid spacings - this assumes a cartesian grid */ - origin = (CCTK_REAL *) malloc (2 * num_dims * sizeof (CCTK_REAL)); - delta = origin + num_dims; - data = (const CCTK_REAL *const *) coord_arrays; - for (dim = 0; dim < num_dims; dim++) + /* check for invalid arguments */ + if (num_dims < 0 || num_interp_points < 0 || + num_input_arrays < 0 || num_output_arrays < 0) { - origin[dim] = data[dim][0]; - delta[dim] = data[dim][1] - data[dim][0]; + return (UTIL_ERROR_BAD_INPUT); } - /* sort the individual interpolation coordinate arrays into a single one */ - coords = (CCTK_REAL *) malloc (num_dims * num_points * sizeof (CCTK_REAL)); - data = (const CCTK_REAL *const *) interp_coord_arrays; - for (point = 0; point < num_points; point++) + /* this interpolation operator computes one output array per input array */ + if (num_input_arrays != num_output_arrays) { - for (dim = 0; dim < num_dims; dim++) - { - *coords++ = data[dim][point]; - } + CCTK_WARN (1, "Number of input arrays must match number of output arrays"); + return (UTIL_ERROR_BAD_INPUT); } - coords -= num_dims * num_points; - /* call the interpolator function */ - retval = LocalInterp_Interpolate (order, - num_points, num_dims, num_out_arrays, - coord_dims, coords, origin, delta, - in_array_types, in_arrays, - out_array_types, out_arrays); - - /* free allocated resources */ - free (coords); - free (origin); - - return (retval); -} - -/**************************************************************************/ -/* local routines */ -/**************************************************************************/ - -/*@@ - @routine LocalInterp_CheckArguments - @date Thu 25 Jan 2001 - @author Thomas Radke - @desc - Checks the interpolation arguments passed in via - the flesh's general interpolation calling interface - - This routine also verifies that the parameters meet - the limitations of LocalInterp's interpolation operators. - @enddesc - - @var GH - @vdesc Pointer to CCTK grid hierarchy - @vtype cGH * - @vio in - @endvar - @var num_dims - @vdesc dimensionality of the underlying grid - @vtype int - @vio in - @endvar - @var num_points - @vdesc number of points to interpolate at - @vtype int - @vio in - @endvar - @var num_in_arrays - @vdesc number of passed input arrays - @vtype int - @vio in - @endvar - @var num_out_arrays - @vdesc number of passed input arrays - @vtype int - @vio in - @endvar - @var interp_coord_array_types - @vdesc types of passed coordinates to interpolate at - @vtype int [num_dims] - @vio in - @endvar - - @returntype int - @returndesc - +1 for success - 0 for success but nothing to do - -1 for failure (wrong parameters passed or limitations not met) - @endreturndesc - @@*/ -static int LocalInterp_CheckArguments (cGH *GH, - int num_dims, - int num_points, - int num_in_arrays, - int num_out_arrays, - const int interp_coord_array_types[]) -{ - int i; - - - /* check for invalid arguments */ - if (num_dims < 0 || num_points < 0 || num_in_arrays < 0 || num_out_arrays < 0) + if (interp_coords_type_code != CCTK_VARIABLE_REAL) { - return (-1); + CCTK_WARN (1, "Interpolation coordinates must be of type CCTK_REAL"); + return (UTIL_ERROR_BAD_INPUT); } /* check if there's anything to do at all */ - /* NOTE: num_points can be 0 in a collective call */ - if (num_dims == 0 || num_in_arrays == 0 || num_out_arrays == 0) + if (num_dims == 0 || num_input_arrays == 0) { return (0); } - /* for now we can only deal with coordinates of type CCTK_REAL */ - for (i = 0; i < num_dims; i++) + + /* get the interpolation order from the user-supplied parameter table */ + order = 1; + if (table >= 0) { - if (interp_coord_array_types[i] != CCTK_VARIABLE_REAL) + /* loop through all table options */ + iterator = Util_TableItCreate (table); + for (iterator = Util_TableItCreate (table); + Util_TableItQueryIsNonNull (iterator) > 0 && + Util_TableItQueryKeyValueInfo (iterator, sizeof (key), key, &type, + &nelems) > 0; + Util_TableItAdvance (iterator)) { - CCTK_WARN (1, "Interpolation coordinates must be of type CCTK_REAL"); - return (-1); + if (CCTK_Equals (key, "order")) + { + if (type == CCTK_VARIABLE_INT && nelems == 1) + { + Util_TableGetInt (table, &order, "order"); + } + else + { + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Invalid value for option 'order' in interpolation " + "parameter options table (must be CCTK_INT scalar value)", + key); + } + } + else if (CCTK_Equals (key, "N_boundary_points_to_omit") || + CCTK_Equals (key, "boundary_off_centering_tolerance") || + CCTK_Equals (key, "boundary_extrapolation_tolerance")) + { + /* warn about unsupported options */ + CCTK_VWarn (4, __LINE__, __FILE__, CCTK_THORNSTRING, + "Option with key '%s' in interpolation parameter options " + "table is not not supported (will be ignored)", key); + } + else + { + /* warn about other options */ + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "Found option with unrecognized key '%s' in interpolation " + "parameter options table (will be ignored)", key); + } } + Util_TableItDestroy (iterator); } - /* LocalInterp's interpolation operators compute one output array - per input array */ - if (num_in_arrays != num_out_arrays) - { - CCTK_WARN (1, "Number of input arrays must match number of output arrays"); - return (-1); - } + /* call the interpolator function */ + retval = LocalInterp_Interpolate (order, num_interp_points, num_dims, + num_output_arrays, input_array_dims, + (const CCTK_REAL *const *) interp_coords, + coord_origin, coord_delta, + input_array_type_codes, input_arrays, + output_array_type_codes, output_arrays); - return (1); + return (retval); } |