Skip to content

Commit

Permalink
Initial commit
Browse files Browse the repository at this point in the history
  • Loading branch information
@certik
certik committed Feb 21, 2013
0 parents commit 030cbf4
Show file tree
Hide file tree
Showing 55 changed files with 3,617 additions and 0 deletions.
96 changes: 96 additions & 0 deletions Makefile
Original file line number Diff line number Diff line change
@@ -0,0 +1,96 @@
# Several changes made by HCP so this would build without trouble
# on a Linux/g77 system.
# (1) changed step to build library to use ar instead of update
# (update must mean something different on someone elses Unix.)
# (2) Added make clean step
# (3) In test step, changed a.out to ./a.out for cautious folk who don't
# have "." in their PATH.
# (4) Change FFLAGS from -O to -O2 -funroll-loops
# (5) Specify FC=gcc in case /usr/bin/f77 is not a link to g77
# (as it won't be if you have f77reorder installed)
# (6) Added targets shared and installshared to make and install a shared
# version of the library. You need /usr/local/lib in /etc/ld.so.conf
# for this to work
# (7) Modified names for dble prec version
LIB=dfftpack

# Use these lines for Linux/g77
FC=g77
FFLAGS=-O2 -funroll-loops -fexpensive-optimizations

# Use these lines for Solaris
#FC=f77
#FFLAGS=-fast -O5

OBJ=\
zfftb.o\
cfftb1.o\
zfftf.o\
cfftf1.o\
zffti.o\
cffti1.o\
dcosqb.o\
cosqb1.o\
dcosqf.o\
cosqf1.o\
dcosqi.o\
dcost.o\
dcosti.o\
ezfft1.o\
dzfftb.o\
dzfftf.o\
dzffti.o\
passb.o\
passb2.o\
passb3.o\
passb4.o\
passb5.o\
passf.o\
passf2.o\
passf3.o\
passf4.o\
passf5.o\
radb2.o\
radb3.o\
radb4.o\
radb5.o\
radbg.o\
radf2.o\
radf3.o\
radf4.o\
radf5.o\
radfg.o\
dfftb.o\
rfftb1.o\
dfftf.o\
rfftf1.o\
dffti.o\
rffti1.o\
dsinqb.o\
dsinqf.o\
dsinqi.o\
dsint.o\
sint1.o\
dsinti.o

lib$(LIB).a: $(OBJ)
ar -rcs lib$(LIB).a $(OBJ)

shared:$(OBJ)
$(FC) -shared -o lib$(LIB).so $(OBJ)

install: lib$(LIB).a
mv lib$(LIB).a /usr/local/lib
rm *.o

installshared:lib$(LIB).so
mv lib$(LIB).so /usr/local/lib
rm *.o
ldconfig

test: test.o
$(FC) test.o -L./ -l$(LIB)
time ./a.out

clean:
rm -f -r *.o *.a *.so
31 changes: 31 additions & 0 deletions README
Original file line number Diff line number Diff line change
@@ -0,0 +1,31 @@
DFFTPACK V1.0
*****************************************************************
A Double precision clone by Hugh C. Pumphrey of:
FFTPACK
version 4 april 1985

The gzipped tar file dp.tgz contains a complete copy of the FORTRAN
sources of fftpack, with everything converted to double precision. If
you do

gunzip dp.tgz
tar xvf dfftpack.tar

You will get a directory called dfftpack with all the source code in
it. There is also:

(*) a Makefile which I have tweaked to work on modern Linux and Solaris
systems. The comments in this file document the changes made.

(*) a file doc which was supplied with fftpack and which has been
altered to reflect the changes made in the change to double precision.

(*) A file doc.double which details the changes I made to the source code

Please send any comments or bug reports to [email protected] . Please
also report if you get dfftpack to build successfully on any system
other than Linux or Solaris.

The original FFTPACK was public domain, so dfftpack is public domain
too. It is released in the hope it will be useful to someone. There is
no warranty of any sort covering this software.
62 changes: 62 additions & 0 deletions cfftb1.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,62 @@
SUBROUTINE CFFTB1 (N,C,CH,WA,IFAC)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION CH(*) ,C(*) ,WA(*) ,IFAC(*)
NF = IFAC(2)
NA = 0
L1 = 1
IW = 1
DO 116 K1=1,NF
IP = IFAC(K1+2)
L2 = IP*L1
IDO = N/L2
IDOT = IDO+IDO
IDL1 = IDOT*L1
IF (IP .NE. 4) GO TO 103
IX2 = IW+IDOT
IX3 = IX2+IDOT
IF (NA .NE. 0) GO TO 101
CALL PASSB4 (IDOT,L1,C,CH,WA(IW),WA(IX2),WA(IX3))
GO TO 102
101 CALL PASSB4 (IDOT,L1,CH,C,WA(IW),WA(IX2),WA(IX3))
102 NA = 1-NA
GO TO 115
103 IF (IP .NE. 2) GO TO 106
IF (NA .NE. 0) GO TO 104
CALL PASSB2 (IDOT,L1,C,CH,WA(IW))
GO TO 105
104 CALL PASSB2 (IDOT,L1,CH,C,WA(IW))
105 NA = 1-NA
GO TO 115
106 IF (IP .NE. 3) GO TO 109
IX2 = IW+IDOT
IF (NA .NE. 0) GO TO 107
CALL PASSB3 (IDOT,L1,C,CH,WA(IW),WA(IX2))
GO TO 108
107 CALL PASSB3 (IDOT,L1,CH,C,WA(IW),WA(IX2))
108 NA = 1-NA
GO TO 115
109 IF (IP .NE. 5) GO TO 112
IX2 = IW+IDOT
IX3 = IX2+IDOT
IX4 = IX3+IDOT
IF (NA .NE. 0) GO TO 110
CALL PASSB5 (IDOT,L1,C,CH,WA(IW),WA(IX2),WA(IX3),WA(IX4))
GO TO 111
110 CALL PASSB5 (IDOT,L1,CH,C,WA(IW),WA(IX2),WA(IX3),WA(IX4))
111 NA = 1-NA
GO TO 115
112 IF (NA .NE. 0) GO TO 113
CALL PASSB (NAC,IDOT,IP,L1,IDL1,C,C,C,CH,CH,WA(IW))
GO TO 114
113 CALL PASSB (NAC,IDOT,IP,L1,IDL1,CH,CH,CH,C,C,WA(IW))
114 IF (NAC .NE. 0) NA = 1-NA
115 L1 = L2
IW = IW+(IP-1)*IDOT
116 CONTINUE
IF (NA .EQ. 0) RETURN
N2 = N+N
DO 117 I=1,N2
C(I) = CH(I)
117 CONTINUE
RETURN
END
62 changes: 62 additions & 0 deletions cfftf1.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,62 @@
SUBROUTINE CFFTF1 (N,C,CH,WA,IFAC)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION CH(*) ,C(*) ,WA(*) ,IFAC(*)
NF = IFAC(2)
NA = 0
L1 = 1
IW = 1
DO 116 K1=1,NF
IP = IFAC(K1+2)
L2 = IP*L1
IDO = N/L2
IDOT = IDO+IDO
IDL1 = IDOT*L1
IF (IP .NE. 4) GO TO 103
IX2 = IW+IDOT
IX3 = IX2+IDOT
IF (NA .NE. 0) GO TO 101
CALL PASSF4 (IDOT,L1,C,CH,WA(IW),WA(IX2),WA(IX3))
GO TO 102
101 CALL PASSF4 (IDOT,L1,CH,C,WA(IW),WA(IX2),WA(IX3))
102 NA = 1-NA
GO TO 115
103 IF (IP .NE. 2) GO TO 106
IF (NA .NE. 0) GO TO 104
CALL PASSF2 (IDOT,L1,C,CH,WA(IW))
GO TO 105
104 CALL PASSF2 (IDOT,L1,CH,C,WA(IW))
105 NA = 1-NA
GO TO 115
106 IF (IP .NE. 3) GO TO 109
IX2 = IW+IDOT
IF (NA .NE. 0) GO TO 107
CALL PASSF3 (IDOT,L1,C,CH,WA(IW),WA(IX2))
GO TO 108
107 CALL PASSF3 (IDOT,L1,CH,C,WA(IW),WA(IX2))
108 NA = 1-NA
GO TO 115
109 IF (IP .NE. 5) GO TO 112
IX2 = IW+IDOT
IX3 = IX2+IDOT
IX4 = IX3+IDOT
IF (NA .NE. 0) GO TO 110
CALL PASSF5 (IDOT,L1,C,CH,WA(IW),WA(IX2),WA(IX3),WA(IX4))
GO TO 111
110 CALL PASSF5 (IDOT,L1,CH,C,WA(IW),WA(IX2),WA(IX3),WA(IX4))
111 NA = 1-NA
GO TO 115
112 IF (NA .NE. 0) GO TO 113
CALL PASSF (NAC,IDOT,IP,L1,IDL1,C,C,C,CH,CH,WA(IW))
GO TO 114
113 CALL PASSF (NAC,IDOT,IP,L1,IDL1,CH,CH,CH,C,C,WA(IW))
114 IF (NAC .NE. 0) NA = 1-NA
115 L1 = L2
IW = IW+(IP-1)*IDOT
116 CONTINUE
IF (NA .EQ. 0) RETURN
N2 = N+N
DO 117 I=1,N2
C(I) = CH(I)
117 CONTINUE
RETURN
END
61 changes: 61 additions & 0 deletions cffti1.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,61 @@
SUBROUTINE CFFTI1 (N,WA,IFAC)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION WA(*) ,IFAC(*) ,NTRYH(4)
DATA NTRYH(1),NTRYH(2),NTRYH(3),NTRYH(4)/3,4,2,5/
NL = N
NF = 0
J = 0
101 J = J+1
IF (J-4) 102,102,103
102 NTRY = NTRYH(J)
GO TO 104
103 NTRY = NTRY+2
104 NQ = NL/NTRY
NR = NL-NTRY*NQ
IF (NR) 101,105,101
105 NF = NF+1
IFAC(NF+2) = NTRY
NL = NQ
IF (NTRY .NE. 2) GO TO 107
IF (NF .EQ. 1) GO TO 107
DO 106 I=2,NF
IB = NF-I+2
IFAC(IB+2) = IFAC(IB+1)
106 CONTINUE
IFAC(3) = 2
107 IF (NL .NE. 1) GO TO 104
IFAC(1) = N
IFAC(2) = NF
TPI = 6.28318530717958647692D0
ARGH = TPI/FLOAT(N)
I = 2
L1 = 1
DO 110 K1=1,NF
IP = IFAC(K1+2)
LD = 0
L2 = L1*IP
IDO = N/L2
IDOT = IDO+IDO+2
IPM = IP-1
DO 109 J=1,IPM
I1 = I
WA(I-1) = 1.0D0
WA(I) = 0.0D0
LD = LD+L1
FI = 0.0D0
ARGLD = FLOAT(LD)*ARGH
DO 108 II=4,IDOT,2
I = I+2
FI = FI+1.D0
ARG = FI*ARGLD
WA(I-1) = COS(ARG)
WA(I) = SIN(ARG)
108 CONTINUE
IF (IP .LE. 5) GO TO 109
WA(I1-1) = WA(I-1)
WA(I1) = WA(I)
109 CONTINUE
L1 = L2
110 CONTINUE
RETURN
END
28 changes: 28 additions & 0 deletions cosqb1.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,28 @@
SUBROUTINE COSQB1 (N,X,W,XH)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION X(1) ,W(1) ,XH(1)
NS2 = (N+1)/2
NP2 = N+2
DO 101 I=3,N,2
XIM1 = X(I-1)+X(I)
X(I) = X(I)-X(I-1)
X(I-1) = XIM1
101 CONTINUE
X(1) = X(1)+X(1)
MODN = MOD(N,2)
IF (MODN .EQ. 0) X(N) = X(N)+X(N)
CALL DFFTB (N,X,XH)
DO 102 K=2,NS2
KC = NP2-K
XH(K) = W(K-1)*X(KC)+W(KC-1)*X(K)
XH(KC) = W(K-1)*X(K)-W(KC-1)*X(KC)
102 CONTINUE
IF (MODN .EQ. 0) X(NS2+1) = W(NS2)*(X(NS2+1)+X(NS2+1))
DO 103 K=2,NS2
KC = NP2-K
X(K) = XH(K)+XH(KC)
X(KC) = XH(K)-XH(KC)
103 CONTINUE
X(1) = X(1)+X(1)
RETURN
END
26 changes: 26 additions & 0 deletions cosqf1.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,26 @@
SUBROUTINE COSQF1 (N,X,W,XH)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION X(1) ,W(1) ,XH(1)
NS2 = (N+1)/2
NP2 = N+2
DO 101 K=2,NS2
KC = NP2-K
XH(K) = X(K)+X(KC)
XH(KC) = X(K)-X(KC)
101 CONTINUE
MODN = MOD(N,2)
IF (MODN .EQ. 0) XH(NS2+1) = X(NS2+1)+X(NS2+1)
DO 102 K=2,NS2
KC = NP2-K
X(K) = W(K-1)*XH(KC)+W(KC-1)*XH(K)
X(KC) = W(K-1)*XH(K)-W(KC-1)*XH(KC)
102 CONTINUE
IF (MODN .EQ. 0) X(NS2+1) = W(NS2)*XH(NS2+1)
CALL DFFTF (N,X,XH)
DO 103 I=3,N,2
XIM1 = X(I-1)-X(I)
X(I) = X(I-1)+X(I)
X(I-1) = XIM1
103 CONTINUE
RETURN
END
14 changes: 14 additions & 0 deletions dcosqb.f
Original file line number Diff line number Diff line change
@@ -0,0 +1,14 @@
SUBROUTINE DCOSQB (N,X,WSAVE)
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
DIMENSION X(*) ,WSAVE(*)
DATA TSQRT2 /2.82842712474619009760D0/
IF (N-2) 101,102,103
101 X(1) = 4.0D0*X(1)
RETURN
102 X1 = 4.0D0*(X(1)+X(2))
X(2) = TSQRT2*(X(1)-X(2))
X(1) = X1
RETURN
103 CALL COSQB1 (N,X,WSAVE,WSAVE(N+1))
RETURN
END
Loading

0 comments on commit 030cbf4

Please sign in to comment.