welcome to the third part of the multivariate polynomials in c++ using templates series. this part will explain how to compute partial derivatives using partial template specialization tricks.

the idea is very simple. given , written as with , we have

and

this allows us to implement partial derivatives iteratively. we assume that the number of the indeterminate by which one should differentiate is already known during compile time. the implementation will be done by a templated function

1template<int N, int n, class T> 2inline poly<n, T> derivative(const poly<n, T> & p);

the purpose of the template parameter `N`

is to say that should be computed. given a polynomial `p`

, the function can be called simply by `derivative<N>(p)`

. to implement `derivative<N>(p)`

, we would like to again use partial template specialization to distinguish between the cases `N == 0`

and `N > 0`

. unfortunately, partial template specialization is not available for template functions. hence, we need to introduce a helper template class:

1template<int N, int n, class T> 2class derivative_computer;

it should contain one static function, `static inline void computeDerivative(const poly<n, T> & src, poly<n, T> & dest);`

which allows to compute the partial derivative of the polynomial in `src`

and store it into `dest`

. then `derviative<N>(p)`

can simply be implemented by

1template<int N, int n, class T> 2inline poly<n, T, std::allocator<T> > derivative(const poly<n, T, std::allocator<T> > & p) 3{ 4 poly<n, T> res; 5 derivative_computer<N, n, T>::computeDerivative(p, res); 6 return res; 7}

now let us focus on the class template `derivative_computer<N, n, T>`

. the general implementation, for `N > 0`

, uses the second formula from above and iteratively calls `derivative_computer<N-1, n-1, T>::computeDerivative()`

:

1template<int N, int n, class T> 2class derivative_computer 3{ 4public: 5 static inline void computeDerivative(const poly<n, T> & src, poly<n, T> & dest) 6 { 7 dest.d_value.resize(src.degree() + 1); 8 for (int i = 0; i <= src.degree(); ++i) 9 derivative_computer<N - 1, n - 1, T>::computeDerivative(src[i], dest[i]); 10 dest.normalize(); 11 } 12};

now we are left with the case that `N = 0`

. the implementation here is also simple:

1template<int n, class T> 2class derivative_computer<0, n, T> 3{ 4public: 5 static inline void computeDerivative(const poly<n, T> & src, poly<n, T> & dest) 6 { 7 dest.d_value.resize(src.degree() >= 0 ? src.degree() : 0); 8 for (int i = 1; i <= src.degree(); ++i) 9 dest[i - 1] = src[i] * (T)i; 10 dest.normalize(); 11 } 12};

here, we cast `i`

to an object of type `T`

before multiplication. from this point on, no further iteration is needed.

we are left with the case that some (malicious or just confused) user calls `derivative<N>(p)`

for a polynomial `p`

of type `poly<n, T>`

with `n <= N`

. in that case, the compiler tries to instanciate `derivative_computer<k, 0, T>::computeDerivative()`

and tries to treat `poly<0, T>`

as a polynomial. this will result in strange error messages, and it would be nice if there is something more readable. for that reason, i added another two specializations:

1template<class T> 2class derivative_computer<0, 0, T> 3{ 4public: 5 class ERROR_N_must_be_strictly_less_than_n_in_derivative_template; 6 7 static inline void computeDerivative(const poly<0, T> &, poly<0, T> &) 8 { 9 ERROR_N_must_be_strictly_less_than_n_in_derivative_template(); 10 } 11}; 12 13template<int N, class T> 14class derivative_computer<N, 0, T> 15{ 16public: 17 class ERROR_N_must_be_strictly_less_than_n_in_derivative_template; 18 19 static inline void computeDerivative(const poly<0, T> &, poly<0, T> &) 20 { 21 ERROR_N_must_be_strictly_less_than_n_in_derivative_template(); 22 } 23};

these two definitions compile fine, but as soon as one of them is instanciated by the compiler, it will try to find out something about the class `ERROR_N_must_be_strictly_less_than_n_in_derivative_template`

and not succeed, and spit out some error message hopefully containing the string `ERROR_N_must_be_strictly_less_than_n_in_derivative_template`

. for example, if i replace the line computing the partial derivative of the intersection equation in `poly-test.cpp`

by `std::cout << derivative<2>(l) << "\n";`

, `g++`

will complain as follows:

1In file included from poly-test.cpp:2:0: 2polynomial.hpp: In static member function ‘static void derivative_computer<N, 0, T>::computeDerivative(const poly<0, T>&, poly<0, T>&) [with int N = 1, T = double]’: 3polynomial.hpp:1200:13: instantiated from ‘static void derivative_computer<N, n, T>::computeDerivative(const poly<n, T>&, poly<n, T>&) [with int N = 2, int n = 1, T = double]’ 4polynomial.hpp:1247:5: instantiated from ‘poly<n, T> derivative(const poly<n, T>&) [with int N = 2, int n = 1, T = double]’ 5poly-test.cpp:19:33: instantiated from here 6polynomial.hpp:1227:9: error: invalid use of incomplete type ‘struct derivative_computer<1, 0, double>::ERROR_N_must_be_strictly_less_than_n_in_derivative_template’ 7polynomial.hpp:1223:11: error: declaration of ‘struct derivative_computer<1, 0, double>::ERROR_N_must_be_strictly_less_than_n_in_derivative_template’ 8make: *** [poly-test.o] Error 1

this is not exactly readable, but the uppercase `ERROR`

quickly lets you look at the important part of that message. (note that i removed the allocators from the error messages, so that it is in sync with the above writing. if you try it out yourself with the source code, you'll see a slightly more complicated error message.)

this completes part three of the series. in the next part, i will discuss how to evaluate polynomials efficiently using templates.