Line data Source code
1 : block
2 : !Local variables-------------------------------
3 : !scalars
4 : integer,parameter :: cplex0=0,option3=3,npwin0=0,ndat1=1
5 : integer :: nx,ny,nz,ldx,ldy,ldz,fftalg,fftalga,fftalgc,fftcache
6 : integer :: idat,rspad,gspad,ix,iy,iz,ifft,ipw
7 : real(dp),parameter :: weight1_i=one,weight1_r=one
8 : !arrays
9 : integer :: kg_kin0(3,0)
10 : real(dp) :: dum_denpot(0,0,0),dum_fofgin(0,0)
11 780045 : real(dp),allocatable :: fofgout(:,:),fofr(:,:,:,:),fg_box(:,:,:,:)
12 : ! *************************************************************************
13 :
14 780045 : fftalg=ngfft(7); fftalga=fftalg/100; fftalgc=MOD(fftalg,10); fftcache=ngfft(8)
15 780045 : nx=ngfft(1); ny=ngfft(2); nz=ngfft(3)
16 :
17 0 : select case (fftalga)
18 : case (FFT_FFTW3)
19 0 : ldx=nx; ldy=ny; ldz=nz ! No augmentation, the caller does not support it.
20 : call fftw3_fftur(fftalg, fftcache, npw_k, nx, ny, nz, ldx, ldy, ldz, nspinor*ndat, istwf_k, mgfft, &
21 0 : kg_k, gbound_k, ur, ug)
22 :
23 : case (FFT_DFTI)
24 780004 : ldx=nx; ldy=ny; ldz=nz ! No augmentation, the caller does not support it.
25 : call dfti_fftur(fftalg, fftcache, npw_k, nx, ny, nz, ldx, ldy, ldz, nspinor*ndat, istwf_k, mgfft, &
26 780004 : kg_k, gbound_k, ur, ug)
27 :
28 : case (FFT_SG, FFT_SG2002)
29 : ! Goedecker routines.
30 : ! 1) sg_fftrisc does not support istwkf >= 2 for the u(r) --> u(G) transform.
31 41 : ldx=ngfft(4); ldy=ngfft(5); ldz=ngfft(6)
32 123 : ABI_MALLOC(fofgout,(2,npw_k))
33 205 : ABI_MALLOC(fofr,(2,ldx,ldy,ldz))
34 :
35 41 : if (istwf_k/=1) then
36 128 : ABI_MALLOC(fg_box,(2,ldx,ldy,ldz))
37 : end if
38 :
39 136 : do idat=1,nspinor*ndat
40 95 : gspad = (idat-1)*npw_k
41 95 : rspad = (idat-1)*nfft
42 : !
43 : ! Fill fofr from input ur array.
44 9545 : do iz=1,nz
45 950045 : do iy=1,ny
46 94594950 : do ix=1,nx
47 93645000 : ifft = ix + (iy-1)*nx + (iz-1)*nx*ny + rspad
48 93645000 : fofr(1,ix,iy,iz) = REAL( ur(ifft))
49 94585500 : fofr(2,ix,iy,iz) = AIMAG(ur(ifft))
50 : end do
51 : end do
52 : end do
53 :
54 95 : if (istwf_k==1) then
55 : !
56 : ! option=3 --> real space to reciprocal space.
57 : ! NOTE that in this case, istwf_k=1 must be used.
58 : ! this is not thread-safe! use ndat1
59 : call sg_fftrisc(cplex0,dum_denpot,dum_fofgin,fofgout,fofr,&
60 : gbound_k,gbound_k,istwf_k,kg_kin0,kg_k,mgfft,ndat1,ngfft,npwin0,npw_k,ldx,ldy,ldz,&
61 15 : option3,weight1_r,weight1_i)
62 :
63 : else
64 : ! Zero-padded on the Box, then transfer data from box to sphere. ndat=1 here.
65 80 : call sg_fftpad(fftcache,mgfft,nx,ny,nz,ldx,ldy,ldz,ndat1,gbound_k,-1,fofr,fg_box)
66 :
67 80 : call cg_box2gsph(nx,ny,nz,ldx,ldy,ldz,ndat1,npw_k,kg_k,fg_box,fofgout,rscal=one/(nx*ny*nz))
68 : end if
69 :
70 3310966 : do ipw=1,npw_k ! Have to convert from REAL to CMPLX
71 3310925 : ug(ipw+gspad) = DCMPLX(fofgout(1,ipw), fofgout(2,ipw))
72 : end do
73 : !
74 : end do ! idat
75 :
76 41 : ABI_FREE(fofgout)
77 41 : ABI_FREE(fofr)
78 41 : ABI_SFREE(fg_box)
79 :
80 : case default
81 780045 : ABI_ERROR(sjoin("Wrong fftalga:", itoa(fftalga)))
82 : end select
83 :
84 : end block
|