OilfieldToolsNet/opm-common
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

synchronize the evaluation code of Spline and Tabulated1DFunction

Browse Source 
in particular out-of-range conditions are now no longer done using
assertations but always trigger an Opm::NumericalProblem to be thrown.
This commit is contained in:
Andreas Lauser 2017-03-09 17:04:51 +01:00
parent d23801b098
commit 0ab39c5f7a
2 changed files with 63 additions and 28 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

46
opm/material/common/Spline.hpp
View File

@ -30,6 +30,7 @@
#include <opm/material/common/TridiagonalMatrix.hpp>
#include <opm/material/common/PolynomialUtils.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/Exceptions.hpp>
#include <opm/common/Unused.hpp>
#include <ostream>
@ -798,7 +799,11 @@ public:
template <class Evaluation>
Evaluation eval(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
typedef Opm::MathToolbox<Evaluation> Toolbox;
if (!extrapolate && !applies(x))
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate a spline outside of its range");
// handle extrapolation
if (extrapolate) {
@ -815,7 +820,7 @@ public:
}
}
return eval_(x, segmentIdx_(x));
return eval_(x, segmentIdx_(Toolbox::scalarValue(x)));
}
/*!
@ -833,7 +838,13 @@ public:
template <class Evaluation>
Evaluation evalDerivative(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
typedef Opm::MathToolbox<Evaluation> Toolbox;
if (!extrapolate && !applies(x))
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate the derivative of a spline outside of its range");
// handle extrapolation
if (extrapolate) {
if (x < xAt(0))
return evalDerivative_(xAt(0), /*segmentIdx=*/0);
@ -841,7 +852,7 @@ public:
return evalDerivative_(xAt(numSamples() - 1), /*segmentIdx=*/numSamples() - 2);
}
return evalDerivative_(x, segmentIdx_(x));
return evalDerivative_(x, segmentIdx_(Toolbox::scalarValue(x)));
}
/*!
@ -859,11 +870,15 @@ public:
template <class Evaluation>
Evaluation evalSecondDerivative(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
if (extrapolate)
typedef Opm::MathToolbox<Evaluation> Toolbox;
if (!extrapolate && !applies(x))
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate the second derivative of a spline outside of its range");
else if (extrapolate)
return 0.0;
return evalDerivative2_(x, segmentIdx_(x));
return evalDerivative2_(x, segmentIdx_(Toolbox::scalarValue(x)));
}
/*!
@ -881,11 +896,15 @@ public:
template <class Evaluation>
Evaluation evalThirdDerivative(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
if (extrapolate)
typedef Opm::MathToolbox<Evaluation> Toolbox;
if (!extrapolate && !applies(x))
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate the third derivative of a spline outside of its range");
else if (extrapolate)
return 0.0;
return evalDerivative3_(x, segmentIdx_(x));
return evalDerivative3_(x, segmentIdx_(Toolbox::scalarValue(x)));
}
/*!
@ -1730,13 +1749,8 @@ protected:
}
// find the segment index for a given x coordinate
template <class Evaluation>
size_t segmentIdx_(const Evaluation& xEval) const
size_t segmentIdx_(Scalar x) const
{
typedef Opm::MathToolbox<Evaluation> Toolbox;
Scalar x = Toolbox::scalarValue(xEval);
// bisection
size_t iLow = 0;
size_t iHigh = numSamples() - 1;

45
opm/material/common/Tabulated1DFunction.hpp
View File

@ -27,10 +27,10 @@
#ifndef OPM_TABULATED_1D_FUNCTION_HPP
#define OPM_TABULATED_1D_FUNCTION_HPP
#include <opm/material/densead/Math.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/Exceptions.hpp>
#include <opm/common/Unused.hpp>
#include <opm/material/densead/Math.hpp>
#include <algorithm>
#include <cassert>
@ -235,7 +235,8 @@ public:
/*!
* \brief Return true iff the given x is in range [x1, xn].
*/
bool applies(Scalar x) const
template <class Evaluation>
bool applies(const Evaluation& x) const
{ return xValues_[0] <= x && x <= xValues_[numSamples() - 1]; }
/*!
@ -252,15 +253,18 @@ public:
{
typedef Opm::MathToolbox<Evaluation> Toolbox;
if (!extrapolate && !applies(x))
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate a tabulated function outside of its range");
size_t segIdx;
if (extrapolate && x < xValues_.front())
segIdx = 0;
else if (extrapolate && x > xValues_.back())
segIdx = numSamples() - 2;
else {
assert(xValues_.front() <= x && x <= xValues_.back());
else
segIdx = findSegmentIndex_(Toolbox::scalarValue(x));
}
Scalar x0 = xValues_[segIdx];
Scalar x1 = xValues_[segIdx + 1];
@ -283,11 +287,18 @@ public:
* cause a failed assertation.
*/
template <class Evaluation>
Evaluation evalDerivative(const Evaluation& x, bool extrapolate OPM_UNUSED = false) const
Evaluation evalDerivative(const Evaluation& x, bool extrapolate = false) const
{
unsigned segIdx = findSegmentIndex_(x);
typedef Opm::MathToolbox<Evaluation> Toolbox;
return evalDerivative(x, segIdx);
if (!extrapolate && !applies(x)) {
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate a derivative of a tabulated"
" function outside of its range");
}
unsigned segIdx = findSegmentIndex_(Toolbox::scalarValue(x));
return evalDerivative_(x, segIdx);
}
/*!
@ -305,9 +316,14 @@ public:
* cause a failed assertation.
*/
template <class Evaluation>
Evaluation evalSecondDerivative(const Evaluation OPM_OPTIM_UNUSED& x, bool extrapolate OPM_OPTIM_UNUSED = false) const
Evaluation evalSecondDerivative(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
if (!extrapolate && !applies(x)) {
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate a second derivative of a tabulated "
" function outside of its range");
}
return 0.0;
}
@ -326,9 +342,14 @@ public:
* cause a failed assertation.
*/
template <class Evaluation>
Evaluation evalThirdDerivative(const Evaluation OPM_OPTIM_UNUSED& x, bool extrapolate OPM_OPTIM_UNUSED = false) const
Evaluation evalThirdDerivative(const Evaluation& x, bool extrapolate = false) const
{
assert(extrapolate || applies(x));
if (!extrapolate && !applies(x)) {
OPM_THROW(Opm::NumericalProblem,
"Tried to evaluate a third derivative of a tabulated "
" function outside of its range");
}
return 0.0;
}