TABLE OF CONTENTS
- ABINIT/m_wvl_psi
- ABINIT/wvl_hpsitopsi
- ABINIT/wvl_nl_gradient
- ABINIT/wvl_psitohpsi
- ABINIT/wvl_tail_corrections
ABINIT/m_wvl_psi [ Modules ]
NAME
m_wvl_psi
FUNCTION
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DC, 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 .
SOURCE
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 #include "abi_common.h" 20 21 module m_wvl_psi 22 23 use defs_basis 24 use defs_wvltypes 25 use m_errors 26 use m_xmpi 27 use m_abicore 28 use m_dtset 29 30 use defs_datatypes, only : pseudopotential_type 31 use defs_abitypes, only : MPI_type 32 use m_energies, only : energies_type 33 use m_pawcprj, only : pawcprj_type, pawcprj_alloc 34 use m_abi2big, only : wvl_vxc_abi2big, wvl_vtrial_abi2big 35 36 implicit none 37 38 private
ABINIT/wvl_hpsitopsi [ Functions ]
NAME
wvl_hpsitopsi
FUNCTION
Heart of the wavelet resolution, compute new wavefunctions mixed witf previous by computing the gradient of the wavefunctions knowing the external potential.
INPUTS
dtset <type(dataset_type)>=input variables. istep=id of the current iteration (first is 1). mpi_enreg=information about MPI parallelization proj <type(wvl_projector_type)>=projectors information for wavelets. vtrial(dtset%nfft)=external potential. xcart(3,natom)=cartesian atomic coordinates
OUTPUT
SIDE EFFECTS
energies <type(energies_type)>=storage for energies computed here : | e_kinetic(OUT)=kinetic energy part of total energy | e_localpsp(OUT)=local pseudopotential part of total energy | e_nlpsp_vfock(OUT)=nonlocal psp + potential Fock ACE part of total energy residm=max value for gradient in the minimisation process. rhor(dtset%nfft)=electron density in r space wfs <type(wvl_projector_type)>=wavefunctions information for wavelets.
SOURCE
80 subroutine wvl_hpsitopsi(cprj,dtset,energies,istep,mcprj,mpi_enreg,residm,wvl,xcart) 81 82 #if defined HAVE_BIGDFT 83 use BigDFT_API, only : hpsitopsi, calculate_energy_and_gradient 84 #endif 85 implicit none 86 87 !Arguments ------------------------------- 88 type(dataset_type), intent(in) :: dtset 89 type(energies_type), intent(inout) :: energies 90 integer, intent(in) :: istep,mcprj 91 type(MPI_type), intent(in) :: mpi_enreg 92 real(dp), intent(inout) :: residm 93 type(wvl_data), intent(inout) :: wvl 94 !arrays 95 real(dp), intent(in) :: xcart(3, dtset%natom) 96 type(pawcprj_type),dimension(dtset%natom,mcprj),intent(inout)::cprj 97 98 !Local variables------------------------------- 99 #if defined HAVE_BIGDFT 100 integer :: iatom,icprj 101 character(len = 500) :: message 102 real(dp), save :: etotal_local 103 integer, save :: ids 104 real(dp) :: gnrm_zero 105 integer :: comm,me,nproc 106 integer :: nlmn(dtset%natom) 107 #endif 108 109 110 ! ********************************************************************* 111 112 DBG_ENTER("COLL") 113 114 #if defined HAVE_BIGDFT 115 116 if(wvl%wfs%ks%orthpar%methOrtho .ne. 0) then 117 write(message,'(2a)') ch10,& 118 & 'wvl_hpsitopsi: the only orthogonalization method supported for PAW+WVL is Cholesky' 119 ABI_ERROR(message) 120 end if 121 122 write(message, '(a,a)' ) ch10,& 123 & ' wvl_hpsitopsi: compute the new wavefunction from the trial potential.' 124 call wrtout(std_out,message,'COLL') 125 126 comm=mpi_enreg%comm_wvl 127 me=xmpi_comm_rank(comm) 128 nproc=xmpi_comm_size(comm) 129 130 !Initialisation of mixing parameter 131 if (istep == 1) then 132 etotal_local = real(1.d100, dp) 133 ids = dtset%nwfshist 134 end if 135 136 !WARNING! e_hartree is taken from the previous iteration as e_xc 137 !Update physical values 138 139 !Precondition, minimise (DIIS or steepest descent) and ortho. 140 !Compute also the norm of the gradient. 141 if(dtset%usepaw==1) then 142 call calculate_energy_and_gradient(istep, me, nproc, wvl%wfs%GPU, dtset%wvl_nprccg, & 143 & dtset%iscf, wvl%e%energs, wvl%wfs%ks, residm, gnrm_zero,wvl%descr%paw) 144 else 145 call calculate_energy_and_gradient(istep, me, nproc, wvl%wfs%GPU, dtset%wvl_nprccg, & 146 & dtset%iscf, wvl%e%energs, wvl%wfs%ks, residm, gnrm_zero) 147 end if 148 etotal_local = wvl%wfs%ks%diis%energy 149 150 if(dtset%usepaw==1) then 151 call hpsitopsi(me, nproc, istep, ids, wvl%wfs%ks,& 152 & wvl%descr%atoms,wvl%projectors%nlpsp,& 153 & wvl%descr%paw,xcart,energies%e_nlpsp_vfock,wvl%projectors%G) 154 else 155 call hpsitopsi(me, nproc, istep, ids, wvl%wfs%ks,& 156 & wvl%descr%atoms,wvl%projectors%nlpsp) 157 end if 158 159 if(dtset%usepaw==1) then 160 ! PENDING : cprj should not be copied' 161 162 ! Cannot use pawcprj_copy because cprj and paw%cprj are not the same objects 163 ! Get nlmn from bigdft cprj, and allocate our copy 164 do iatom=1,dtset%natom 165 nlmn(iatom) = wvl%descr%paw%cprj(iatom,1)%nlmn 166 end do 167 call pawcprj_alloc(cprj,mcprj,nlmn) 168 169 do iatom=1,dtset%natom 170 do icprj=1,mcprj 171 cprj(iatom,icprj)%cp(:,:)= wvl%descr%paw%cprj(iatom,icprj)%cp(:,:) 172 end do 173 end do 174 175 176 end if 177 178 #else 179 BIGDFT_NOTENABLED_ERROR() 180 if (.false.) write(std_out,*) dtset%nstep,energies%e_ewald,istep,mcprj,mpi_enreg%nproc,residm,& 181 & wvl%wfs%ks,xcart(1,1),cprj(1,1)%nlmn 182 #endif 183 184 DBG_EXIT("COLL") 185 186 end subroutine wvl_hpsitopsi
ABINIT/wvl_nl_gradient [ Functions ]
NAME
wvl_nl_gradient
FUNCTION
Compute the non local part of the wavefunction gradient.
INPUTS
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
328 subroutine wvl_nl_gradient(grnl, mpi_enreg, natom, rprimd, wvl, xcart) 329 330 #if defined HAVE_BIGDFT 331 use BigDFT_API, only: nonlocal_forces 332 #endif 333 implicit none 334 335 !Arguments ------------------------------------ 336 !scalars 337 integer, intent(in) :: natom 338 type(MPI_type),intent(in) :: mpi_enreg 339 type(wvl_data),intent(inout) :: wvl 340 !arrays 341 real(dp),intent(in) :: xcart(3,natom),rprimd(3,3) 342 real(dp),intent(inout) :: grnl(3,natom) 343 344 !Local variables------------------------------- 345 #if defined HAVE_BIGDFT 346 !scalars 347 integer :: ia,ierr,igeo,me,nproc,spaceComm 348 character(len=500) :: message 349 !arrays 350 real(dp),allocatable :: gxyz(:,:) 351 real(dp)::strtens(6,4) 352 #endif 353 354 ! ************************************************************************* 355 356 #if defined HAVE_BIGDFT 357 358 !Compute forces 359 write(message, '(a,a)' ) ' wvl_nl_gradient(): compute non-local part to gradient.' 360 call wrtout(std_out,message,'COLL') 361 362 !Nullify output arrays. 363 grnl(:, :) = zero 364 strtens(:,:)=zero 365 366 ABI_MALLOC(gxyz,(3, natom)) 367 gxyz(:,:) = zero 368 369 !Add the nonlocal part of the forces to grtn (BigDFT routine) 370 spaceComm=mpi_enreg%comm_wvl 371 me=xmpi_comm_rank(spaceComm) 372 nproc=xmpi_comm_size(spaceComm) 373 call nonlocal_forces(wvl%descr%Glr, & 374 & wvl%descr%h(1), wvl%descr%h(2), wvl%descr%h(3), wvl%descr%atoms, & 375 & xcart, wvl%wfs%ks%orbs, wvl%projectors%nlpsp, wvl%wfs%ks%Lzd%Glr%wfd, & 376 & wvl%wfs%ks%psi, gxyz, .true.,strtens(1,2), & 377 & proj_G=wvl%projectors%G,paw=wvl%descr%paw) 378 379 if (nproc > 1) then 380 call xmpi_sum(gxyz, spaceComm, ierr) 381 end if 382 383 !Forces should be in reduced coordinates. 384 do ia = 1, natom, 1 385 do igeo = 1, 3, 1 386 grnl(igeo, ia) = - rprimd(1, igeo) * gxyz(1, ia) - & 387 & rprimd(2, igeo) * gxyz(2, ia) - & 388 & rprimd(3, igeo) * gxyz(3, ia) 389 end do 390 end do 391 ABI_FREE(gxyz) 392 393 #else 394 BIGDFT_NOTENABLED_ERROR() 395 if (.false.) write(std_out,*) natom,mpi_enreg%nproc,wvl%wfs%ks,xcart(1,1),rprimd(1,1),grnl(1,1) 396 #endif 397 398 end subroutine wvl_nl_gradient
ABINIT/wvl_psitohpsi [ Functions ]
NAME
wvl_psitohpsi
FUNCTION
Compute new trial potential and calculate the hamiltionian application into hpsi.
INPUTS
mpi_enreg=information about MPI parallelization
OUTPUT
vxc(nfft,nspden)=exchange-correlation potential (hartree) vtrial(nfft,nspden)=new potential
NOTES
SOURCE
208 subroutine wvl_psitohpsi(alphamix,eexctX, eexcu, ehart, ekin_sum, epot_sum, eproj_sum, eSIC_DC, & 209 & itrp, iter, iscf, me, natom, nfft, nproc, nspden, rpnrm, scf, & 210 & vexcu, wvl, wvlbigdft, xcart, xcstr,vtrial,vxc) 211 212 #if defined HAVE_BIGDFT 213 use BigDFT_API, only: psitohpsi, KS_POTENTIAL, total_energies 214 #endif 215 implicit none 216 217 !Arguments------------------------------- 218 !scalars 219 integer, intent(in) :: me, nproc, itrp, iter, iscf, natom, nfft, nspden 220 real(dp), intent(in) :: alphamix 221 real(dp), intent(out) :: rpnrm 222 logical, intent(in) :: scf 223 logical, intent(in) :: wvlbigdft 224 type(wvl_data), intent(inout) :: wvl 225 real(dp), intent(inout) :: eexctX,eSIC_DC,ehart,eexcu,vexcu, ekin_sum, epot_sum, eproj_sum 226 real(dp), dimension(6), intent(out) :: xcstr 227 real(dp), intent(inout) :: xcart(3, natom) 228 !arrays 229 real(dp),intent(out), optional :: vxc(nfft,nspden) 230 real(dp),intent(out), optional :: vtrial(nfft,nspden) 231 232 !Local variables------------------------------- 233 !scalars 234 #if defined HAVE_BIGDFT 235 character(len=500) :: message 236 integer :: linflag = 0 237 character(len=3), parameter :: unblock_comms = "OFF" 238 #endif 239 240 ! ************************************************************************* 241 242 DBG_ENTER("COLL") 243 244 #if defined HAVE_BIGDFT 245 246 if(wvl%descr%atoms%npspcode(1)==7) then 247 call psitohpsi(me,nproc,wvl%descr%atoms,scf,wvl%den%denspot, & 248 & itrp, iter, iscf, alphamix,& 249 & wvl%projectors%nlpsp,xcart,linflag,unblock_comms, & 250 & wvl%wfs%GPU,wvl%wfs%ks,wvl%e%energs,rpnrm,xcstr,& 251 & wvl%projectors%G,wvl%descr%paw) 252 else 253 call psitohpsi(me,nproc,wvl%descr%atoms,scf,wvl%den%denspot, & 254 & itrp, iter, iscf, alphamix,& 255 & wvl%projectors%nlpsp,xcart,linflag,unblock_comms, & 256 & wvl%wfs%GPU,wvl%wfs%ks,wvl%e%energs,rpnrm,xcstr) 257 end if 258 259 if(scf) then 260 ehart = wvl%e%energs%eh 261 eexcu = wvl%e%energs%exc 262 vexcu = wvl%e%energs%evxc 263 end if 264 eexctX = wvl%e%energs%eexctX 265 eSIC_DC = wvl%e%energs%evsic 266 ekin_sum = wvl%e%energs%ekin 267 eproj_sum = wvl%e%energs%eproj 268 epot_sum = wvl%e%energs%epot 269 270 !Correct local potential, since in BigDFT 271 !this variable contains more terms 272 !Do the following only if sumpion==.true. in psolver_rhohxc. 273 !For the moment it is set to false. 274 275 epot_sum=epot_sum-real(2,dp)*wvl%e%energs%eh 276 epot_sum=epot_sum-wvl%e%energs%evxc 277 278 if(wvlbigdft) then 279 call total_energies(wvl%e%energs, iter, me) 280 end if 281 282 !Note: if evxc is not rested here, 283 !we have to rest this from etotal in prtene, afterscfcv and etotfor. 284 !check ABINIT-6.15.1. 285 286 if(scf) then 287 if (present(vxc)) then 288 write(message, '(a,a,a,a)' ) ch10, ' wvl_psitohpsi : but why are you copying vxc :..o(' 289 call wrtout(std_out,message,'COLL') 290 call wvl_vxc_abi2big(2,vxc,wvl%den) 291 end if 292 if (wvl%den%denspot%rhov_is == KS_POTENTIAL .and. present(vtrial)) then 293 write(message, '(a,a,a,a)' ) ch10, ' wvl_psitohpsi : but why are you copying vtrial :..o(' 294 call wrtout(std_out,message,'COLL') 295 call wvl_vtrial_abi2big(2,vtrial,wvl%den) 296 end if 297 end if 298 299 #else 300 BIGDFT_NOTENABLED_ERROR() 301 if (.false.) write(std_out,*) me,nproc,itrp,iter,iscf,natom,nfft,nspden,alphamix,rpnrm,scf,& 302 & wvlbigdft,wvl%wfs%ks,eexctX,eSIC_DC,ehart,eexcu,vexcu,ekin_sum,& 303 & epot_sum,eproj_sum,xcstr(1),xcart(1,1),vxc(1,1),vtrial(1,1) 304 #endif 305 306 DBG_EXIT("COLL") 307 308 end subroutine wvl_psitohpsi
ABINIT/wvl_tail_corrections [ Functions ]
NAME
wvl_tail_corrections
FUNCTION
Perform a minimization on the wavefunctions (especially the treatment of the kinetic operator) with exponentialy decreasing functions on boundaries.
INPUTS
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
420 subroutine wvl_tail_corrections(dtset, energies, etotal, mpi_enreg, psps, wvl, xcart) 421 422 #if defined HAVE_BIGDFT 423 use BigDFT_API, only: CalculateTailCorrection 424 #endif 425 implicit none 426 427 !Arguments ------------------------------------ 428 !scalars 429 real(dp),intent(out) :: etotal 430 type(MPI_type),intent(in) :: mpi_enreg 431 type(dataset_type),intent(in) :: dtset 432 type(energies_type),intent(inout) :: energies 433 type(pseudopotential_type),intent(in) :: psps 434 type(wvl_data),intent(inout) :: wvl 435 !arrays 436 real(dp),intent(in) :: xcart(3,dtset%natom) 437 438 !Local variables------------------------------- 439 #if defined HAVE_BIGDFT 440 !scalars 441 integer :: ierr,me,nbuf,nproc,nsize,spaceComm 442 real(dp) :: ekin_sum,epot_sum,eproj_sum 443 logical :: parallel 444 character(len=500) :: message 445 !arrays 446 integer :: ntails(3) 447 real(dp) :: atails(3) 448 #endif 449 450 ! ************************************************************************* 451 452 #if defined HAVE_BIGDFT 453 454 spaceComm=mpi_enreg%comm_wvl 455 me=xmpi_comm_rank(spaceComm) 456 nproc=xmpi_comm_size(spaceComm) 457 parallel = (nproc > 1) 458 459 !Write a message with the total energy before tail corrections. 460 etotal = energies%e_kinetic + energies%e_hartree + energies%e_xc + & 461 & energies%e_localpsp + energies%e_corepsp + energies%e_fock+& 462 & energies%e_entropy + energies%e_elecfield + energies%e_magfield+& 463 & energies%e_ewald + energies%e_chempot + energies%e_vdw_dftd 464 if (dtset%usepaw==0) etotal = etotal + energies%e_nlpsp_vfock 465 if (dtset%usepaw/=0) etotal = etotal + energies%e_paw 466 write(message,'(a,2x,e19.12)') ' Total energy before tail correction', etotal 467 call wrtout(std_out, message, 'COLL') 468 469 !Calculate kinetic energy correction due to boundary conditions 470 nbuf = nint(dtset%tl_radius / dtset%wvl_hgrid) 471 ntails = (/ wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3 /) + 2 * nbuf 472 atails = real(ntails, dp) * dtset%wvl_hgrid 473 write(message,'(a,a,i6,a,A,A,3F12.6,A,A,3I12,A)') ch10,& 474 & ' Tail requires ',nbuf,' additional grid points around cell.', ch10, & 475 & ' | new acell:', atails, ch10, & 476 & ' | new box size for wavelets:', ntails, ch10 477 call wrtout(std_out,message,'COLL') 478 call wrtout(ab_out,message,'COLL') 479 480 481 !Calculate energy correction due to finite size effects 482 !---reformat potential 483 nsize = wvl%descr%Glr%d%n1i * wvl%descr%Glr%d%n2i 484 ABI_MALLOC(wvl%den%denspot%pot_work, (nsize * wvl%descr%Glr%d%n3i * dtset%nsppol)) 485 486 if (parallel) then 487 call xmpi_allgatherv(wvl%den%denspot%rhov, & 488 & nsize * wvl%den%denspot%dpbox%nscatterarr(me, 2), & 489 & wvl%den%denspot%pot_work, nsize * wvl%den%denspot%dpbox%ngatherarr(:,2), & 490 & nsize * wvl%den%denspot%dpbox%ngatherarr(:,3),spaceComm,ierr) 491 else 492 call dcopy(wvl%descr%Glr%d%n1i * wvl%descr%Glr%d%n2i * & 493 & wvl%descr%Glr%d%n3i * dtset%nsppol,wvl%den%denspot%rhov,1,wvl%den%denspot%pot_work,1) 494 end if 495 496 if(dtset%usepaw==1) then 497 call CalculateTailCorrection(me, nproc, wvl%descr%atoms, dtset%tl_radius, & 498 & wvl%wfs%ks%orbs, wvl%wfs%ks%lzd%Glr, wvl%projectors%nlpsp, dtset%tl_nprccg, & 499 & wvl%den%denspot%pot_work, dtset%wvl_hgrid, xcart, psps%gth_params%radii_cf, & 500 & dtset%wvl_crmult, dtset%wvl_frmult, dtset%nsppol, & 501 & wvl%wfs%ks%psi, .false., ekin_sum, epot_sum, eproj_sum,& 502 & wvl%projectors%G,wvl%descr%paw) 503 else 504 call CalculateTailCorrection(me, nproc, wvl%descr%atoms, dtset%tl_radius, & 505 & wvl%wfs%ks%orbs, wvl%wfs%ks%lzd%Glr, wvl%projectors%nlpsp, dtset%tl_nprccg, & 506 & wvl%den%denspot%pot_work, dtset%wvl_hgrid, xcart, psps%gth_params%radii_cf, & 507 & dtset%wvl_crmult, dtset%wvl_frmult, dtset%nsppol, & 508 & wvl%wfs%ks%psi, .false., ekin_sum, epot_sum, eproj_sum) 509 end if 510 511 ABI_FREE(wvl%den%denspot%pot_work) 512 513 energies%e_kinetic = ekin_sum 514 energies%e_localpsp = epot_sum - two * energies%e_hartree 515 energies%e_nlpsp_vfock = eproj_sum 516 energies%e_corepsp = zero 517 energies%e_chempot = zero 518 #if defined HAVE_BIGDFT 519 energies%e_localpsp = energies%e_localpsp - wvl%e%energs%evxc 520 #endif 521 522 write(message,'(a,3(1x,e18.11))') ' ekin_sum,epot_sum,eproj_sum', & 523 ekin_sum,epot_sum,eproj_sum 524 call wrtout(std_out, message, 'COLL') 525 write(message,'(a,2(1x,e18.11))') ' ehart,eexcu', & 526 & energies%e_hartree,energies%e_xc 527 call wrtout(std_out, message, 'COLL') 528 529 etotal = energies%e_kinetic + energies%e_hartree + energies%e_xc + & 530 & energies%e_localpsp + energies%e_corepsp + energies%e_fock+& 531 & energies%e_entropy + energies%e_elecfield + energies%e_magfield+& 532 & energies%e_ewald + energies%e_vdw_dftd 533 if (dtset%usepaw==0) etotal = etotal + energies%e_nlpsp_vfock 534 if (dtset%usepaw/=0) etotal = etotal + energies%e_paw 535 536 write(message,'(a,2x,e19.12)') ' Total energy with tail correction', etotal 537 call wrtout(std_out, message, 'COLL') 538 539 !--- End if of tail calculation 540 541 #else 542 BIGDFT_NOTENABLED_ERROR() 543 if (.false.) write(std_out,*) etotal,mpi_enreg%nproc,dtset%nstep,energies%e_ewald,psps%npsp,& 544 & wvl%wfs%ks,xcart(1,1) 545 #endif 546 547 end subroutine wvl_tail_corrections