C++ constexpr for polynomial evaluation with Horner's method

Question

I have want to be able go evaluate the derivative of a polynomial using Horner's method and use the result as a constexpr. This seems incredibly mundane, but I am missing something obvious because the compiler says that I am exceeding maximum recursion depth. The core recursion happens here:

template<size_t d, size_t i, typename C, typename X>
constexpr X evalImpl(const C &c, const X &x) {
    return i >= (C::SizeAtCompileTime - 1 - d) ? 1 : evalImpl<d, i + 1, C, X>(c, x);
}

You should not need to even know Horner's method to know what it is going on here, so I have stripped the code as much as possible, even removing how x is used, because it does not seem relevant to the problem I am having.

The idea is that when an index i equals the degree of the polynomial Degree<C>::value minus the order of the derivative d, then the recursion should stop. Otherwise, it should increment the index i and try again.

I am calling the above recursion with a call of the form

 eval<derivative, 0>(c, x)

where c is an Eigen Matrix of type Eigen::Matrix<double,1,7>, and x is a double. The idea is to start at 0 and count up to essentially the degree of the polynomial.

The compiler error message is of the form

In file included from /mnt/c/proj/src/main.cpp:11:0:
/mnt/c/proj/src/polynomial.h: In instantiation of 'constexpr X {anonymous}::evalImpl(const C&, const X&) [with long unsigned int d = 1ul; long unsigned int i = 898ul; C = Eigen::Matrix<double, 1, 7>; X = double]':
/mnt/c/proj/src/polynomial.h:74:108:   recursively required from 'constexpr X {anonymous}::evalImpl(const C&, const X&) [with long unsigned int d = 1ul; long unsigned int i = 1ul; C = Eigen::Matrix<double, 1, 7>; X = double]'
/mnt/c/proj/src/polynomial.h:74:108:   required from 'constexpr X {anonymous}::evalImpl(const C&, const X&) [with long unsigned int d = 1ul; long unsigned int i = 0ul; C = Eigen::Matrix<double, 1, 7>; X = double]'
/mnt/c/proj/src/polynomial.h:109:39:   required from 'constexpr X Polynomial::eval(const C&, const X&) [with long unsigned int d = 1ul; C = Eigen::Matrix<double, 1, 7>; X = double]'
/mnt/c/proj/src/main.cpp:306:66:   required from here
/mnt/c/proj/src/polynomial.h:74:108: fatal error: template instantiation depth exceeds maximum of 900 (use -ftemplate-depth= to increase the maximum)
         return i >= (C::SizeAtCompileTime - 1 - d) ? 1 : evalImpl<d, i + 1, C, X>(c, x);

Answer 1

The condition here:

return i >= (C::SizeAtCompileTime - 1 - d) ? 1 : evalImpl<d, i + 1, C, X>(c, x);

is not a if constexpr. Thus no matter whether i >= (C::SizeAtCompileTime - 1 - d) is true or false, the remaining one will always be instantiated. Thus recursion won't stop as you want.

Change to this:

if constexpr (i >= (C::SizeAtCompileTime - 1 - d)) {
    return 1;
} else { 
    return evalImpl<d, i + 1, C, X>(c, x);  
}

EDIT:

If you don't have access to C++17, use tag dispatch:

template<size_t d, size_t i, typename C, typename X>
constexpr X evalImpl_impl(const C &c, const X &x, std::true_type) {
    return 1;
}

template<size_t d, size_t i, typename C, typename X>
constexpr X evalImpl_impl(const C &c, const X &x, std::false_type) {
    return evalImpl_impl<d, i + 1, C, X>(c, x, std::integral_constant<bool, C::SizeAtCompileTime-1-d <= i+1>{});
}

template<size_t d, size_t i, typename C, typename X>
constexpr X evalImpl(const C &c, const X &x) {
    return evalImpl_impl(c, x, std::integral_constant<bool, C::SizeAtCompileTime-1-d <= i>{});
}