JSim RealNData Interpolation

Implementation Specification

The RealNData.realVal(double[] xvals) method implements N-dimensional interpolation based on grid values (not grid array indexes). Various interpolation algorithms may be implemented, subject to the following guidelines.

The method may be called concurrently by several threads, (even on the same RealNData), and so must be reentrant.

Interpolation must efficiently support both regular and irregular grids. Currently all model data is defined on regular grids, but this may change in the future. Data curves in project datasets is often defined on irregular grids. GridData.nearestInx(double xval) should be used to find the surrounding grid indexes for an xval in an efficient manner. Note that IrregularGridData.nearestInx() currently (Sep 06) uses an inefficient low-to-high search strategy. This should be replaced by an efficient binary search at some point in the future.

Interpolation algorithms may use as many or as few nearby grid points as is appropriate. A linear algorithm may use as few as 2, or as many as 2^N nearest neighbors. A spline or other smooth algorithm may use even more by drawing on non-nearest neighbors. This specification is silent on the number of neighbors used for interpolation. However, in does require robustness when NaN's are encountered in one or more neighbors. All model data arrays (RealNData) in initialized to NaNs at the start of a run. Values are calculated in an unpredictable order, depending upon the planner/compiler. When interpolation is requested, data values for some neighbor grid points may be NaNs. Currently (due to compiler limitations) all JSim model loops move from low index to high index, but this should not be assumed by interpolation algorithms as other calculation orders will, no doubt, be implemented in the future.

Linear algorithms assign weight factors for each of the 2^N nearest neighbors based on distance to the query xvals[], with the returned value being a summation of data values times weights. Algorithms should not add zero weights to summation in case the associated data value is NaN (IEEE spec says 0*NaN gives NaN). Furthermore, weights that are very nearly zero should be set to zero. Here "very nearly zero" mean abs(x) < THRESH where THRESH is a static parameter currently recommended to have value 1e-7.

Note: The current THRESH is 1e-7 times the minimum grid delta, which is unit incompatible with the dimensionless weights. Therefore, the current threshold is non-sensical.

How to extend this NaN-protection policy to non-linear algorithms that sample non-nearest neighbors has yet to be defined.

It may be assumed that the GridData's overwhich the queried RealNData assumes values have non-NaN values.

xvals[] query outside grid bounds but within one grid delta should return values at grid bounds. When the query is more than one delta outside grid bounds, this the value returned may be either NaN or the value at grid bounds. Future versions of this spec will probably pick one or the other but, as of this writing, it is not completely clear which is best.

Verification specification

There are 3 types of verification required: accuracy, robustness and performance.

Accuracy verification

Needed: model/project files that demonstrate smooth and correct interpolation behaviour for NaN-less RealNData objects. These files will be incorporated into the JSIMSRC/VERIFY/interp suite.

Robustness

Needed: model files that check for undesirable NaNs generated at run-time. These files will be incorporated into the JSIMSRC/VERIFY/interp suite. Any additional NaN problem models discovered by others during testing are appreciated.

Performance

A JSim/data/InterpTest.java for assessing algorithm performance will be provided. Run via $JSIMSRC/build/jsinterpperf