Skip to content

Commit c17acdc

Browse files
authored
Merge pull request #7 from fortran-fans/update_wsave_intent
Add `dzffti/dzfftf/dzfftb` interfaces
2 parents 0c4ecc4 + 3151a4f commit c17acdc

File tree

5 files changed

+330
-21
lines changed

5 files changed

+330
-21
lines changed

doc/doc

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -230,7 +230,7 @@ but slower version of dfftf.
230230
input parameters
231231

232232
n the length of the array r to be transformed. the method
233-
is must efficient when n is the product of small primes.
233+
is most efficient when n is the product of small primes.
234234

235235
r a real array of length n which contains the sequence
236236
to be transformed. r is not destroyed.

doc/specs/fftpack.md

Lines changed: 247 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,7 @@ Experimental.
2323

2424
Pure function.
2525

26-
#### Snytax
26+
#### Syntax
2727

2828
`call [[fftpack(module):zffti(interface)]](n, wsave)`
2929

@@ -75,7 +75,7 @@ Experimental.
7575

7676
Pure function.
7777

78-
#### Snytax
78+
#### Syntax
7979

8080
`call [[fftpack(module):zfftf(interface)]](n, c, wsave)`
8181

@@ -99,7 +99,7 @@ for j=1,...,n
9999
```
100100

101101
`wsave`: Shall be a `real` array.
102-
This argument is `intent(inout)`.
102+
This argument is `intent(in)`.
103103
A `real` work array which must be dimensioned at least `4n+15` in the program that calls `zfftf`.
104104
The wsave array must be initialized by calling subroutine `zffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`.
105105
This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first.
@@ -145,7 +145,7 @@ Experimental.
145145

146146
Pure function.
147147

148-
#### Snytax
148+
#### Syntax
149149

150150
`call [[fftpack(module):zfftb(interface)]](n, c, wsave)`
151151

@@ -169,7 +169,7 @@ for j=1,...,n
169169
```
170170

171171
`wsave`: Shall be a `real` array.
172-
This argument is `intent(inout)`.
172+
This argument is `intent(in)`.
173173
A `real` work array which must be dimensioned at least `4n+15` in the program that calls `zfftf`. The `wsave` array must be initialized by calling subroutine `zffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`. This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first. The same `wsave` array can be used by `zfftf` and `zfftb`.
174174
Contains initialization calculations which must not be destroyed between calls of subroutine `zfftf` or `zfftb`.
175175

@@ -205,7 +205,7 @@ Experimental.
205205

206206
Pure function.
207207

208-
#### Snytax
208+
#### Syntax
209209

210210
`result = [[fftpack(module):fft(interface)]](x [, n])`
211211

@@ -253,7 +253,7 @@ Experimental.
253253

254254
Pure function.
255255

256-
#### Snytax
256+
#### Syntax
257257

258258
`result = [[fftpack(module):ifft(interface)]](x [, n])`
259259

@@ -298,7 +298,7 @@ Experimental.
298298

299299
Pure function.
300300

301-
#### Snytax
301+
#### Syntax
302302

303303
`call [[fftpack(module):dffti(interface)]](n, wsave)`
304304

@@ -348,7 +348,7 @@ Experimental.
348348

349349
Pure function.
350350

351-
#### Snytax
351+
#### Syntax
352352

353353
`call [[fftpack(module):dfftf(interface)]](n, r, wsave)`
354354

@@ -385,7 +385,7 @@ if n is even
385385
```
386386

387387
`wsave`: Shall be a `real` array.
388-
This argument is `intent(inout)`.
388+
This argument is `intent(in)`.
389389
A `real` work array which must be dimensioned at least `4n+15` in the program that calls `dfftf`.
390390
The wsave array must be initialized by calling subroutine `dffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`.
391391
This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first.
@@ -430,7 +430,7 @@ Experimental.
430430

431431
Pure function.
432432

433-
#### Snytax
433+
#### Syntax
434434

435435
`call [[fftpack(module):dfftb(interface)]](n, r, wsave)`
436436

@@ -464,7 +464,7 @@ for n odd and for i = 1,...,n
464464
```
465465

466466
`wsave`: Shall be a `real` array.
467-
This argument is `intent(inout)`.
467+
This argument is `intent(in)`.
468468
A `real` work array which must be dimensioned at least `2n+15` in the program that calls `dfftf`. The `wsave` array must be initialized by calling subroutine `dffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`. This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first. The same `wsave` array can be used by `dfftf` and `dfftb`.
469469
Contains initialization calculations which must not be destroyed between calls of subroutine `dfftf` or `dfftb`.
470470

@@ -499,7 +499,7 @@ Experimental.
499499

500500
Pure function.
501501

502-
#### Snytax
502+
#### Syntax
503503

504504
`result = [[fftpack(module):rfft(interface)]](x [, n])`
505505

@@ -547,7 +547,7 @@ Experimental.
547547

548548
Pure function.
549549

550-
#### Snytax
550+
#### Syntax
551551

552552
`result = [[fftpack(module):irfft(interface)]](x [, n])`
553553

@@ -576,6 +576,237 @@ program demo_irfft
576576
end program demo_irfft
577577
```
578578

579+
## Simplified fourier transform of double real periodic sequences
580+
581+
### `dzffti`
582+
583+
#### Description
584+
585+
Initializes the array `wsave` which is used in both `dzfftf` and `dzfftb`.
586+
The prime factorization of `n` together with a tabulation of the trigonometric functions are computed and stored in `wsave`.
587+
588+
#### Status
589+
590+
Experimental
591+
592+
#### Class
593+
594+
Prue function.
595+
596+
#### Syntax
597+
598+
`call [[fftpack(module):dzffti(interface)]](n, wsave)`
599+
600+
#### Arguments
601+
602+
`n`: Shall be an `integer` scalar.
603+
This argument is `intent(in)`.
604+
The length of the sequence to be transformed.
605+
606+
`wsave`: Shall be a `real` and rank-1 array.
607+
This argument is `intent(out)`.
608+
A work array which must be dimensioned at least `3*n+15`.
609+
The same work array can be used for both `dzfftf` and `dzfftb` as long as n remains unchanged.
610+
Different `wsave` arrays are required for different values of `n`.
611+
612+
##### Warning
613+
614+
The contents of `wsave` must not be changed between calls of `dzfftf` or `dzfftb`.
615+
616+
#### Example
617+
618+
```fortran
619+
program demo_dzffti
620+
use fftpack, only: dzffti
621+
real(kind=8) :: x(4) = [1, 2, 3, 4]
622+
real(kind=8) :: w(3*4 + 15)
623+
call dzffti(4, w) !! Initializes the array `w` which is used in both `dzfftf` and `dzfftb`.
624+
end program demo_dzffti
625+
```
626+
627+
### `dzfftf`
628+
629+
#### Description
630+
631+
Computes the fourier coefficients of a `real` perodic sequence (fourier analysis).
632+
The transform is defined below at output parameters `azero`, `a` and `b`.
633+
`dzfftf` is a simplified but **slower version** of `dfftf`.
634+
635+
#### Status
636+
637+
Experimental
638+
639+
#### Class
640+
641+
Pure function.
642+
643+
#### Syntax
644+
645+
`call [[fftpack(module):dzfftf(interface)]](n, r, azero, a, b, wsave)`
646+
647+
#### Arguments
648+
649+
`n`: Shall be an `integer` scalar.
650+
This argument is `intent(in)`.
651+
The length of the array `r` to be transformed.
652+
The method is most efficient when `n` is the product of small primes.
653+
654+
`r`: Shall be a `real` and rank-1 array.
655+
This argument is `intent(in)`.
656+
A `real` array of length `n` which contains the sequence to be transformed. `r` is not destroyed.
657+
658+
`azero`: Shall be a `real` scalar.
659+
This argument is `intent(out)`.
660+
The sum from `i=1` to `i=n` of `r(i)/n`.
661+
662+
`a`, `b`: Shall be a `real` and rank-1 array.
663+
This argument is `intent(out)`.
664+
```
665+
for n even b(n/2)=0. and a(n/2) is the sum from i=1 to i=n of (-1)**(i-1)*r(i)/n
666+
667+
for n even define kmax=n/2-1
668+
for n odd define kmax=(n-1)/2
669+
670+
then for k=1,...,kmax
671+
672+
a(k) equals the sum from i=1 to i=n of
673+
674+
2./n*r(i)*cos(k*(i-1)*2*pi/n)
675+
676+
b(k) equals the sum from i=1 to i=n of
677+
678+
2./n*r(i)*sin(k*(i-1)*2*pi/n)
679+
```
680+
681+
`wsave`: Shall be a `real` and rank-1 array.
682+
This argument is `intent(in)`.
683+
A work array which must be dimensioned at least `3*n+15`.
684+
In the program that calls `dzfftf`. The `wsave` array must be initialized by calling subroutine `dzffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`.
685+
This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first.
686+
The same `wsave` array can be used by `dzfftf` and `dzfftb`.
687+
688+
#### Example
689+
690+
```fortran
691+
program demo_dzfftf
692+
use fftpack, only: dzffti, dzfftf
693+
real(kind=8) :: x(4) = [1, 2, 3, 4]
694+
real(kind=8) :: w(3*4 + 15)
695+
real(kind=8) :: azero, a(4/2), b(4/2)
696+
call dzffti(4, w)
697+
call dzfftf(4, x, azero, a, b, w) !! `azero`: 2.5; `a`: [-1.0, -0.5]; `b`: [-1.0, -0.0]
698+
end program demo_dzfftf
699+
```
700+
701+
### `dzfftb`
702+
703+
#### Description
704+
705+
Computes a `real` perodic sequence from its fourier coefficients (fourier synthesis).
706+
The transform is defined below at output parameter `r`.
707+
`dzfftb` is a simplified but **slower version** of `dfftb`.
708+
709+
#### Status
710+
711+
Experimental
712+
713+
#### Class
714+
715+
Pure function.
716+
717+
#### Syntax
718+
719+
`call [[fftpack(module):dzfftb(interface)]](n, r, azero, a, b, wsave)`
720+
721+
#### Arguments
722+
723+
`n`: Shall be an `integer` scalar.
724+
This argument is `intent(in)`.
725+
The length of the output array `r`.
726+
The method is most efficient when `n` is the product of small primes.
727+
728+
`r`: Shall be a `real` and rank-1 array.
729+
This argument is `intent(out)`.
730+
```
731+
if n is even define kmax=n/2
732+
if n is odd define kmax=(n-1)/2
733+
734+
then for i=1,...,n
735+
736+
r(i)=azero plus the sum from k=1 to k=kmax of
737+
738+
a(k)*cos(k*(i-1)*2*pi/n)+b(k)*sin(k*(i-1)*2*pi/n)
739+
```
740+
Complex notation:
741+
```
742+
for j=1,...,n
743+
744+
r(j) equals the sum from k=-kmax to k=kmax of
745+
746+
c(k)*exp(i*k*(j-1)*2*pi/n)
747+
748+
where
749+
750+
c(k) = .5*cmplx(a(k),-b(k)) for k=1,...,kmax
751+
752+
c(-k) = conjg(c(k))
753+
754+
c(0) = azero
755+
756+
and i=sqrt(-1)
757+
```
758+
Amplitude - phase notation:
759+
```
760+
for i=1,...,n
761+
762+
r(i) equals azero plus the sum from k=1 to k=kmax of
763+
764+
alpha(k)*cos(k*(i-1)*2*pi/n+beta(k))
765+
766+
where
767+
768+
alpha(k) = sqrt(a(k)*a(k)+b(k)*b(k))
769+
770+
cos(beta(k))=a(k)/alpha(k)
771+
772+
sin(beta(k))=-b(k)/alpha(k)
773+
```
774+
775+
`azero`: Shall be a `real` scalar.
776+
This argument is `intent(in)`.
777+
The constant fourier coefficient.
778+
779+
`a`, `b`: Shall be a `real` and rank-1 array.
780+
This argument is `intent(in)`.
781+
Arrays which contain the remaining fourier coefficients these arrays are not destroyed.
782+
The length of these arrays depends on whether `n` is even or odd.
783+
```
784+
if n is even n/2 locations are required
785+
if n is odd (n-1)/2 locations are required
786+
```
787+
788+
`wsave`: Shall be a `real` and rank-1 array.
789+
This argument is `intent(in)`.
790+
A work array which must be dimensioned at least `3*n+15`.
791+
In the program that calls `dzfftf`. The `wsave` array must be initialized by calling subroutine `dzffti(n,wsave)` and a different `wsave` array must be used for each different value of `n`.
792+
This initialization does not have to be repeated so long as `n` remains unchanged thus subsequent transforms can be obtained faster than the first.
793+
The same `wsave` array can be used by `dzfftf` and `dzfftb`.
794+
795+
#### Example
796+
797+
```fortran
798+
program demo_dzfftb
799+
use fftpack, only: dzffti, dzfftf, dzfftb
800+
real(kind=8) :: x(4) = [1, 2, 3, 4]
801+
real(kind=8) :: w(3*4 + 15)
802+
real(kind=8) :: azero, a(4/2), b(4/2)
803+
call dzffti(4, w)
804+
call dzfftf(4, x, azero, a, b, w) !! `azero`: 2.5; `a`: [-1.0, -0.5]; `b`: [-1.0, -0.0]
805+
x = 0.0
806+
call dzfftb(4, x, azero, a, b, w) !! `x`: [1.0, 2.0, 3.0, 4.0]
807+
end program demo_dzfftb
808+
```
809+
579810
## Utility functions
580811

581812
### `fftshift`
@@ -592,7 +823,7 @@ Experimental.
592823

593824
Pure function.
594825

595-
#### Snytax
826+
#### Syntax
596827

597828
`result = [[fftpack(module):fftshift(interface)]](x)`
598829

@@ -633,7 +864,7 @@ Experimental.
633864

634865
Pure function.
635866

636-
#### Snytax
867+
#### Syntax
637868

638869
`result = [[fftpack(module):ifftshift(interface)]](x)`
639870

0 commit comments

Comments
 (0)