TABLE OF CONTENTS
- ABINIT/m_pawfgr
- m_paw_ij/pawfgr_nullify
- m_pawfgr/indgrid
- m_pawfgr/pawfgr_destroy
- m_pawfgr/pawfgr_init
- m_pawfgr/pawfgr_type
ABINIT/m_pawfgr [ Modules ]
NAME
m_pawfgr
FUNCTION
This module contains the definition of the pawfgr_type structured datatype, as well as related functions and methods. pawfgr_type variables define Fine rectangular GRid parameters and related data.
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (MT, FJ) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
OUTPUT
NOTES
* Routines tagged with "@type_name" are strongly connected to the definition of the data type. Strongly connected means that the proper functioning of the implementation relies on the assumption that the tagged procedure is consistent with the type declaration. Every time a developer changes the structure "type_name" adding new entries, he/she has to make sure that all the strongly connected routines are changed accordingly to accommodate the modification of the data type Typical examples of strongly connected routines are creation, destruction or reset methods.
SOURCE
30 #if defined HAVE_CONFIG_H 31 #include "config.h" 32 #endif 33 34 #include "abi_common.h" 35 36 MODULE m_pawfgr 37 38 use defs_basis 39 use m_errors 40 use m_abicore 41 use m_xmpi 42 use m_dtset 43 44 use m_kg, only : getcut 45 46 implicit none 47 48 private 49 50 !public procedures. 51 public :: pawfgr_init 52 public :: pawfgr_destroy 53 public :: pawfgr_nullify 54 public :: indgrid
m_paw_ij/pawfgr_nullify [ Functions ]
[ Top ] [ m_paw_ij ] [ Functions ]
NAME
pawfgr_nullify
FUNCTION
Nullify pointers and flags in a pawfgr structure
SIDE EFFECTS
Pawfgr<type(pawfgr_type)>=Fine GRid parameters and related data. Nullified in output
SOURCE
324 subroutine pawfgr_nullify(Pawfgr) 325 326 !Arguments ------------------------------------ 327 !arrays 328 type(Pawfgr_type),intent(inout) :: Pawfgr 329 330 !Local variables------------------------------- 331 332 ! ************************************************************************* 333 334 !@Pawfgr_type 335 336 nullify(Pawfgr%coatofin) 337 nullify(Pawfgr%fintocoa) 338 339 Pawfgr%usefinegrid=0 340 341 end subroutine pawfgr_nullify
m_pawfgr/indgrid [ Functions ]
[ Top ] [ m_pawfgr ] [ Functions ]
NAME
indgrid
FUNCTION
Calculate the correspondance between the coarse grid and the fine grid
INPUTS
nfftc=total number of FFt grid=n1*n2*n3 for the coarse grid nfftf=total number of FFt grid=n1*n2*n3 for the fine grid ngfftc(18)=contain all needed information about 3D FFT, for the coarse grid, see ~abinit/doc/variables/vargs.htm#ngfft ngfftf(18)=contain all needed information about 3D FFT, for the fine grid, see ~abinit/doc/variables/vargs.htm#ngfft
OUTPUT
coatofin(nfftc)= index of the points of the coarse grid on the fine grid fintocoa(nfftf)=index of the points of the fine grid on the coarse grid (=0 if the point of the fine grid does not belong to the coarse grid).
SOURCE
370 subroutine indgrid(coatofin,fintocoa,nfftc,nfftf,ngfftc,ngfftf) 371 372 !Arguments ------------------------------------ 373 !scalars 374 integer,intent(in) :: nfftc,nfftf 375 !arrays 376 integer,intent(in) :: ngfftc(18),ngfftf(18) 377 integer,intent(out) :: coatofin(nfftc),fintocoa(nfftf) 378 379 !Local variables------------------------------- 380 !scalars 381 integer :: i1,i2,i3,if1,if2,if3,ii,ing,n1c,n1f,n2c,n2f,n3c,n3f,narg1,narg2 382 !arrays 383 integer :: id(3) 384 integer,allocatable :: gc(:,:),gf(:,:) 385 character(len=500) :: msg 386 387 ! ************************************************************************* 388 389 DBG_ENTER("COLL") 390 391 n1c=ngfftc(1);n2c=ngfftc(2);n3c=ngfftc(3) 392 n1f=ngfftf(1);n2f=ngfftf(2);n3f=ngfftf(3) 393 394 ABI_MALLOC(gc,(3,max(n1c,n2c,n3c))) 395 do ii=1,3 396 id(ii)=ngfftc(ii)/2+2 397 do ing=1,ngfftc(ii) 398 gc(ii,ing)=ing-(ing/id(ii))*ngfftc(ii)-1 399 end do 400 end do 401 402 ABI_MALLOC(gf,(3,max(n1f,n2f,n3f))) 403 do ii=1,3 404 id(ii)=ngfftf(ii)/2+2 405 do ing=1,ngfftf(ii) 406 gf(ii,ing)=ing-(ing/id(ii))*ngfftf(ii)-1 407 end do 408 end do 409 410 coatofin=0;fintocoa=0 411 do i1=1,n1c 412 do if1=1,n1f 413 if(gc(1,i1)==gf(1,if1)) then 414 do i2=1,n2c 415 do if2=1,n2f 416 if(gc(2,i2)==gf(2,if2)) then 417 do i3=1,n3c 418 do if3=1,n3f 419 if(gc(3,i3)==gf(3,if3)) then 420 narg1=i1+n1c*(i2-1+n2c*(i3-1)) 421 narg2=if1+n1f*(if2-1+n2f*(if3-1)) 422 coatofin(narg1)=narg2 423 fintocoa(narg2)=narg1 424 exit 425 end if 426 end do 427 end do 428 exit ! To avoid N_fine * N_coarse scaling 429 end if 430 end do 431 end do 432 exit ! To avoid N_fine * N_coarse scaling 433 end if 434 end do 435 end do 436 437 !Check coatofin to make sure there are no zeros! 438 do ii=1,ubound(coatofin,1) 439 if (coatofin(ii)==0) then 440 msg = 'A zero was found in coatofin. Check that the fine FFT mesh is finer in each dimension than the coarse FFT mesh.' 441 ABI_ERROR(msg) 442 end if 443 end do 444 445 ABI_FREE(gf) 446 ABI_FREE(gc) 447 448 DBG_EXIT("COLL") 449 450 end subroutine indgrid
m_pawfgr/pawfgr_destroy [ Functions ]
[ Top ] [ m_pawfgr ] [ Functions ]
NAME
pawfgr_destroy
FUNCTION
Deallocate pointers and nullify flags in a pawfgr structure
SIDE EFFECTS
pawfgr<type(pawfgr_type)>= Fine GRid parameters and related data
SOURCE
282 subroutine pawfgr_destroy(Pawfgr) 283 284 !Arguments ------------------------------------ 285 !arrays 286 type(Pawfgr_type),intent(inout) :: Pawfgr 287 288 !Local variables------------------------------- 289 290 ! ************************************************************************* 291 292 DBG_ENTER("COLL") 293 294 !@Pawfgr_type 295 296 if (associated(Pawfgr%coatofin)) then 297 ABI_FREE(Pawfgr%coatofin) 298 end if 299 if (associated(Pawfgr%fintocoa)) then 300 ABI_FREE(Pawfgr%fintocoa) 301 end if 302 303 Pawfgr%usefinegrid=0 304 305 DBG_EXIT("COLL") 306 307 end subroutine pawfgr_destroy
m_pawfgr/pawfgr_init [ Functions ]
[ Top ] [ m_pawfgr ] [ Functions ]
NAME
pawfgr_init
FUNCTION
Initialize a pawfgr_type datatype, reporting also info on the mesh according to the method used (norm-conserving PSP or PAW)
INPUTS
k0(3)=input k vector for k+G sphere Dtset <type(dataset_type)>=all input variables for this dataset %dilatmx %usepaw %usewvl %natom %ngfft %ngfftdg %nfft %mgfft %mgfftdg %dilatmx %pawecutdg %ecut gmet(3,3)=reciprocal space metric (bohr^-2)
OUTPUT
ecut_eff=effective energy cutoff (hartree) for coarse planewave basis sphere ecutdg_eff=effective energy cutoff (hartree) for dense planewave basis sphere gsqcutc_eff=(PAW) Fourier cutoff on G^2 for "large sphere" of radius double for the coarse FFT grid nfftf=(effective) number of FFT grid points (for this proc), for dense FFT mesh mgfftf=maximum size of 1D FFTs, for dense FFT mesh ngfftc(18),ngfftf(18)=contain all needed information about 3D FFT, for coarse and dense FFT mesh, resp. see ~abinit/doc/variables/vargs.htm#ngfft Pawfgr<pawfgr_type>=For PAW, Fine rectangular GRid parameters and related data
SOURCE
158 subroutine pawfgr_init(Pawfgr,Dtset,mgfftf,nfftf,ecut_eff,ecutdg_eff,ngfftc,ngfftf,& 159 & gsqcutc_eff,gsqcutf_eff,gmet,k0) ! optional 160 161 !Arguments ------------------------------------ 162 !scalars 163 integer,intent(out) :: nfftf,mgfftf 164 real(dp),intent(out) :: ecut_eff,ecutdg_eff 165 real(dp),intent(out),optional :: gsqcutf_eff,gsqcutc_eff 166 type(dataset_type),intent(in) :: Dtset 167 type(Pawfgr_type),intent(out) :: Pawfgr 168 !arrays 169 real(dp),intent(in),optional :: gmet(3,3) 170 integer,intent(out) :: ngfftc(18),ngfftf(18) 171 real(dp),intent(in),optional :: k0(3) 172 173 !Local variables------------------------------- 174 integer :: ii,nfftc_tot,nfftf_tot 175 real(dp) :: boxcut,boxcutc 176 character(len=500) :: msg 177 178 !************************************************************************ 179 180 DBG_ENTER("COLL") 181 182 !@Pawfgr_type 183 if ((present(gsqcutc_eff).or.present(gsqcutf_eff)).and.& 184 ((.not.present(gmet)).or.(.not.present(k0)))) then 185 ABI_BUG('To compute gsqcut[c,f]_eff, both k0 and gmet must be present as argument !') 186 end if 187 188 ngfftc(:)=Dtset%ngfft(:) 189 190 SELECT CASE (Dtset%usepaw) 191 192 CASE (0) 193 ! === Norm-conserving pseudopotentials === 194 nfftf=Dtset%nfft ; mgfftf=Dtset%mgfft ; ngfftf(:)=Dtset%ngfft(:) 195 Pawfgr%usefinegrid=0 196 ABI_MALLOC(Pawfgr%coatofin,(0)) 197 ABI_MALLOC(Pawfgr%fintocoa,(0)) 198 ecut_eff =Dtset%ecut*Dtset%dilatmx**2 199 ecutdg_eff=ecut_eff 200 201 CASE (1) 202 ! == PAW calculation === 203 if (any(Dtset%ngfftdg(1:3)/=Dtset%ngfft(1:3)) .and. Dtset%usewvl==0) then 204 ! Use fine FFT grid generated according to pawecutdg. 205 nfftf=Dtset%nfftdg ; mgfftf=Dtset%mgfftdg ; ngfftf(:)=Dtset%ngfftdg(:) 206 nfftc_tot =ngfftc(1)*ngfftc(2)*ngfftc(3) 207 nfftf_tot =ngfftf(1)*ngfftf(2)*ngfftf(3) 208 Pawfgr%usefinegrid=1 209 ABI_MALLOC(Pawfgr%coatofin,(nfftc_tot)) 210 ABI_MALLOC(Pawfgr%fintocoa,(nfftf_tot)) 211 call indgrid(Pawfgr%coatofin,Pawfgr%fintocoa,nfftc_tot,nfftf_tot,ngfftc,ngfftf) 212 213 else 214 ! Do not use fine FFT mesh. Simple transfer that can be done in parallel with only local info. 215 nfftf=Dtset%nfft ; mgfftf=Dtset%mgfft ; ngfftf(:)=Dtset%ngfft(:) 216 Pawfgr%usefinegrid=0 217 ABI_MALLOC(Pawfgr%coatofin,(Dtset%nfft)) 218 ABI_MALLOC(Pawfgr%fintocoa,(Dtset%nfft)) 219 do ii=1,Dtset%nfft 220 Pawfgr%coatofin(ii)=ii ; Pawfgr%fintocoa(ii)=ii 221 end do 222 end if 223 ecutdg_eff=Dtset%pawecutdg*Dtset%dilatmx**2 224 ecut_eff =Dtset%ecut*Dtset%dilatmx**2 225 226 CASE DEFAULT 227 write(msg,'(a,i4)')' Wrong value of usepaw: ',Dtset%usepaw 228 ABI_BUG(msg) 229 END SELECT 230 231 ! Store useful dimensions in Pawfgr 232 Pawfgr%nfftc=Dtset%nfft ; Pawfgr%mgfftc=Dtset%mgfft ; Pawfgr%ngfftc(:)=Dtset%ngfft(:) 233 Pawfgr%nfft=nfftf ; Pawfgr%mgfft=mgfftf ; Pawfgr%ngfft (:)=ngfftf(:) 234 235 ! Get boxcut for given gmet, ngfft, and ecut (center at k0) === 236 ! boxcut=ratio of basis sphere diameter to fft box side 237 boxcut=-one 238 if (Dtset%usepaw==1) then 239 if (present(gsqcutc_eff)) then 240 write(msg,'(2a)')ch10,' Coarse grid specifications (used for wave-functions):' 241 call wrtout(std_out,msg) 242 call getcut(boxcutc,ecut_eff,gmet,gsqcutc_eff,Dtset%iboxcut,std_out,k0,ngfftc) 243 end if 244 if (present(gsqcutf_eff)) then 245 write(msg,'(2a)')ch10,' Fine grid specifications (used for densities):' 246 call wrtout(std_out,msg) 247 call getcut(boxcut,ecutdg_eff,gmet,gsqcutf_eff,Dtset%iboxcut,std_out,k0,ngfftf) 248 end if 249 else if (present(gsqcutc_eff)) then 250 call getcut(boxcut,ecut_eff,gmet,gsqcutc_eff,Dtset%iboxcut,std_out,k0,ngfftc) 251 gsqcutf_eff=gsqcutc_eff 252 end if 253 254 ! Check that boxcut>=2 if intxc=1; otherwise intxc must be set=0. 255 if (boxcut>=zero .and. boxcut<two .and. Dtset%intxc==1) then 256 write(msg,'(a,es12.4,5a)')& 257 ' boxcut=',boxcut,' is < 2.0 => intxc must be 0;',ch10,& 258 ' Need larger ngfft to use intxc=1.',ch10,& 259 ' Action: you could increase ngfft, or decrease ecut, or put intxc=0.' 260 ABI_ERROR(msg) 261 end if 262 263 DBG_EXIT("COLL") 264 265 end subroutine pawfgr_init
m_pawfgr/pawfgr_type [ Types ]
[ Top ] [ m_pawfgr ] [ Types ]
NAME
pawfgr_type
FUNCTION
For PAW, Fine GRid parameters and related data
SOURCE
69 type,public :: pawfgr_type 70 71 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 72 ! declared in another part of ABINIT, that might need to take into account your modification. 73 74 !Integer scalars 75 76 integer :: mgfft, nfft 77 ! Values of mffft and nfft for the fine grid: 78 ! mgfft= max(ngfft(i)) [max. size of 1D FFT grid] 79 ! nfft=ngfft1*ngfft2*ngfft3 [number of pts in the FFT box] 80 81 integer :: mgfftc, nfftc 82 ! Values of mffft and nfft for the COARSE grid: 83 ! mgfftc= max(ngfftc(i)) [max. size of 1D FFT grid] 84 ! nfftc=ngfftc1*ngfftc2*ngfftc3 [number of pts in the FFT box] 85 86 integer :: usefinegrid 87 ! Flag: =1 if a double-grid is used to convert spherical data 88 ! to Fourier grid. =0 otherwise 89 90 !Integer arrays 91 92 ! MG TODO: Replace with allocatable 93 integer, pointer :: coatofin(:) 94 ! coatofin(nfftc) 95 ! Index of the points of the coarse grid on the fine grid 96 97 integer, pointer :: fintocoa(:) 98 ! fintocoa(nfft) 99 ! Index of the points of the fine grid on the coarse grid 100 ! (=0 if the point of the fine grid does not belong to the coarse grid) 101 102 integer :: ngfft(18) 103 ! ngfft(1:18)=integer array with FFT box dimensions and other 104 ! information on FFTs, for the fine rectangular grid 105 106 integer :: ngfftc(18) 107 ! ngfft(1:18)=integer array with FFT box dimensions and other 108 ! information on FFTs, for the COARSE rectangular grid 109 110 end type pawfgr_type