TABLE OF CONTENTS
ABINIT/green_kernel [ Functions ]
NAME
green_kernel
FUNCTION
this routine compute the fourrier transform of the Green kernel for the recursion method
COPYRIGHT
Copyright (C) 2008-2018 ABINIT group ( ). 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
inf_rmet=define the infinitesimal metric : rprimd*(transpose(rprimd)) divided by the number of discretisation point inf_ucvol=volume of infinitesimal cell mult=variance of the Gaussian (=rtrotter/beta) mpi_enreg=information about MPI parallelization ngfft=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft nfft=total number of fft grid points debug_rec=debugging variable
OUTPUT
ZT_p=fourier transforme of the Green kernel
PARENTS
first_rec
CHILDREN
fourdp,timab,wrtout
NOTES
at this time : - need a rectangular box
SOURCE
41 #if defined HAVE_CONFIG_H 42 #include "config.h" 43 #endif 44 45 #include "abi_common.h" 46 47 subroutine green_kernel(ZT_p,inf_rmet,inf_ucvol,mult,mpi_enreg,ngfft,nfft) 48 49 use defs_basis 50 use defs_abitypes 51 use m_profiling_abi 52 53 !This section has been created automatically by the script Abilint (TD). 54 !Do not modify the following lines by hand. 55 #undef ABI_FUNC 56 #define ABI_FUNC 'green_kernel' 57 use interfaces_14_hidewrite 58 use interfaces_18_timing 59 use interfaces_53_ffts 60 !End of the abilint section 61 62 implicit none 63 64 !Arguments ------------------------------- 65 !scalars 66 integer,intent(in) :: nfft 67 real(dp),intent(in) :: inf_ucvol,mult 68 type(MPI_type),intent(in) :: mpi_enreg 69 !arrays 70 integer,intent(in) :: ngfft(18) 71 real(dp),intent(in) :: inf_rmet(3,3) 72 real(dp),intent(out) :: ZT_p(1:2,0:nfft-1) 73 74 !Local variables------------------------------- 75 !scalars 76 integer,parameter :: n_green_max=5 77 integer :: ii,isign,jj,kk,n_green,xx,yy,zz 78 real(dp) :: acc, norme 79 character(len=500) :: msg 80 !arrays 81 real(dp) :: tsec(2) 82 real(dp),allocatable :: T_p(:) 83 84 ! ************************************************************************* 85 86 call timab(603,1,tsec) 87 88 norme = (mult/pi)**(onehalf) 89 90 ABI_ALLOCATE(T_p,(0:nfft-1)) 91 92 !n_green should be better chosen for non rectangular cell 93 do xx=1, n_green_max 94 n_green = xx 95 if(exp(-mult*dsq_green(xx*ngfft(1),0,0,inf_rmet))<tol14 & 96 & .and. exp(-mult*dsq_green(0,xx*ngfft(2),0,inf_rmet))<tol14 & 97 & .and. exp(-mult*dsq_green(0,0,xx*ngfft(3),inf_rmet))<tol14 ) exit 98 end do 99 100 acc = zero 101 T_p = zero 102 do kk = 0,ngfft(3)-1 103 do jj = 0,ngfft(2)-1 104 do ii = 0,ngfft(1)-1 105 106 do xx=-n_green,n_green-1 107 do yy=-n_green,n_green-1 108 do zz=-n_green,n_green-1 109 110 T_p(ii+ngfft(1)*jj+ngfft(1)*ngfft(2)*kk) = T_p(ii+ngfft(1)*jj+ngfft(1)*ngfft(2)*kk)+ & 111 & exp(-mult*dsq_green(ii+xx*ngfft(1),jj+yy*ngfft(2),kk+zz*ngfft(3),inf_rmet)) 112 113 end do 114 end do 115 end do 116 117 T_p(ii+ngfft(1)*jj+ngfft(1)*ngfft(2)*kk) = norme*T_p(ii+ngfft(1)*jj+ngfft(1)*ngfft(2)*kk) 118 acc = acc + inf_ucvol* T_p(ii+ngfft(1)*jj+ngfft(1)*ngfft(2)*kk) 119 120 end do 121 end do 122 end do 123 124 T_p(:)= (one/acc)*T_p(:) 125 126 !if(debug_rec)then 127 write(msg,'(a,d12.3,2(2a,i8),2(2a,3d12.3),2a,d16.6)')& 128 & ' on the boundary ', exp(-mult*dsq_green(ngfft(1),0,0,inf_rmet)),ch10, & 129 & ' no zero ', count(T_p>tol14),ch10, & 130 & ' n_green ', n_green,ch10, & 131 & ' erreur_n_green ', exp(-mult*dsq_green(n_green*ngfft(1),0,0,inf_rmet)), & 132 & exp(-mult*dsq_green(0,n_green*ngfft(2),0,inf_rmet)), & 133 & exp(-mult*dsq_green(0,0,n_green*ngfft(3),inf_rmet)),ch10,& 134 & ' erreur_troncat ', T_p(ngfft(1)/2), & 135 & T_p(ngfft(1)*(ngfft(2)/2)), & 136 & T_P(ngfft(1)*ngfft(2)*(ngfft(3)/2)),ch10, & 137 & ' erreurT_p ',abs(acc-1.d0) 138 call wrtout(std_out,msg,'COLL') 139 !endif 140 141 142 isign = -1 143 call fourdp(1,ZT_p,T_p,isign,mpi_enreg,nfft,ngfft,1,0) 144 145 ABI_DEALLOCATE(T_p) 146 147 ZT_p(:,:) = real(nfft,dp)*ZT_p 148 149 150 call timab(603,2,tsec) 151 152 contains 153 154 function dsq_green(ii,jj,kk,inf_rmet) 155 156 157 !This section has been created automatically by the script Abilint (TD). 158 !Do not modify the following lines by hand. 159 #undef ABI_FUNC 160 #define ABI_FUNC 'dsq_green' 161 !End of the abilint section 162 163 real(dp) :: dsq_green 164 integer,intent(in) :: ii,jj,kk 165 real(dp),intent(in) :: inf_rmet(3,3) 166 dsq_green= inf_rmet(1,1)*dble(ii**2)& 167 & +inf_rmet(2,2)*dble(jj**2)& 168 & +inf_rmet(3,3)*dble(kk**2)& 169 & +two*(inf_rmet(1,2)*dble(ii*jj)& 170 & +inf_rmet(2,3)*dble(jj*kk)& 171 & +inf_rmet(3,1)*dble(kk*ii)) 172 end function dsq_green 173 174 end subroutine green_kernel