TABLE OF CONTENTS
ABINIT/fftpac [ Functions ]
NAME
fftpac
FUNCTION
Allow for data copying to modify the stride (dimensioning) of a three- dimensional array, for more efficient three dimensional fft. NOTE that the arrays are in REAL space. Note that arrays aa and bb may be the same array (start at the same address). The array aa also incorporate a spin variable. MG FIXME: THIS IS **VERY BAD** AS FORTRAN DOES NOT ALLOW FOR ALIASING
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, MF, XG, GMR). This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
ispden=actual spin-density of interest nspden=number of spin-density components n1,n2,n3=actual data dimensions, dimensions of complex array a nd1,nd2,nd3=array dimensions of (larger) array b ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft option= see description of side effects
OUTPUT
(see side effects)
SIDE EFFECTS
aa & bb arrays are treated as input or output depending on option: option=1 aa(n1*n2*n3,ispden) <-- bb(nd1,nd2,nd3) real case option=2 aa(n1*n2*n3,ispden) --> bb(nd1,nd2,nd3) real case option=10 aa(n1*n2*n3,ispden) <-- bb(nd1,nd2,nd3) complex case like option 1 real part option=11 aa(n1*n2*n3,ispden) <-- bb(nd1,nd2,nd3) complex case like option 1 imag part
PARENTS
dfpt_mkrho,dfpt_nstpaw,dfpt_rhofermi,dfpt_vtorho,dfptnl_resp,energy fock_getghc,getgh1c,gwls_hamiltonian,ks_ddiago,m_epjdos,m_io_kss,mkrho suscep_stat,vtorho
CHILDREN
ptabs_fourdp
SOURCE
50 #if defined HAVE_CONFIG_H 51 #include "config.h" 52 #endif 53 54 #include "abi_common.h" 55 56 subroutine fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,nd1,nd2,nd3,ngfft,aa,bb,option) 57 58 use defs_basis 59 use m_errors 60 61 use defs_abitypes, only : MPI_type 62 use m_mpinfo, only : ptabs_fourdp 63 64 !This section has been created automatically by the script Abilint (TD). 65 !Do not modify the following lines by hand. 66 #undef ABI_FUNC 67 #define ABI_FUNC 'fftpac' 68 !End of the abilint section 69 70 implicit none 71 72 !Arguments ------------------------------------ 73 !scalars 74 integer,intent(in) :: ispden,n1,n2,n3,nd1,nd2,nd3,nspden,option 75 type(mpi_type),intent(in) :: mpi_enreg 76 !arrays 77 integer,intent(in) :: ngfft(18) 78 real(dp),intent(inout) :: aa(n1*n2*n3/ngfft(10),nspden),bb(nd1,nd2,nd3) 79 80 !Local variables------------------------------- 81 !scalars 82 integer :: i1,i2,i3,index,me_fft,nproc_fft 83 character(len=500) :: message 84 !arrays 85 integer, ABI_CONTIGUOUS pointer :: fftn2_distrib(:),ffti2_local(:) 86 integer, ABI_CONTIGUOUS pointer :: fftn3_distrib(:),ffti3_local(:) 87 88 ! ************************************************************************* 89 90 me_fft=ngfft(11); nproc_fft=ngfft(10) 91 92 if (option==1.or.option==2) then 93 if (nd1<n1.or.nd2<n2.or.nd3<n3) then 94 write(message,'(a,3i0,2a,3i0,a)')& 95 & 'Each of nd1,nd2,nd3=',nd1,nd2,nd3,ch10,& 96 & 'must be >= n1, n2, n3 =',n1,n2,n3,'.' 97 MSG_BUG(message) 98 end if 99 else 100 if (2*nd1<n1.or.nd2<n2.or.nd3<n3) then 101 write(message,'(a,3i0,2a,3i0,a)')& 102 & 'Each of 2*nd1,nd2,nd3=',2*nd1,nd2,nd3,ch10,& 103 & 'must be >= (n1, n2, n3) =',n1,n2,n3,'.' 104 MSG_BUG(message) 105 end if 106 end if 107 108 ! Get the distrib associated with this fft_grid 109 call ptabs_fourdp(mpi_enreg,n2,n3,fftn2_distrib,ffti2_local,fftn3_distrib,ffti3_local) 110 111 if (option==1) then 112 do i3=1,n3 113 if (me_fft==fftn3_distrib(i3)) then 114 do i2=1,n2 115 do i1=1,n1 116 aa(i1+n1*(i2-1+n2*(ffti3_local(i3)-1)),ispden)=bb(i1,i2,i3) 117 end do 118 end do 119 end if 120 end do 121 122 else if (option==2) then 123 ! Here we avoid corrupting the data in a while writing to b in the 124 ! case in which a and b are same array. 125 ! Also: replace "trash" data with 0 s to avoid floating point 126 ! exceptions when this data is actually manipulated in fft. 127 do i3=nd3,n3+1,-1 128 do i2=nd2,1,-1 129 do i1=nd1,1,-1 130 bb(i1,i2,i3)=0.d0 131 end do 132 end do 133 end do 134 do i3=n3,1,-1 135 if (me_fft==fftn3_distrib(i3)) then 136 do i2=nd2,n2+1,-1 137 do i1=nd1,1,-1 138 bb(i1,i2,i3)=0.d0 139 end do 140 end do 141 do i2=n2,1,-1 142 do i1=nd1,n1+1,-1 143 bb(i1,i2,i3)=0.d0 144 end do 145 do i1=n1,1,-1 146 bb(i1,i2,i3)=aa(i1+n1*(i2-1+n2*(ffti3_local(i3) - 1)),ispden) 147 end do 148 end do 149 end if 150 end do 151 ! MF 152 else if (option==10 .or. option==11) then 153 index=1 154 if(option==11) index=2 155 do i3=1,n3 156 do i2=1,n2 157 do i1=1,n1/2 158 aa(index,ispden)=bb(i1,i2,i3) 159 index=index+2 160 end do 161 end do 162 end do 163 ! MF 164 else 165 write(message,'(a,i0,a)')' Bad option =',option,'.' 166 MSG_BUG(message) 167 end if 168 169 end subroutine fftpac