[SciPy-User] Chebeshev Polynomial implimentation python

Eric Carlson ecarlson at eng.ua.edu
Mon Dec 20 22:55:09 EST 2010


On 12/18/2010 12:35 PM, Pramod wrote:
> Matlab imple mentation :
> for N = 1:Nmax;
> [D,x] = cheb(N);
>
> How to impliment above (written in matlab ) chebshev  polynomial in
> python

cheb is not a standard matlab function, but if this is it:

function [D,x]=cheb(N)

if N==0, D=0; x=1; return; end
x=cos(pi*(0:N)/N)';
c=[2; ones(N-1,1); 2].*(-1).^(0:N)';
X=repmat(x,1,N+1);
dX=X-X';
D=(c*(1./c)')./(dX+eye(N+1)); % off diagonals
D=D-diag(sum(D')); % diagonals


Then a python version could be given by:

from numpy import cos,pi,linspace,array,matrix,ones,hstack,eye,diag,tile

def cheb(N):

     if N==0:
         D=0.0
         x=1.0
     else:
         x=cos(pi*linspace(0,N,N+1)/N)
         c=matrix(hstack(([2.],ones(N-1),[2.]))*(-1)**linspace(0,N,N+1)).T
         X=tile(x,[N+1,1])
         dX=(X-X.T).T
         D=array(c*(1/c).T)/(dX+eye(N+1)); # off diagonals
         D=D-diag(sum(D,axis=1)); # diagonals

     return D,x

(D,x)=cheb(3)
print D
(D,x)=cheb(4)
print D


I almost literally translated the matlab code, so I do not know how 
efficient it all is, and have given zero thought to better ways to do 
it. You need to be very careful with matrix and array types. Unless you 
KNOW you want a matrix, you probably want an array. One of the biggest 
transitions from matlab to python is learning to not worry about whether 
something is a column or row array, except when you are doing matrix 
multiplication.

HTH








More information about the SciPy-User mailing list