Hide
Easily highlight source code for your blog with our Syntax Highlighter. Join Siafoo Now or Learn More

Recursive template for the Legendre Polynomials Atom Feed 0

In Brief A template implementation of the Bonnet recursion formula for generating of Legendre polynomials. I wrote this as a sort of template meta-programming exercise but have verified that the first ~10 generated polynomials produce identical output to that of the legendre_p function from the boost Math Toolkit.... more
# 's
 1/*
2 * Legendre.hpp
3 *
4 * Created on: Oct 27, 2009
5 * Author: stou@icapsid.net
6 * License: Public Domain (send me improvements though)
7 */
8
9namespace Legendre
10{
11
12template<class T, unsigned int L> struct P
13{
14 T operator()(T x)
15 {
16 return ((2 * L - 1) / T(L)) * x * P<T, L - 1> ()(x) - ((L - 1) / T(L)) * P<T, L - 2> ()(x);
17 }
18};
19
20template<class T> struct P<T, 0>
21{
22 T operator()(T x)
23 {
24 return T(1);
25 }
26};
27
28template<class T> struct P<T, 1>
29{
30 T operator()(T x)
31 {
32 return x;
33 }
34};
35
36// Sum of the first N legendre polynomials
37template<class T, unsigned int N> struct P_sum
38{
39 T operator()(T x)
40 {
41 return P_sum<T, N - 1> ()(x) + P<T, N - 1> ()(x);
42 }
43};
44
45template<class T> struct P_sum<T, 0>
46{
47 T operator()(T x)
48 {
49 return T(0);
50 }
51};
52
53}

A template implementation of the Bonnet recursion formula for generating of Legendre polynomials. I wrote this as a sort of template meta-programming exercise but have verified that the first ~10 generated polynomials produce identical output to that of the legendre_p function from the boost Math Toolkit.

Use the template like this:

# 's
1P<float, 5>()(x)
2
3// OR
4P<float, 10> p();
5p_10(x)

It is also possible to sum the first N polynomials using:

# 's
1P_sum<float, 5>()(x)

You are much better off using the boost Math Toolkit when summing the polynomials because:

Warning

The polynomial sums are VERY VERY slow above > 10, much slower than simply looping over legendre_p from the boost Math Toolkit.

Also this isn't exactly the Bonnet recursion formula, it is a rearranged version: