Boiled it down a bit more to only include code that actually takes time.
First time around I found the other variant more instructive because it
shows the discrepancy between the DCT and the loop but might be
confusing. Thus here the bare minimum that correctly calculates the
coefficients of the first derivative from the coefficients of the
Chebyshev polynomials.
Cheers
Robert
On 25.07.2011 00:43, Robert Elsner wrote:
> Hey Everybody,
>
> I am approximating the derivative of nonperiodic functions on [-1,1]
> using Chebyshev polynomials. The implementation is straightforward and
> works well but is painfully slow. The routine wastes most of its time on
> a trivial operation (see comment in the code)
> Unfortunately the spectral coefficients are calculated recursively and
> thus I haven't found a way to speed up the code. I am aware of list
> comprehensions, the speed advantage of calling native numpy functions
> etc. But no luck yet finding an alternate formulation that speeds up the
> calculation. I attached a stripped down variant of the code. I am very
> grateful for tips or directions where to look for a solution (well I
> know writing a small C extension might be the way to go but Id love to
> have speedy Python)
>
> Note: The DCT takes almost no time compared to the statement in the loop.
>
> Cheers
> Robert
>
>
>
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> [email protected]
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
import scipy.fftpack
import numpy as np
import timeit
def chebyshev_derivative_gc():
""" Calculate the first derivative of a function
"""
N = 64
a = np.random.random( ( 6, N ) )
b = np.zeros( a.shape )
b[:, N - 2] = 2.0 * ( N - 1 ) * a[:, N - 1]
for j in xrange( N - 2, 0, -1 ):
# The next line is the offender - try by commenting out.
b[:, j - 1] = ( 2.0 * j ) * a[:, j] + b[:, j + 1]
pass
# Timing code
t = timeit.Timer( chebyshev_derivative_gc )
print t.repeat( 3, 1000 )
_______________________________________________
NumPy-Discussion mailing list
[email protected]
http://mail.scipy.org/mailman/listinfo/numpy-discussion
