diff --git a/MathLib/ODE/CVodeSolver.cpp b/MathLib/ODE/CVodeSolver.cpp index 37eee55cb2d5482d9fb2022490fb3829ecb572c0..0261d3b1f22cfdef4122d0fbd24d5907753529bf 100644 --- a/MathLib/ODE/CVodeSolver.cpp +++ b/MathLib/ODE/CVodeSolver.cpp @@ -90,7 +90,7 @@ public: ~CVodeSolverImpl(); private: - void setFunction(std::unique_ptr<FunctionHandles>&& f); + void setFunction(std::unique_ptr<detail::FunctionHandles>&& f); void preSolve(); void solve(const double t_end); @@ -113,7 +113,7 @@ private: unsigned _num_equations; void* _cvode_mem; - std::unique_ptr<FunctionHandles> _f; + std::unique_ptr<detail::FunctionHandles> _f; int _linear_multistep_method = CV_ADAMS; int _nonlinear_solver_iteration = CV_FUNCTIONAL; @@ -178,8 +178,9 @@ CVodeSolverImpl::CVodeSolverImpl(const BaseLib::ConfigTree& config, auto f_wrapped = [](const realtype t, const N_Vector y, N_Vector ydot, void* function_handles) -> int { - bool successful = static_cast<FunctionHandles*>(function_handles) - ->call(t, NV_DATA_S(y), NV_DATA_S(ydot)); + bool successful = + static_cast<detail::FunctionHandles*>(function_handles) + ->call(t, NV_DATA_S(y), NV_DATA_S(ydot)); return successful ? 0 : 1; }; @@ -206,7 +207,7 @@ void CVodeSolverImpl::setTolerance(const double abstol, const double reltol) _reltol = reltol; } -void CVodeSolverImpl::setFunction(std::unique_ptr<FunctionHandles>&& f) +void CVodeSolverImpl::setFunction(std::unique_ptr<detail::FunctionHandles>&& f) { _f = std::move(f); assert(_num_equations == _f->getNumEquations()); @@ -251,7 +252,7 @@ void CVodeSolverImpl::preSolve() // Caution: by calling the DENSE_COL() macro we assume that matrices // are stored contiguously in memory! bool successful = - static_cast<FunctionHandles*>(function_handles) + static_cast<detail::FunctionHandles*>(function_handles) ->callJacobian(t, NV_DATA_S(y), NV_DATA_S(ydot), DENSE_COL(jac, 0)); return successful ? 0 : 1; @@ -309,7 +310,7 @@ void CVodeSolver::setTolerance(const double abstol, const double reltol) _impl->setTolerance(abstol, reltol); } -void CVodeSolver::setFunction(std::unique_ptr<FunctionHandles>&& f) +void CVodeSolver::setFunction(std::unique_ptr<detail::FunctionHandles>&& f) { _impl->setFunction(std::move(f)); } diff --git a/MathLib/ODE/CVodeSolver.h b/MathLib/ODE/CVodeSolver.h index 196492b2fedf731fd635741c8b4d577b87453b9d..40fc1e6456beb3d73b75ca470b48f1386244b305 100644 --- a/MathLib/ODE/CVodeSolver.h +++ b/MathLib/ODE/CVodeSolver.h @@ -13,6 +13,7 @@ #include "BaseLib/ConfigTree.h" #include "OdeSolverTypes.h" +#include "FunctionHandles.h" namespace MathLib { @@ -33,7 +34,7 @@ protected: void setTolerance(double const* const abstol, const double reltol); void setTolerance(const double abstol, const double reltol); - void setFunction(std::unique_ptr<FunctionHandles>&& f); + void setFunction(std::unique_ptr<detail::FunctionHandles>&& f); void setIC(const double t0, double const* const y0); diff --git a/MathLib/ODE/ConcreteOdeSolver.h b/MathLib/ODE/ConcreteOdeSolver.h index 92e0d0d34543d55761d8f4720f6bf5170f49c828..c6ccce9c3c5c673319a9a47fd7612abe9aafd18b 100644 --- a/MathLib/ODE/ConcreteOdeSolver.h +++ b/MathLib/ODE/ConcreteOdeSolver.h @@ -15,7 +15,7 @@ #include "BaseLib/ConfigTree.h" #include "OdeSolver.h" -#include "Handles.h" +#include "FunctionHandles.h" #ifdef CVODE_FOUND #include "CVodeSolver.h" @@ -53,9 +53,8 @@ public: void setFunction(Function f, JacobianFunction df) override { Implementation::setFunction( - std::unique_ptr< // TODO unique_ptr not needed - detail::Handles<NumEquations>>{ - new detail::Handles<NumEquations>{f, df}}); + std::unique_ptr<detail::FunctionHandlesImpl<NumEquations>>{ + new detail::FunctionHandlesImpl<NumEquations>{f, df}}); } void setTolerance(const std::array<double, NumEquations>& abstol, diff --git a/MathLib/ODE/Handles.h b/MathLib/ODE/FunctionHandles.h similarity index 68% rename from MathLib/ODE/Handles.h rename to MathLib/ODE/FunctionHandles.h index 3036c803834bf910f79fc2fd30649e75382314e3..c268f21b68f043a66b3245cbf50fd16e1aa83fb8 100644 --- a/MathLib/ODE/Handles.h +++ b/MathLib/ODE/FunctionHandles.h @@ -16,14 +16,31 @@ namespace MathLib { namespace detail { +class FunctionHandles +{ +public: + virtual bool call(const double t, double const* const y, + double* const ydot) = 0; + virtual bool callJacobian(const double t, + double const* const y, + double* const ydot, + double* const jac) = 0; + + virtual bool hasJacobian() const = 0; + + virtual unsigned getNumEquations() const = 0; + + virtual ~FunctionHandles() = default; +}; + /// Function handles for N equations. template <unsigned N> -struct Handles : public MathLib::FunctionHandles +struct FunctionHandlesImpl : FunctionHandles { using Function = MathLib::Function<N>; using JacobianFunction = MathLib::JacobianFunction<N>; - Handles(Function& f, JacobianFunction& df) : f(f), df(df) {} + FunctionHandlesImpl(Function& f, JacobianFunction& df) : f(f), df(df) {} bool call(const double t, const double* const y, double* const ydot) override { diff --git a/MathLib/ODE/OdeSolverTypes.h b/MathLib/ODE/OdeSolverTypes.h index ba7978d6a6c76164e5a2d7a5015cc5ca8321062b..5d8affc7876e402e8eba01aedb6901521ab865db 100644 --- a/MathLib/ODE/OdeSolverTypes.h +++ b/MathLib/ODE/OdeSolverTypes.h @@ -38,23 +38,6 @@ using JacobianFunction = std::function<bool(const double t, MappedConstVector<N> ydot, MappedMatrix<N, N> jac)>; -// This is an internal detail -class FunctionHandles -{ -public: - virtual bool call(const double t, double const* const y, - double* const ydot) = 0; - virtual bool callJacobian(const double t, - double const* const y, - double* const ydot, - double* const jac) = 0; - - virtual bool hasJacobian() const = 0; - - virtual unsigned getNumEquations() const = 0; - - virtual ~FunctionHandles() = default; -}; -} +} // namespace MathLib #endif // MATHLIB_ODE_ODESOLVERTYPES_H