Class implementing cubic splines. More...
#include <Spline.hpp>
Classes | |
struct | ComparatorX_ |
Helper class needed to sort the input sampling points. More... | |
Public Types | |
enum | SplineType { Full, Natural, Periodic, Monotonic } |
The type of the spline to be created. More... | |
Public Member Functions | |
Spline () | |
Default constructor for a spline. | |
Spline (Scalar x0, Scalar x1, Scalar y0, Scalar y1, Scalar m0, Scalar m1) | |
Convenience constructor for a full spline with just two sampling points. | |
template<class ScalarArrayX , class ScalarArrayY > | |
Spline (size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, SplineType splineType=Natural, bool sortInputs=true) | |
Convenience constructor for a natural or a periodic spline. | |
template<class PointArray > | |
Spline (size_t nSamples, const PointArray &points, SplineType splineType=Natural, bool sortInputs=true) | |
Convenience constructor for a natural or a periodic spline. | |
template<class ScalarContainer > | |
Spline (const ScalarContainer &x, const ScalarContainer &y, SplineType splineType=Natural, bool sortInputs=true) | |
Convenience constructor for a natural or a periodic spline. | |
template<class PointContainer > | |
Spline (const PointContainer &points, SplineType splineType=Natural, bool sortInputs=true) | |
Convenience constructor for a natural or a periodic spline. | |
template<class ScalarArray > | |
Spline (size_t nSamples, const ScalarArray &x, const ScalarArray &y, Scalar m0, Scalar m1, bool sortInputs=true) | |
Convenience constructor for a full spline. | |
template<class PointArray > | |
Spline (size_t nSamples, const PointArray &points, Scalar m0, Scalar m1, bool sortInputs=true) | |
Convenience constructor for a full spline. | |
template<class ScalarContainerX , class ScalarContainerY > | |
Spline (const ScalarContainerX &x, const ScalarContainerY &y, Scalar m0, Scalar m1, bool sortInputs=true) | |
Convenience constructor for a full spline. | |
template<class PointContainer > | |
Spline (const PointContainer &points, Scalar m0, Scalar m1, bool sortInputs=true) | |
Convenience constructor for a full spline. | |
void | set (Scalar x0, Scalar x1, Scalar y0, Scalar y1, Scalar m0, Scalar m1) |
Set the sampling points and the boundary slopes of the spline with two sampling points. | |
template<class ScalarArrayX , class ScalarArrayY > | |
void | setXYArrays (size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, Scalar m0, Scalar m1, bool sortInputs=true) |
Set the sampling points and the boundary slopes of a full spline using C-style arrays. | |
template<class ScalarContainerX , class ScalarContainerY > | |
void | setXYContainers (const ScalarContainerX &x, const ScalarContainerY &y, Scalar m0, Scalar m1, bool sortInputs=true) |
Set the sampling points and the boundary slopes of a full spline using STL-compatible containers. | |
template<class PointArray > | |
void | setArrayOfPoints (size_t nSamples, const PointArray &points, Scalar m0, Scalar m1, bool sortInputs=true) |
Set the sampling points and the boundary slopes of a full spline using a C-style array. | |
template<class XYContainer > | |
void | setContainerOfPoints (const XYContainer &points, Scalar m0, Scalar m1, bool sortInputs=true) |
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of array-like objects. | |
template<class XYContainer > | |
void | setContainerOfTuples (const XYContainer &points, Scalar m0, Scalar m1, bool sortInputs=true) |
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of tuple-like objects. | |
template<class ScalarArrayX , class ScalarArrayY > | |
void | setXYArrays (size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, SplineType splineType=Natural, bool sortInputs=true) |
Set the sampling points natural spline using C-style arrays. | |
template<class ScalarContainerX , class ScalarContainerY > | |
void | setXYContainers (const ScalarContainerX &x, const ScalarContainerY &y, SplineType splineType=Natural, bool sortInputs=true) |
Set the sampling points of a natural spline using STL-compatible containers. | |
template<class PointArray > | |
void | setArrayOfPoints (size_t nSamples, const PointArray &points, SplineType splineType=Natural, bool sortInputs=true) |
Set the sampling points of a natural spline using a C-style array. | |
template<class XYContainer > | |
void | setContainerOfPoints (const XYContainer &points, SplineType splineType=Natural, bool sortInputs=true) |
Set the sampling points of a natural spline using a STL-compatible container of array-like objects. | |
template<class XYContainer > | |
void | setContainerOfTuples (const XYContainer &points, SplineType splineType=Natural, bool sortInputs=true) |
Set the sampling points of a natural spline using a STL-compatible container of tuple-like objects. | |
template<class Evaluation > | |
bool | applies (const Evaluation &x) const |
Return true iff the given x is in range [x1, xn]. | |
size_t | numSamples () const |
Return the number of (x, y) values. | |
Scalar | xAt (size_t sampleIdx) const |
Return the x value of a given sampling point. | |
Scalar | valueAt (size_t sampleIdx) const |
Return the x value of a given sampling point. | |
void | printCSV (Scalar xi0, Scalar xi1, size_t k, std::ostream &os=std::cout) const |
Prints k tuples of the format (x, y, dx/dy, isMonotonic) to stdout. | |
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 spline's second derivative at a given position. | |
template<class Evaluation > | |
Evaluation | evalThirdDerivative (const Evaluation &x, bool extrapolate=false) const |
Evaluate the spline's third derivative at a given position. | |
template<class Evaluation > | |
Evaluation | intersect (const Evaluation &a, const Evaluation &b, const Evaluation &c, const Evaluation &d) const |
Find the intersections of the spline with a cubic polynomial in the whole interval, throws Opm::MathError exception if there is more or less than one solution. | |
template<class Evaluation > | |
Evaluation | intersectInterval (Scalar x0, Scalar x1, const Evaluation &a, const Evaluation &b, const Evaluation &c, const Evaluation &d) const |
Find the intersections of the spline with a cubic polynomial in a sub-interval of the spline, throws Opm::MathError exception if there is more or less than one solution. | |
int | monotonic (Scalar x0, Scalar x1, bool extrapolate OPM_OPTIM_UNUSED=false) const |
Returns 1 if the spline is monotonically increasing, -1 if the spline is mononously decreasing and 0 if the spline is not monotonous in the interval (x0, x1). | |
int | monotonic () const |
Same as monotonic(x0, x1), but with the entire range of the spline as interval. | |
Protected Member Functions | |
void | sortInput_ () |
Sort the sample points in ascending order of their x value. | |
void | reverseSamplingPoints_ () |
Reverse order of the elements in the arrays which contain the sampling points. | |
void | setNumSamples_ (size_t nSamples) |
Resizes the internal vectors to store the sample points. | |
void | makeFullSpline_ (Scalar m0, Scalar m1) |
Create a natural spline from the already set sampling points. | |
void | makeNaturalSpline_ () |
Create a natural spline from the already set sampling points. | |
void | makePeriodicSpline_ () |
Create a periodic spline from the already set sampling points. | |
template<class DestVector , class SourceVector > | |
void | assignSamplingPoints_ (DestVector &destX, DestVector &destY, const SourceVector &srcX, const SourceVector &srcY, unsigned nSamples) |
Set the sampling point vectors. | |
template<class DestVector , class ListIterator > | |
void | assignFromArrayList_ (DestVector &destX, DestVector &destY, const ListIterator &srcBegin, const ListIterator &srcEnd, unsigned nSamples) |
template<class DestVector , class ListIterator > | |
void | assignFromTupleList_ (DestVector &destX, DestVector &destY, ListIterator srcBegin, ListIterator srcEnd, unsigned nSamples) |
Set the sampling points. | |
template<class Vector , class Matrix > | |
void | makeFullSystem_ (Matrix &M, Vector &d, Scalar m0, Scalar m1) |
Make the linear system of equations Mx = d which results in the moments of the full spline. | |
template<class Vector , class Matrix > | |
void | makeNaturalSystem_ (Matrix &M, Vector &d) |
Make the linear system of equations Mx = d which results in the moments of the natural spline. | |
template<class Matrix , class Vector > | |
void | makePeriodicSystem_ (Matrix &M, Vector &d) |
Make the linear system of equations Mx = d which results in the moments of the periodic spline. | |
template<class Vector > | |
void | makeMonotonicSpline_ (Vector &slopes) |
Create a monotonic spline from the already set sampling points. | |
template<class MomentsVector , class SlopeVector > | |
void | setSlopesFromMoments_ (SlopeVector &slopes, const MomentsVector &moments) |
Convert the moments at the sample points to slopes. | |
template<class Evaluation > | |
Evaluation | eval_ (const Evaluation &x, size_t i) const |
template<class Evaluation > | |
Evaluation | evalDerivative_ (const Evaluation &x, size_t i) const |
template<class Evaluation > | |
Evaluation | evalDerivative2_ (const Evaluation &x, size_t i) const |
template<class Evaluation > | |
Evaluation | evalDerivative3_ (const Evaluation &x, size_t i) const |
template<class Evaluation > | |
Evaluation | h00_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h10_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h01_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h11_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h00_prime_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h10_prime_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h01_prime_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h11_prime_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h00_prime2_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h10_prime2_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h01_prime2_ (const Evaluation &t) const |
template<class Evaluation > | |
Evaluation | h11_prime2_ (const Evaluation &t) const |
template<class Evaluation > | |
Scalar | h00_prime3_ (const Evaluation &t OPM_UNUSED) const |
template<class Evaluation > | |
Scalar | h10_prime3_ (const Evaluation &t OPM_UNUSED) const |
template<class Evaluation > | |
Scalar | h01_prime3_ (const Evaluation &t OPM_UNUSED) const |
template<class Evaluation > | |
Scalar | h11_prime3_ (const Evaluation &t OPM_UNUSED) const |
int | monotonic_ (size_t i, Scalar x0, Scalar x1, int &r) const |
template<class Evaluation > | |
size_t | intersectSegment_ (Evaluation *sol, size_t segIdx, const Evaluation &a, const Evaluation &b, const Evaluation &c, const Evaluation &d, Scalar x0=-1e30, Scalar x1=1e30) const |
Find all the intersections of a segment of the spline with a cubic polynomial within a specified interval. | |
size_t | segmentIdx_ (Scalar x) const |
Scalar | h_ (size_t i) const |
Returns x[i] - x[i - 1]. | |
Scalar | x_ (size_t i) const |
Returns the y coordinate of the i-th sampling point. | |
Scalar | y_ (size_t i) const |
Returns the y coordinate of the i-th sampling point. | |
Scalar | slope_ (size_t i) const |
Returns the slope (i.e. | |
Scalar | a_ (size_t i) const |
Scalar | b_ (size_t i) const |
Scalar | c_ (size_t i) const |
Scalar | d_ (size_t i) const |
Protected Attributes | |
Vector | xPos_ |
Vector | yPos_ |
Vector | slopeVec_ |
Class implementing cubic splines.
This class supports full, natural, periodic and monotonic cubic splines.
Full a splines are splines which, given
sampling points
, fulfill the following conditions
for any given boundary slopes and
.
Natural splines which are defined by
For periodic splines of splines the slopes at the boundaries are identical:
Finally, there are monotonic splines which guarantee that the curve is confined by its sampling points, i.e.,
For more information on monotonic splines, see http://en.wikipedia.org/wiki/Monotone_cubic_interpolation
Full, natural and periodic splines are continuous in their first and second derivatives, i.e.,
holds for such splines. Monotonic splines are only continuous up to their first derivative, i.e., for these only
is true.
enum Opm::Spline::SplineType |
The type of the spline to be created.
To specfiy the acutal curve, use one of the set() methods.
Opm::Spline< Scalar >::Spline | ( | ) | [inline] |
Default constructor for a spline.
To specfiy the acutal curve, use one of the set() methods.
Opm::Spline< Scalar >::Spline | ( | Scalar | x0, | |
Scalar | x1, | |||
Scalar | y0, | |||
Scalar | y1, | |||
Scalar | m0, | |||
Scalar | m1 | |||
) | [inline] |
Convenience constructor for a full spline with just two sampling points.
x0 | The ![]() | |
x1 | The ![]() | |
y0 | The ![]() | |
y1 | The ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() |
Opm::Spline< Scalar >::Spline | ( | size_t | nSamples, | |
const ScalarArrayX & | x, | |||
const ScalarArrayY & | y, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a natural or a periodic spline.
nSamples | The number of sampling points (must be > 2) | |
x | An array containing the ![]() | |
y | An array containing the ![]() | |
periodic | Indicates whether a natural or a periodic spline should be created |
Opm::Spline< Scalar >::Spline | ( | size_t | nSamples, | |
const PointArray & | points, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a natural or a periodic spline.
nSamples | The number of sampling points (must be > 2) | |
points | An array of ![]() | |
periodic | Indicates whether a natural or a periodic spline should be created |
Opm::Spline< Scalar >::Spline | ( | const ScalarContainer & | x, | |
const ScalarContainer & | y, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a natural or a periodic spline.
x | An array containing the ![]() | |
y | An array containing the ![]() | |
periodic | Indicates whether a natural or a periodic spline should be created |
Opm::Spline< Scalar >::Spline | ( | const PointContainer & | points, | |
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a natural or a periodic spline.
points | An array of ![]() | |
periodic | Indicates whether a natural or a periodic spline should be created |
Opm::Spline< Scalar >::Spline | ( | size_t | nSamples, | |
const ScalarArray & | x, | |||
const ScalarArray & | y, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a full spline.
nSamples | The number of sampling points (must be >= 2) | |
x | An array containing the ![]() | |
y | An array containing the ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() | |
sortInputs | Indicates whether the sample points should be sorted (this is not necessary if they are already sorted in ascending or descending order) |
Opm::Spline< Scalar >::Spline | ( | size_t | nSamples, | |
const PointArray & | points, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a full spline.
nSamples | The number of sampling points (must be >= 2) | |
points | An array containing the ![]() ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() | |
sortInputs | Indicates whether the sample points should be sorted (this is not necessary if they are already sorted in ascending or descending order) |
Opm::Spline< Scalar >::Spline | ( | const ScalarContainerX & | x, | |
const ScalarContainerY & | y, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a full spline.
x | An array containing the ![]() | |
y | An array containing the ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() | |
sortInputs | Indicates whether the sample points should be sorted (this is not necessary if they are already sorted in ascending or descending order) |
Opm::Spline< Scalar >::Spline | ( | const PointContainer & | points, | |
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Convenience constructor for a full spline.
points | An array of ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() | |
sortInputs | Indicates whether the sample points should be sorted (this is not necessary if they are already sorted in ascending or descending order) |
void Opm::Spline< Scalar >::assignFromTupleList_ | ( | DestVector & | destX, | |
DestVector & | destY, | |||
ListIterator | srcBegin, | |||
ListIterator | srcEnd, | |||
unsigned | nSamples | |||
) | [inline, protected] |
Set the sampling points.
Here we assume that the elements of the source vector have an [] operator where v[0] is the x value and v[1] is the y value if the sampling point.
void Opm::Spline< Scalar >::assignSamplingPoints_ | ( | DestVector & | destX, | |
DestVector & | destY, | |||
const SourceVector & | srcX, | |||
const SourceVector & | srcY, | |||
unsigned | nSamples | |||
) | [inline, protected] |
Set the sampling point vectors.
This takes care that the order of the x-values is ascending, although the input must be ordered already!
Evaluation Opm::Spline< 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 spline 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::Spline< 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 spline'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::Spline< Scalar >::evalSecondDerivative | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the spline's second derivative at a given position.
x | The value on the abscissa where the spline'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::Spline< Scalar >::evalThirdDerivative | ( | const Evaluation & | x, | |
bool | extrapolate = false | |||
) | const [inline] |
Evaluate the spline's third derivative at a given position.
x | The value on the abscissa where the spline'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 ![]() |
void Opm::Spline< Scalar >::makeFullSpline_ | ( | Scalar | m0, | |
Scalar | m1 | |||
) | [inline, protected] |
Create a natural spline from the already set sampling points.
This creates a temporary matrix and right hand side vector.
void Opm::Spline< Scalar >::makeMonotonicSpline_ | ( | Vector & | slopes | ) | [inline, protected] |
Create a monotonic spline from the already set sampling points.
This code is inspired by opm-core's "MonotCubicInterpolator" class and also uses the approach by Fritsch and Carlson, see
void Opm::Spline< Scalar >::makeNaturalSpline_ | ( | ) | [inline, protected] |
Create a natural spline from the already set sampling points.
This creates a temporary matrix and right hand side vector.
void Opm::Spline< Scalar >::makePeriodicSpline_ | ( | ) | [inline, protected] |
Create a periodic spline from the already set sampling points.
This creates a temporary matrix and right hand side vector.
int Opm::Spline< Scalar >::monotonic | ( | Scalar | x0, | |
Scalar | x1, | |||
bool extrapolate | OPM_OPTIM_UNUSED = false | |||
) | const [inline] |
Returns 1 if the spline is monotonically increasing, -1 if the spline is mononously decreasing and 0 if the spline is not monotonous in the interval (x0, x1).
In the corner case that the spline is constant within the given interval, this method returns 3.
void Opm::Spline< Scalar >::printCSV | ( | Scalar | xi0, | |
Scalar | xi1, | |||
size_t | k, | |||
std::ostream & | os = std::cout | |||
) | const [inline] |
Prints k tuples of the format (x, y, dx/dy, isMonotonic) to stdout.
If the spline 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 > spline.csv gnuplot
gnuplot> plot "spline.csv" using 1:2 w l ti "Curve", \ "spline.csv" using 1:3 w l ti "Derivative", \ "spline.csv" using 1:4 w p ti "Monotonic" ----------- snap -----------
void Opm::Spline< Scalar >::set | ( | Scalar | x0, | |
Scalar | x1, | |||
Scalar | y0, | |||
Scalar | y1, | |||
Scalar | m0, | |||
Scalar | m1 | |||
) | [inline] |
Set the sampling points and the boundary slopes of the spline with two sampling points.
x0 | The ![]() | |
x1 | The ![]() | |
y0 | The ![]() | |
y1 | The ![]() | |
m0 | The slope of the spline at ![]() | |
m1 | The slope of the spline at ![]() |
void Opm::Spline< Scalar >::setArrayOfPoints | ( | size_t | nSamples, | |
const PointArray & | points, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points of a natural spline using a C-style array.
This method uses a single array of sampling points, which are seen as an array-like object which provides access to the X and Y coordinates. In this context 'array-like' means that an access to the members is provided via the [] operator. (e.g. C arrays, std::vector, std::array, etc.) The array containing the sampling points must be of size 'nSamples' at least. Also, the number of sampling points must be larger than 1.
void Opm::Spline< Scalar >::setArrayOfPoints | ( | size_t | nSamples, | |
const PointArray & | points, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points and the boundary slopes of a full spline using a C-style array.
This method uses a single array of sampling points, which are seen as an array-like object which provides access to the X and Y coordinates. In this context 'array-like' means that an access to the members is provided via the [] operator. (e.g. C arrays, std::vector, std::array, etc.) The array containing the sampling points must be of size 'nSamples' at least. Also, the number of sampling points must be larger than 1.
void Opm::Spline< Scalar >::setContainerOfPoints | ( | const XYContainer & | points, | |
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points of a natural spline using a STL-compatible container of array-like objects.
This method uses a single STL-compatible container of sampling points, which are assumed to be array-like objects storing the X and Y coordinates. "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::Spline< Scalar >::setContainerOfPoints | ( | const XYContainer & | points, | |
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of array-like objects.
This method uses a single STL-compatible container of sampling points, which are assumed to be array-like objects storing the X and Y coordinates. "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::Spline< Scalar >::setContainerOfTuples | ( | const XYContainer & | points, | |
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points of a natural spline 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::Spline< Scalar >::setContainerOfTuples | ( | const XYContainer & | points, | |
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points and the boundary slopes of a full spline 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::Spline< Scalar >::setSlopesFromMoments_ | ( | SlopeVector & | slopes, | |
const MomentsVector & | moments | |||
) | [inline, protected] |
Convert the moments at the sample points to slopes.
This requires to use cubic Hermite interpolation, but it is required because for monotonic splines the second derivative is not continuous.
void Opm::Spline< Scalar >::setXYArrays | ( | size_t | nSamples, | |
const ScalarArrayX & | x, | |||
const ScalarArrayY & | y, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points natural spline using C-style arrays.
This method uses separate array-like objects for the values of the X and Y coordinates. In this context 'array-like' means that an access to the members is provided via the [] operator. (e.g. C arrays, std::vector, std::array, etc.) Each array must be of size 'nSamples' at least. Also, the number of sampling points must be larger than 1.
void Opm::Spline< Scalar >::setXYArrays | ( | size_t | nSamples, | |
const ScalarArrayX & | x, | |||
const ScalarArrayY & | y, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points and the boundary slopes of a full spline using C-style arrays.
This method uses separate array-like objects for the values of the X and Y coordinates. In this context 'array-like' means that an access to the members is provided via the [] operator. (e.g. C arrays, std::vector, std::array, etc.) Each array must be of size 'nSamples' at least. Also, the number of sampling points must be larger than 1.
void Opm::Spline< Scalar >::setXYContainers | ( | const ScalarContainerX & | x, | |
const ScalarContainerY & | y, | |||
SplineType | splineType = Natural , |
|||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points of a natural spline using STL-compatible containers.
This method uses separate STL-compatible containers for the values of the sampling points' X and Y coordinates. "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::Spline< Scalar >::setXYContainers | ( | const ScalarContainerX & | x, | |
const ScalarContainerY & | y, | |||
Scalar | m0, | |||
Scalar | m1, | |||
bool | sortInputs = true | |||
) | [inline] |
Set the sampling points and the boundary slopes of a full spline using STL-compatible containers.
This method uses separate STL-compatible containers for the values of the sampling points' X and Y coordinates. "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.
Scalar Opm::Spline< Scalar >::slope_ | ( | size_t | i | ) | const [inline, protected] |
Returns the slope (i.e.
first derivative) of the spline at the i-th sampling point.