TABLE OF CONTENTS


ABINIT/m_dfti [ Modules ]

[ Top ] [ Modules ]

NAME

 m_dfti

FUNCTION

  This module provides wrappers for the MKL DFTI routines: in-place and out-of-place version.

COPYRIGHT

 Copyright (C) 2009-2018 ABINIT group (MG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

NOTES

  1) MPI parallelism is not supported
  2) For better performance the FFT divisions should contain small factors  (/2, 3, 5, 7, 11, 13/)
     see http://software.intel.com/sites/products/documentation/doclib/mkl_sa/11/mkl_userguide_lnx/index.htm

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 #ifdef HAVE_FFT_DFTI
29 
30 ! Include and generate MKL_DFTI module
31 #include "mkl_dfti.f90"
32 
33 ! Macros for template files.
34 #define FFTLIB "DFTI"
35 #define FFT_PREF(name) CONCAT(dfti_,name)
36 #define SPAWN_THREADS_HERE(ndat, nthreads) dfti_spawn_threads_here(ndat, nthreads)
37 
38 #endif
39 
40 MODULE m_dfti
41 
42  use defs_basis
43  use m_abicore
44  use m_errors
45  use m_xomp
46  use m_cgtools
47  use m_cplxtools
48  use m_fftcore
49  use iso_c_binding
50 #ifdef HAVE_FFT_DFTI
51  use MKL_DFTI
52 #endif
53 
54  use m_fstrings,  only : basename, strcat, int2char10
55  use m_hide_blas, only : xcopy
56  use m_fft_mesh,  only : zpad_t, zpad_init, zpad_free
57 
58  implicit none
59 
60  private
61 
62 ! Entry points for client code
63  public :: dfti_seqfourdp         ! 3D FFT of lengths nx, ny, nz. Mainly used for densities or potentials.
64  public :: dfti_seqfourwf      ! FFT transform of wavefunctions (high-level interface).
65  public :: dfti_fftrisc
66  public :: dfti_fftug          ! G-->R, 3D zero-padded FFT of lengths nx, ny, nz. Mainly used for wavefunctions
67  public :: dfti_fftur          ! R-->G, 3D zero-padded FFT of lengths nx, ny, nz. Mainly used for wavefunctions
68 
69 ! Low-level routines.
70  public :: dfti_r2c_op         ! Real to complex transform (out-of-place version).
71  public :: dfti_c2r_op         ! Complex to real transform (out-of-place version).
72  public :: dfti_c2c_op         ! complex to complex transform (out-of-place version).
73  public :: dfti_c2c_ip         ! complex to complex transform (in-place version).
74  public :: dfti_many_dft_op    ! Driver routine for many out-of-place 3D complex-to-complex FFTs.
75  public :: dfti_many_dft_ip    ! Driver routine for many in-place 3D complex-to-complex FFTs.
76  public :: dfti_fftpad         ! Driver routines for zero-padded FFT of wavefunctions.
77 
78  !FIXME I don't know why gcc does not recognize this one
79  ! Perhaps I have to provide an interfaces for fofr(:,:,:,:)
80  public :: dfti_fftpad_dp      ! Driver routines for zero-padded FFT of wavefunctions.
81  public :: dfti_fftug_dp       ! Driver routines for zero-padded FFT of wavefunctions.
82  public :: dfti_use_lib_threads

m_dfti/dfti_alloc_complex_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_alloc_complex_dpc

FUNCTION

INPUTS

OUTPUT

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2546 #ifdef HAVE_FFT_DFTI
2547 
2548 subroutine dfti_alloc_complex_dpc(size,cptr,fptr)
2549 
2550 
2551 !This section has been created automatically by the script Abilint (TD).
2552 !Do not modify the following lines by hand.
2553 #undef ABI_FUNC
2554 #define ABI_FUNC 'dfti_alloc_complex_dpc'
2555 !End of the abilint section
2556 
2557  implicit none
2558 
2559 !Arguments ------------------------------------
2560 !scalars
2561  integer,intent(in) :: size
2562  complex(dpc),ABI_CONTIGUOUS pointer :: fptr(:)
2563  type(C_PTR),intent(out) :: cptr
2564 !arrays
2565 
2566 ! *************************************************************************
2567 
2568  cptr = mkl_malloc( INT(2*size*C_DOUBLE, KIND=C_SIZE_T), DFTI_DEFAULT_ALIGNMENT_DP)
2569  if (.not. C_ASSOCIATED(cptr)) then
2570    MSG_ERROR("mkl_malloc returned NULL!")
2571  end if
2572 
2573  call c_f_pointer(cptr, fptr, [size])
2574 
2575 end subroutine dfti_alloc_complex_dpc

m_dfti/dfti_alloc_complex_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_alloc_complex_spc

FUNCTION

INPUTS

OUTPUT

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2493 #ifdef HAVE_FFT_DFTI
2494 
2495 subroutine dfti_alloc_complex_spc(size,cptr,fptr)
2496 
2497 
2498 !This section has been created automatically by the script Abilint (TD).
2499 !Do not modify the following lines by hand.
2500 #undef ABI_FUNC
2501 #define ABI_FUNC 'dfti_alloc_complex_spc'
2502 !End of the abilint section
2503 
2504  implicit none
2505 
2506 !Arguments ------------------------------------
2507 !scalars
2508  integer,intent(in) :: size
2509  complex(spc),ABI_CONTIGUOUS pointer :: fptr(:)
2510  type(C_PTR),intent(out) :: cptr
2511 !arrays
2512 
2513 ! *************************************************************************
2514 
2515  cptr = mkl_malloc( INT(2*size*C_FLOAT, KIND=C_SIZE_T), DFTI_DEFAULT_ALIGNMENT_SP)
2516  if (.not. C_ASSOCIATED(cptr)) then
2517    MSG_ERROR("mkl_malloc returned NULL!")
2518  end if
2519 
2520  call c_f_pointer(cptr, fptr, [size])
2521 
2522 end subroutine dfti_alloc_complex_spc

m_dfti/dfti_alloc_real_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_alloc_real_dp

FUNCTION

INPUTS

OUTPUT

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2440 #ifdef HAVE_FFT_DFTI
2441 
2442 subroutine dfti_alloc_real_dp(size,cptr,fptr)
2443 
2444 
2445 !This section has been created automatically by the script Abilint (TD).
2446 !Do not modify the following lines by hand.
2447 #undef ABI_FUNC
2448 #define ABI_FUNC 'dfti_alloc_real_dp'
2449 !End of the abilint section
2450 
2451  implicit none
2452 
2453 !Arguments ------------------------------------
2454 !scalars
2455  integer,intent(in) :: size
2456  real(dp),ABI_CONTIGUOUS pointer :: fptr(:)
2457  type(C_PTR),intent(out) :: cptr
2458 !arrays
2459 
2460 ! *************************************************************************
2461 
2462  cptr = mkl_malloc( INT(size*C_DOUBLE, KIND=C_SIZE_T), DFTI_DEFAULT_ALIGNMENT_DP)
2463  if (.not. C_ASSOCIATED(cptr)) then
2464    MSG_ERROR("mkl_malloc returned NULL!")
2465  end if
2466 
2467  call c_f_pointer(cptr, fptr, [size])
2468 
2469 end subroutine dfti_alloc_real_dp

m_dfti/dfti_c2c_ip_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2c_ip_dpc

FUNCTION

 Driver routine for in-place 3D complex-complex FFT. TARGET: DPC arrays

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the array.
 ndat=Number of FFTs to be done.
 isign= +1 : ff(G) => ff(R); -1 : ff(R) => ff(G)

SIDE EFFECTS

  ff(ldx*ldy*ldz*ndat)=
    In input: the complex array to be transformed.
    In output: the Fourier transformed in the space specified by isign.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1369 subroutine dfti_c2c_ip_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,isign,ff)
1370 
1371 
1372 !This section has been created automatically by the script Abilint (TD).
1373 !Do not modify the following lines by hand.
1374 #undef ABI_FUNC
1375 #define ABI_FUNC 'dfti_c2c_ip_dpc'
1376 !End of the abilint section
1377 
1378  implicit none
1379 
1380 !Arguments ------------------------------------
1381 !scalars
1382  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,isign
1383 !arrays
1384  complex(dpc),intent(inout) :: ff(ldx*ldy*ldz*ndat)
1385 
1386 ! *************************************************************************
1387 
1388 ! Include Fortran template
1389 #undef DEV_DFTI_PRECISION
1390 #define DEV_DFTI_PRECISION DFTI_DOUBLE
1391 #include "dfti_c2c_ip.finc"
1392 
1393 end subroutine dfti_c2c_ip_dpc

m_dfti/dfti_c2c_ip_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2c_ip_spc

FUNCTION

 Driver routine for in-place 3D complex-complex FFT. TARGET: SPC arrays

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the array.
 ndat=Number of FFTs to be done.
 isign= +1 : ff(G) => ff(R); -1 : ff(R) => ff(G)

SIDE EFFECTS

  ff(ldx*ldy*ldz*ndat)=
    In input: the complex array to be transformed.
    In output: the Fourier transform in the space specified by isign.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1315 subroutine dfti_c2c_ip_spc(nx,ny,nz,ldx,ldy,ldz,ndat,isign,ff)
1316 
1317 
1318 !This section has been created automatically by the script Abilint (TD).
1319 !Do not modify the following lines by hand.
1320 #undef ABI_FUNC
1321 #define ABI_FUNC 'dfti_c2c_ip_spc'
1322 !End of the abilint section
1323 
1324  implicit none
1325 
1326 !Arguments ------------------------------------
1327 !scalars
1328  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,isign
1329 !arrays
1330  complex(spc),intent(inout) :: ff(ldx*ldy*ldz*ndat)
1331 
1332 ! *************************************************************************
1333 
1334 ! Include Fortran template
1335 #undef DEV_DFTI_PRECISION
1336 #define DEV_DFTI_PRECISION DFTI_SINGLE
1337 #include "dfti_c2c_ip.finc"
1338 
1339 end subroutine dfti_c2c_ip_spc

m_dfti/dfti_c2c_op_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2c_op_dpc

FUNCTION

 Driver routine for out-of-place 3D complex-complex FFT of lengths nx, ny, nz.
 TARGET: DPC arrays

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the array.
 ndat=Number of FFTs to be done.
 isign= +1 : ff(G) => gg(R); -1 : ff(R) => gg(G)
 ff(ldx*ldy*ldz*ndat)=The array to be transformed.

OUTPUT

   gg(ldx*ldy*ldz*ndat)=The FFT of ff.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1478 subroutine dfti_c2c_op_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,isign,ff,gg)
1479 
1480 
1481 !This section has been created automatically by the script Abilint (TD).
1482 !Do not modify the following lines by hand.
1483 #undef ABI_FUNC
1484 #define ABI_FUNC 'dfti_c2c_op_dpc'
1485 !End of the abilint section
1486 
1487  implicit none
1488 
1489 !Arguments ------------------------------------
1490 !scalars
1491  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,isign,ndat
1492 !arrays
1493  complex(dpc),intent(in) :: ff(ldx*ldy*ldz*ndat)
1494  complex(dpc),intent(out) :: gg(ldx*ldy*ldz*ndat)
1495 
1496 ! *************************************************************************
1497 
1498 ! Include Fortran template
1499 #undef DEV_DFTI_PRECISION
1500 #define DEV_DFTI_PRECISION DFTI_DOUBLE
1501 #include "dfti_c2c_op.finc"
1502 
1503 end subroutine dfti_c2c_op_dpc

m_dfti/dfti_c2c_op_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2c_op_spc

FUNCTION

 Driver routine for out-of-place 3D complex-complex FFT of lengths nx, ny, nz.
 TARGET: spc arrays

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the array.
 ndat=Number of FFTs to be done.
 isign= +1 : ff(G) => gg(R); -1 : ff(R) => gg(G)
 ff(ldx*ldy*ldz*ndat)=The array to be transformed.

OUTPUT

   gg(ldx*ldy*ldz*ndat)=The FFT of ff.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1423 subroutine dfti_c2c_op_spc(nx,ny,nz,ldx,ldy,ldz,ndat,isign,ff,gg)
1424 
1425 
1426 !This section has been created automatically by the script Abilint (TD).
1427 !Do not modify the following lines by hand.
1428 #undef ABI_FUNC
1429 #define ABI_FUNC 'dfti_c2c_op_spc'
1430 !End of the abilint section
1431 
1432  implicit none
1433 
1434 !Arguments ------------------------------------
1435 !scalars
1436  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,isign,ndat
1437 !arrays
1438  complex(spc),intent(in) :: ff(ldx*ldy*ldz*ndat)
1439  complex(spc),intent(out) :: gg(ldx*ldy*ldz*ndat)
1440 
1441 ! *************************************************************************
1442 
1443 ! Include Fortran template
1444 #undef DEV_DFTI_PRECISION
1445 #define DEV_DFTI_PRECISION DFTI_SINGLE
1446 #include "dfti_c2c_op.finc"
1447 
1448 end subroutine dfti_c2c_op_spc

m_dfti/dfti_c2r_op_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2r_op_dp

FUNCTION

 Driver routine for out-of-place 3D complex-to-real FFT of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Physical dimension of the f array (to avoid cache conflicts).
 ndat=Number of FFTs to be done.
 ff(2,ldx*ldy*ldz*ndat)=The complex array to be transformed.

OUTPUT

 gg(ldx*ldy*ldz*ndat)=The backwards real FFT of ff.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2217 subroutine dfti_c2r_op_dp(nx,ny,nz,ldx,ldy,ldz,ndat,ff,gg)
2218 
2219 
2220 !This section has been created automatically by the script Abilint (TD).
2221 !Do not modify the following lines by hand.
2222 #undef ABI_FUNC
2223 #define ABI_FUNC 'dfti_c2r_op_dp'
2224 !End of the abilint section
2225 
2226  implicit none
2227 
2228 !Arguments ------------------------------------
2229 !scalars
2230  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
2231 !arrays
2232  real(dp),target,intent(in) :: ff(2*ldx*ldy*ldz*ndat)
2233  real(dp),intent(inout) :: gg(ldx*ldy*ldz*ndat)    !vz_i
2234 
2235 #ifdef HAVE_FFT_DFTI
2236 !Local variables-------------------------------
2237 !scalars
2238  type(C_ptr) :: ff_cptr
2239 !arrays
2240  complex(dpc),ABI_CONTIGUOUS pointer :: ff_fptr(:)
2241 
2242 ! *************************************************************************
2243 
2244  ff_cptr = C_loc(ff)
2245  call C_F_pointer(ff_cptr,ff_fptr, shape=(/ldx*ldy*ldz*ndat/))
2246 
2247  call dfti_c2r_op_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,ff_fptr,gg)
2248 
2249 #else
2250  MSG_ERROR("FFT_DFTI support not activated")
2251  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz/))
2252  ABI_UNUSED((/ff(1),gg(1)/))
2253 #endif
2254 
2255 end subroutine dfti_c2r_op_dp

m_dfti/dfti_c2r_op_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_c2r_op_dpc

FUNCTION

 Driver routine for out-of-place 3D complex-to-real FFT of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Physical dimension of the f array (to avoid cache conflicts).
 ndat=Number of FFTs to be done.
 ff(2,ldx*ldy*ldz*ndat)=The complex array to be transformed.

OUTPUT

 gg(2,ldx*ldy*ldz*ndat)=The backwards real FFT of ff.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

2106 subroutine dfti_c2r_op_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,ff,gg)
2107 
2108 
2109 !This section has been created automatically by the script Abilint (TD).
2110 !Do not modify the following lines by hand.
2111 #undef ABI_FUNC
2112 #define ABI_FUNC 'dfti_c2r_op_dpc'
2113 !End of the abilint section
2114 
2115  implicit none
2116 
2117 !Arguments ------------------------------------
2118 !scalars
2119  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
2120 !arrays
2121  complex(dpc),intent(in) :: ff(ldx*ldy*ldz*ndat)
2122  real(dp),intent(out) :: gg(ldx*ldy*ldz*ndat)
2123 
2124 #ifdef HAVE_FFT_DFTI
2125 !Local variables-------------------------------
2126 !scalars
2127  integer :: status,nhp,padx,i2,i3,igp,igf,idat,padatf,padatp,ii
2128  type(DFTI_DESCRIPTOR),pointer :: Desc
2129  type(C_PTR) :: cptr
2130 !arrays
2131  complex(dpc),ABI_CONTIGUOUS pointer :: ff_hp(:)
2132 
2133 ! *************************************************************************
2134 
2135  !stride  = 1
2136  !idist   = nhp
2137  !odist   = nx*ny*nz
2138  !n       = (/nx,ny,nz/)
2139  !inembed = (/(nx/2+1),ny,nz/)
2140  !onembed = (/nx,ny,nz/) ! check this
2141  !my_plan = retrieve_plan3(n,ndat,inembed,stride,idist,onembed,stride,odist,FFTW_BACKWARD,my_flags,Saved_plans)
2142 
2143  ! Fill the Hermitian part: Hermitian redundancy: out[i] is the conjugate of out[n-i]
2144  padx    = (nx/2+1)
2145  nhp     = padx*ny*nz
2146 
2147  call dfti_alloc_complex(nhp*ndat,cptr,ff_hp)
2148 
2149 !!$OMP PARALLEL DO PRIVATE(padatf,padatp,igf,igp)
2150  do idat=1,ndat
2151    padatf=(idat-1)*ldx*ldy*ldz
2152    padatp=(idat-1)*padx*ny*nz
2153    do i3=1,nz
2154      do i2=1,ny
2155        igf = (i3-1)*ldx*ldy + (i2-1)*ldx   + padatf
2156        igp = (i3-1)*padx*ny + (i2-1)*padx  + padatp
2157        ff_hp(igp+1:igp+padx) = ff(igf+1:igf+padx)
2158      end do
2159    end do
2160  end do
2161 
2162  status = DftiCreateDescriptor(Desc, DFTI_DOUBLE, DFTI_REAL, 3, (/nx,ny,nz/) )
2163  DFTI_CHECK(status)
2164 
2165  status = DftiSetValue(Desc, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX)
2166  status = DftiSetValue(Desc, DFTI_PLACEMENT, DFTI_NOT_INPLACE )
2167  status = DftiSetValue(Desc, DFTI_INPUT_STRIDES, (/0, 1, padx, padx*ny/))
2168  status = DftiSetValue(Desc, DFTI_INPUT_DISTANCE, nhp)
2169  status = DftiSetValue(Desc, DFTI_OUTPUT_STRIDES, (/0, 1, ldx, ldx*ldy/))
2170  status = DftiSetValue(Desc, DFTI_OUTPUT_DISTANCE, ldx*ldy*ldz)
2171  status = DftiSetValue(Desc, DFTI_NUMBER_OF_TRANSFORMS, ndat)
2172  DFTI_CHECK(status)
2173 
2174  DFTI_CHECK( DftiCommitDescriptor(Desc) )
2175 
2176  DFTI_CHECK( DftiComputeBackward(Desc, ff_hp, gg) )
2177 
2178  DFTI_CHECK( DftiFreeDescriptor(Desc) )
2179 
2180  call dfti_free(cptr)
2181 
2182 #else
2183  MSG_ERROR("FFT_DFTI support not activated")
2184  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz/))
2185  ABI_UNUSED(ff(1))
2186  ABI_UNUSED(gg(1))
2187 #endif
2188 
2189 end subroutine dfti_c2r_op_dpc

m_dfti/dfti_check_status [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_check_status

FUNCTION

  Error handler for DFTI wrappers. Print error message and abort.

INPUTS

  status = status error reported by dfti.
  [file] = file name
  [line] = line number

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2279 #ifdef HAVE_FFT_DFTI
2280 
2281 subroutine dfti_check_status(status,file,line)
2282 
2283 
2284 !This section has been created automatically by the script Abilint (TD).
2285 !Do not modify the following lines by hand.
2286 #undef ABI_FUNC
2287 #define ABI_FUNC 'dfti_check_status'
2288 !End of the abilint section
2289 
2290  implicit none
2291 
2292 !Arguments ------------------------------------
2293 !scalars
2294  integer,intent(in) :: status
2295  integer,optional,intent(in) :: line
2296  character(len=*),optional,intent(in) :: file
2297 
2298 !Local variables-------------------------------
2299  integer :: f90line
2300  character(len=10) :: lnum
2301  character(len=500) :: f90name
2302  character(len=500) :: my_msg
2303  character(len=DFTI_MAX_MESSAGE_LENGTH+500) :: err_msg
2304 
2305 ! *************************************************************************
2306 
2307  if (PRESENT(line)) then
2308    f90line=line
2309  else
2310    f90line=0
2311  end if
2312  call int2char10(f90line,lnum)
2313 
2314  if (PRESENT(file)) then
2315    f90name = basename(file)
2316  else
2317    f90name='Subroutine Unknown'
2318  end if
2319 
2320  my_msg=strcat(f90name,":",lnum,":")
2321 
2322  if (status /= 0) then
2323    if (.not. DftiErrorClass(status, DFTI_NO_ERROR)) then
2324      err_msg = strcat(my_msg," Error: ",DftiErrorMessage(status))
2325      MSG_ERROR(err_msg)
2326    end if
2327  end if
2328 
2329 end subroutine dfti_check_status

m_dfti/dfti_fftpad_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_fftpad_dp

FUNCTION

  This routine transforms wavefunctions using 3D zero-padded FFTs with DFTI.
  The 3D ffts are computed only on lines and planes which have non zero elements (see zpad_init)
  FFT transform is in-place.

INPUTS

   nx,ny,nz=Logical dimensions of the FFT mesh.
   ldx,ldy,ldz=Physical dimension of the f array (to avoid cache conflicts).
   ndat=Numer of FFTs
   mgfft=MAX(nx,ny,nz), only used to dimension gbound
   isign=The sign of the transform.
   gbound(2*mgfft+8,2)= The boundaries of the basis sphere of G vectors at a given k-point.
     See sphereboundary for more info.

SIDE EFFECTS

   ff(2*ldx*ldy*ldz*ndat)=
     input: The array with the data to be transformed.
     output: The results of the FFT.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

1686 subroutine dfti_fftpad_dp(ff,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign,gbound)
1687 
1688 
1689 !This section has been created automatically by the script Abilint (TD).
1690 !Do not modify the following lines by hand.
1691 #undef ABI_FUNC
1692 #define ABI_FUNC 'dfti_fftpad_dp'
1693 !End of the abilint section
1694 
1695  implicit none
1696 
1697 !Arguments ------------------------------------
1698 !scalars
1699  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign
1700 !arrays
1701  integer,intent(in) :: gbound(2*mgfft+8,2)
1702  real(dp),target,intent(inout) :: ff(2*ldx*ldy*ldz*ndat)
1703 
1704 !Local variables-------------------------------
1705 #ifdef HAVE_FFT_DFTI
1706 !scalars
1707  type(C_ptr) :: cptr
1708 !arrays
1709  complex(dpc),ABI_CONTIGUOUS pointer :: fptr(:)
1710 
1711 ! *************************************************************************
1712 
1713  ! Associate complex fptr with real ff via the C pointer
1714  cptr = C_loc(ff)
1715  call C_F_pointer(cptr,fptr, SHAPE=(/ldx*ldy*ldz*ndat/))
1716 
1717  ! Call complex version --> a lot of boilerplate code avoided
1718  call dfti_fftpad_dpc(fptr,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign,gbound)
1719 
1720 #else
1721  MSG_ERROR("FFT_DFTI support not activated")
1722  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign/))
1723  ABI_UNUSED(gbound(1,1))
1724  ABI_UNUSED(ff(1))
1725 #endif
1726 
1727 end subroutine dfti_fftpad_dp

m_dfti/dfti_fftpad_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_fftpad_dpc

FUNCTION

  This routine transforms wavefunctions using 3D zero-padded FFTs with DFTI.
  The 3D ffts are computed only on lines and planes which have non zero elements (see zpad_init)
  FFT transform is in-place. Target: complex DPC arrays.

INPUTS

   nx,ny,nz=Logical dimensions of the FFT mesh.
   ldx,ldy,ldz=Physical dimension of the f array (to avoid cache conflicts).
   ndat=Number of FFTs.
   mgfft=MAX(nx,ny,nz), only used to dimension gbound
   isign=The sign of the transform.
   gbound(2*mgfft+8,2)= The boundaries of the basis sphere of G vectors at a given k-point.
     See sphereboundary for more info.

SIDE EFFECTS

   ff(ldx*ldy*ldz*ndat)=
     input: The array with the data to be transformed.
     output: The results of the FFT.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

1763 subroutine dfti_fftpad_dpc(ff,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign,gbound)
1764 
1765 
1766 !This section has been created automatically by the script Abilint (TD).
1767 !Do not modify the following lines by hand.
1768 #undef ABI_FUNC
1769 #define ABI_FUNC 'dfti_fftpad_dpc'
1770 !End of the abilint section
1771 
1772  implicit none
1773 
1774 !Arguments ------------------------------------
1775 !scalars
1776  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign
1777 !arrays
1778  integer,intent(in) :: gbound(2*mgfft+8,2)
1779  complex(dpc),intent(inout) :: ff(ldx*ldy*ldz*ndat)
1780 
1781 !Local variables-------------------------------
1782 !scalars
1783 #ifdef HAVE_FFT_DFTI
1784 
1785 ! *************************************************************************
1786 
1787 ! Include Fortran template
1788 #undef DEV_DFTI_PRECISION
1789 #define DEV_DFTI_PRECISION DFTI_DOUBLE
1790 
1791 #include "dfti_fftpad.finc"
1792 
1793 #else
1794  MSG_ERROR("FFT_DFTI support not activated")
1795  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign/))
1796  ABI_UNUSED(gbound(1,1))
1797  ABI_UNUSED(ff(1))
1798 #endif
1799 
1800 end subroutine dfti_fftpad_dpc

m_dfti/dfti_fftpad_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_fftpad_spc

FUNCTION

  This routine transforms wavefunctions using 3D zero-padded FFTs with DFTI.
  The 3D ffts are computed only on lines and planes which have non zero elements (see zpad_init)
  FFT transform is in-place. Target: complex SPC arrays.

INPUTS

   nx,ny,nz=Logical dimensions of the FFT mesh.
   ldx,ldy,ldz=Physical dimension of the f array (to avoid cache conflicts).
   ndat=Number of FFTs.
   mgfft=MAX(nx,ny,nz), only used to dimension gbound
   isign=The sign of the transform.
   gbound(2*mgfft+8,2)= The boundaries of the basis sphere of G vectors at a given k-point.
     See sphereboundary for more info.

SIDE EFFECTS

   ff(ldx*ldy*ldz*ndat)=
     input: The array with the data to be transformed.
     output: The results of the FFT.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1835 subroutine dfti_fftpad_spc(ff,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign,gbound)
1836 
1837 
1838 !This section has been created automatically by the script Abilint (TD).
1839 !Do not modify the following lines by hand.
1840 #undef ABI_FUNC
1841 #define ABI_FUNC 'dfti_fftpad_spc'
1842 !End of the abilint section
1843 
1844  implicit none
1845 
1846 !Arguments ------------------------------------
1847 !scalars
1848  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign
1849 !arrays
1850  integer,intent(in) :: gbound(2*mgfft+8,2)
1851  complex(spc),intent(inout) :: ff(ldx*ldy*ldz*ndat)
1852 
1853 !Local variables-------------------------------
1854 !scalars
1855 #ifdef HAVE_FFT_DFTI
1856 
1857 ! *************************************************************************
1858 
1859 ! Include Fortran template
1860 #undef DEV_DFTI_PRECISION
1861 #define DEV_DFTI_PRECISION DFTI_SINGLE
1862 
1863 #include "dfti_fftpad.finc"
1864 
1865 #else
1866  MSG_ERROR("FFT_DFTI support not activated")
1867  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,isign/))
1868  ABI_UNUSED(gbound(1,1))
1869  ABI_UNUSED(ff(1))
1870 #endif
1871 
1872 end subroutine dfti_fftpad_spc

m_dfti/dfti_fftrisc_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftrisc_dp

FUNCTION

 Carry out Fourier transforms between real and reciprocal (G) space,
 for wavefunctions, contained in a sphere in reciprocal space,
 in both directions. Also accomplish some post-processing.

NOTES

 Specifically uses rather sophisticated algorithms, based on S Goedecker
 routines, specialized for superscalar RISC architecture.
 Zero padding : saves 7/12 execution time
 Bi-dimensional data locality in most of the routine : cache reuse
 For k-point (0 0 0) : takes advantage of symmetry of data.
 Note however that no blocking is used, in both 1D z-transform
 or subsequent 2D transform. This should be improved.

INPUTS

  cplex= if 1 , denpot is real, if 2 , denpot is complex
     (cplex=2 only allowed for option=2 when istwf_k=1)
     one can also use cplex=0 if option=0 or option=3
  fofgin(2,npwin)=holds input wavefunction in G vector basis sphere.
  gboundin(2*mgfft+8,2)=sphere boundary info for reciprocal to real space
  gboundout(2*mgfft+8,2)=sphere boundary info for real to reciprocal space
  istwf_k=option parameter that describes the storage of wfs
  kg_kin(3,npwin)=reduced planewave coordinates, input
  kg_kout(3,npwout)=reduced planewave coordinates, output
  mgfft=maximum size of 1D FFTs
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  npwin=number of elements in fofgin array (for option 0, 1 and 2)
  npwout=number of elements in fofgout array (for option 2 and 3)
  ldx,ldy,ldz=ngfft(4),ngfft(5),ngfft(6), dimensions of fofr.
  option= if 0: do direct FFT
          if 1: do direct FFT, then sum the density
          if 2: do direct FFT, multiply by the potential, then do reverse FFT
          if 3: do reverse FFT only
  weight=weight to be used for the accumulation of the density in real space
          (needed only when option=1)

OUTPUT

  (see side effects)

OPTIONS

  The different options are:
  - reciprocal to real space and output the result (when option=0),
  - reciprocal to real space and accumulate the density (when option=1) or
  - reciprocal to real space, apply the local potential to the wavefunction
    in real space and produce the result in reciprocal space (when option=2)
  - real space to reciprocal space (when option=3).
  option=0 IS NOT ALLOWED when istwf_k>2
  option=3 IS NOT ALLOWED when istwf_k>=2

SIDE EFFECTS

  for option==0, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 fofr(2,ldx,ldy,ldz) contains the Fourier Transform of fofgin;
                 no use of denpot, fofgout and npwout.
  for option==1, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 denpot(cplex*ldx,ldy,ldz) contains the input density at input,
                 and the updated density at output;
                 no use of fofgout and npwout.
  for option==2, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 denpot(cplex*ldx,ldy,ldz) contains the input local potential;
                 fofgout(2,npwout) contains the output function;
  for option==3, fofr(2,ldx,ldy,ldz) contains the real space wavefunction;
                 fofgout(2,npwout) contains its Fourier transform;
                 no use of fofgin and npwin.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

753 subroutine dfti_fftrisc_dp(cplex,denpot,fofgin,fofgout,fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
754 & mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
755 
756 
757 !This section has been created automatically by the script Abilint (TD).
758 !Do not modify the following lines by hand.
759 #undef ABI_FUNC
760 #define ABI_FUNC 'dfti_fftrisc_dp'
761 !End of the abilint section
762 
763  implicit none
764 
765 !Arguments ------------------------------------
766 !scalars
767  integer,intent(in) :: cplex,istwf_k,mgfft,ldx,ldy,ldz,npwin,npwout,option
768  real(dp),intent(in) :: weight_i,weight_r
769 !arrays
770  integer,intent(in) :: gboundin(2*mgfft+8,2),gboundout(2*mgfft+8,2)
771  integer,intent(in) :: kg_kin(3,npwin),kg_kout(3,npwout),ngfft(18)
772  real(dp),intent(in) :: fofgin(2,npwin)
773  real(dp),intent(inout) :: denpot(cplex*ldx,ldy,ldz)
774  real(dp),intent(inout) :: fofr(2,ldx*ldy*ldz)
775  real(dp),intent(inout) :: fofgout(2,npwout)    !vz_i
776 
777 ! *************************************************************************
778 
779 #ifdef HAVE_FFT_DFTI
780 
781 #undef  FFT_PRECISION
782 #undef  MYKIND
783 #undef  MYCZERO
784 #undef  MYCMPLX
785 #undef  MYCONJG
786 
787 #define FFT_PRECISION DFTI_DOUBLE
788 #define MYKIND DPC
789 #define MYCZERO (0._dp,0._dp)
790 #define MYCMPLX  DCMPLX
791 #define MYCONJG  DCONJG
792 
793 #include "dfti_fftrisc.finc"
794 
795 #else
796  MSG_ERROR("DFTI support not activated")
797  ABI_UNUSED((/cplex,gboundin(1,1),gboundout(1,1),istwf_k,kg_kin(1,1),kg_kout(1,1)/))
798  ABI_UNUSED((/mgfft,ngfft(1),npwin,npwout,ldx,ldy,ldz,option/))
799  ABI_UNUSED((/denpot(1,1,1),fofgin(1,1),fofgout(1,1),fofr(1,1),weight_r,weight_i/))
800 #endif
801 
802 end subroutine dfti_fftrisc_dp

m_dfti/dfti_fftrisc_sp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftrisc_sp

FUNCTION

 Carry out Fourier transforms between real and reciprocal (G) space,
 for wavefunctions, contained in a sphere in reciprocal space,
 in both directions. Also accomplish some post-processing.

NOTES

 Specifically uses rather sophisticated algorithms, based on S Goedecker
 routines, specialized for superscalar RISC architecture.
 Zero padding : saves 7/12 execution time
 Bi-dimensional data locality in most of the routine : cache reuse
 For k-point (0 0 0) : takes advantage of symmetry of data.
 Note however that no blocking is used, in both 1D z-transform
 or subsequent 2D transform. This should be improved.

INPUTS

  cplex= if 1 , denpot is real, if 2 , denpot is complex
     (cplex=2 only allowed for option=2 when istwf_k=1)
     one can also use cplex=0 if option=0 or option=3
  fofgin(2,npwin)=holds input wavefunction in G vector basis sphere.
  gboundin(2*mgfft+8,2)=sphere boundary info for reciprocal to real space
  gboundout(2*mgfft+8,2)=sphere boundary info for real to reciprocal space
  istwf_k=option parameter that describes the storage of wfs
  kg_kin(3,npwin)=reduced planewave coordinates, input
  kg_kout(3,npwout)=reduced planewave coordinates, output
  mgfft=maximum size of 1D FFTs
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  npwin=number of elements in fofgin array (for option 0, 1 and 2)
  npwout=number of elements in fofgout array (for option 2 and 3)
  ldx,ldy,ldz=ngfft(4),ngfft(5),ngfft(6), dimensions of fofr.
  option= if 0: do direct FFT
          if 1: do direct FFT, then sum the density
          if 2: do direct FFT, multiply by the potential, then do reverse FFT
          if 3: do reverse FFT only
  weight=weight to be used for the accumulation of the density in real space
          (needed only when option=1)

OUTPUT

  (see side effects)

OPTIONS

  The different options are:
  - reciprocal to real space and output the result (when option=0),
  - reciprocal to real space and accumulate the density (when option=1) or
  - reciprocal to real space, apply the local potential to the wavefunction
    in real space and produce the result in reciprocal space (when option=2)
  - real space to reciprocal space (when option=3).
  option=0 IS NOT ALLOWED when istwf_k>2
  option=3 IS NOT ALLOWED when istwf_k>=2

SIDE EFFECTS

  for option==0, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 fofr(2,ldx,ldy,ldz) contains the Fourier Transform of fofgin;
                 no use of denpot, fofgout and npwout.
  for option==1, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 denpot(cplex*ldx,ldy,ldz) contains the input density at input,
                 and the updated density at output;
                 no use of fofgout and npwout.
  for option==2, fofgin(2,npwin)=holds input wavefunction in G sphere;
                 denpot(cplex*ldx,ldy,ldz) contains the input local potential;
                 fofgout(2,npwout) contains the output function;
  for option==3, fofr(2,ldx,ldy,ldz) contains the real space wavefunction;
                 fofgout(2,npwout) contains its Fourier transform;
                 no use of fofgin and npwin.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

622 subroutine dfti_fftrisc_sp(cplex,denpot,fofgin,fofgout,fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
623 & mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
624 
625 
626 !This section has been created automatically by the script Abilint (TD).
627 !Do not modify the following lines by hand.
628 #undef ABI_FUNC
629 #define ABI_FUNC 'dfti_fftrisc_sp'
630 !End of the abilint section
631 
632  implicit none
633 
634 !Arguments ------------------------------------
635 !scalars
636  integer,intent(in) :: cplex,istwf_k,mgfft,ldx,ldy,ldz,npwin,npwout,option
637  real(dp),intent(in) :: weight_i,weight_r
638 !arrays
639  integer,intent(in) :: gboundin(2*mgfft+8,2),gboundout(2*mgfft+8,2)
640  integer,intent(in) :: kg_kin(3,npwin),kg_kout(3,npwout),ngfft(18)
641  real(sp),intent(in) :: fofgin(2,npwin)
642  real(dp),intent(inout) :: denpot(cplex*ldx,ldy,ldz)
643  real(sp),intent(inout) :: fofr(2,ldx*ldy*ldz)
644  real(sp),intent(inout) :: fofgout(2,npwout)    !vz_i
645 
646 ! *************************************************************************
647 
648 #ifdef HAVE_FFT_DFTI
649 
650 #undef FFT_PRECISION
651 #undef MYKIND
652 #undef MYCZERO
653 #undef MYCMPLX
654 #undef MYCONJG
655 
656 #define FFT_PRECISION DFTI_SINGLE
657 #define MYKIND SPC
658 #define MYCZERO (0._sp,0._sp)
659 #define MYCMPLX  CMPLX
660 #define MYCONJG  CONJG
661 
662 #include "dfti_fftrisc.finc"
663 
664 #else
665  MSG_ERROR("DFTI support not activated")
666  ABI_UNUSED((/cplex,gboundin(1,1),gboundout(1,1),istwf_k,kg_kin(1,1),kg_kout(1,1)/))
667  ABI_UNUSED((/mgfft,ngfft(1),npwin,npwout,ldx,ldy,ldz,option/))
668  ABI_UNUSED((/denpot(1,1,1),weight_r,weight_i/))
669  ABI_UNUSED((/fofgin(1,1),fofgout(1,1),fofr(1,1)/))
670 #endif
671 
672 end subroutine dfti_fftrisc_sp

m_dfti/dfti_fftug_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftug_dp

FUNCTION

 Compute ndat zero-padded FFTs from G to R space.
 Mainly used for the transform of wavefunctions.
 TARGET: dp arrays with real and imaginary part

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.
  ug(npw_k*ndat)=wavefunctions in reciprocal space.

OUTPUT

  ur(ldx*ldy*ldz*ndat)=wavefunctions in real space.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

840 subroutine dfti_fftug_dp(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ug,ur)
841 
842 
843 !This section has been created automatically by the script Abilint (TD).
844 !Do not modify the following lines by hand.
845 #undef ABI_FUNC
846 #define ABI_FUNC 'dfti_fftug_dp'
847 !End of the abilint section
848 
849  implicit none
850 
851 !Arguments ------------------------------------
852 !scalars
853  integer,intent(in) :: fftalg,fftcache
854  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
855 !arrays
856  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
857  real(dp),target,intent(in) :: ug(2*npw_k*ndat)
858  real(dp),target,intent(inout) :: ur(2*ldx*ldy*ldz*ndat)
859 
860 #ifdef HAVE_FFT_DFTI
861 !Local variables-------------------------------
862 !scalars
863  integer,parameter :: dist=2
864  real(dp) :: fofgout(2,0)
865  real(dp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
866 
867 ! *************************************************************************
868 
869 #undef TK_PREF
870 #define TK_PREF(name) CONCAT(cg_,name)
871 
872 #include "fftug.finc"
873 
874 #else
875  ! Silence compiler warning
876  MSG_ERROR("FFT_DFTI support not activated")
877  ABI_UNUSED((/fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
878  ABI_UNUSED((/ug(1),ur(1)/))
879 #endif
880 
881 end subroutine dfti_fftug_dp

m_dfti/dfti_fftug_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftug_dpc

FUNCTION

 Compute ndat zero-padded FFTs from G ro R.
 Mainly used for the transform of wavefunctions.
 TARGET: DPC arrays

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.
  ug(npw_k*ndat)=wavefunctions in reciprocal space

OUTPUT

  ur(ldx*ldy*ldz*ndat)=wavefunctions in real space.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

 996 subroutine dfti_fftug_dpc(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ug,ur)
 997 
 998 
 999 !This section has been created automatically by the script Abilint (TD).
1000 !Do not modify the following lines by hand.
1001 #undef ABI_FUNC
1002 #define ABI_FUNC 'dfti_fftug_dpc'
1003 !End of the abilint section
1004 
1005  implicit none
1006 
1007 !Arguments ------------------------------------
1008 !scalars
1009  integer,intent(in) :: fftalg,fftcache
1010  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
1011 !arrays
1012  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
1013  complex(dpc),target,intent(in) :: ug(npw_k*ndat)
1014  complex(dpc),target,intent(inout) :: ur(ldx*ldy*ldz*ndat)    !vz_i
1015 
1016 #ifdef HAVE_FFT_DFTI
1017 !Local variables-------------------------------
1018 !scalars
1019  integer,parameter :: dist=1
1020 !arrays
1021  real(dp) :: fofgout(2,0)
1022  real(dp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
1023 
1024 ! *************************************************************************
1025 
1026 #undef TK_PREF
1027 #define TK_PREF(name) CONCAT(cplx_,name)
1028 
1029 #include "fftug.finc"
1030 
1031 #else
1032  ! Silence compiler warning
1033  MSG_ERROR("FFT_DFTI support not activated")
1034  ABI_UNUSED((/fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
1035  ABI_UNUSED((/ug(1),ur(1)/))
1036 #endif
1037 
1038 end subroutine dfti_fftug_dpc

m_dfti/dfti_fftug_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftug_spc

FUNCTION

 Compute ndat zero-padded FFTs from G- to R-space .
 Mainly used for the transform of wavefunctions.
 TARGET: spc arrays

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.
  ug(npw_k*ndat)=wavefunctions in reciprocal space.

OUTPUT

  ur(ldx*ldy*ldz*ndat)=wavefunctions in real space.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

918 subroutine dfti_fftug_spc(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ug,ur)
919 
920 
921 !This section has been created automatically by the script Abilint (TD).
922 !Do not modify the following lines by hand.
923 #undef ABI_FUNC
924 #define ABI_FUNC 'dfti_fftug_spc'
925 !End of the abilint section
926 
927  implicit none
928 
929 !Arguments ------------------------------------
930 !scalars
931  integer,intent(in) :: fftalg,fftcache
932  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
933 !arrays
934  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
935  complex(spc),target,intent(in) :: ug(npw_k*ndat)
936  complex(spc),target,intent(inout) :: ur(ldx*ldy*ldz*ndat)    !vz_i
937 
938 #ifdef HAVE_FFT_DFTI
939 !Local variables-------------------------------
940 !scalars
941  integer,parameter :: dist=1
942  real(sp) :: fofgout(2,0)
943  real(sp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
944 
945 ! *************************************************************************
946 
947 #undef TK_PREF
948 #define TK_PREF(name) CONCAT(cplx_,name)
949 
950 #include "fftug.finc"
951 
952 #else
953  ! Silence compiler warning
954  MSG_ERROR("FFT_DFTI support not activated")
955  ABI_UNUSED((/fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
956  ABI_UNUSED((/ug(1),ur(1)/))
957 #endif
958 
959 end subroutine dfti_fftug_spc

m_dfti/dfti_fftur_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftur_dp

FUNCTION

 Compute ndat zero-padded FFTs from R- to G-space .
 Mainly used for the transform of wavefunctions.
 TARGET: dp arrays with real and imaginary part.

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.

 SIDE EFFECT
 ur(2,ldx*ldy*ldz*ndat)= In input: wavefunctions in real space.
                         Destroyed in output. Do not use ur anymore!

OUTPUT

 ug(2,npw_k*ndat)=wavefunctions in reciprocal space.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1077 subroutine dfti_fftur_dp(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ur,ug)
1078 
1079 
1080 !This section has been created automatically by the script Abilint (TD).
1081 !Do not modify the following lines by hand.
1082 #undef ABI_FUNC
1083 #define ABI_FUNC 'dfti_fftur_dp'
1084 !End of the abilint section
1085 
1086  implicit none
1087 
1088 !Arguments ------------------------------------
1089 !scalars
1090  integer,intent(in) :: fftalg,fftcache
1091  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
1092 !arrays
1093  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
1094  real(dp),target,intent(inout) :: ur(2*ldx*ldy*ldz*ndat)
1095  real(dp),target,intent(inout) :: ug(2*npw_k*ndat)    !vz_i
1096 
1097 #ifdef HAVE_FFT_DFTI
1098 !Local variables-------------------------------
1099 !scalars
1100  integer,parameter :: dist=2
1101 !arrays
1102  real(dp) :: dum_ugin(2,0)
1103  real(dp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
1104 
1105 ! *************************************************************************
1106 
1107 #undef TK_PREF
1108 #define TK_PREF(name) CONCAT(cg_,name)
1109 
1110 #include "fftur.finc"
1111 
1112 #else
1113  ! Silence compiler warning
1114  MSG_ERROR("FFT_DFTI support not activated")
1115  ABI_UNUSED((/fftalg,fftcache/))
1116  ABI_UNUSED((/npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
1117  ABI_UNUSED((/ug(1),ur(1)/))
1118 #endif
1119 
1120 end subroutine dfti_fftur_dp

m_dfti/dfti_fftur_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftur_dpc

FUNCTION

 Compute ndat zero-padded FFTs from R ro G.
 Mainly used for the transform of wavefunctions.
 TARGET: DPC arrays

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.

 SIDE EFFECT
 ur(ldx*ldy*ldz*ndat)= In input: wavefunctions in real space.
                       Destroyed in output. Do not use ur anymore!

OUTPUT

 ug(npw_k*ndat)=wavefunctions in reciprocal space

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1242 subroutine dfti_fftur_dpc(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ur,ug)
1243 
1244 
1245 !This section has been created automatically by the script Abilint (TD).
1246 !Do not modify the following lines by hand.
1247 #undef ABI_FUNC
1248 #define ABI_FUNC 'dfti_fftur_dpc'
1249 !End of the abilint section
1250 
1251  implicit none
1252 
1253 !Arguments ------------------------------------
1254 !scalars
1255  integer,intent(in) :: fftalg,fftcache
1256  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
1257 !arrays
1258  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
1259  complex(dpc),target,intent(inout) :: ur(ldx*ldy*ldz*ndat)
1260  complex(dpc),target,intent(inout) :: ug(npw_k*ndat)    !vz_i
1261 
1262 #ifdef HAVE_FFT_DFTI
1263 !Local variables-------------------------------
1264 !scalars
1265  integer,parameter :: dist=1
1266 !arrays
1267  real(dp) :: dum_ugin(2,0)
1268  real(dp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
1269 
1270 ! *************************************************************************
1271 
1272 #undef TK_PREF
1273 #define TK_PREF(name) CONCAT(cplx_,name)
1274 
1275 #include "fftur.finc"
1276 
1277 #else
1278  ! Silence compiler warning
1279  MSG_ERROR("FFT_DFTI support not activated")
1280  ABI_UNUSED((/fftalg,fftcache/))
1281  ABI_UNUSED((/npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
1282  ABI_UNUSED((/ug(1),ur(1)/))
1283 #endif
1284 
1285 end subroutine dfti_fftur_dpc

m_dfti/dfti_fftur_spc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_fftur_spc

FUNCTION

 Compute ndat zero-padded FFTs from R- to G-space .
 Mainly used for the transform of wavefunctions.
 TARGET: spc arrays

INPUTS

 fftalg=FFT algorith (see input variable)
 fftcache=size of the cache (kB)
 npw_k=number of plane waves for this k-point.
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat=Number of transforms
 istwf_k=Option describing the storage of the wavefunction.
 mgfft=Max number of FFT divisions (used to dimension gbound)
 kg_k(3,npw_k)=G-vectors in reduced coordinates
 gbound(2*mgfft+8,2)=Table for padded-FFT. See sphereboundary.

 SIDE EFFECT
 ur(ldx*ldy*ldz*ndat)= In input: wavefunctions in real space.
                       Destroyed in output. Do not use ur anymore!

OUTPUT

 ug(npw_k*ndat)=wavefunctions in reciprocal space.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1160 subroutine dfti_fftur_spc(fftalg,fftcache,npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k,gbound,ur,ug)
1161 
1162 
1163 !This section has been created automatically by the script Abilint (TD).
1164 !Do not modify the following lines by hand.
1165 #undef ABI_FUNC
1166 #define ABI_FUNC 'dfti_fftur_spc'
1167 !End of the abilint section
1168 
1169  implicit none
1170 
1171 !Arguments ------------------------------------
1172 !scalars
1173  integer,intent(in) :: fftalg,fftcache
1174  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft
1175 !arrays
1176  integer,intent(in) :: gbound(2*mgfft+8,2),kg_k(3,npw_k)
1177  complex(spc),target,intent(inout) :: ur(ldx*ldy*ldz*ndat)
1178  complex(spc),target,intent(inout) :: ug(npw_k*ndat)    !vz_i
1179 
1180 #ifdef HAVE_FFT_DFTI
1181 !Local variables-------------------------------
1182 !scalars
1183  integer,parameter :: dist=1
1184 !arrays
1185  real(sp) :: dum_ugin(2,0)
1186  real(sp),ABI_CONTIGUOUS pointer :: real_ug(:,:),real_ur(:,:)
1187 
1188 ! *************************************************************************
1189 
1190 #undef TK_PREF
1191 #define TK_PREF(name) CONCAT(cplx_,name)
1192 
1193 #include "fftur.finc"
1194 
1195 #else
1196  ! Silence compiler warning
1197  MSG_ERROR("FFT_DFTI support not activated")
1198  ABI_UNUSED((/fftalg,fftcache/))
1199  ABI_UNUSED((/npw_k,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_k(1,1),gbound(1,1)/))
1200  ABI_UNUSED((/ug(1),ur(1)/))
1201 #endif
1202 
1203 end subroutine dfti_fftur_spc

m_dfti/dfti_many_dft_ip [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_many_dft_ip

FUNCTION

 Driver routine for many in-place 3D complex-to-complex FFTs of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimension of the finout array (to avoid cache conflicts).
 ndat=Number of FFTs to be done.
 isign=sign of Fourier transform exponent: current convention uses
   +1 for transforming from G to r,
   -1 for transforming from r to G.

OUTPUT

 finout(2,ldx*ldy*ldz*ndat)=
   In input: The complex array to be transformed.
   In output: The FFT results.

PARENTS

CHILDREN

      c_f_pointer

SOURCE

1611 subroutine dfti_many_dft_ip(nx,ny,nz,ldx,ldy,ldz,ndat,isign,finout)
1612 
1613 
1614 !This section has been created automatically by the script Abilint (TD).
1615 !Do not modify the following lines by hand.
1616 #undef ABI_FUNC
1617 #define ABI_FUNC 'dfti_many_dft_ip'
1618 !End of the abilint section
1619 
1620  implicit none
1621 
1622 !Arguments ------------------------------------
1623 !scalars
1624  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,isign
1625 !arrays
1626  real(dp),target,intent(inout) :: finout(2*ldx*ldy*ldz*ndat)
1627 
1628 #ifdef HAVE_FFT_DFTI
1629 !Local variables-------------------------------
1630 !scalars
1631  type(C_ptr) :: finout_cptr
1632 !arrays
1633  complex(dpc),ABI_CONTIGUOUS pointer :: finout_fptr(:)
1634 
1635 ! *************************************************************************
1636 
1637  ! Associate complex finout_fptr with real ffinout via the C pointer
1638  finout_cptr = C_loc(finout)
1639  call C_F_pointer(finout_cptr,finout_fptr, shape=(/ldx*ldy*ldz*ndat/))
1640 
1641  ! Call complex version --> a lot of boilerplate code avoided
1642  call dfti_c2c_ip(nx,ny,nz,ldx,ldy,ldz,ndat,isign,finout_fptr)
1643 
1644 #else
1645  MSG_ERROR("FFT_DFTI support not activated")
1646  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz,ndat,isign/))
1647  ABI_UNUSED(finout(1))
1648 #endif
1649 
1650 end subroutine dfti_many_dft_ip

m_dfti/dfti_many_dft_op [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_many_dft_op

FUNCTION

 Driver routine for many out-of-place 3D complex-to-complex FFTs of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimension of the fin and fout arrays (to avoid cache conflicts).
 ndat=Number of FFTs to be done.
 fin(2*ldx*ldy*ldz*ndat)=The complex array to be transformed.
 isign=sign of Fourier transform exponent: current convention uses
   +1 for transforming from G to r,
   -1 for transforming from r to G.

OUTPUT

 fout(2,ldx*ldy*ldz*ndat)=The Fourier transform of fin.

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

1535 subroutine dfti_many_dft_op(nx,ny,nz,ldx,ldy,ldz,ndat,isign,fin,fout)
1536 
1537 
1538 !This section has been created automatically by the script Abilint (TD).
1539 !Do not modify the following lines by hand.
1540 #undef ABI_FUNC
1541 #define ABI_FUNC 'dfti_many_dft_op'
1542 !End of the abilint section
1543 
1544  implicit none
1545 
1546 !Arguments ------------------------------------
1547 !scalars
1548  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat,isign
1549 !arrays
1550  real(dp),target,intent(in) :: fin(2*ldx*ldy*ldz*ndat)
1551  real(dp),target,intent(out) :: fout(2*ldx*ldy*ldz*ndat)
1552 
1553 #ifdef HAVE_FFT_DFTI
1554 !Local variables-------------------------------
1555 !scalars
1556  type(C_ptr) :: fin_cptr, fout_cptr
1557 !arrays
1558  complex(dpc),ABI_CONTIGUOUS pointer :: fin_fptr(:),fout_fptr(:)
1559 
1560 ! *************************************************************************
1561 
1562  ! Associate complex pointers with real inputs via the C pointers
1563  fin_cptr = C_loc(fin)
1564  call C_F_pointer(fin_cptr,fin_fptr, shape=(/ldx*ldy*ldz*ndat/))
1565 
1566  fout_cptr = C_loc(fout)
1567  call C_F_pointer(fout_cptr,fout_fptr, shape=(/ldx*ldy*ldz*ndat/))
1568 
1569  ! Call complex version --> a lot of boilerplate code avoided
1570  call dfti_c2c_op(nx,ny,nz,ldx,ldy,ldz,ndat,isign,fin_fptr,fout_fptr)
1571 
1572 #else
1573  MSG_ERROR("FFT_DFTI support not activated")
1574  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz,ndat,isign/))
1575  ABI_UNUSED(fin(1))
1576  ABI_UNUSED(fout(1))
1577 #endif
1578 
1579 end subroutine dfti_many_dft_op

m_dfti/dfti_r2c_op_dp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_r2c_op_dp

FUNCTION

 Driver routine for out-of-place 3D real-to-complex FFT of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the f array (to avoid cache conflicts).
 ndat=Number of FFTs to be done.
 ff(ldx*ldy*ldz*ndat)=The real array to be transformed.

OUTPUT

 gg(2*ldx*ldy*ldz*ndat)=The forward FFT of ff (real valued)

PARENTS

CHILDREN

      c_f_pointer

SOURCE

2038 subroutine dfti_r2c_op_dp(nx,ny,nz,ldx,ldy,ldz,ndat,ff,gg)
2039 
2040 
2041 !This section has been created automatically by the script Abilint (TD).
2042 !Do not modify the following lines by hand.
2043 #undef ABI_FUNC
2044 #define ABI_FUNC 'dfti_r2c_op_dp'
2045 !End of the abilint section
2046 
2047  implicit none
2048 
2049 !Arguments ------------------------------------
2050 !scalars
2051  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
2052 !arrays
2053  real(dp),intent(in) :: ff(ldx*ldy*ldz*ndat)
2054  real(dp),target,intent(out) :: gg(2*ldx*ldy*ldz*ndat)
2055 
2056 #ifdef HAVE_FFT_DFTI
2057 !Local variables-------------------------------
2058 !scalars
2059  type(C_ptr) :: gg_cptr
2060 !arrays
2061  complex(dpc),ABI_CONTIGUOUS pointer :: gg_fptr(:)
2062 
2063 ! *************************************************************************
2064 
2065  gg_cptr = C_loc(gg)
2066  call C_F_pointer(gg_cptr,gg_fptr, shape=(/ldx*ldy*ldz*ndat/))
2067 
2068  call dfti_r2c_op_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,ff,gg_fptr)
2069 
2070 #else
2071  MSG_ERROR("FFT_DFTI support not activated")
2072  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz/))
2073  ABI_UNUSED(ff)
2074  ABI_UNUSED(gg(1))
2075 #endif
2076 
2077 end subroutine dfti_r2c_op_dp

m_dfti/dfti_r2c_op_dpc [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_r2c_op_dpc

FUNCTION

 Driver routine for out-of-place 3D real-to-complex FFT of lengths nx, ny, nz.

INPUTS

 nx,ny,nz=Number of points along the three directions.
 ldx,ldy,ldz=Physical dimensions of the f array (to avoid cache conflicts).
 ff(ldx*ldy*ldz*ndat)=The real array to be transformed.
 ndat=Number of FFTs to be done.

OUTPUT

 gg(ldx*ldy*ldz*ndat)=The forward FFT of ff (complex valued)

PARENTS

      m_dfti

CHILDREN

      c_f_pointer

SOURCE

1901 subroutine dfti_r2c_op_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,ff,gg)
1902 
1903 
1904 !This section has been created automatically by the script Abilint (TD).
1905 !Do not modify the following lines by hand.
1906 #undef ABI_FUNC
1907 #define ABI_FUNC 'dfti_r2c_op_dpc'
1908 !End of the abilint section
1909 
1910  implicit none
1911 
1912 !Arguments ------------------------------------
1913 !scalars
1914  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
1915 !arrays
1916  real(dp),intent(in) :: ff(ldx*ldy*ldz*ndat)
1917  complex(dpc),intent(out) :: gg(ldx*ldy*ldz*ndat)
1918 
1919 #ifdef HAVE_FFT_DFTI
1920 !Local variables-------------------------------
1921 !scalars
1922  integer :: status,nhp,padx,i1,i2,i3,igp,igf,imgf,ii
1923  integer :: i1inv,i2inv,i3inv,idat,padatf
1924  type(DFTI_DESCRIPTOR),pointer :: Desc
1925  type(C_PTR) :: cptr
1926 !arrays
1927  integer,allocatable :: i1inver(:),i2inver(:),i3inver(:)
1928  complex(dpc),ABI_CONTIGUOUS pointer :: gg_hp(:)
1929 
1930 ! *************************************************************************
1931 
1932  padx = (nx/2+1)
1933  nhp = (nx/2+1)*ny*nz
1934 
1935  call dfti_alloc_complex(nhp*ndat,cptr,gg_hp)
1936 
1937  status = DftiCreateDescriptor(Desc, DFTI_DOUBLE, DFTI_REAL, 3, (/nx,ny,nz/) )
1938  DFTI_CHECK(status)
1939 
1940  status = DftiSetValue(Desc, DFTI_CONJUGATE_EVEN_STORAGE, DFTI_COMPLEX_COMPLEX)
1941  status = DftiSetValue(Desc, DFTI_PLACEMENT, DFTI_NOT_INPLACE )
1942  status = DftiSetValue(Desc, DFTI_INPUT_STRIDES,  (/0, 1, ldx,  ldx*ldy/))
1943  status = DftiSetValue(Desc, DFTI_INPUT_DISTANCE, ldx*ldy*ldz)
1944  status = DftiSetValue(Desc, DFTI_OUTPUT_STRIDES, (/0, 1, padx, padx*ny/))
1945  status = DftiSetValue(Desc, DFTI_OUTPUT_DISTANCE, nhp)
1946  status = DftiSetValue(Desc, DFTI_NUMBER_OF_TRANSFORMS, ndat)
1947  status = DftiSetValue(Desc, DFTI_FORWARD_SCALE, one / DBLE(nx*ny*nz) )
1948  DFTI_CHECK(status)
1949 
1950  DFTI_CHECK( DftiCommitDescriptor(Desc) )
1951 
1952  DFTI_CHECK( DftiComputeForward(Desc, ff, gg_hp) )
1953 
1954  DFTI_CHECK( DftiFreeDescriptor(Desc) )
1955 
1956  ! Reconstruct full FFT: Hermitian redundancy: out[i] is the conjugate of out[n-i]
1957  ABI_MALLOC(i1inver,(padx))
1958  ABI_MALLOC(i2inver,(ny))
1959  ABI_MALLOC(i3inver,(nz))
1960 
1961  i1inver(1)=1
1962  do i1=2,padx
1963    i1inver(i1)=nx+2-i1
1964  end do
1965 
1966  i2inver(1)=1
1967  do i2=2,ny
1968    i2inver(i2)=ny+2-i2
1969  end do
1970 
1971  i3inver(1)=1
1972  do i3=2,nz
1973    i3inver(i3)=nz+2-i3
1974  end do
1975 
1976  igp=0
1977  do idat=1,ndat
1978    padatf=(idat-1)*ldx*ldy*ldz
1979    do i3=1,nz
1980      i3inv = i3inver(i3)
1981      do i2=1,ny
1982        i2inv = i2inver(i2)
1983        do i1=1,padx
1984          igp=igp+1
1985          igf = i1 + (i3-1)*ldx*ldy + (i2-1)*ldx + padatf
1986          gg(igf) =  gg_hp(igp)
1987          i1inv = i1inver(i1)
1988          if (i1inv/=i1) then
1989            imgf = i1inv + (i3inv-1)*ldx*ldy + (i2inv-1)*ldx + padatf
1990            gg(imgf) = DCONJG(gg_hp(igp))
1991          end if
1992        end do
1993      end do
1994    end do
1995  end do
1996 
1997  ABI_FREE(i1inver)
1998  ABI_FREE(i2inver)
1999  ABI_FREE(i3inver)
2000 
2001  call dfti_free(cptr)
2002 
2003 #else
2004  MSG_ERROR("FFT_DFTI support not activated")
2005  ABI_UNUSED((/nx,ny,nz,ldx,ldy,ldz/))
2006  ABI_UNUSED(ff)
2007  ABI_UNUSED(gg(1))
2008 #endif
2009 
2010 end subroutine dfti_r2c_op_dpc

m_dfti/dfti_seqfourdp [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

  dfti_seqfourdp

FUNCTION

 Driver routine for 3D FFT of lengths nx, ny, nz. Mainly used for densities or potentials.
 FFT Transform is out-of-place

INPUTS

 cplex=1 if fofr is real, 2 if fofr is complex
 nx,ny,nz=Number of point along the three directions.
 ldx,ldy,ldz=Leading dimensions of the array.
 ndat = Number of FFTS
 isign= +1 : fofg(G) => fofr(R);
        -1 : fofr(R) => fofg(G)
 fofg(2,ldx*ldy*ldz*ndat)=The array to be transformed.

OUTPUT

 fofr(cplex,ldx*ldy*ldz*ndat)=The FFT of fofg

PARENTS

      fourdp,m_fft

CHILDREN

      c_f_pointer

SOURCE

206 subroutine dfti_seqfourdp(cplex,nx,ny,nz,ldx,ldy,ldz,ndat,isign,fofg,fofr)
207 
208 
209 !This section has been created automatically by the script Abilint (TD).
210 !Do not modify the following lines by hand.
211 #undef ABI_FUNC
212 #define ABI_FUNC 'dfti_seqfourdp'
213 !End of the abilint section
214 
215  implicit none
216 
217 !Arguments ------------------------------------
218 !scalars
219  integer,intent(in) :: cplex,nx,ny,nz,ldx,ldy,ldz,ndat,isign
220 !arrays
221  real(dp),intent(inout) :: fofg(2*ldx*ldy*ldz*ndat)
222  real(dp),intent(inout) :: fofr(cplex*ldx*ldy*ldz*ndat)
223 
224 ! *************************************************************************
225 
226  select case (cplex)
227  case (2) ! Complex to Complex.
228 
229    select case (isign)
230    case (+1)
231      call dfti_many_dft_op(nx,ny,nz,ldx,ldy,ldz,ndat,isign,fofg,fofr)
232 
233    case (-1) ! -1
234      call dfti_many_dft_op(nx,ny,nz,ldx,ldy,ldz,ndat,isign,fofr,fofg)
235 
236    case default
237      MSG_BUG("Wrong isign")
238    end select
239 
240  case (1) ! Real case.
241 
242    select case (isign)
243    case (+1) ! G --> R
244      call dfti_c2r_op(nx,ny,nz,ldx,ldy,ldz,ndat,fofg,fofr)
245 
246    case (-1) ! R --> G
247      call dfti_r2c_op(nx,ny,nz,ldx,ldy,ldz,ndat,fofr,fofg)
248 
249    case default
250      MSG_BUG("Wrong isign")
251    end select
252 
253  case default
254    MSG_BUG(" Wrong value for cplex")
255  end select
256 
257 end subroutine dfti_seqfourdp

m_dfti/dfti_seqfourwf [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_seqfourwf

FUNCTION

 Carry out composite Fourier transforms between real and reciprocal (G) space.
 Wavefunctions, contained in a sphere in reciprocal space,
 can be FFT to real space. They can also be FFT from real space
 to a sphere. Also, the density maybe accumulated, and a local potential can be applied.

 The different options are :
 - option=0 --> reciprocal to real space and output the result.
 - option=1 --> reciprocal to real space and accumulate the density.
 - option=2 --> reciprocal to real space, apply the local potential to the wavefunction
                in real space and produce the result in reciprocal space.
 - option=3 --> real space to reciprocal space.
                NOTE that in this case, fftalg=1x1 MUST be used. This may be changed in the future.

INPUTS

 cplex= if 1 , denpot is real, if 2 , denpot is complex
    (cplex=2 only allowed for option=2, and istwf_k=1)
    not relevant if option=0 or option=3, so cplex=0 can be used to minimize memory
 fofgin(2,npwin)=holds input wavefunction in G vector basis sphere.
                 (intent(in) but the routine sphere can modify it for another iflag)
 gboundin(2*mgfft+8,2)=sphere boundary info for reciprocal to real space
 gboundout(2*mgfft+8,2)=sphere boundary info for real to reciprocal space
 istwf_k=option parameter that describes the storage of wfs
 kg_kin(3,npwin)=reduced planewave coordinates, input
 kg_kout(3,npwout)=reduced planewave coordinates, output
 mgfft=maximum size of 1D FFTs
 ndat=number of FFT to do in //
 ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
 npwin=number of elements in fofgin array (for option 0, 1 and 2)
 npwout=number of elements in fofgout array (for option 2 and 3)
 ldx,ldy,ldz=ngfft(4),ngfft(5),ngfft(6), dimensions of fofr.
 option= if 0: do direct FFT
         if 1: do direct FFT, then sum the density
         if 2: do direct FFT, multiply by the potential, then do reverse FFT
         if 3: do reverse FFT only
 weight_r=weight to be used for the accumulation of the density in real space
         (needed only when option=1)

OUTPUT

  (see side effects)

SIDE EFFECTS

 Input/Output
 for option==0, fofgin(2,npwin*ndat)=holds input wavefunction in G sphere;
                fofr(2,ldx,ldy,ldz) contains the output Fourier Transform of fofgin;
                no use of denpot, fofgout and npwout.
 for option==1, fofgin(2,npwin*ndat)=holds input wavefunction in G sphere;
                denpot(cplex*ldx,ldy,ldz) contains the input density at input,
                and the updated density at output (accumulated);
                no use of fofgout and npwout.
 for option==2, fofgin(2,npwin*ndat)=holds input wavefunction in G sphere;
                denpot(cplex*ldx,ldy,ldz) contains the input local potential;
                fofgout(2,npwout*ndat) contains the output function;
 for option==3, fofr(2,ldx,ldy,ldz*ndat) contains the input real space wavefunction;
                fofgout(2,npwout*ndat) contains its output Fourier transform;
                no use of fofgin and npwin.

PARENTS

      fourwf

CHILDREN

      c_f_pointer

SOURCE

331 subroutine dfti_seqfourwf(cplex,denpot,fofgin,fofgout,fofr,gboundin,gboundout,istwf_k,&
332 &  kg_kin,kg_kout,mgfft,ndat,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
333 
334 
335 !This section has been created automatically by the script Abilint (TD).
336 !Do not modify the following lines by hand.
337 #undef ABI_FUNC
338 #define ABI_FUNC 'dfti_seqfourwf'
339 !End of the abilint section
340 
341  implicit none
342 
343 !Arguments ------------------------------------
344 !scalars
345  integer,intent(in) :: cplex,istwf_k,ldx,ldy,ldz,ndat,npwin,npwout,option,mgfft
346  real(dp),intent(in) :: weight_i,weight_r
347 !arrays
348  integer,intent(in) :: gboundin(2*mgfft+8,2),gboundout(2*mgfft+8,2)
349  integer,intent(in) :: kg_kin(3,npwin),kg_kout(3,npwout),ngfft(18)
350  real(dp),intent(inout) :: denpot(cplex*ldx,ldy,ldz),fofgin(2,npwin*ndat)
351  real(dp),intent(inout) :: fofr(2,ldx*ldy*ldz*ndat)
352  real(dp),intent(out) :: fofgout(2,npwout*ndat)
353 
354 !Local variables-------------------------------
355 !scalars
356  integer,parameter :: ndat1=1
357  integer :: nx,ny,nz,fftalg,fftalga,fftalgc,fftcache,dat,ptg,ptr,ptgin,ptgout,nthreads
358  logical :: use_fftrisc
359  character(len=500) :: msg
360 
361 ! *************************************************************************
362 
363  if (option==1 .and. cplex/=1) then
364    write(msg,'(a,i0)')' With the option number 1, cplex must be 1 but it is cplex=',cplex
365    MSG_ERROR(msg)
366  end if
367 
368  if (option==2 .and. (cplex/=1 .and. cplex/=2)) then
369    write(msg,'(a,i0)')' With the option number 2, cplex must be 1 or 2, but it is cplex=',cplex
370    MSG_ERROR(msg)
371  end if
372 
373  nx=ngfft(1); ny=ngfft(2); nz=ngfft(3)
374  fftalg=ngfft(7); fftalga=fftalg/100; fftalgc=MOD(fftalg,10)
375  fftcache=ngfft(8)
376 
377  use_fftrisc = (fftalgc==2)
378  if (istwf_k==2.and.option==3) use_fftrisc = .FALSE.
379  if (istwf_k>2.and.ANY(option==(/0,3/))) use_fftrisc = .FALSE.
380 
381  nthreads = xomp_get_num_threads(open_parallel=.TRUE.)
382 
383  if (use_fftrisc) then
384    !call wrtout(std_out,strcat(ABI_FUNC," calls dfti_fftrisc"),"COLL")
385    if (ndat==1) then
386      call dfti_fftrisc_dp(cplex,denpot,fofgin,fofgout,fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
387 &      mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
388 
389    else
390      ! All this boilerplate code is needed because the caller might pass zero-sized arrays
391      ! for the arguments that are not referenced and we don't want to have problems at run-time.
392      ! Moreover option 1 requires a special treatment when threads are started at this level.
393 
394      SELECT CASE (option)
395      CASE (0)
396        !
397        ! fofgin -> fofr, no use of denpot, fofgout and npwout.
398        if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
399          do dat=1,ndat
400            ptg = 1 + (dat-1)*npwin
401            ptr = 1 + (dat-1)*ldx*ldy*ldz
402            call dfti_fftrisc_dp(cplex,denpot,fofgin(1,ptg),fofgout,fofr(1,ptr),gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
403 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
404          end do
405        else
406 !$OMP PARALLEL DO PRIVATE(ptg,ptr)
407          do dat=1,ndat
408            ptg = 1 + (dat-1)*npwin
409            ptr = 1 + (dat-1)*ldx*ldy*ldz
410            call dfti_fftrisc_dp(cplex,denpot,fofgin(1,ptg),fofgout,fofr(1,ptr),gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
411 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
412          end do
413        end if
414 
415      CASE (1)
416        !fofgin -> local ur and accumulate density in denpot
417        ! TODO this is delicate part to do in parallel, as one should OMP reduce denpot.
418        do dat=1,ndat
419          ptg = 1 + (dat-1)*npwin
420          ptr = 1 + (dat-1)*ldx*ldy*ldz
421          call dfti_fftrisc_dp(cplex,denpot,fofgin(1,ptg),fofgout,fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
422 &          mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
423        end do
424 
425      CASE (2)
426        ! <G|vloc(r)|fofgin(r)> in fofgout
427        if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
428          do dat=1,ndat
429            ptgin  = 1 + (dat-1)*npwin
430            ptgout = 1 + (dat-1)*npwout
431            call dfti_fftrisc_dp(cplex,denpot,fofgin(1,ptgin),fofgout(1,ptgout),fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
432 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
433          end do
434        else
435 !$OMP PARALLEL DO PRIVATE(ptgin,ptgout)
436          do dat=1,ndat
437            ptgin  = 1 + (dat-1)*npwin
438            ptgout = 1 + (dat-1)*npwout
439            call dfti_fftrisc_dp(cplex,denpot,fofgin(1,ptgin),fofgout(1,ptgout),fofr,gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
440 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
441          end do
442        end if
443 
444      CASE (3)
445        !fofr -> fofgout
446        if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
447          do dat=1,ndat
448            ptr    = 1 + (dat-1)*ldx*ldy*ldz
449            ptgout = 1 + (dat-1)*npwout
450            call dfti_fftrisc_dp(cplex,denpot,fofgin,fofgout(1,ptgout),fofr(1,ptr),gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
451 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
452          end do
453        else
454 !$OMP PARALLEL DO PRIVATE(ptr,ptgout)
455          do dat=1,ndat
456            ptr    = 1 + (dat-1)*ldx*ldy*ldz
457            ptgout = 1 + (dat-1)*npwout
458            call dfti_fftrisc_dp(cplex,denpot,fofgin,fofgout(1,ptgout),fofr(1,ptr),gboundin,gboundout,istwf_k,kg_kin,kg_kout,&
459 &            mgfft,ngfft,npwin,npwout,ldx,ldy,ldz,option,weight_r,weight_i)
460          end do
461        end if
462 
463      CASE DEFAULT
464        write(msg,'(a,i0,a)')'Option',option,' is not allowed. Only option=0, 1, 2 or 3 are allowed presently.'
465        MSG_ERROR(msg)
466      END SELECT
467 
468    end if
469 
470  else
471    SELECT CASE (option)
472    CASE (0)
473      !
474      ! FFT u(g) --> u(r)
475      if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
476        call dfti_fftug_dp(fftalg,fftcache,npwin,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_kin,gboundin,fofgin,fofr)
477      else
478 !$OMP PARALLEL DO PRIVATE(ptg, ptr)
479        do dat=1,ndat
480          ptg = 1 + (dat-1)*npwin
481          ptr = 1 + (dat-1)*ldx*ldy*ldz
482          call dfti_fftug_dp(fftalg,fftcache,npwin,nx,ny,nz,ldx,ldy,ldz,ndat1,&
483 &          istwf_k,mgfft,kg_kin,gboundin,fofgin(1,ptg),fofr(1,ptr))
484        end do
485      end if
486 
487    CASE (1)
488      ! TODO this is delicate part to do in parallel, as one should OMP reduce denpot.
489      call dfti_fftug_dp(fftalg,fftcache,npwin,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_kin,gboundin,fofgin,fofr)
490      call cg_addtorho(nx,ny,nz,ldx,ldy,ldz,ndat,weight_r,weight_i,fofr,denpot)
491 
492    CASE (2)
493 
494      if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
495        call dfti_fftug_dp(fftalg,fftcache,npwin,nx,ny,nz,ldx,ldy,ldz,ndat,istwf_k,mgfft,kg_kin,gboundin,fofgin,fofr)
496        call cg_vlocpsi(nx,ny,nz,ldx,ldy,ldz,ndat,cplex,denpot,fofr)
497 
498        !  The data for option==2 is now in fofr.
499        call dfti_fftpad_dp(fofr,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,-1,gboundout)
500 
501        call cg_box2gsph(nx,ny,nz,ldx,ldy,ldz,ndat,npwout,kg_kout,fofr,fofgout)
502      else
503 
504 !$OMP PARALLEL DO PRIVATE(ptg, ptr)
505        do dat=1,ndat
506          ptg = 1 + (dat-1)*npwin
507          ptr = 1 + (dat-1)*ldx*ldy*ldz
508          call dfti_fftug_dp(fftalg,fftcache,npwin,nx,ny,nz,ldx,ldy,ldz,ndat1,&
509 &          istwf_k,mgfft,kg_kin,gboundin,fofgin(1,ptg),fofr(1,ptr))
510 
511          call cg_vlocpsi(nx,ny,nz,ldx,ldy,ldz,ndat1,cplex,denpot,fofr(1,ptr))
512 
513          !  The data for option==2 is now in fofr.
514          call dfti_fftpad_dp(fofr(1,ptr),nx,ny,nz,ldx,ldy,ldz,ndat1,mgfft,-1,gboundout)
515 
516          ptg = 1 + (dat-1)*npwout
517          call cg_box2gsph(nx,ny,nz,ldx,ldy,ldz,ndat1,npwout,kg_kout,fofr(1,ptr),fofgout(1,ptg))
518        end do
519      end if
520 
521    CASE (3)
522      !  The data for option==3 is already in fofr.
523      if (.not.dfti_spawn_threads_here(ndat,nthreads)) then
524        call dfti_fftpad_dp(fofr,nx,ny,nz,ldx,ldy,ldz,ndat,mgfft,-1,gboundout)
525        call cg_box2gsph(nx,ny,nz,ldx,ldy,ldz,ndat,npwout,kg_kout,fofr,fofgout)
526      else
527 !$OMP PARALLEL DO PRIVATE(ptg, ptr)
528        do dat=1,ndat
529          ptg = 1 + (dat-1)*npwout
530          ptr = 1 + (dat-1)*ldx*ldy*ldz
531          call dfti_fftpad_dp(fofr(1,ptr),nx,ny,nz,ldx,ldy,ldz,ndat1,mgfft,-1,gboundout)
532          call cg_box2gsph(nx,ny,nz,ldx,ldy,ldz,ndat1,npwout,kg_kout,fofr(1,ptr),fofgout(1,ptg))
533        end do
534      end if
535 
536    CASE DEFAULT
537      write(msg,'(a,i0,a)')'Option',option,' is not allowed. Only option=0, 1, 2 or 3 are allowed presently.'
538      MSG_ERROR(msg)
539    END SELECT
540  end if
541 
542 end subroutine dfti_seqfourwf

m_dfti/dfti_spawn_threads_here [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_spawn_threads_here

FUNCTION

  Helper function that returns true if FFT calls should be OMP
  parallelized in the client code.

INPUTS

  ndat=Number of FFT transforms to do
  nthreads = Number of threads available

PARENTS

SOURCE

2352 function dfti_spawn_threads_here(ndat,nthreads) result(ans)
2353 
2354 
2355 !This section has been created automatically by the script Abilint (TD).
2356 !Do not modify the following lines by hand.
2357 #undef ABI_FUNC
2358 #define ABI_FUNC 'dfti_spawn_threads_here'
2359 !End of the abilint section
2360 
2361  implicit none
2362 
2363 !Arguments ------------------------------------
2364 !scalars
2365  integer,intent(in) :: ndat,nthreads
2366  logical :: ans
2367 !arrays
2368 
2369 ! *************************************************************************
2370 
2371  ans = .FALSE.
2372 #ifdef HAVE_OPENMP
2373  ans = (nthreads > 1 .and. MOD(ndat,nthreads) == 0 .and. .not. USE_LIB_THREADS)
2374 #else
2375  ABI_UNUSED((/ndat,nthreads/))
2376 #endif
2377 
2378 end function dfti_spawn_threads_here

m_dfti/dfti_use_lib_threads [ Functions ]

[ Top ] [ m_dfti ] [ Functions ]

NAME

 dfti_use_lib_threads

FUNCTION

INPUTS

PARENTS

      m_fft

CHILDREN

      c_f_pointer

SOURCE

2399 subroutine dfti_use_lib_threads(logvar)
2400 
2401 
2402 !This section has been created automatically by the script Abilint (TD).
2403 !Do not modify the following lines by hand.
2404 #undef ABI_FUNC
2405 #define ABI_FUNC 'dfti_use_lib_threads'
2406 !End of the abilint section
2407 
2408  implicit none
2409 
2410 !Arguments ------------------------------------
2411 !scalars
2412  logical,intent(in) :: logvar
2413 
2414 ! *************************************************************************
2415 
2416  USE_LIB_THREADS = logvar
2417 
2418 end subroutine dfti_use_lib_threads