TABLE OF CONTENTS
ABINIT/laplacian [ Functions ]
NAME
laplacian
FUNCTION
compute the laplacian of a function defined in real space the code is written in the way of /3xc/xcden.F90
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, MT) 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/Infos/contributors .
INPUTS
gprimd(3,3)=dimensional reciprocal space primitive translations mpi_enreg=informations about MPI parallelization nfft=number of points of the fft grid nfunc=number of functions on the grid for which the laplacian is to be calculated ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft paral_kgb=flag controlling (k,g,bands) parallelization (optional) rdfuncr(nfft,nfunc)=real(dp) discretized functions in real space rdfuncg_in TO BE DESCRIBED SB 090901 laplacerdfuncg_in TO BE DESCRIBED SB 090901 (optional) g2cart_in(nfft) = G**2 on the grid
OUTPUT
(optional) laplacerdfuncr = laplacian in real space of the functions in rdfuncr (optional) rdfuncg = real(dp) discretized functions in fourier space (optional) laplacerdfuncg = real(dp) discretized laplacian of the functions in fourier space (optional) g2cart_out(nfft) = G**2 on the grid rdfuncg_out TO BE DESCRIBED SB 090901 laplacerdfuncg_out TO BE DESCRIBED SB 090901
SIDE EFFECTS
WARNINGS
NOTES
PARENTS
frskerker1,frskerker2,moddiel_csrb,prcrskerker1,prcrskerker2
CHILDREN
fourdp,ptabs_fourdp
SOURCE
51 #if defined HAVE_CONFIG_H 52 #include "config.h" 53 #endif 54 55 #include "abi_common.h" 56 57 subroutine laplacian(gprimd,mpi_enreg,nfft,nfunc,ngfft,paral_kgb,rdfuncr,& 58 & laplacerdfuncr,rdfuncg_out,laplacerdfuncg_out,g2cart_out,rdfuncg_in,g2cart_in) 59 60 use defs_basis 61 use defs_abitypes 62 use m_profiling_abi 63 use m_errors 64 65 use m_mpinfo, only : ptabs_fourdp 66 67 !This section has been created automatically by the script Abilint (TD). 68 !Do not modify the following lines by hand. 69 #undef ABI_FUNC 70 #define ABI_FUNC 'laplacian' 71 use interfaces_53_ffts 72 !End of the abilint section 73 74 implicit none 75 76 !Arguments ------------------------------------ 77 !scalars 78 integer,intent(in) :: nfft,nfunc,paral_kgb 79 type(MPI_type),intent(in) :: mpi_enreg 80 !arrays 81 integer,intent(in) :: ngfft(18) 82 real(dp),intent(in) :: gprimd(3,3) 83 real(dp),intent(inout),optional :: laplacerdfuncr(nfft,nfunc) 84 real(dp),intent(inout),optional,target :: rdfuncr(nfft,nfunc) 85 real(dp),intent(in),optional,target :: g2cart_in(nfft) !vz_i 86 real(dp),intent(out),optional,target :: g2cart_out(nfft) !vz_i 87 real(dp),intent(out),optional,target :: laplacerdfuncg_out(2,nfft,nfunc) 88 real(dp),intent(in),optional,target :: rdfuncg_in(2,nfft,nfunc) !vz_i 89 real(dp),intent(out),optional,target :: rdfuncg_out(2,nfft,nfunc) 90 91 !Local variables------------------------------- 92 !scalars 93 integer :: count,i1,i2,i3,id1,id2,id3,ifft,ifunc,ig1,ig2,ig3,ii1,n1,n2 94 integer :: n3 95 real(dp) :: b11,b12,b13,b21,b22,b23,b31,b32,b33 96 !arrays 97 integer, ABI_CONTIGUOUS pointer :: fftn2_distrib(:),ffti2_local(:) 98 integer, ABI_CONTIGUOUS pointer :: fftn3_distrib(:),ffti3_local(:) 99 real(dp),ABI_CONTIGUOUS pointer :: g2cart(:),laplacerdfuncg(:,:,:),rdfuncg(:,:,:) 100 101 ! ************************************************************************* 102 103 !Keep local copy of fft dimensions 104 n1=ngfft(1) 105 n2=ngfft(2) 106 n3=ngfft(3) 107 108 if(present(laplacerdfuncg_out)) then 109 laplacerdfuncg => laplacerdfuncg_out 110 else 111 ABI_ALLOCATE(laplacerdfuncg,(2,nfft,nfunc)) 112 end if 113 114 ! Get the distrib associated with this fft_grid 115 call ptabs_fourdp(mpi_enreg,n2,n3,fftn2_distrib,ffti2_local,fftn3_distrib,ffti3_local) 116 117 !change the real density rdfuncr on real space on the real density 118 !rdfuncg in reciprocal space 119 if(.not.present(rdfuncg_in)) then 120 if(present(rdfuncg_out)) then 121 rdfuncg => rdfuncg_out 122 else 123 ABI_ALLOCATE(rdfuncg,(2,nfft,nfunc)) 124 end if 125 if(present(rdfuncr)) then 126 do ifunc=1,nfunc 127 call fourdp(1,rdfuncg(:,:,ifunc),rdfuncr(:,ifunc),-1,mpi_enreg,nfft,ngfft,paral_kgb,0) 128 end do 129 end if 130 else 131 rdfuncg => rdfuncg_in 132 end if 133 134 135 136 !apply the laplacian on laplacerdfuncr 137 !code from /3xc/xcden.F90 138 !see also src/5common/hatre.F90 and src/5common/moddiel.F90 139 !Keep local copy of fft dimensions 140 !Initialize computation of G^2 in cartesian coordinates 141 if(.not.present(g2cart_in)) then 142 if(present(g2cart_out)) then 143 g2cart => g2cart_out 144 else 145 ABI_ALLOCATE(g2cart,(nfft)) 146 end if 147 id1=int(n1/2)+2 148 id2=int(n2/2)+2 149 id3=int(n3/2)+2 150 count=0 151 do i3=1,n3 152 ifft=(i3-1)*n1*(n2/mpi_enreg%nproc_fft) 153 ig3=i3-int(i3/id3)*n3-1 154 do i2=1,n2 155 if (fftn2_distrib(i2)==mpi_enreg%me_fft) then 156 ig2=i2-int(i2/id2)*n2-1 157 158 ii1=1 159 do i1=ii1,n1 160 ig1=i1-int(i1/id1)*n1-1 161 ifft=ifft+1 162 163 b11=gprimd(1,1)*real(ig1,dp) 164 b21=gprimd(2,1)*real(ig1,dp) 165 b31=gprimd(3,1)*real(ig1,dp) 166 b12=gprimd(1,2)*real(ig2,dp) 167 b22=gprimd(2,2)*real(ig2,dp) 168 b32=gprimd(3,2)*real(ig2,dp) 169 b13=gprimd(1,3)*real(ig3,dp) 170 b23=gprimd(2,3)*real(ig3,dp) 171 b33=gprimd(3,3)*real(ig3,dp) 172 173 g2cart(ifft)=( & 174 & (b11+b12+b13)**2& 175 & +(b21+b22+b23)**2& 176 & +(b31+b32+b33)**2& 177 & ) 178 do ifunc=1,nfunc 179 ! compute the laplacien in Fourier space that is * (i x 2pi x G)**2 180 laplacerdfuncg(1,ifft,ifunc) = -rdfuncg(1,ifft,ifunc)*g2cart(ifft)*two_pi*two_pi 181 laplacerdfuncg(2,ifft,ifunc) = -rdfuncg(2,ifft,ifunc)*g2cart(ifft)*two_pi*two_pi 182 end do 183 end do 184 end if 185 end do 186 end do 187 if(.not.present(g2cart_out)) then 188 ABI_DEALLOCATE(g2cart) 189 end if 190 else 191 g2cart => g2cart_in 192 do ifunc=1,nfunc 193 do ifft=1,nfft 194 ! compute the laplacien in Fourier space that is * (i x 2pi x G)**2 195 laplacerdfuncg(1,ifft,ifunc) = -rdfuncg(1,ifft,ifunc)*g2cart(ifft)*two_pi*two_pi 196 laplacerdfuncg(2,ifft,ifunc) = -rdfuncg(2,ifft,ifunc)*g2cart(ifft)*two_pi*two_pi 197 end do 198 end do 199 end if 200 201 !get the result back into real space 202 if(present(laplacerdfuncr)) then 203 do ifunc=1,nfunc 204 call fourdp(1,laplacerdfuncg(:,:,ifunc),laplacerdfuncr(:,ifunc),1,mpi_enreg,nfft,ngfft,paral_kgb,0) 205 end do 206 end if 207 208 !deallocate pointers 209 if((.not.present(rdfuncg_in)).and.(.not.present(rdfuncg_in))) then 210 ABI_DEALLOCATE(rdfuncg) 211 end if 212 if(.not.present(laplacerdfuncg_out)) then 213 ABI_DEALLOCATE(laplacerdfuncg) 214 end if 215 216 end subroutine laplacian