C++ array as template parameter

Question

I want to write a simple polynomial class that can take an array of coefficients and expand it into a function a compile time so I don't need to loop over the coefficients at run time. I want to do something like this:

template <PARAM_TYPE, PARAMS>
class P {
public:
    PARAM_TYPE eval(PARAM_TYPE p){
        //Does PARAMS[0] * pow(p, PARAMS.length() -1) + ... + PARAMS[N-1] 
    }
}

Sample call

P<double,{2,4,3}> quadratic;
quadratic.eval(5); //returns 73

I don't want to be doing the loop since that will take time. Ideally I want to be able to form the expression above at compile time. Is this possible? Thanks

something like this maybe? http://stackoverflow.com/questions/7877293/non-type-variadic-function-templates-in-c11 not an exact answer, but gives the idea.. basically you get a parameter pack of integers. You can operate over them using the ... notation but you might need an http://en.cppreference.com/w/cpp/utility/integer_sequence as well. but I don't actually know how to do what you want... — xaxxon, Mar 17 '16 at 23:35
Re "I don't want to be doing the loop since that will take time.", there will be a loop in the underlying machine code no matter what. — Cheers and hth. - Alf, Mar 17 '16 at 23:49
Do you want `quadratic.eval(5)` to be calculated at runtime or compile time? — Kurt Stutsman, Mar 18 '16 at 01:05
I want the compiler to spit out the unrolled version of the loop at compile time and substitute in the cas at runtime — GXR, Mar 18 '16 at 03:35

score 1 · Accepted Answer · edited May 23 '17 at 12:23

Here is an example of doing what you want. The compiler is finicky about whether or not it optimizes away all the code into constants depending on the usage I noticed and which compiler you use.

test here

#include <type_traits>

template<class T, unsigned Exponent>
inline constexpr typename std::enable_if<Exponent == 0, T>::type
pow2(const T base)
{
    return 1;
}

template<class T, unsigned Exponent>
inline constexpr typename std::enable_if<Exponent % 2 != 0, T>::type
pow2(const T base)
{
      return base * pow2<T, (Exponent-1)/2>(base) * pow2<T, (Exponent-1)/2>(base);
}

template<class T, unsigned Exponent>
inline constexpr typename std::enable_if<Exponent != 0 && Exponent % 2 == 0, T>::type
pow2(const T base)
{
    return pow2<T, Exponent / 2>(base) * pow2<T, Exponent / 2>(base);
}

template<typename ParamType>
inline constexpr ParamType polynomial(const ParamType&, const ParamType& c0)
{
  return c0;
}

template<typename ParamType, typename Coeff0, typename ...Coeffs>
inline constexpr ParamType polynomial(const ParamType& x, const Coeff0& c0, const Coeffs& ...cs)
{
    return (static_cast<ParamType>(c0) * pow2<ParamType, sizeof...(cs)>(x)) + polynomial(x, static_cast<ParamType>(cs)...);
}

unsigned run(unsigned x)
{
   return polynomial(x, 2, 4, 3);
}

double run(double x)
{
  return polynomial(x, 2, 4, 3);
}

unsigned const_unsigned()
{
    static const unsigned value = polynomial(5, 2, 4, 3);
    return value;
}

double const_double()
{
    static const double value = polynomial(5, 2, 4, 3);
    return value;
}

EDIT: I have updated the code to a use a tweaked version of pow2<>() that aggressively performs calculations at compile time. This version optimizes so well at -O2 that it actually surprised me. You can see the generated assembly for the full program using the button above the code. If all arguments are constant, the compiler will generate the entire constant value at compile time. If the first argument is runtime-dependent, it generates very tight code for it still.

(Thanks to @dyp on this question for the inspiration to pow)

So if I have `polynomial(x, 0.5, 4, 3, 2, 7, 9)` the compiler will generate the unrolled version of `0.5*x^5 + ... + 9`? — GXR, Mar 18 '16 at 04:16
I have updated the code in the answer to more aggressively inline in the context of `x` not being a compile-time constant. This code produces surprisingly optimal code and removes all the function calls the previous version would generate when `x` was not known at compile time. — Kurt Stutsman, Mar 18 '16 at 14:59
For your question about `pow`, the standard versions will convert both arguments to `double` and use a floating point power function when it's not necessary. Also this version allows for better inlining. — Kurt Stutsman, Mar 18 '16 at 15:00

score 0 · Answer 2 · answered Mar 18 '16 at 16:57

To evaluate a polynom, a good algorithm is Horner (See https://en.wikipedia.org/wiki/Horner%27s_method). The main idea is to compute the polynom recursively. Let's a polynom P of order n with coefficient ai. It is easy to see that the sequence Pk = Pk-1*x0 + an-k with P0 = an, that P(x0) = Pn.

So let's implement this algorithm using constexpr function:

template<class T>
constexpr double horner(double x, T an) { return an; }

template<class... T, class U = T>
constexpr double horner(double x, U an, T... a) { return horner(x, a...) * x + an; }

std::cout << horner(5., 1, 2, 1) << std::endl;

//test if the invocation of the constexpr function takes the constant expression branch
std::cout << noexcept(horner(5., 1, 2, 1)) << std::endl;

As you see, it is really easy to implement the evaluation of a polynom with constexpr functions using the recursive formula.

C++ array as template parameter

2 Answers2