Gromacs  2025.0-dev-20241011-013a99c
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
List of all members | Public Member Functions | Static Public Attributes
gmx::CubicSplineTable Class Reference

#include <gromacs/tables/cubicsplinetable.h>

Description

Cubic spline interpolation table.

This class interpolates a function specified either as an analytical expression or from user-provided table data.

At initialization, you provide the reference function of vectors as a list of tuples that contain a brief name, the function, and derivative for each function to tabulate. To create a table with two functions this initializer list can for instance look like

{ {"LJ6", lj6Func, lj6Der}, {"LJ12", lj12Func, lj12Der} }

The names are only used so exceptions during initialization can be traced to a specific table.

When interpolating, there are methods to interpolate either 1, 2, or 3 functions in one go. By default these interpolation routines will operate on tables with the same number of functions as specified in the interpolation method (debug builds check that this is consistent with the table). However, it is also possible to use optional template parameters that specify the total number of functions in a table, and what function index to interpolate. For instance, to interpolate the derivative of the second function (i.e., index 1) in a multi-function-table with three functions in total, you can write

table.evaluateDerivative<3,1>(x,&der);

Here too, debug builds will check that the template parameters are consistent with the table.

This class interpolates a function specified either as an analytical expression or from user-provided table data. The coefficients for each table point are precalculated such that we simply evaluate

\begin{eqnarray*} V(x) = Y + F \epsilon + G \epsilon^2 + H \epsilon^3 V'(x) = (F + 2 G \epsilon + 3 H \epsilon^2)/h \end{eqnarray*}

Where h is the spacing and epsilon the fractional offset from table point.

While it is possible to create tables only from function values (i.e., no derivatives), it is recommended to provide derivatives for higher accuracy and to avoid issues with numerical differentiation. Note that the table input should be smooth, i.e. it should not contain noise e.g. from an (iterative) Boltzmann inversion procedure - you have been warned.

Note
This class is responsible for fundamental interpolation of any function, which might or might not correspond to a potential. For this reason both input and output derivatives are proper function derivatives, and we do not swap any signs to get forces directly from the table.
There will be a small additional accuracy loss from the internal operation where we calculate the epsilon offset from the nearest table point, since the integer part we subtract can get large in those cases. While this is technically possible to solve with extended precision arithmetics, that would introduce extra instructions in some highly performance-sensitive code parts. For typical GROMACS interaction functions the derivatives will decay faster than the potential, which means it will never play any role. For other functions it will only cause a small increase in the relative error for arguments where the magnitude of the function or derivative is very small. Since we typically sum several results in GROMACS, this should never show up in any real cases, and for this reason we choose not to do the extended precision arithmetics.
These routines are not suitable for table ranges starting far away from zero, since we allocate memory and calculate indices starting from range zero for efficiency reasons.

Public Member Functions

 CubicSplineTable (std::initializer_list< AnalyticalSplineTableInput > analyticalInputList, const std::pair< real, real > &range, real tolerance=defaultTolerance)
 Initialize table data from function. More...
 
 CubicSplineTable (std::initializer_list< NumericalSplineTableInput > numericalInputList, const std::pair< real, real > &range, real tolerance=defaultTolerance)
 Initialize table data from tabulated values and derivatives. More...
 
template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void evaluateFunctionAndDerivative (T r, T *functionValue, T *derivativeValue) const
 Evaluate both function and derivative, single table function. More...
 
template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void evaluateFunction (T r, T *functionValue) const
 Evaluate function value only, single table function. More...
 
template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void evaluateDerivative (T r, T *derivativeValue) const
 Evaluate function derivative only, single table function. More...
 
template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void evaluateFunctionAndDerivative (T r, T *functionValue0, T *derivativeValue0, T *functionValue1, T *derivativeValue1) const
 Evaluate both function and derivative, two table functions. More...
 
template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void evaluateFunction (T r, T *functionValue0, T *functionValue1) const
 Evaluate function value only, two table functions. More...
 
template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void evaluateDerivative (T r, T *derivativeValue0, T *derivativeValue1) const
 Evaluate function derivative only, two table functions. More...
 
template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void evaluateFunctionAndDerivative (T r, T *functionValue0, T *derivativeValue0, T *functionValue1, T *derivativeValue1, T *functionValue2, T *derivativeValue2) const
 Evaluate both function and derivative, three table functions. More...
 
template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void evaluateFunction (T r, T *functionValue0, T *functionValue1, T *functionValue2) const
 Evaluate function value only, three table functions. More...
 
template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void evaluateDerivative (T r, T *derivativeValue0, T *derivativeValue1, T *derivativeValue2) const
 Evaluate function derivative only, three table functions. More...
 
real tableSpacing () const
 Return the table spacing (distance between points) More...
 

Static Public Attributes

static LIBGROMACS_EXPORT const real defaultTolerance = 10.0 * GMX_FLOAT_EPS
 Default tolerance for cubic spline tables. More...
 

Constructor & Destructor Documentation

gmx::CubicSplineTable::CubicSplineTable ( std::initializer_list< AnalyticalSplineTableInput analyticalInputList,
const std::pair< real, real > &  range,
real  tolerance = defaultTolerance 
)

Initialize table data from function.

Parameters
analyticalInputListInitializer list with one or more functions to tabulate, specified as elements with a string description and the function as well as derivative. The function will also be called for values smaller than the lower limit of the range, but we avoid calling it for 0.0 if that value is not included in the range. Constructor will throw gmx::APIError for negative values. Due to the way the numerical derivative evaluation depends on machine precision internally, this range must be at least 0.001, or the constructor throws gmx::APIError.
rangeRange over which the function will be tabulated. Constructor will throw gmx::APIError for negative values, or if the value/derivative vector does not cover the range.
toleranceRequested accuracy of the table. This will be used to calculate the required internal spacing. If this cannot be achieved (for instance because the table would require too much memory) the constructor will throw gmx::ToleranceError.
Note
The functions are always defined in double precision to avoid losing accuracy when constructing tables.
Since we fill the table for values below range.first, you can achieve a smaller table by using a smaller range where the tolerance has to be met, and accept that a few function calls below range.first do not quite reach the tolerance.
Warning
For efficiency reasons (since this code is used in some inner (kernels), we always allocate memory and calculate table indices for the complete interval [0,range.second], although the data will not be valid outside the definition range to avoid calling the function there. This means you should not use this class to tabulate functions for small ranges very far away from zero, since you would both waste a huge amount of memory and incur truncation errors when calculating the index.
Exceptions
gmx::ToleranceErrorif the requested tolerance cannot be achieved, and gmx::APIError for other incorrect input.
gmx::CubicSplineTable::CubicSplineTable ( std::initializer_list< NumericalSplineTableInput numericalInputList,
const std::pair< real, real > &  range,
real  tolerance = defaultTolerance 
)

Initialize table data from tabulated values and derivatives.

Parameters
numericalInputListInitializer list with one or more functions to tabulate, specified as a string description, vectors with function and derivative values, and the input spacing. Data points are separated by the spacing parameter, starting from 0. Values below the lower limit of the range will be used to attempt defining the table, but we avoid using index 0 unless 0.0 is included in the range. Some extra points beyond range.second are required to re-interpolate values, so add some margin. The constructor will throw gmx::APIError if the input vectors are too short to cover the requested range (and they must always be at least five points).
rangeRange over which the function will be tabulated. Constructor will throw gmx::APIError for negative values, or if the value/derivative vector does not cover the range.
toleranceRequested accuracy of the table. This will be used to calculate the required internal spacing and possibly re-interpolate. The constructor will throw gmx::ToleranceError if the input spacing is too coarse to achieve this accuracy.
Note
The input data vectors are always double precision to avoid losing accuracy when constructing tables.
Since we fill the table for values below range.first, you can achieve a smaller table by using a smaller range where the tolerance has to be met, and accept that a few function calls below range.first do not quite reach the tolerance.
Warning
For efficiency reasons (since this code is used in some inner (kernels), we always allocate memory and calculate table indices for the complete interval [0,range.second], although the data will not be valid outside the definition range to avoid calling the function there. This means you should not use this class to tabulate functions for small ranges very far away from zero, since you would both waste a huge amount of memory and incur truncation errors when calculating the index.

Member Function Documentation

template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void gmx::CubicSplineTable::evaluateDerivative ( r,
T *  derivativeValue 
) const
inline

Evaluate function derivative only, single table function.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 1
funcIndexIndex of function to evaluate in table, default is 0
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function derivative
[out]derivativeValueFunction derivative

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void gmx::CubicSplineTable::evaluateDerivative ( r,
T *  derivativeValue0,
T *  derivativeValue1 
) const
inline

Evaluate function derivative only, two table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 2
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function derivative
[out]derivativeValue0Interpolated derivative for first function
[out]derivativeValue1Interpolated derivative for second function

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void gmx::CubicSplineTable::evaluateDerivative ( r,
T *  derivativeValue0,
T *  derivativeValue1,
T *  derivativeValue2 
) const
inline

Evaluate function derivative only, three table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 3
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
funcIndex2Index of 3rd function to evaluate in table, default is 2
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function derivative
[out]derivativeValue0Interpolated derivative for first function
[out]derivativeValue1Interpolated derivative for second function
[out]derivativeValue2Interpolated derivative for third function

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void gmx::CubicSplineTable::evaluateFunction ( r,
T *  functionValue 
) const
inline

Evaluate function value only, single table function.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 1
funcIndexIndex of function to evaluate in table, default is 0
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function value
[out]functionValueFunction value

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void gmx::CubicSplineTable::evaluateFunction ( r,
T *  functionValue0,
T *  functionValue1 
) const
inline

Evaluate function value only, two table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 2
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function value
[out]functionValue0Interpolated value for first function
[out]functionValue1Interpolated value for second function

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void gmx::CubicSplineTable::evaluateFunction ( r,
T *  functionValue0,
T *  functionValue1,
T *  functionValue2 
) const
inline

Evaluate function value only, three table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 3
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
funcIndex2Index of 3rd function to evaluate in table, default is 2
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function value
[out]functionValue0Interpolated value for first function
[out]functionValue1Interpolated value for second function
[out]functionValue2Interpolated value for third function

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 1, int funcIndex = 0, typename T >
void gmx::CubicSplineTable::evaluateFunctionAndDerivative ( r,
T *  functionValue,
T *  derivativeValue 
) const
inline

Evaluate both function and derivative, single table function.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 1
funcIndexIndex of function to evaluate in table, default is 0
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function and derivative
[out]functionValueFunction value
[out]derivativeValueFunction derivative

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 2, int funcIndex0 = 0, int funcIndex1 = 1, typename T >
void gmx::CubicSplineTable::evaluateFunctionAndDerivative ( r,
T *  functionValue0,
T *  derivativeValue0,
T *  functionValue1,
T *  derivativeValue1 
) const
inline

Evaluate both function and derivative, two table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 2
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function and derivative
[out]functionValue0Interpolated value for first function
[out]derivativeValue0Interpolated derivative for first function
[out]functionValue1Interpolated value for second function
[out]derivativeValue1Interpolated derivative for second function

For debug builds we assert that the input values fall in the range specified when constructing the table.

template<int numFuncInTable = 3, int funcIndex0 = 0, int funcIndex1 = 1, int funcIndex2 = 2, typename T >
void gmx::CubicSplineTable::evaluateFunctionAndDerivative ( r,
T *  functionValue0,
T *  derivativeValue0,
T *  functionValue1,
T *  derivativeValue1,
T *  functionValue2,
T *  derivativeValue2 
) const
inline

Evaluate both function and derivative, three table functions.

This is a templated method where the template can be either real or SimdReal.

Template Parameters
numFuncInTableNumber of separate functions in table, default is 3
funcIndex0Index of 1st function to evaluate in table, default is 0
funcIndex1Index of 2nd function to evaluate in table, default is 1
funcIndex2Index of 3rd function to evaluate in table, default is 2
TType (SimdReal or real) of lookup and result
Parameters
rPoints for which to evaluate function and derivative
[out]functionValue0Interpolated value for first function
[out]derivativeValue0Interpolated derivative for first function
[out]functionValue1Interpolated value for second function
[out]derivativeValue1Interpolated derivative for second function
[out]functionValue2Interpolated value for third function
[out]derivativeValue2Interpolated derivative for third function

For debug builds we assert that the input values fall in the range specified when constructing the table.

real gmx::CubicSplineTable::tableSpacing ( ) const
inline

Return the table spacing (distance between points)

You should never have to use this for normal code, but due to the way tables are constructed internally we need this in the unit tests to check relative tolerances over each interval.

Returns
table spacing.

Member Data Documentation

const real gmx::CubicSplineTable::defaultTolerance = 10.0 * GMX_FLOAT_EPS
static

Default tolerance for cubic spline tables.

This is 10*GMX_FLOAT_EPS in single precision, and 1e-10 for double precision. It might not be worth setting this tolerance lower than 1e-10 in double precision, both because you will end up with very large tables, and because functions like r^-12 become so large for small values of r the table generation code will lead to some precision loss even in double precision.


The documentation for this class was generated from the following files: