On Nov 21, 2:52�pm, Jerry Avins <j...@ieee.org> wrote:
> c...@claysturner.com wrote:
> > On Nov 21, 11:42 am, "sasidhar" <sasidha...@gmail.com> wrote:
> >> Hi all,
>
> >> I am trying to calculate derivatives of chebyshev coefficients using FFTW
> >> for cosine(x). However I am not getting back the right values when I
> >> perform inverse chebyshev transform. Can anyone check the code below and
> >> tell me the possible errors. Thank you.
>
> >> � � � program test
> >> � � � implicit none
> >> #include "fftw3.f"
>
> >> � � � INTEGER, PARAMETER �:: �N=4
> >> � � � INTEGER(8) :: �plan
> >> � � � REAL(8), DIMENSION(N) :: �in, out, out1, in1
> >> � � � REAL(8), DIMENSION(N) :: dout
> >> � � � REAL :: �pi, y
> >> � � � INTEGER :: j, p
>
> >> � � � pi = 3.14
> >> � � � write(*,*) 'Input array:'
> >> � � � do j = 1,N,1
> >> � � � � y =-1*cos(pi*(j)/N)
> >> � � � � in(j) = cos(2*pi*y/N)
> >> � � � enddo
>
> >> � � � call dfftw_plan_r2r_1d ( plan, N, in, out,
> >> � � �^ � � � � � � � � � � � � FFTW_REDFT00, FFTW_ESTIMATE )
> >> � � � call dfftw_execute ( plan )
> >> � � � call dfftw_destroy_plan ( plan )
> >> � � � out = out/6
> >> � � � do j = 1,N,1
> >> � � � � write(*,*) ' � �out(',j,') = ',out(j)
> >> � � � enddo
>
> >> � � � dout(N) = 0
> >> � � � dout(N-1) = 2*(N)*out(N)
>
> >> � � � do j = N-2,1,-1
> >> � � � � dout(j) = dout(j+2) + 2*(j+1)*out(j+1)
> >> � � � enddo
>
> >> � � � call dfftw_plan_r2r_1d(plan, N+1, dout, out1,
> >> � � �^ � � � � � � � � � � � � FFTW_REDFT00, FFTW_ESTIMATE)
> >> � � � call dfftw_execute ( plan )
> >> � � � call dfftw_destroy_plan ( plan )
>
> >> � � � write(*,*) 'Output array after inverse FFT:'
> >> � � � do j = 1,N,1
> >> � � � � write(*,*) ' � �out1(',j,') = ',out1(j)
> >> � � � enddo
>
> >> � � � do j = 1,N,1
> >> � � � � y =-1*cos(pi*(j)/N)
> >> � � � � in1(j) = -1*sin(2*pi*y/N)
> >> � � � � write(*,*) ' � �in1(',j,') = ',in1(j)
> >> � � � enddo
> >> � � � end
>
> >> �Regards,
> >> �Sasidhar
>
> > Hello Sasidhar,
>
> > I don't know about your code, but derivatives of coefficients
> > (constants?) are zero. Or are you trying to find derivatives of
> > Chebyshev polynomials? These do have well defined forms in terms
> > Chebyshev polynomials (although both 1st and 2nd kind are needed
> > though). Or you can just stay with polynomials.
>
> > Maybe if you can give a clearer description of what you are really
> > trying to do, then we may be better able to help.
>
> > Are you trying to do differentiation via FFTs? You know aliasing can
> > be an issue here!
>
> Clay, What is a Chebychev transform?
>
> Jerry
> --
> Engineering is the art of making what you want from things you can get.
> �����������������������������������������������������������������������- Hide quoted text -
>
> - Show quoted text -
Hey Jerry,
Basically it is a form of expansion of a function in terms of
Chebyshev polys where they are sampled on a nonuniformly spaced set of
nodes. The nodes are found by cos(k*pi/N) where k goes from
0,1,2,3,.., N. Using a set from -1 to +1. Scaling and translation may
be used for other domains. The transform is also self inverting.
On the "chebyshev nodes", the discrete chebyshev functions (formed by
sampling the continuous ones) are orthogonal with respect to the
chebyshev weighting function.
Chebys have a bunch of neat properties such as a composition of chebys
is also a chebys! But by doing a cosine transformation on the domain,
the cheby transform can be found via a FFT. And since it is self
inverting, the OP wants to take a sample function and go around the
loop and see if all of his steps are working properly.
He needs to do a discrete cheby transform by the direct method so as
to compare that halfway step to what he gets from the FFT method. Then
he can try the inverse to see if he gets back to where he started as a
sanity check.
IHTH,
Clay
Clay