fft3d.f90 2.57 KB
Newer Older
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
      MODULE m_fft3d
      CONTAINS
      SUBROUTINE fft3d(&
     &                 afft,bfft,fg3,&
     &                 stars,isn,scaled)

!************************************************************
!*                                                          *
!* interface for fg3(star) -- FFT --> (a,b)fft (r) (isn=+1) *
!*         or (a,b)fft (r) -- FFT --> fg3(star)    (isn=-1) *
!*                                                          *
!* dimension of (a,b)fft is (3*k1d x 3*k2d x 3*k3d)         *
!* afft and bfft contain the real/imaginary part of the FFT *
!* igfft(i,1) is the pointer from the G-sphere to stars     *
!* igfft(i,2) is the pointer from the G-sphere to fft-grid  *
!* pgfft(i)   contains the phases of the G-vectors of sph.  *
!*                                                          *
!************************************************************
      USE m_cfft
      USE m_types
      IMPLICIT NONE
    
      INTEGER, INTENT (IN) :: isn
      TYPE(t_stars),INTENT(IN):: stars
      REAL,    INTENT (INOUT) :: afft(0:27*stars%k1d*stars%k2d*stars%k3d-1)
      REAL,    INTENT (INOUT) :: bfft(0:27*stars%k1d*stars%k2d*stars%k3d-1)
      COMPLEX                 :: fg3(stars%ng3)
      LOGICAL,INTENT(IN),OPTIONAL :: scaled !< determines if coefficients are scaled by stars%nstr

      INTEGER i,ifftd
      REAL scale
32
      COMPLEX ctmp
33 34 35 36 37 38 39 40 41 42

      ifftd=27*stars%k1d*stars%k2d*stars%k3d
     
      IF (isn.GT.0) THEN
!
!  ---> put stars onto the fft-grid 
!
        afft=0.0
        bfft=0.0
        DO i=0,stars%kimax
43 44 45
          ctmp = fg3(stars%igfft(i,1))*stars%pgfft(i)
          afft(stars%igfft(i,2))=real(ctmp)
          bfft(stars%igfft(i,2))=aimag(ctmp)
46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62
        ENDDO
      ENDIF

!---> now do the fft (isn=+1 : G -> r ; isn=-1 : r -> G)

      CALL cfft(afft,bfft,ifftd,3*stars%k1d,3*stars%k1d,isn)
      CALL cfft(afft,bfft,ifftd,3*stars%k2d,9*stars%k1d*stars%k2d,isn)
      CALL cfft(afft,bfft,ifftd,3*stars%k3d,ifftd,isn)

      IF (isn.LT.0) THEN
!
!  ---> collect stars from the fft-grid
!
        DO i=1,stars%ng3
          fg3(i) = cmplx(0.0,0.0)
        ENDDO
        DO i=0,stars%kimax
63 64
          fg3(stars%igfft(i,1)) = fg3(stars%igfft(i,1)) + CONJG( stars%pgfft(i) ) * &
     &                CMPLX(afft(stars%igfft(i,2)),bfft(stars%igfft(i,2)))
65 66 67 68 69 70 71 72 73 74 75 76 77 78 79
        ENDDO
        scale=1.0/ifftd
        IF (PRESENT(scaled)) THEN
           IF (scaled) THEN
               fg3=scale*fg3/stars%nstr
           ELSE
               fg3=scale*fg3    
           ENDIF
        ELSE
           fg3=scale*fg3/stars%nstr
        ENDIF
      ENDIF

      END SUBROUTINE fft3d
      END MODULE m_fft3d