c c Numerical Analysis: c The Mathematics of Scientific Computing c D.R. Kincaid & E.W. Cheney c Brooks/Cole Publ., 1990 c c Section 6.12 c c Example of Fast Fourier Transform method c c c file: fft.f c parameter (m=4,n3=15) complex z(0:n3),c(0:n3),d(0:n3),w,f,u,v c print * print *,' Fast Fourier Transform example' print *,' Section 6.12, Kincaid-Cheney' print * c nn=2**m m2=nn/2 pi=4.0*atan(1.0) y=2.0*pi/float(nn) w=cmplx(0.0 , -y) w=cexp(w) z(0)=cmplx(1.0,0.0) c(0)=f(0.0) do 2 k=1,nn-1 z(k)=w*z(k-1) c(k)=f(y*float(k)) 2 continue print *,' Initial c-vector' print 5,c print * do 4 n=0,m-1 print *,' n = ',n print 5,c print * mn=2**(m-n-1) n2=2**(n-1) do 3 k=0,mn-1 do 3 j=0,n2 u=c(2*n2*k+j) v=z(j*mn)*c(2*n2*k+j+m2) d(4*n2*k+j)=(u+v)*cmplx(0.5,0.0) d(4*n2*k+j+2*n2)=(u-v)*cmplx(0.5,0.0) 3 continue do 4 j=0,nn-1 c(j)=d(j) 4 continue print *,' Final c-vector' print 5,c c 5 format (2x,16( '(',e13.6,',',e13.6,')',4x)) stop end c function f(x) complex f,z,d u=cos(x) v=sin(x) z=cmplx(u,v) d=cmplx(16.0,0.0) do 2 k=1,15 y=float(k) d=d*z + cmplx(16.0-y,0.0) 2 continue f=d c return end