[1403] | 1 | ! |
---|
| 2 | ! $Id: mod_fft_fftw.F90 1403 2010-07-01 09:02:53Z evignon $ |
---|
| 3 | ! |
---|
| 4 | |
---|
[986] | 5 | MODULE mod_fft_fftw |
---|
| 6 | |
---|
| 7 | #ifdef FFT_FFTW |
---|
| 8 | |
---|
[1403] | 9 | REAL, SAVE :: scale_factor |
---|
| 10 | INTEGER, SAVE :: vsize |
---|
| 11 | INTEGER, PARAMETER :: inc=1 |
---|
[986] | 12 | |
---|
[1403] | 13 | INTEGER*8, ALLOCATABLE, DIMENSION(:), SAVE :: plan_forward |
---|
| 14 | INTEGER*8, ALLOCATABLE, DIMENSION(:), SAVE :: plan_backward |
---|
[986] | 15 | |
---|
| 16 | CONTAINS |
---|
| 17 | |
---|
[1403] | 18 | SUBROUTINE Init_fft(iim,nvectmax) |
---|
[986] | 19 | IMPLICIT NONE |
---|
[1403] | 20 | #include <fftw3.f> |
---|
[986] | 21 | INTEGER :: iim |
---|
[1403] | 22 | INTEGER :: nvectmax |
---|
| 23 | |
---|
| 24 | INTEGER :: itmp |
---|
| 25 | |
---|
| 26 | INTEGER :: rank |
---|
| 27 | INTEGER :: howmany |
---|
| 28 | INTEGER :: istride, idist |
---|
| 29 | INTEGER :: ostride, odist |
---|
| 30 | INTEGER, DIMENSION(1) :: n_array, inembed, onembed |
---|
| 31 | |
---|
| 32 | REAL, DIMENSION(iim+1,nvectmax) :: dbidon |
---|
| 33 | COMPLEX, DIMENSION(iim/2+1,nvectmax) :: cbidon |
---|
| 34 | |
---|
| 35 | vsize = iim |
---|
| 36 | scale_factor = 1./SQRT(1.*vsize) |
---|
| 37 | |
---|
| 38 | dbidon = 0 |
---|
| 39 | cbidon = 0 |
---|
| 40 | |
---|
| 41 | ALLOCATE(plan_forward(nvectmax)) |
---|
| 42 | ALLOCATE(plan_backward(nvectmax)) |
---|
[986] | 43 | |
---|
[1403] | 44 | WRITE(*,*)"!---------------------!" |
---|
| 45 | WRITE(*,*)"! !" |
---|
| 46 | WRITE(*,*)"! INITIALISATION FFTW !" |
---|
| 47 | WRITE(*,*)"! !" |
---|
| 48 | WRITE(*,*)"!---------------------!" |
---|
[986] | 49 | |
---|
[1403] | 50 | ! On initialise tous les plans |
---|
| 51 | DO itmp = 1, nvectmax |
---|
| 52 | rank = 1 |
---|
| 53 | n_array(1) = iim |
---|
| 54 | howmany = itmp |
---|
| 55 | inembed(1) = iim + 1 ; onembed(1) = iim/2 + 1 |
---|
| 56 | istride = 1 ; ostride = 1 |
---|
| 57 | idist = iim + 1 ; odist = iim/2 + 1 |
---|
[986] | 58 | |
---|
[1403] | 59 | CALL dfftw_plan_many_dft_r2c(plan_forward(itmp), rank, n_array, howmany, & |
---|
| 60 | & dbidon, inembed, istride, idist, & |
---|
| 61 | & cbidon, onembed, ostride, odist, & |
---|
| 62 | & FFTW_ESTIMATE) |
---|
| 63 | |
---|
| 64 | rank = 1 |
---|
| 65 | n_array(1) = iim |
---|
| 66 | howmany = itmp |
---|
| 67 | inembed(1) = iim/2 + 1 ; onembed(1) = iim + 1 |
---|
| 68 | istride = 1 ; ostride = 1 |
---|
| 69 | idist = iim/2 + 1 ; odist = iim + 1 |
---|
| 70 | CALL dfftw_plan_many_dft_c2r(plan_backward(itmp), rank, n_array, howmany, & |
---|
| 71 | & cbidon, inembed, istride, idist, & |
---|
| 72 | & dbidon, onembed, ostride, odist, & |
---|
| 73 | & FFTW_ESTIMATE) |
---|
| 74 | |
---|
| 75 | ENDDO |
---|
| 76 | |
---|
| 77 | WRITE(*,*)"!-------------------------!" |
---|
| 78 | WRITE(*,*)"! !" |
---|
| 79 | WRITE(*,*)"! FIN INITIALISATION FFTW !" |
---|
| 80 | WRITE(*,*)"! !" |
---|
| 81 | WRITE(*,*)"!-------------------------!" |
---|
| 82 | |
---|
[986] | 83 | END SUBROUTINE Init_fft |
---|
| 84 | |
---|
| 85 | |
---|
| 86 | SUBROUTINE fft_forward(vect,TF_vect,nb_vect) |
---|
| 87 | IMPLICIT NONE |
---|
[1403] | 88 | #include <fftw3.f> |
---|
| 89 | INTEGER,INTENT(IN) :: nb_vect |
---|
| 90 | REAL,INTENT(IN) :: vect(vsize+inc,nb_vect) |
---|
| 91 | COMPLEX,INTENT(OUT) :: TF_vect(vsize/2+1,nb_vect) |
---|
| 92 | |
---|
| 93 | CALL dfftw_execute_dft_r2c(plan_forward(nb_vect),vect,TF_vect) |
---|
| 94 | |
---|
| 95 | TF_vect = scale_factor * TF_vect |
---|
| 96 | |
---|
[986] | 97 | END SUBROUTINE fft_forward |
---|
| 98 | |
---|
| 99 | SUBROUTINE fft_backward(TF_vect,vect,nb_vect) |
---|
| 100 | IMPLICIT NONE |
---|
[1403] | 101 | #include <fftw3.f> |
---|
| 102 | INTEGER,INTENT(IN) :: nb_vect |
---|
| 103 | REAL,INTENT(OUT) :: vect(vsize+inc,nb_vect) |
---|
| 104 | COMPLEX,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect) |
---|
[986] | 105 | |
---|
[1403] | 106 | CALL dfftw_execute_dft_c2r(plan_backward(nb_vect),TF_vect,vect) |
---|
| 107 | |
---|
| 108 | vect = scale_factor * vect |
---|
| 109 | |
---|
[986] | 110 | END SUBROUTINE fft_backward |
---|
| 111 | |
---|
| 112 | #endif |
---|
| 113 | |
---|
| 114 | END MODULE mod_fft_fftw |
---|