1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
|
*DECK RADB5
SUBROUTINE RADB5 (IDO, L1, CC, CH, WA1, WA2, WA3, WA4)
C***BEGIN PROLOGUE RADB5
C***SUBSIDIARY
C***PURPOSE Calculate the fast Fourier transform of subvectors of
C length five.
C***LIBRARY SLATEC (FFTPACK)
C***TYPE SINGLE PRECISION (RADB5-S)
C***AUTHOR Swarztrauber, P. N., (NCAR)
C***ROUTINES CALLED (NONE)
C***REVISION HISTORY (YYMMDD)
C 790601 DATE WRITTEN
C 830401 Modified to use SLATEC library source file format.
C 860115 Modified by Ron Boisvert to adhere to Fortran 77 by
C (a) changing dummy array size declarations (1) to (*),
C (b) changing definition of variables PI, TI11, TI12,
C TR11, TR12 by using FORTRAN intrinsic functions ATAN
C and SIN instead of DATA statements.
C 881128 Modified by Dick Valent to meet prologue standards.
C 890831 Modified array declarations. (WRB)
C 891214 Prologue converted to Version 4.0 format. (BAB)
C 900402 Added TYPE section. (WRB)
C***END PROLOGUE RADB5
DIMENSION CC(IDO,5,*), CH(IDO,L1,5), WA1(*), WA2(*), WA3(*),
+ WA4(*)
C***FIRST EXECUTABLE STATEMENT RADB5
PI = 4.*ATAN(1.)
TR11 = SIN(.1*PI)
TI11 = SIN(.4*PI)
TR12 = -SIN(.3*PI)
TI12 = SIN(.2*PI)
DO 101 K=1,L1
TI5 = CC(1,3,K)+CC(1,3,K)
TI4 = CC(1,5,K)+CC(1,5,K)
TR2 = CC(IDO,2,K)+CC(IDO,2,K)
TR3 = CC(IDO,4,K)+CC(IDO,4,K)
CH(1,K,1) = CC(1,1,K)+TR2+TR3
CR2 = CC(1,1,K)+TR11*TR2+TR12*TR3
CR3 = CC(1,1,K)+TR12*TR2+TR11*TR3
CI5 = TI11*TI5+TI12*TI4
CI4 = TI12*TI5-TI11*TI4
CH(1,K,2) = CR2-CI5
CH(1,K,3) = CR3-CI4
CH(1,K,4) = CR3+CI4
CH(1,K,5) = CR2+CI5
101 CONTINUE
IF (IDO .EQ. 1) RETURN
IDP2 = IDO+2
IF((IDO-1)/2.LT.L1) GO TO 104
DO 103 K=1,L1
CDIR$ IVDEP
DO 102 I=3,IDO,2
IC = IDP2-I
TI5 = CC(I,3,K)+CC(IC,2,K)
TI2 = CC(I,3,K)-CC(IC,2,K)
TI4 = CC(I,5,K)+CC(IC,4,K)
TI3 = CC(I,5,K)-CC(IC,4,K)
TR5 = CC(I-1,3,K)-CC(IC-1,2,K)
TR2 = CC(I-1,3,K)+CC(IC-1,2,K)
TR4 = CC(I-1,5,K)-CC(IC-1,4,K)
TR3 = CC(I-1,5,K)+CC(IC-1,4,K)
CH(I-1,K,1) = CC(I-1,1,K)+TR2+TR3
CH(I,K,1) = CC(I,1,K)+TI2+TI3
CR2 = CC(I-1,1,K)+TR11*TR2+TR12*TR3
CI2 = CC(I,1,K)+TR11*TI2+TR12*TI3
CR3 = CC(I-1,1,K)+TR12*TR2+TR11*TR3
CI3 = CC(I,1,K)+TR12*TI2+TR11*TI3
CR5 = TI11*TR5+TI12*TR4
CI5 = TI11*TI5+TI12*TI4
CR4 = TI12*TR5-TI11*TR4
CI4 = TI12*TI5-TI11*TI4
DR3 = CR3-CI4
DR4 = CR3+CI4
DI3 = CI3+CR4
DI4 = CI3-CR4
DR5 = CR2+CI5
DR2 = CR2-CI5
DI5 = CI2-CR5
DI2 = CI2+CR5
CH(I-1,K,2) = WA1(I-2)*DR2-WA1(I-1)*DI2
CH(I,K,2) = WA1(I-2)*DI2+WA1(I-1)*DR2
CH(I-1,K,3) = WA2(I-2)*DR3-WA2(I-1)*DI3
CH(I,K,3) = WA2(I-2)*DI3+WA2(I-1)*DR3
CH(I-1,K,4) = WA3(I-2)*DR4-WA3(I-1)*DI4
CH(I,K,4) = WA3(I-2)*DI4+WA3(I-1)*DR4
CH(I-1,K,5) = WA4(I-2)*DR5-WA4(I-1)*DI5
CH(I,K,5) = WA4(I-2)*DI5+WA4(I-1)*DR5
102 CONTINUE
103 CONTINUE
RETURN
104 DO 106 I=3,IDO,2
IC = IDP2-I
CDIR$ IVDEP
DO 105 K=1,L1
TI5 = CC(I,3,K)+CC(IC,2,K)
TI2 = CC(I,3,K)-CC(IC,2,K)
TI4 = CC(I,5,K)+CC(IC,4,K)
TI3 = CC(I,5,K)-CC(IC,4,K)
TR5 = CC(I-1,3,K)-CC(IC-1,2,K)
TR2 = CC(I-1,3,K)+CC(IC-1,2,K)
TR4 = CC(I-1,5,K)-CC(IC-1,4,K)
TR3 = CC(I-1,5,K)+CC(IC-1,4,K)
CH(I-1,K,1) = CC(I-1,1,K)+TR2+TR3
CH(I,K,1) = CC(I,1,K)+TI2+TI3
CR2 = CC(I-1,1,K)+TR11*TR2+TR12*TR3
CI2 = CC(I,1,K)+TR11*TI2+TR12*TI3
CR3 = CC(I-1,1,K)+TR12*TR2+TR11*TR3
CI3 = CC(I,1,K)+TR12*TI2+TR11*TI3
CR5 = TI11*TR5+TI12*TR4
CI5 = TI11*TI5+TI12*TI4
CR4 = TI12*TR5-TI11*TR4
CI4 = TI12*TI5-TI11*TI4
DR3 = CR3-CI4
DR4 = CR3+CI4
DI3 = CI3+CR4
DI4 = CI3-CR4
DR5 = CR2+CI5
DR2 = CR2-CI5
DI5 = CI2-CR5
DI2 = CI2+CR5
CH(I-1,K,2) = WA1(I-2)*DR2-WA1(I-1)*DI2
CH(I,K,2) = WA1(I-2)*DI2+WA1(I-1)*DR2
CH(I-1,K,3) = WA2(I-2)*DR3-WA2(I-1)*DI3
CH(I,K,3) = WA2(I-2)*DI3+WA2(I-1)*DR3
CH(I-1,K,4) = WA3(I-2)*DR4-WA3(I-1)*DI4
CH(I,K,4) = WA3(I-2)*DI4+WA3(I-1)*DR4
CH(I-1,K,5) = WA4(I-2)*DR5-WA4(I-1)*DI5
CH(I,K,5) = WA4(I-2)*DI5+WA4(I-1)*DR5
105 CONTINUE
106 CONTINUE
RETURN
END
|
