mattdeboard
mattdeboard

Reputation: 700

Iterative polynomial multiplication -- Chebyshev polynomials in Python

My question is: What is the best approach to iterative polynomial multiplication in Python?

I thought an interesting project would be to write a function in Python to generate the coefficients and exponents of each term for a Chebyshev polynomial of a given degree. The recursive function to generate such a polynomial (represented by Tn(x)) is:

With:

T0(x) = 1

and

T1(x) = x:

Tn(x) = 2xTn-1(x) - Tn-2(x)

What I have so far isn't very useful, but I am having trouble kind of wrapping my brain around how to get this going. What I want to happen is the following:

>> chebyshev(4)
[[8,4], [8,2], [1,0]]

This list represents the Chebyshev polynomial of the 4th degree: T4(x) = 8x4 - 8x2 + 1

import sys
def chebyshev(n, a=[1,0], b=[1,1]):
    z = [2,1]
    result = []
    if n == 0:
        return a
    if n == 1:
        return b
    print >> sys.stderr, ([z[0]*b[0], 
                           z[1]+b[1]],
                          a) # This displays the proper result for n = 2
    return result

The one solution I found on the web didn't work, so I am hoping someone can shed some light.

p.s. More information on Chebyshev polynomials: CSU Fullteron, Wikipedia - Chebyshev polynomials. They are very cool/useful, and tie together some really interesting trig functions/properties; worth a read.

Upvotes: 0

Views: 2662

Answers (3)

Nico Schlömer
Nico Schlömer

Reputation: 58841

orthopy (a project of mine) also supports computation of Chebyshev polynomials. With

import orthopy

# from sympy.abc import x
x = 0.5

normalization = "normal"   # or "classical", "monic"
evaluator = orthopy.c1.chebyshev1.Eval(x, normalization)
for _ in range(10):
    print(next(evaluator))
0.5641895835477564
0.39894228040143276
-0.39894228040143265
...

you get the values of the polynomials with increasing degree at x = 0.5. You can use a list/vector of multiple values, or even sympy symbolics.

Computation happens with recurrence relations of course. If you're interested in the coefficients, check out

rc = orthopy.c1.chebyshev1.RecurrenceCoefficients("monic", symbolic=True)

Upvotes: 0

bobobobo
bobobobo

Reputation: 67276

The best implementation for Chebyshev is :

// Computes T_n(x), with -1 <= x <= 1
real T( int n, real x )
{
  return cos( n*acos(x) ) ;
}

If you test this against other implementations, including explicit polynomial evaluation and iteratively computing the recurrence relation, this is actually just as fast. Try it yourself..

Generally:

  • Explicit polynomial evaluation is the worst (for large n)
  • Recursive evaluation is a little better
  • cosine evaluation is the best

Upvotes: 0

Joshua Smith
Joshua Smith

Reputation: 6631

SciPy has an implementation for Chebyshev

http://www.scipy.org/doc/api_docs/SciPy.special.orthogonal.html

I would suggest looking at their code.

Upvotes: 2

Related Questions