25 #ifndef OPM_SPLINE_HPP
26 #define OPM_SPLINE_HPP
30 #include <opm/common/ErrorMacros.hpp>
88 template<
class Scalar>
92 typedef std::vector<Scalar> Vector;
126 Scalar y0, Scalar y1,
127 Scalar m0, Scalar m1)
128 {
set(x0, x1, y0, y1, m0, m1); }
138 template <
class ScalarArrayX,
class ScalarArrayY>
140 const ScalarArrayX &x,
141 const ScalarArrayY &y,
143 bool sortInputs =
false)
144 { this->
setXYArrays(nSamples, x, y, splineType, sortInputs); }
153 template <
class Po
intArray>
155 const PointArray &points,
157 bool sortInputs =
false)
167 template <
class ScalarContainer>
169 const ScalarContainer &y,
171 bool sortInputs =
false)
180 template <
class Po
intContainer>
183 bool sortInputs =
false)
196 template <
class ScalarArray>
198 const ScalarArray &x,
199 const ScalarArray &y,
202 bool sortInputs =
false)
203 { this->
setXYArrays(nSamples, x, y, m0, m1, sortInputs); }
214 template <
class Po
intArray>
216 const PointArray &points,
219 bool sortInputs =
false)
231 template <
class ScalarContainerX,
class ScalarContainerY>
233 const ScalarContainerY &y,
236 bool sortInputs =
false)
247 template <
class Po
intContainer>
251 bool sortInputs =
false)
258 {
return xPos_.size(); }
271 void set(Scalar x0, Scalar x1,
272 Scalar y0, Scalar y1,
273 Scalar m0, Scalar m1)
326 template <
class ScalarArrayX,
class ScalarArrayY>
328 const ScalarArrayX &x,
329 const ScalarArrayY &y,
330 Scalar m0, Scalar m1,
331 bool sortInputs =
false)
333 assert(nSamples > 1);
336 for (
size_t i = 0; i < nSamples; ++i) {
360 template <
class ScalarContainerX,
class ScalarContainerY>
362 const ScalarContainerY &y,
363 Scalar m0, Scalar m1,
364 bool sortInputs =
false)
366 assert(x.size() == y.size());
367 assert(x.size() > 1);
371 std::copy(x.begin(), x.end(),
xPos_.begin());
372 std::copy(y.begin(), y.end(),
yPos_.begin());
394 template <
class Po
intArray>
396 const PointArray &points,
399 bool sortInputs =
false)
403 assert(nSamples > 1);
406 for (
size_t i = 0; i < nSamples; ++i) {
407 xPos_[i] = points[i][0];
408 yPos_[i] = points[i][1];
431 template <
class XYContainer>
435 bool sortInputs =
false)
439 assert(points.size() > 1);
442 typename XYContainer::const_iterator it = points.begin();
443 typename XYContainer::const_iterator endIt = points.end();
444 for (
size_t i = 0; it != endIt; ++i, ++it) {
473 template <
class XYContainer>
477 bool sortInputs =
false)
481 typename XYContainer::const_iterator it = points.begin();
482 typename XYContainer::const_iterator endIt = points.end();
483 for (
int i = 0; it != endIt; ++i, ++it) {
484 xPos_[i] = std::get<0>(*it);
485 yPos_[i] = std::get<1>(*it);
514 template <
class ScalarArrayX,
class ScalarArrayY>
516 const ScalarArrayX &x,
517 const ScalarArrayY &y,
519 bool sortInputs =
false)
521 assert(nSamples > 1);
524 for (
size_t i = 0; i < nSamples; ++i) {
536 else if (splineType ==
Natural)
541 OPM_THROW(std::runtime_error,
"Spline type " << splineType <<
" not supported at this place");
555 template <
class ScalarContainerX,
class ScalarContainerY>
557 const ScalarContainerY &y,
559 bool sortInputs =
false)
561 assert(x.size() == y.size());
562 assert(x.size() > 1);
565 std::copy(x.begin(), x.end(),
xPos_.begin());
566 std::copy(y.begin(), y.end(),
yPos_.begin());
575 else if (splineType ==
Natural)
580 OPM_THROW(std::runtime_error,
"Spline type " << splineType <<
" not supported at this place");
595 template <
class Po
intArray>
597 const PointArray &points,
599 bool sortInputs =
false)
603 assert(nSamples > 1);
606 for (
size_t i = 0; i < nSamples; ++i) {
607 xPos_[i] = points[i][0];
608 yPos_[i] = points[i][1];
618 else if (splineType ==
Natural)
623 OPM_THROW(std::runtime_error,
"Spline type " << splineType <<
" not supported at this place");
637 template <
class XYContainer>
640 bool sortInputs =
false)
644 assert(points.size() > 1);
647 typename XYContainer::const_iterator it = points.begin();
648 typename XYContainer::const_iterator endIt = points.end();
649 for (
size_t i = 0; it != endIt; ++ i, ++it) {
661 else if (splineType ==
Natural)
666 OPM_THROW(std::runtime_error,
"Spline type " << splineType <<
" not supported at this place");
683 template <
class XYContainer>
686 bool sortInputs =
false)
690 typename XYContainer::const_iterator it = points.begin();
691 typename XYContainer::const_iterator endIt = points.end();
692 for (
unsigned i = 0; it != endIt; ++i, ++it) {
693 xPos_[i] = std::get<0>(*it);
694 yPos_[i] = std::get<1>(*it);
704 else if (splineType ==
Natural)
709 OPM_THROW(std::runtime_error,
"Spline type " << splineType <<
" not supported at this place");
761 void printCSV(Scalar xi0, Scalar xi1,
size_t k, std::ostream &os = std::cout)
const
766 for (
size_t i = 0; i <= k; ++i) {
767 double x = i*(x1 - x0)/k + x0;
768 double x_p1 = x + (x1 - x0)/k;
775 y = (x -
x_(0))*dy_dx +
y_(0);
776 mono = (dy_dx>0)?1:-1;
778 else if (x >
x_(n)) {
780 y = (x -
x_(n))*dy_dx +
y_(n);
781 mono = (dy_dx>0)?1:-1;
784 OPM_THROW(std::runtime_error,
785 "The sampling points given to a spline must be sorted by their x value!");
794 os << x <<
" " << y <<
" " << dy_dx <<
" " << mono <<
"\n";
809 Scalar
eval(Scalar x,
bool extrapolate=
false)
const
811 assert(extrapolate ||
applies(x));
818 return y0 + m*(x -
xMin());
820 else if (x >
xMax()) {
823 return y0 + m*(x -
xMax());
841 template <
class Evaluation>
842 Evaluation
eval(
const Evaluation& x,
bool extrapolate=
false)
const
844 assert(extrapolate ||
applies(x.value));
848 if (x.value <
xMin()) {
851 return Evaluation::createConstant(y0 + m*(x.value -
xMin()));
853 else if (x >
xMax()) {
856 return Evaluation::createConstant(y0 + m*(x.value -
xMax()));
877 assert(extrapolate ||
applies(x));
902 assert(extrapolate ||
applies(x));
923 assert(extrapolate ||
applies(x));
936 template <
class Evaluation>
940 const Evaluation& d)
const
951 template <
class Evaluation>
956 const Evaluation& d)
const
960 Evaluation tmpSol[3], sol = 0;
964 for (
size_t i = iFirst; i <= iLast; ++i)
972 OPM_THROW(std::runtime_error,
973 "Spline has more than one intersection");
978 OPM_THROW(std::runtime_error,
979 "Spline has no intersection");
1004 assert(extrapolate);
1012 if (
x_(i + 1) >= x1) {
1027 for (; i < iEnd - 1; ++i) {
1037 assert(extrapolate);
1041 return (r < 0 || r==3)?-1:0;
1043 return (r > 0 || r==3)?1:0;
1072 {
return x_.at(idxA) <
x_.at(idxB); }
1074 const std::vector<Scalar> &
x_;
1085 std::vector<unsigned> idxVector(static_cast<unsigned>(n));
1086 for (
unsigned i = 0; i < n; ++i)
1087 idxVector[static_cast<unsigned>(i)] = i;
1092 std::sort(idxVector.begin(), idxVector.end(), cmp);
1095 std::vector<Scalar> tmpX(n), tmpY(n);
1096 for (
size_t i = 0; i < idxVector.size(); ++ i) {
1097 tmpX[i] =
xPos_[idxVector[i]];
1098 tmpY[i] =
yPos_[idxVector[i]];
1112 for (
unsigned i = 0; i <= (n - 1)/2; ++i) {
1123 xPos_.resize(nSamples);
1124 yPos_.resize(nSamples);
1143 M.
solve(moments, d);
1164 M.
solve(moments, d);
1188 M.
solve(moments, d);
1191 for (
int i = static_cast<int>(
numSamples()) - 2; i >= 0; --i)
1192 moments[static_cast<unsigned>(i+1)] = moments[
static_cast<unsigned>(i)];
1205 template <
class DestVector,
class SourceVector>
1208 const SourceVector &srcX,
1209 const SourceVector &srcY,
1212 assert(numSamples >= 2);
1218 if (srcX[0] > srcX[numSamples - 1])
1219 idx = numSamples - i - 1;
1220 destX[i] = srcX[idx];
1221 destY[i] = srcY[idx];
1225 template <
class DestVector,
class ListIterator>
1228 const ListIterator &srcBegin,
1229 const ListIterator &srcEnd,
1232 assert(numSamples >= 2);
1235 ListIterator it = srcBegin;
1237 bool reverse =
false;
1238 if ((*srcBegin)[0] > (*it)[0])
1243 for (
int i = 0; it != srcEnd; ++i, ++it) {
1246 idx = numSamples - i - 1;
1247 destX[i] = (*it)[0];
1248 destY[i] = (*it)[1];
1259 template <
class DestVector,
class ListIterator>
1262 ListIterator srcBegin,
1263 ListIterator srcEnd,
1266 assert(numSamples >= 2);
1272 ListIterator it = srcBegin;
1274 bool reverse =
false;
1275 if (std::get<0>(*srcBegin) > std::get<0>(*it))
1280 for (
int i = 0; it != srcEnd; ++i, ++it) {
1283 idx = numSamples - i - 1;
1284 destX[i] = std::get<0>(*it);
1285 destY[i] = std::get<1>(*it);
1293 template <
class Vector,
class Matrix>
1301 d[0] = 6/
h_(1) * ( (
y_(1) -
y_(0))/
h_(1) - m0);
1310 (m1 - (
y_(n) -
y_(n - 1))/
h_(n));
1317 template <
class Vector,
class Matrix>
1327 for (
size_t i = 1; i < n; ++i) {
1328 Scalar lambda_i =
h_(i + 1) / (
h_(i) +
h_(i + 1));
1329 Scalar mu_i = 1 - lambda_i;
1331 6 / (
h_(i) +
h_(i + 1))
1333 ( (
y_(i + 1) -
y_(i))/
h_(i + 1) - (
y_(i) -
y_(i - 1))/
h_(i));
1337 M[i][i + 1] = lambda_i;
1342 Scalar lambda_0 = 0;
1363 template <
class Matrix,
class Vector>
1372 assert(M.
rows() == n);
1375 for (
size_t i = 2; i < n; ++i) {
1376 Scalar lambda_i =
h_(i + 1) / (
h_(i) +
h_(i + 1));
1377 Scalar mu_i = 1 - lambda_i;
1379 6 / (
h_(i) +
h_(i + 1))
1381 ( (
y_(i + 1) -
y_(i))/
h_(i + 1) - (
y_(i) -
y_(i - 1))/
h_(i));
1385 M[i-1][i] = lambda_i;
1389 Scalar lambda_n =
h_(1) / (
h_(n) +
h_(1));
1390 Scalar lambda_1 =
h_(2) / (
h_(1) +
h_(2));;
1391 Scalar mu_1 = 1 - lambda_1;
1392 Scalar mu_n = 1 - lambda_n;
1411 M[n-1][0] = lambda_n;
1425 template <
class Vector>
1431 std::vector<Scalar> delta(n);
1432 for (
size_t k = 0; k < n - 1; ++k)
1433 delta[k] = (
y_(k + 1) -
y_(k))/(
x_(k + 1) -
x_(k));
1436 for (
size_t k = 1; k < n - 1; ++k)
1437 slopes[k] = (delta[k - 1] + delta[k])/2;
1438 slopes[0] = delta[0];
1439 slopes[n - 1] = delta[n - 2];
1442 for (
size_t k = 0; k < n - 1; ++k) {
1451 Scalar alpha = slopes[k] / delta[k];
1452 Scalar beta = slopes[k + 1] / delta[k];
1454 if (alpha < 0 || (k > 0 && slopes[k] / delta[k - 1] < 0)) {
1458 else if (alpha*alpha + beta*beta > 3*3) {
1459 Scalar tau = 3.0/
std::sqrt(alpha*alpha + beta*beta);
1460 slopes[k] = tau*alpha*delta[k];
1461 slopes[k + 1] = tau*beta*delta[k];
1474 template <
class MomentsVector,
class SlopeVector>
1485 Scalar h = this->
h_(n - 1);
1490 (
y_(n - 1) -
y_(n - 2))/h
1492 h/6*(moments[n-1] - moments[n - 2]);
1497 moments[n - 1] * x*x / (2 * h)
1503 for (
size_t i = 0; i < n - 1; ++ i) {
1506 Scalar h_i = this->
h_(i + 1);
1511 (
y_(i+1) -
y_(i))/h_i
1513 h_i/6*(moments[i+1] - moments[i]);
1516 - moments[i] * x_i1*x_i1 / (2 * h_i)
1523 slopes[n - 1] = mRight;
1532 Scalar delta =
h_(i + 1);
1533 Scalar t = (x -
x_(i))/delta;
1544 template <
class Evaluation>
1545 Evaluation
eval_(
const Evaluation& x,
size_t i)
const
1548 Scalar delta =
h_(i + 1);
1549 Scalar t = (x.value -
x_(i))/delta;
1559 for (
unsigned varIdx = 0; varIdx < result.derivatives.size(); ++ varIdx)
1560 result.derivatives[varIdx] = df_dg*x.derivatives[varIdx];
1570 Scalar delta =
h_(i + 1);
1571 Scalar t = (x -
x_(i))/delta;
1572 Scalar alpha = 1 / delta;
1587 Scalar delta =
h_(i + 1);
1588 Scalar t = (x -
x_(i))/delta;
1589 Scalar alpha = 1 / delta;
1604 Scalar delta =
h_(i + 1);
1605 Scalar t = (x -
x_(i))/delta;
1606 Scalar alpha = 1 / delta;
1618 {
return (2*t - 3)*t*t + 1; }
1621 {
return ((t - 2)*t + 1)*t; }
1624 {
return (-2*t + 3)*t*t; }
1627 {
return (t - 1)*t*t; }
1631 {
return (3*2*t - 2*3)*t; }
1634 {
return (3*t - 2*2)*t + 1; }
1637 {
return (-3*2*t + 2*3)*t; }
1640 {
return (3*t - 2)*t; }
1644 {
return 2*3*2*t - 2*3; }
1647 {
return 2*3*t - 2*2; }
1650 {
return -2*3*2*t + 2*3; }
1653 {
return 2*3*t - 2; }
1686 Scalar disc = b*b - 4*a*c;
1690 if (x0*(x0*a + b) + c > 0) {
1691 r = (r==3 || r == 1)?1:0;
1695 r = (r==3 || r == -1)?-1:0;
1700 Scalar xE1 = (-b + disc)/(2*a);
1701 Scalar xE2 = (-b - disc)/(2*a);
1710 if (x0*(x0*a + b) + c > 0) {
1711 r = (r==3 || r == 1)?1:0;
1715 r = (r==3 || r == -1)?-1:0;
1719 if ((x0 < xE1 && xE1 < x1) ||
1720 (x0 < xE2 && xE2 < x1))
1729 if (x0*(x0*a + b) + c > 0) {
1730 r = (r==3 || r == 1)?1:0;
1734 r = (r==3 || r == -1)?-1:0;
1743 template <
class Evaluation>
1746 const Evaluation& a,
1747 const Evaluation& b,
1748 const Evaluation& c,
1749 const Evaluation& d,
1750 Scalar x0 = -1e100, Scalar x1 = 1e100)
const
1762 for (
int j = 0; j < n; ++j) {
1763 if (x0 <= sol[j] && sol[j] <= x1) {
1778 while (iLow + 1 < iHigh) {
1779 size_t i = (iLow + iHigh) / 2;
1791 Scalar
h_(
size_t i)
const
1793 assert(
x_(i) >
x_(i-1));
1795 return x_(i) -
x_(i - 1);
1801 Scalar
x_(
size_t i)
const
1802 {
return xPos_[i]; }
1807 Scalar
y_(
size_t i)
const
1808 {
return yPos_[i]; }
1819 Scalar
a_(
size_t i)
const
1824 Scalar
b_(
size_t i)
const
1829 Scalar
c_(
size_t i)
const
1834 Scalar
d_(
size_t i)
const
1835 {
return eval_(Scalar(0.0), i); }
int monotonic() const
Same as monotonic(x0, x1), but with the entire range of the spline as interval.
Definition: Spline.hpp:1058
void solve(XVector &x, const BVector &b) const
Calculate the solution for a linear system of equations.
Definition: TridiagonalMatrix.hpp:712
Scalar h00_(Scalar t) const
Definition: Spline.hpp:1617
Scalar h10_(Scalar t) const
Definition: Spline.hpp:1620
Provides a tridiagonal matrix that also supports non-zero entries in the upper right and lower left...
Definition: TridiagonalMatrix.hpp:48
void assignSamplingPoints_(DestVector &destX, DestVector &destY, const SourceVector &srcX, const SourceVector &srcY, int numSamples)
Set the sampling point vectors.
Definition: Spline.hpp:1206
Definition: Spline.hpp:104
Spline()
Default constructor for a spline.
Definition: Spline.hpp:112
Scalar evalDerivative2_(Scalar x, size_t i) const
Definition: Spline.hpp:1584
int monotonic(Scalar x0, Scalar x1, OPM_OPTIM_UNUSED bool extrapolate=false) const
Returns 1 if the spline is monotonically increasing, -1 if the spline is mononously decreasing and 0 ...
Definition: Spline.hpp:992
Scalar yLast() const
Return the y value of the rightmost sampling point.
Definition: Spline.hpp:741
void setXYArrays(size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using C-style arrays.
Definition: Spline.hpp:327
Scalar d_(size_t i) const
Definition: Spline.hpp:1834
Vector slopeVec_
Definition: Spline.hpp:1839
Definition: Air_Mesitylene.hpp:31
void assignFromArrayList_(DestVector &destX, DestVector &destY, const ListIterator &srcBegin, const ListIterator &srcEnd, int numSamples)
Definition: Spline.hpp:1226
ComparatorX_(const std::vector< Scalar > &x)
Definition: Spline.hpp:1067
void setContainerOfPoints(const XYContainer &points, SplineType splineType=Natural, bool sortInputs=false)
Set the sampling points of a natural spline using a STL-compatible container of array-like objects...
Definition: Spline.hpp:638
Scalar y_(size_t i) const
Returns the y coordinate of the i-th sampling point.
Definition: Spline.hpp:1807
bool applies(Scalar x) const
Return true iff the given x is in range [x1, xn].
Definition: Spline.hpp:715
Evaluation< Scalar, VarSetTag, numVars > max(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:114
Scalar h11_prime_(Scalar t) const
Definition: Spline.hpp:1639
Evaluation< Scalar, VarSetTag, numVars > sqrt(const Evaluation< Scalar, VarSetTag, numVars > &x)
Definition: Math.hpp:278
void setContainerOfTuples(const XYContainer &points, SplineType splineType=Natural, bool sortInputs=false)
Set the sampling points of a natural spline using a STL-compatible container of tuple-like objects...
Definition: Spline.hpp:684
void setContainerOfTuples(const XYContainer &points, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of ...
Definition: Spline.hpp:474
void setContainerOfPoints(const XYContainer &points, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using a STL-compatible container of ...
Definition: Spline.hpp:432
Spline(size_t nSamples, const ScalarArray &x, const ScalarArray &y, Scalar m0, Scalar m1, bool sortInputs=false)
Convenience constructor for a full spline.
Definition: Spline.hpp:197
SplineType
The type of the spline to be created.
Definition: Spline.hpp:100
Spline(size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, SplineType splineType=Natural, bool sortInputs=false)
Convenience constructor for a natural or a periodic spline.
Definition: Spline.hpp:139
Spline(const PointContainer &points, Scalar m0, Scalar m1, bool sortInputs=false)
Convenience constructor for a full spline.
Definition: Spline.hpp:248
Scalar evalSecondDerivative(Scalar x, bool extrapolate=false) const
Evaluate the spline's second derivative at a given position.
Definition: Spline.hpp:900
void makeNaturalSystem_(Matrix &M, Vector &d)
Make the linear system of equations Mx = d which results in the moments of the natural spline...
Definition: Spline.hpp:1318
Scalar evalDerivative3_(Scalar x, size_t i) const
Definition: Spline.hpp:1601
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.
Definition: Spline.hpp:271
Evaluation eval(const Evaluation &x, bool extrapolate=false) const
Evaluate the spline for a given function evaluation.
Definition: Spline.hpp:842
Vector xPos_
Definition: Spline.hpp:1837
int invertCubicPolynomial(SolContainer *sol, Scalar a, Scalar b, Scalar c, Scalar d)
Invert a cubic polynomial analytically.
Definition: PolynomialUtils.hpp:151
void makePeriodicSpline_()
Create a periodic spline from the already set sampling points.
Definition: Spline.hpp:1175
Spline(size_t nSamples, const PointArray &points, Scalar m0, Scalar m1, bool sortInputs=false)
Convenience constructor for a full spline.
Definition: Spline.hpp:215
void makeMonotonicSpline_(Vector &slopes)
Create a monotonic spline from the already set sampling points.
Definition: Spline.hpp:1426
Scalar h11_(Scalar t) const
Definition: Spline.hpp:1626
Scalar xMax() const
Return the x value of the rightmost sampling point.
Definition: Spline.hpp:729
Scalar evalDerivative_(Scalar x, size_t i) const
Definition: Spline.hpp:1567
Scalar eval_(Scalar x, size_t i) const
Definition: Spline.hpp:1529
void assignFromTupleList_(DestVector &destX, DestVector &destY, ListIterator srcBegin, ListIterator srcEnd, int numSamples)
Set the sampling points.
Definition: Spline.hpp:1260
Scalar xMin() const
Return the x value of the leftmost sampling point.
Definition: Spline.hpp:723
void setSlopesFromMoments_(SlopeVector &slopes, const MomentsVector &moments)
Convert the moments at the sample points to slopes.
Definition: Spline.hpp:1475
Scalar h01_prime3_(Scalar) const
Definition: Spline.hpp:1662
void setArrayOfPoints(size_t nSamples, const PointArray &points, SplineType splineType=Natural, bool sortInputs=false)
Set the sampling points of a natural spline using a C-style array.
Definition: Spline.hpp:596
void setXYContainers(const ScalarContainerX &x, const ScalarContainerY &y, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using STL-compatible containers...
Definition: Spline.hpp:361
Scalar yFirst() const
Return the y value of the leftmost sampling point.
Definition: Spline.hpp:735
Definition: Spline.hpp:102
Scalar h00_prime_(Scalar t) const
Definition: Spline.hpp:1630
void sortInput_()
Sort the sample points in ascending order of their x value.
Definition: Spline.hpp:1080
Scalar h00_prime2_(Scalar t) const
Definition: Spline.hpp:1643
Vector yPos_
Definition: Spline.hpp:1838
Provides a tridiagonal matrix that also supports non-zero entries in the upper right and lower left...
Scalar h00_prime3_(Scalar) const
Definition: Spline.hpp:1656
size_t segmentIdx_(Scalar x) const
Definition: Spline.hpp:1772
Evaluation< Scalar, VarSetTag, numVars > abs(const Evaluation< Scalar, VarSetTag, numVars > &)
Definition: Math.hpp:41
Scalar h01_(Scalar t) const
Definition: Spline.hpp:1623
Definition: Spline.hpp:103
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...
Definition: Spline.hpp:952
Class implementing cubic splines.
Definition: Spline.hpp:89
Provides free functions to invert polynomials of degree 1, 2 and 3.
Scalar h01_prime2_(Scalar t) const
Definition: Spline.hpp:1649
Evaluation< Scalar, VarSetTag, numVars > min(const Evaluation< Scalar, VarSetTag, numVars > &x1, const Evaluation< Scalar, VarSetTag, numVars > &x2)
Definition: Math.hpp:61
Scalar evalDerivative(Scalar x, bool extrapolate=false) const
Evaluate the spline's derivative at a given position.
Definition: Spline.hpp:875
Definition: Spline.hpp:101
Scalar h_(size_t i) const
Returns x[i] - x[i - 1].
Definition: Spline.hpp:1791
void setArrayOfPoints(size_t nSamples, const PointArray &points, Scalar m0, Scalar m1, bool sortInputs=false)
Set the sampling points and the boundary slopes of a full spline using a C-style array.
Definition: Spline.hpp:395
#define OPM_OPTIM_UNUSED
Definition: Unused.hpp:42
Spline(const PointContainer &points, SplineType splineType=Natural, bool sortInputs=false)
Convenience constructor for a natural or a periodic spline.
Definition: Spline.hpp:181
const std::vector< Scalar > & x_
Definition: Spline.hpp:1074
void reverseSamplingPoints_()
Reverse order of the elements in the arrays which contain the sampling points.
Definition: Spline.hpp:1108
Scalar eval(Scalar x, bool extrapolate=false) const
Evaluate the spline at a given position.
Definition: Spline.hpp:809
bool operator()(unsigned idxA, unsigned idxB) const
Definition: Spline.hpp:1071
Scalar a_(size_t i) const
Definition: Spline.hpp:1819
void makeNaturalSpline_()
Create a natural spline from the already set sampling points.
Definition: Spline.hpp:1154
Scalar h11_prime3_(Scalar) const
Definition: Spline.hpp:1665
void makeFullSpline_(Scalar m0, Scalar m1)
Create a natural spline from the already set sampling points.
Definition: Spline.hpp:1133
Spline(const ScalarContainer &x, const ScalarContainer &y, SplineType splineType=Natural, bool sortInputs=false)
Convenience constructor for a natural or a periodic spline.
Definition: Spline.hpp:168
Scalar x_(size_t i) const
Returns the y coordinate of the i-th sampling point.
Definition: Spline.hpp:1801
Scalar h10_prime3_(Scalar) const
Definition: Spline.hpp:1659
int monotonic_(size_t i, Scalar x0, Scalar x1, int &r) const
Definition: Spline.hpp:1676
Scalar evalThirdDerivative(Scalar x, bool extrapolate=false) const
Evaluate the spline's third derivative at a given position.
Definition: Spline.hpp:921
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.
Definition: Spline.hpp:761
Scalar b_(size_t i) const
Definition: Spline.hpp:1824
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.
Definition: Spline.hpp:937
Scalar h10_prime_(Scalar t) const
Definition: Spline.hpp:1633
Helper class needed to sort the input sampling points.
Definition: Spline.hpp:1065
void setNumSamples_(size_t nSamples)
Resizes the internal vectors to store the sample points.
Definition: Spline.hpp:1121
void makePeriodicSystem_(Matrix &M, Vector &d)
Make the linear system of equations Mx = d which results in the moments of the periodic spline...
Definition: Spline.hpp:1364
Scalar h01_prime_(Scalar t) const
Definition: Spline.hpp:1636
size_t rows() const
Return the number of rows of the matrix.
Definition: TridiagonalMatrix.hpp:138
Spline(size_t nSamples, const PointArray &points, SplineType splineType=Natural, bool sortInputs=false)
Convenience constructor for a natural or a periodic spline.
Definition: Spline.hpp:154
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...
Definition: Spline.hpp:1294
Scalar h10_prime2_(Scalar t) const
Definition: Spline.hpp:1646
Scalar c_(size_t i) const
Definition: Spline.hpp:1829
Provides the OPM_UNUSED macro.
Evaluation eval_(const Evaluation &x, size_t i) const
Definition: Spline.hpp:1545
Spline(const ScalarContainerX &x, const ScalarContainerY &y, Scalar m0, Scalar m1, bool sortInputs=false)
Convenience constructor for a full spline.
Definition: Spline.hpp:232
Scalar slope_(size_t i) const
Returns the slope (i.e. first derivative) of the spline at the i-th sampling point.
Definition: Spline.hpp:1814
Scalar h11_prime2_(Scalar t) const
Definition: Spline.hpp:1652
Spline(Scalar x0, Scalar x1, Scalar y0, Scalar y1, Scalar m0, Scalar m1)
Convenience constructor for a full spline with just two sampling points.
Definition: Spline.hpp:125
size_t intersectSegment_(Evaluation *sol, size_t segIdx, const Evaluation &a, const Evaluation &b, const Evaluation &c, const Evaluation &d, Scalar x0=-1e100, Scalar x1=1e100) const
Find all the intersections of a segment of the spline with a cubic polynomial within a specified inte...
Definition: Spline.hpp:1744
void setXYArrays(size_t nSamples, const ScalarArrayX &x, const ScalarArrayY &y, SplineType splineType=Natural, bool sortInputs=false)
Set the sampling points natural spline using C-style arrays.
Definition: Spline.hpp:515
void setXYContainers(const ScalarContainerX &x, const ScalarContainerY &y, SplineType splineType=Natural, bool sortInputs=false)
Set the sampling points of a natural spline using STL-compatible containers.
Definition: Spline.hpp:556
size_t numSamples() const
Returns the number of sampling points.
Definition: Spline.hpp:257