diff --git a/ecl_geometry/include/ecl/geometry/cubic_spline.hpp b/ecl_geometry/include/ecl/geometry/cubic_spline.hpp index aac50db..e093957 100644 --- a/ecl_geometry/include/ecl/geometry/cubic_spline.hpp +++ b/ecl_geometry/include/ecl/geometry/cubic_spline.hpp @@ -177,6 +177,7 @@ class ECL_PUBLIC CubicSpline : public BluePrintFactory< CubicSpline > { * @exception : StandardException : throws if input x value is outside the spline range [debug mode only]. */ double operator()(const double &x) const; + std::map operator()(int last_index, const double &x); /** * @brief Spline derivative. * @@ -186,6 +187,7 @@ class ECL_PUBLIC CubicSpline : public BluePrintFactory< CubicSpline > { * @exception : StandardException : throws if input x value is outside the spline range [debug mode only]. */ double derivative(double x) const; + std::map derivative(const int &last_index, const double &x) const; /** * @brief Spline second derivative. * @@ -195,7 +197,7 @@ class ECL_PUBLIC CubicSpline : public BluePrintFactory< CubicSpline > { * @exception : StandardException : throws if input x value is outside the spline range [debug mode only]. */ double dderivative(double x) const; - + std::map dderivative(const int &last_index, const double &x) const; /** * @brief The discretised domain for this spline. * diff --git a/ecl_geometry/include/ecl/geometry/smooth_linear_spline.hpp b/ecl_geometry/include/ecl/geometry/smooth_linear_spline.hpp index 9804cc1..bb4473d 100644 --- a/ecl_geometry/include/ecl/geometry/smooth_linear_spline.hpp +++ b/ecl_geometry/include/ecl/geometry/smooth_linear_spline.hpp @@ -88,6 +88,8 @@ class ecl_geometry_PUBLIC SmoothLinearSpline { * @exception : StandardException : throws if x is outside the spline range [debug mode only]. */ double operator()(const double &x) const; + std::map operator()(const int &last_index,const double &x) const; + /** * @brief Spline derivative. * @@ -99,6 +101,8 @@ class ecl_geometry_PUBLIC SmoothLinearSpline { * @exception : StandardException : throws if x is outside the spline range [debug mode only]. */ double derivative(const double &x) const; + std::map derivative(const int &last_index, const double &x) const; + /** * @brief Spline second derivative. * @@ -108,6 +112,7 @@ class ecl_geometry_PUBLIC SmoothLinearSpline { * @exception : StandardException : throws if x is outside the spline range [debug mode only]. */ double dderivative(const double &x) const; + std::map dderivative(const int &last_index, const double &x) const; /** * @brief The discretised domain for this spline. @@ -150,7 +155,7 @@ class ecl_geometry_PUBLIC SmoothLinearSpline { * @brief Streaming output insertion operator for smoothed linear splines. * * Streaming output insertion operator for smoothed linear splines. This - * simply lists the spline segments and corners (linear functions and + * simply lists the spline segments and corners (linear functions andTensionSpline * quintic polynomials) in algebraic form. * * @tparam OutputStream : the type of stream being used. diff --git a/ecl_geometry/include/ecl/geometry/tension_spline.hpp b/ecl_geometry/include/ecl/geometry/tension_spline.hpp index 325565d..09064b7 100644 --- a/ecl_geometry/include/ecl/geometry/tension_spline.hpp +++ b/ecl_geometry/include/ecl/geometry/tension_spline.hpp @@ -25,6 +25,7 @@ #include #include #include "macros.hpp" +#include /***************************************************************************** ** Namespaces @@ -156,6 +157,8 @@ class ecl_geometry_PUBLIC TensionSpline : public BluePrintFactory< TensionSpline * @exception : StandardException : throws if x is outside the spline range [debug mode only]. */ double operator()(const double &x) const; + + std::map operator()(const int &last_index, const double &x) const; /** * @brief Spline derivative. * @@ -167,6 +170,8 @@ class ecl_geometry_PUBLIC TensionSpline : public BluePrintFactory< TensionSpline * @exception : StandardException : throws if x is outside the spline range [debug mode only]. */ double derivative(const double &x) const; + + std::map derivative(const int &last_index, const double &x) const; /** * @brief Spline second derivative. * @@ -177,6 +182,8 @@ class ecl_geometry_PUBLIC TensionSpline : public BluePrintFactory< TensionSpline */ double dderivative(const double &x) const; + std::map dderivative(const int &last_index, const double &x) const; + /** * @brief The discretised domain for this spline. * diff --git a/ecl_geometry/src/lib/cubic_spline.cpp b/ecl_geometry/src/lib/cubic_spline.cpp index 33e66b8..0f76a96 100644 --- a/ecl_geometry/src/lib/cubic_spline.cpp +++ b/ecl_geometry/src/lib/cubic_spline.cpp @@ -35,6 +35,16 @@ double CubicSpline::operator()(const double &x) const { return cubic_polynomials[index](x); } +std::map CubicSpline::operator()(int last_index, const double &x) { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, cubic_polynomials[index](x)}}; + return m; +} + double CubicSpline::derivative(double x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -44,6 +54,16 @@ double CubicSpline::derivative(double x) const { return cubic_polynomials[index].derivative()(x); } +std::map CubicSpline::derivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, cubic_polynomials[index].derivative()(x)}}; + return m; +} + double CubicSpline::dderivative(double x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -53,4 +73,15 @@ double CubicSpline::dderivative(double x) const { return cubic_polynomials[index].derivative().derivative()(x); } +std::map CubicSpline::dderivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, cubic_polynomials[index].derivative().derivative()(x)}}; + return m; +} + + } // namespace ecl diff --git a/ecl_geometry/src/lib/smooth_linear_spline.cpp b/ecl_geometry/src/lib/smooth_linear_spline.cpp index fec1839..9bdc857 100644 --- a/ecl_geometry/src/lib/smooth_linear_spline.cpp +++ b/ecl_geometry/src/lib/smooth_linear_spline.cpp @@ -122,6 +122,21 @@ double SmoothLinearSpline::operator()(const double &x) const { } } +std::map SmoothLinearSpline::operator()(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + if ( index % 2 == 0 ) { // linear + std::map m = {{index,segments[index/2](x)}}; + return m; + } else { // quintic + std::map m = {{index, corners[(index-1)/2](x)}}; + return m; + } +} + double SmoothLinearSpline::derivative(const double &x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -135,6 +150,23 @@ double SmoothLinearSpline::derivative(const double &x) const { } } +std::map SmoothLinearSpline::derivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m; + if ( index % 2 == 0 ) + { // linear + std::map m = {{index,segments[index/2].derivative(x)}}; + return m; + } else { // quintic + std::map m = {{index,corners[(index-1)/2].derivative(x)}}; + return m; + } +} + double SmoothLinearSpline::dderivative(const double &x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -148,6 +180,23 @@ double SmoothLinearSpline::dderivative(const double &x) const { } } +std::map SmoothLinearSpline::dderivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m; + if ( index % 2 == 0 ) + { // linear + std::map m = {{index,segments[index/2].dderivative(x)}}; + return m; + } else { // quintic + std::map m = {{index,corners[(index-1)/2].dderivative(x)}}; + return m; + } +} + } // namespace ecl diff --git a/ecl_geometry/src/lib/tension_spline.cpp b/ecl_geometry/src/lib/tension_spline.cpp index 388872f..6d520d3 100644 --- a/ecl_geometry/src/lib/tension_spline.cpp +++ b/ecl_geometry/src/lib/tension_spline.cpp @@ -35,6 +35,16 @@ double TensionSpline::operator()(const double &x) const { return functions[index](tension,x); } +std::map TensionSpline::operator()(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, functions[index](tension,x)}}; + return m; +} + double TensionSpline::derivative(const double &x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -44,6 +54,16 @@ double TensionSpline::derivative(const double &x) const { return functions[index].derivative(tension,x); } + std::map TensionSpline::derivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, functions[index].derivative(tension,x)}}; + return m; + } + double TensionSpline::dderivative(const double &x) const { ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); int index = 0; @@ -53,4 +73,14 @@ double TensionSpline::dderivative(const double &x) const { return functions[index].dderivative(tension,x); } +std::map TensionSpline::dderivative(const int &last_index, const double &x) const { + ecl_assert_throw( ( ( x >= discretised_domain.front() ) && ( x <= discretised_domain.back() ) ), StandardException(LOC,OutOfRangeError) ); + int index = last_index; + while ( x > discretised_domain[index+1] ) { + ++index; + } + std::map m = {{index, functions[index].dderivative(tension, x)}}; + return m; +} + } // namespace ecl