Implements a linearly interpolated scalar function that depends on one variable. More...
#include <Tabulated1DFunction.hpp>
Classes | |
struct | ComparatorX_ |
Helper class needed to sort the input sampling points. | |
Public Member Functions | |
Tabulated1DFunction () | |
Default constructor for a piecewise linear function. | |
template<class ScalarArrayX , class ScalarArrayY > | |
Tabulated1DFunction (size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, bool sortInputs=true) | |
Convenience constructor for a piecewise linear function. | |
template<class ScalarContainer > | |
Tabulated1DFunction (const ScalarContainer &x, const ScalarContainer &y, bool sortInputs=true) | |
Convenience constructor for a piecewise linear function. | |
template<class PointContainer > | |
Tabulated1DFunction (const PointContainer &points, bool sortInputs=true) | |
Convenience constructor for a piecewise linear function. | |
template<class ScalarArrayX , class ScalarArrayY > | |
void | setXYArrays (size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, bool sortInputs=true) |
Set the sampling points for the piecewise linear function. | |
template<class ScalarContainerX , class ScalarContainerY > | |
void | setXYContainers (const ScalarContainerX &x, const ScalarContainerY &y, bool sortInputs=true) |
Set the sampling points for the piecewise linear function. | |
template<class PointArray > | |
void | setArrayOfPoints (size_t nSamples, const PointArray &points, bool sortInputs=true) |
Set the sampling points for the piecewise linear function. | |
template<class XYContainer > | |
void | setContainerOfTuples (const XYContainer &points, bool sortInputs=true) |
Set the sampling points of the piecewise linear function using a STL-compatible container of tuple-like objects. | |
size_t | numSamples () const |
Returns the number of sampling points. | |
Scalar | xMin () const |
Return the x value of the leftmost sampling point. | |
Scalar | xMax () const |
Return the x value of the rightmost sampling point. | |
Scalar | xAt (size_t i) const |
Return the x value of the a sample point with a given index. | |
Scalar | valueAt (size_t i) const |
Return the value of the a sample point with a given index. | |
template<class Evaluation > | |
bool | applies (const Evaluation &x) const |
Return true iff the given x is in range [x1, xn]. | |
template<class Evaluation > | |
Evaluation | eval (const Evaluation &x, bool extrapolate=false) const |
Evaluate the spline at a given position. | |
template<class Evaluation > | |
Evaluation | evalDerivative (const Evaluation &x, bool extrapolate=false) const |
Evaluate the spline's derivative at a given position. | |
template<class Evaluation > | |
Evaluation | evalSecondDerivative (const Evaluation &x, bool extrapolate=false) const |
Evaluate the function's second derivative at a given position. | |
template<class Evaluation > | |
Evaluation | evalThirdDerivative (const Evaluation &x, bool extrapolate=false) const |
Evaluate the function's third derivative at a given position. | |
int | monotonic (Scalar x0, Scalar x1, bool extrapolate OPM_OPTIM_UNUSED=false) const |
Returns 1 if the function is monotonically increasing, -1 if the function is mononously decreasing and 0 if the function is not monotonous in the interval (x0, x1). | |
int | monotonic () const |
Same as monotonic(x0, x1), but with the entire range of the function as interval. | |
void | printCSV (Scalar xi0, Scalar xi1, unsigned k, std::ostream &os=std::cout) const |
Prints k tuples of the format (x, y, dx/dy, isMonotonic) to stdout. |
Implements a linearly interpolated scalar function that depends on one variable.
Opm::Tabulated1DFunction< Scalar >::Tabulated1DFunction | ( | ) | [inline] |
Default constructor for a piecewise linear function.
To specfiy the acutal curve, use one of the set() methods.
Opm::Tabulated1DFunction< Scalar >::Tabulated1DFunction | ( | size_t | nSamples, | |
const ScalarArrayX & | x, | |||
const ScalarArrayY & | y, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a piecewise linear function.
nSamples | The number of sampling points (must be >= 2) | |
x | An array containing the ![]() | |
y | An array containing the ![]() |
Opm::Tabulated1DFunction< Scalar >::Tabulated1DFunction | ( | const ScalarContainer & | x, | |
const ScalarContainer & | y, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a piecewise linear function.
x | An array containing the ![]() | |
y | An array containing the ![]() |
Opm::Tabulated1DFunction< Scalar >::Tabulated1DFunction | ( | const PointContainer & | points, | |
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a piecewise linear function.
points | A container of ![]() |
Evaluation Opm::Tabulated1DFunction< Scalar >::eval | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the spline at a given position.
x | The value on the abscissa where the function ought to be evaluated | |
extrapolate | If this parameter is set to true, the function will be extended beyond its range by straight lines, if false calling extrapolate for ![]() |
Evaluation Opm::Tabulated1DFunction< Scalar >::evalDerivative | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the spline's derivative at a given position.
x | The value on the abscissa where the function's derivative ought to be evaluated | |
extrapolate | If this parameter is set to true, the spline will be extended beyond its range by straight lines, if false calling extrapolate for ![]() |
Evaluation Opm::Tabulated1DFunction< Scalar >::evalSecondDerivative | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the function's second derivative at a given position.
Since this class represents a piecewise linear function, this method will always return 0.
x | The value on the abscissa where the function's derivative ought to be evaluated | |
extrapolate | If this parameter is set to true, the function will be extended beyond its range by straight lines, if false calling extrapolate for ![]() |
Evaluation Opm::Tabulated1DFunction< Scalar >::evalThirdDerivative | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the function's third derivative at a given position.
Since this class represents a piecewise linear function, this method will always return 0.
x | The value on the abscissa where the function's derivative ought to be evaluated | |
extrapolate | If this parameter is set to true, the function will be extended beyond its range by straight lines, if false calling extrapolate for ![]() |
int Opm::Tabulated1DFunction< Scalar >::monotonic | ( | Scalar | x0, | |
Scalar | x1, | |||
bool extrapolate | OPM_OPTIM_UNUSED = false | |||
) | const [inline] |
Returns 1 if the function is monotonically increasing, -1 if the function is mononously decreasing and 0 if the function is not monotonous in the interval (x0, x1).
In the corner case that the function is constant within the given interval, this method returns 3.
void Opm::Tabulated1DFunction< Scalar >::printCSV | ( | Scalar | xi0, | |
Scalar | xi1, | |||
unsigned | k, | |||
std::ostream & | os = std::cout | |||
) | const [inline] |
Prints k tuples of the format (x, y, dx/dy, isMonotonic) to stdout.
If the function does not apply for parts of [x0, x1] it is extrapolated using a straight line. The result can be inspected using the following commands:
----------- snip ----------- ./yourProgramm > function.csv gnuplot
gnuplot> plot "function.csv" using 1:2 w l ti "Curve", \ "function.csv" using 1:3 w l ti "Derivative" ----------- snap -----------
void Opm::Tabulated1DFunction< Scalar >::setContainerOfTuples | ( | const XYContainer & | points, | |
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points of the piecewise linear function using a STL-compatible container of tuple-like objects.
This method uses a single STL-compatible container of sampling points, which are assumed to be tuple-like objects storing the X and Y coordinates. "tuple-like" means that the objects provide access to the x values via std::get<0>(obj) and to the y value via std::get<1>(obj) (e.g. std::tuple or std::pair). "STL-compatible" means that the container provides access to iterators using the begin(), end() methods and also provides a size() method. Also, the number of entries in the X and the Y containers must be equal and larger than 1.
void Opm::Tabulated1DFunction< Scalar >::setXYArrays | ( | size_t | nSamples, | |
const ScalarArrayX & | x, | |||
const ScalarArrayY & | y, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points for the piecewise linear function.
This method takes C-style arrays (pointers) as arguments.
void Opm::Tabulated1DFunction< Scalar >::setXYContainers | ( | const ScalarContainerX & | x, | |
const ScalarContainerY & | y, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points for the piecewise linear function.
This method takes STL compatible containers as arguments.