Line data    Source code

       1              : block 
       2              : 
       3              : !Local variables ------------------------------
       4              : !scalars
       5              :  integer :: ifft,ix,iy,iz,nx,ny,nz,ldx,ldy,ldz,dat,fftalga,padat
       6              :  real(dp) :: fact
       7              :  character(len=500) :: msg
       8              : !arrays
       9         2958 :  real(dp),allocatable :: arr(:,:,:,:),ftarr(:,:,:,:) 
      10              : 
      11              : ! *************************************************************************
      12              : 
      13         2958 :   nx = plan%dims(1);   ny=plan%dims(2);   nz=plan%dims(3)
      14         2958 :  ldx = plan%embed(1); ldy=plan%embed(2); ldz=plan%embed(3)
      15         2958 :  fftalga = plan%fftalg/100
      16              : 
      17            0 :  select case (fftalga)
      18              :  case (FFT_FFTW3) 
      19            0 :    call fftw3_c2c_op(nx, ny, nz, ldx, ldy, ldz, ndat__, iscale__, isign, ff, gg)
      20              : 
      21              :  case (FFT_DFTI) 
      22         2890 :    call dfti_c2c_op(nx, ny, nz, ldx, ldy, ldz, ndat__, iscale__, isign, ff, gg)
      23              : 
      24              :  case (FFT_SG, FFT_SG2002)
      25              :    ! Fallback to sg_fft_cc 
      26              :    ! 1) we have to change shape and type complex -> double because we have an explicit interface. 
      27              :    ! 2) Cannot use [ZC]copy: this is a template used both for single and double precision.
      28              : 
      29          340 :    ABI_MALLOC(arr,   (2,ldx,ldy,ldz))
      30          272 :    ABI_MALLOC(ftarr, (2,ldx,ldy,ldz))
      31              : 
      32          208 :    do dat=1,ndat__
      33          140 :      padat = (dat-1) * plan%ldxyz
      34              :      !  
      35              :      ! Copy input data in arr
      36         3740 :      do iz=1,nz
      37       198140 :        do iy=1,ny
      38     15166800 :          do ix=1,nx
      39     14968800 :            ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
      40     14968800 :            arr(1,ix,iy,iz) = DBLE (ff(ifft))
      41     15163200 :            arr(2,ix,iy,iz) = AIMAG(ff(ifft)) 
      42              :          end do 
      43              :        end do
      44              :      end do
      45              : 
      46              :      ! c2c with ndat = 1
      47          140 :      call sg_fft_cc(plan%fftcache,nx,ny,nz,ldx,ldy,ldz,1,isign,arr,ftarr)
      48              : 
      49          208 :      if (isign == -1) then  
      50              :        ! Copy and scale the transform
      51           70 :        fact = one/plan%nfft
      52           70 :        if (iscale__ == 0) fact = one
      53         1870 :        do iz=1,nz
      54        99070 :          do iy=1,ny
      55      7583400 :            do ix=1,nx
      56      7484400 :              ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
      57      7581600 :              gg(ifft) = fact * DCMPLX(ftarr(1,ix,iy,iz), ftarr(2,ix,iy,iz))
      58              :            end do 
      59              :          end do
      60              :        end do
      61              :        !
      62              :      else
      63              :        ! Direct copy.
      64         1870 :        do iz=1,nz
      65        99070 :          do iy=1,ny
      66      7583400 :            do ix=1,nx
      67      7484400 :              ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
      68      7581600 :              gg(ifft) = DCMPLX(ftarr(1,ix,iy,iz), ftarr(2,ix,iy,iz))
      69              :            end do 
      70              :          end do
      71              :        end do
      72              :      end if
      73              : 
      74              :    end do
      75              : 
      76           68 :    ABI_FREE(arr)
      77           68 :    ABI_FREE(ftarr)
      78              : 
      79              :  case default
      80            0 :    write(msg,'(a,i0)')"Wrong fftalga ",fftalga
      81         2958 :    ABI_ERROR(msg)
      82              :  end select
      83              : 
      84              : end block