Branch data 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 : 12223 : real(dp),allocatable :: arr(:,:,:,:),ftarr(:,:,:,:)
10 : : ! *************************************************************************
11 : :
12 : 12223 : nx = plan%dims(1); ny=plan%dims(2); nz=plan%dims(3)
13 : 12223 : ldx = plan%embed(1); ldy=plan%embed(2); ldz=plan%embed(3)
14 : 12223 : fftalga = plan%fftalg/100
15 : :
16 : : ! Cpu version
17 : 0 : select case (fftalga)
18 : : case (FFT_FFTW3)
19 : 0 : call fftw3_c2c_ip(nx, ny, nz, ldx, ldy, ldz, ndat__, iscale__, isign, ff)
20 : :
21 : : case (FFT_DFTI)
22 : 12145 : call dfti_c2c_ip(nx, ny, nz, ldx, ldy, ldz, ndat__, iscale__, isign, ff)
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) transform is out-of-place here
28 : : ! 3) Cannot use [ZC]copy: this is a template used both for single and double precision.
29 : :
30 : : ! MG: Well, now one can use, intrinsic :: iso_c_binding for dp version but there's no point in optimizing
31 : : ! this part as fftw3/dfti are much faster.
32 : :
33 [ - + + - : 390 : ABI_MALLOC(arr, (2, ldx, ldy, ldz))
- + + - -
+ + - + -
+ - - + -
+ - + + -
- + + - -
+ + - + -
+ - - + -
+ ]
34 [ - + + - : 312 : ABI_MALLOC(ftarr, (2, ldx, ldy, ldz))
- + + - -
+ + - - +
+ - - + +
- - + + -
- + + - -
+ + - ]
35 : :
36 [ + + + + ]: 228 : do dat=1,ndat__
37 : 150 : padat = (dat-1) * plan%ldxyz
38 : :
39 : : ! Copy input data in arr
40 [ + + + + ]: 4650 : do iz=1,nz
41 [ + + + + ]: 280050 : do iy=1,ny
42 [ + + + + ]: 22538700 : do ix=1,nx
43 : 22258800 : ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
44 : 22258800 : arr(1,ix,iy,iz) = DBLE (ff(ifft))
45 : 22534200 : arr(2,ix,iy,iz) = AIMAG(ff(ifft))
46 : : end do
47 : : end do
48 : : end do
49 : :
50 : : ! c2c with ndat = 1
51 : 150 : call sg_fft_cc(plan%fftcache, nx, ny, nz, ldx, ldy, ldz, 1, isign, arr, ftarr)
52 : : !
53 : : ! Copy results stored in ftarr
54 : : !
55 [ + + + + ]: 228 : if (isign == -1) then
56 : : ! Copy and scale the transform
57 : 80 : fact = one/plan%nfft
58 [ - + + + ]: 80 : if (iscale__ == 0) fact = one
59 [ + + + + ]: 2780 : do iz=1,nz
60 [ + + + + ]: 180980 : do iy=1,ny
61 [ + + + + ]: 14955300 : do ix=1,nx
62 : 14774400 : ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
63 : 14952600 : ff(ifft) = fact * DCMPLX(ftarr(1,ix,iy,iz), ftarr(2,ix,iy,iz))
64 : : end do
65 : : end do
66 : : end do
67 : : !
68 : : else
69 : : ! Direct copy.
70 [ + + + + ]: 1870 : do iz=1,nz
71 [ + + + + ]: 99070 : do iy=1,ny
72 [ + + + + ]: 7583400 : do ix=1,nx
73 : 7484400 : ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy + padat
74 : 7581600 : ff(ifft) = DCMPLX(ftarr(1,ix,iy,iz), ftarr(2,ix,iy,iz))
75 : : end do
76 : : end do
77 : : end do
78 : : end if
79 : :
80 : : end do
81 : :
82 [ - + - + ]: 78 : ABI_FREE(arr)
83 [ - + - + ]: 78 : ABI_FREE(ftarr)
84 : :
85 : : case default
86 : 0 : write(msg,'(a,i0)')"Wrong fftalga: ",fftalga
87 [ - + + - : 12223 : ABI_ERROR(msg)
- + + - ]
88 : : end select
89 : :
90 : : end block
|
