TABLE OF CONTENTS


ABINIT/m_wvl_psi [ Modules ]

[ Top ] [ Modules ]

NAME

  m_wvl_psi

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 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 .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_wvl_psi
27 
28   use defs_basis
29   use defs_abitypes
30   use defs_datatypes
31   use defs_wvltypes
32   use m_errors
33   use m_xmpi
34   use m_abicore
35 
36   use m_energies, only : energies_type
37   use m_pawcprj,  only : pawcprj_type, pawcprj_alloc
38   use m_abi2big,  only : wvl_vxc_abi2big, wvl_vtrial_abi2big
39 
40  implicit none
41 
42  private

ABINIT/wvl_hpsitopsi [ Functions ]

[ Top ] [ 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=informations about MPI parallelization
  proj <type(wvl_projector_type)>=projectors informations 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_nonlocalpsp(OUT)=nonlocal pseudopotential 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 informations for wavelets.

PARENTS

      vtorho

CHILDREN

      calculate_energy_and_gradient,hpsitopsi,pawcprj_alloc,wrtout

SOURCE

 90 subroutine wvl_hpsitopsi(cprj,dtset,energies,istep,mcprj,mpi_enreg,residm,wvl,xcart)
 91 
 92 #if defined HAVE_BIGDFT
 93   use BigDFT_API, only : hpsitopsi, calculate_energy_and_gradient
 94 #endif
 95 
 96 !This section has been created automatically by the script Abilint (TD).
 97 !Do not modify the following lines by hand.
 98 #undef ABI_FUNC
 99 #define ABI_FUNC 'wvl_hpsitopsi'
100 !End of the abilint section
101 
102   implicit none
103 
104 !Arguments -------------------------------
105   type(dataset_type), intent(in)         :: dtset
106   type(energies_type), intent(inout)     :: energies
107   integer, intent(in)                    :: istep,mcprj
108   type(MPI_type), intent(in)             :: mpi_enreg
109   real(dp), intent(inout)                :: residm
110   type(wvl_data), intent(inout)          :: wvl
111   !arrays
112   real(dp), intent(in) :: xcart(3, dtset%natom)
113   type(pawcprj_type),dimension(dtset%natom,mcprj),intent(inout)::cprj
114 
115 !Local variables-------------------------------
116 #if defined HAVE_BIGDFT
117   integer               :: iatom,icprj
118   character(len = 500)  :: message
119   real(dp), save        :: etotal_local
120   integer, save         :: ids
121   real(dp)              :: gnrm_zero
122   integer               :: comm,me,nproc
123   integer :: nlmn(dtset%natom)
124 #endif
125 
126 
127 ! *********************************************************************
128 
129  DBG_ENTER("COLL")
130 
131 #if defined HAVE_BIGDFT
132 
133  if(wvl%wfs%ks%orthpar%methOrtho .ne. 0) then
134    write(message,'(2a)') ch10,&
135 &   'wvl_hpsitopsi: the only orthogonalization method supported for PAW+WVL is Cholesky'
136    MSG_ERROR(message)
137  end if
138 
139  write(message, '(a,a)' ) ch10,&
140 & ' wvl_hpsitopsi: compute the new wavefunction from the trial potential.'
141  call wrtout(std_out,message,'COLL')
142 
143  comm=mpi_enreg%comm_wvl
144  me=xmpi_comm_rank(comm)
145  nproc=xmpi_comm_size(comm)
146 
147 !Initialisation of mixing parameter
148  if (istep == 1) then
149    etotal_local = real(1.d100, dp)
150    ids          = dtset%nwfshist
151  end if
152 
153 !WARNING! e_hartree is taken from the previous iteration as e_xc
154 !Update physical values
155 
156 !Precondition, minimise (DIIS or steepest descent) and ortho.
157 !Compute also the norm of the gradient.
158  if(dtset%usepaw==1) then
159    call calculate_energy_and_gradient(istep, me, nproc, wvl%wfs%GPU, dtset%wvl_nprccg, &
160 &   dtset%iscf, wvl%e%energs, wvl%wfs%ks, residm, gnrm_zero,wvl%descr%paw)
161  else
162    call calculate_energy_and_gradient(istep, me, nproc, wvl%wfs%GPU, dtset%wvl_nprccg, &
163 &   dtset%iscf, wvl%e%energs, wvl%wfs%ks, residm, gnrm_zero)
164  end if
165  etotal_local = wvl%wfs%ks%diis%energy
166 
167  if(dtset%usepaw==1) then
168    call hpsitopsi(me, nproc, istep, ids, wvl%wfs%ks,&
169 &   wvl%descr%atoms,wvl%projectors%nlpsp,&
170 &   wvl%descr%paw,xcart,energies%e_nonlocalpsp,wvl%projectors%G)
171  else
172    call hpsitopsi(me, nproc, istep, ids, wvl%wfs%ks,&
173 &   wvl%descr%atoms,wvl%projectors%nlpsp)
174  end if
175 
176  if(dtset%usepaw==1) then
177 !  PENDING : cprj should not be copied'
178 
179    ! Cannot use pawcprj_copy because cprj and paw%cprj are not the same objects
180    ! Get nlmn from bigdft cprj, and allocate our copy
181    do iatom=1,dtset%natom
182      nlmn(iatom) = wvl%descr%paw%cprj(iatom,1)%nlmn
183    end do
184    call pawcprj_alloc(cprj,mcprj,nlmn)
185 
186    do iatom=1,dtset%natom
187      do icprj=1,mcprj
188        cprj(iatom,icprj)%cp(:,:)= wvl%descr%paw%cprj(iatom,icprj)%cp(:,:)
189      end do
190    end do
191 
192 
193  end if
194 
195 #else
196  BIGDFT_NOTENABLED_ERROR()
197  if (.false.) write(std_out,*) dtset%nstep,energies%e_ewald,istep,mcprj,mpi_enreg%nproc,residm,&
198 & wvl%wfs%ks,xcart(1,1),cprj(1,1)%nlmn
199 #endif
200 
201  DBG_EXIT("COLL")
202 
203 end subroutine wvl_hpsitopsi

ABINIT/wvl_nl_gradient [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_nl_gradient

FUNCTION

 Compute the non local part of the wavefunction gradient.

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      forstr,vtorho

CHILDREN

      nonlocal_forces,wrtout,xmpi_sum

SOURCE

364 subroutine wvl_nl_gradient(grnl, mpi_enreg, natom, rprimd, wvl, xcart)
365 
366 #if defined HAVE_BIGDFT
367  use BigDFT_API, only: nonlocal_forces
368 #endif
369 
370 !This section has been created automatically by the script Abilint (TD).
371 !Do not modify the following lines by hand.
372 #undef ABI_FUNC
373 #define ABI_FUNC 'wvl_nl_gradient'
374 !End of the abilint section
375 
376  implicit none
377 
378 !Arguments ------------------------------------
379 !scalars
380  integer, intent(in) :: natom
381  type(MPI_type),intent(in) :: mpi_enreg
382  type(wvl_data),intent(inout) :: wvl
383 !arrays
384  real(dp),intent(in) :: xcart(3,natom),rprimd(3,3)
385  real(dp),intent(inout) :: grnl(3,natom)
386 
387 !Local variables-------------------------------
388 #if defined HAVE_BIGDFT
389 !scalars
390  integer :: ia,ierr,igeo,me,nproc,spaceComm
391  character(len=500) :: message
392 !arrays
393  real(dp),allocatable :: gxyz(:,:)
394  real(dp)::strtens(6,4)
395 #endif
396 
397 ! *************************************************************************
398 
399 #if defined HAVE_BIGDFT
400 
401 !Compute forces
402  write(message, '(a,a)' ) ' wvl_nl_gradient(): compute non-local part to gradient.'
403  call wrtout(std_out,message,'COLL')
404 
405 !Nullify output arrays.
406  grnl(:, :) = zero
407  strtens(:,:)=zero
408 
409  ABI_ALLOCATE(gxyz,(3, natom))
410  gxyz(:,:) = zero
411 
412 !Add the nonlocal part of the forces to grtn (BigDFT routine)
413  spaceComm=mpi_enreg%comm_wvl
414  me=xmpi_comm_rank(spaceComm)
415  nproc=xmpi_comm_size(spaceComm)
416  call nonlocal_forces(wvl%descr%Glr, &
417 & wvl%descr%h(1), wvl%descr%h(2), wvl%descr%h(3), wvl%descr%atoms, &
418 & xcart, wvl%wfs%ks%orbs, wvl%projectors%nlpsp, wvl%wfs%ks%Lzd%Glr%wfd, &
419 & wvl%wfs%ks%psi, gxyz, .true.,strtens(1,2), &
420 & proj_G=wvl%projectors%G,paw=wvl%descr%paw)
421 
422  if (nproc > 1) then
423    call xmpi_sum(gxyz, spaceComm, ierr)
424  end if
425 
426 !Forces should be in reduced coordinates.
427  do ia = 1, natom, 1
428    do igeo = 1, 3, 1
429      grnl(igeo, ia) = - rprimd(1, igeo) * gxyz(1, ia) - &
430 &     rprimd(2, igeo) * gxyz(2, ia) - &
431 &     rprimd(3, igeo) * gxyz(3, ia)
432    end do
433  end do
434  ABI_DEALLOCATE(gxyz)
435 
436 #else
437  BIGDFT_NOTENABLED_ERROR()
438  if (.false.) write(std_out,*) natom,mpi_enreg%nproc,wvl%wfs%ks,xcart(1,1),rprimd(1,1),grnl(1,1)
439 #endif
440 
441 end subroutine wvl_nl_gradient

ABINIT/wvl_psitohpsi [ Functions ]

[ Top ] [ 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

PARENTS

      afterscfloop,rhotov,setvtr,vtorho

CHILDREN

      psitohpsi,total_energies,wrtout,wvl_vtrial_abi2big,wvl_vxc_abi2big

SOURCE

231 subroutine wvl_psitohpsi(alphamix,eexctX, eexcu, ehart, ekin_sum, epot_sum, eproj_sum, eSIC_DC, &
232      & itrp, iter, iscf, me, natom, nfft, nproc, nspden, rpnrm, scf, &
233      & vexcu, wvl, wvlbigdft, xcart, xcstr,vtrial,vxc)
234 
235 #if defined HAVE_BIGDFT
236  use BigDFT_API, only: psitohpsi, KS_POTENTIAL, total_energies
237 #endif
238 
239 !This section has been created automatically by the script Abilint (TD).
240 !Do not modify the following lines by hand.
241 #undef ABI_FUNC
242 #define ABI_FUNC 'wvl_psitohpsi'
243 !End of the abilint section
244 
245  implicit none
246 
247 !Arguments-------------------------------
248 !scalars
249  integer, intent(in) :: me, nproc, itrp, iter, iscf, natom, nfft, nspden
250  real(dp), intent(in) :: alphamix
251  real(dp), intent(out) :: rpnrm
252  logical, intent(in) :: scf
253  logical, intent(in) :: wvlbigdft
254  type(wvl_data), intent(inout) :: wvl
255  real(dp), intent(inout) :: eexctX,eSIC_DC,ehart,eexcu,vexcu, ekin_sum, epot_sum, eproj_sum
256  real(dp), dimension(6), intent(out) :: xcstr
257  real(dp), intent(inout) :: xcart(3, natom)
258 !arrays
259  real(dp),intent(out), optional :: vxc(nfft,nspden)
260  real(dp),intent(out), optional :: vtrial(nfft,nspden)
261 
262 !Local variables-------------------------------
263 !scalars
264 #if defined HAVE_BIGDFT
265  character(len=500) :: message
266  integer :: linflag = 0
267  character(len=3), parameter :: unblock_comms = "OFF"
268 #endif
269 
270 ! *************************************************************************
271 
272  DBG_ENTER("COLL")
273 
274 #if defined HAVE_BIGDFT
275 
276  if(wvl%descr%atoms%npspcode(1)==7) then
277    call psitohpsi(me,nproc,wvl%descr%atoms,scf,wvl%den%denspot, &
278 &   itrp, iter, iscf, alphamix,&
279 &   wvl%projectors%nlpsp,xcart,linflag,unblock_comms, &
280 &   wvl%wfs%GPU,wvl%wfs%ks,wvl%e%energs,rpnrm,xcstr,&
281 &   wvl%projectors%G,wvl%descr%paw)
282  else
283    call psitohpsi(me,nproc,wvl%descr%atoms,scf,wvl%den%denspot, &
284 &   itrp, iter, iscf, alphamix,&
285 &   wvl%projectors%nlpsp,xcart,linflag,unblock_comms, &
286 &   wvl%wfs%GPU,wvl%wfs%ks,wvl%e%energs,rpnrm,xcstr)
287  end if
288 
289  if(scf) then
290    ehart     = wvl%e%energs%eh
291    eexcu     = wvl%e%energs%exc
292    vexcu     = wvl%e%energs%evxc
293  end if
294  eexctX    = wvl%e%energs%eexctX
295  eSIC_DC   = wvl%e%energs%evsic
296  ekin_sum  = wvl%e%energs%ekin
297  eproj_sum = wvl%e%energs%eproj
298  epot_sum  = wvl%e%energs%epot
299 
300 !Correct local potential, since in BigDFT
301 !this variable contains more terms
302 !Do the following only if sumpion==.true. in psolver_rhohxc.
303 !For the moment it is set to false.
304 
305  epot_sum=epot_sum-real(2,dp)*wvl%e%energs%eh
306  epot_sum=epot_sum-wvl%e%energs%evxc
307 
308  if(wvlbigdft) then
309    call total_energies(wvl%e%energs, iter, me)
310  end if
311 
312 !Note: if evxc is not rested here,
313 !we have to rest this from etotal in prtene, afterscfcv and etotfor.
314 !check ABINIT-6.15.1.
315 
316  if(scf) then
317    if (present(vxc)) then
318      write(message, '(a,a,a,a)' ) ch10, ' wvl_psitohpsi : but why are you copying vxc :..o('
319      call wrtout(std_out,message,'COLL')
320      call wvl_vxc_abi2big(2,vxc,wvl%den)
321    end if
322    if (wvl%den%denspot%rhov_is == KS_POTENTIAL .and. present(vtrial)) then
323      write(message, '(a,a,a,a)' ) ch10, ' wvl_psitohpsi : but why are you copying vtrial :..o('
324      call wrtout(std_out,message,'COLL')
325      call wvl_vtrial_abi2big(2,vtrial,wvl%den)
326    end if
327  end if
328 
329 #else
330  BIGDFT_NOTENABLED_ERROR()
331  if (.false.) write(std_out,*) me,nproc,itrp,iter,iscf,natom,nfft,nspden,alphamix,rpnrm,scf,&
332 & wvlbigdft,wvl%wfs%ks,eexctX,eSIC_DC,ehart,eexcu,vexcu,ekin_sum,&
333 & epot_sum,eproj_sum,xcstr(1),xcart(1,1),vxc(1,1),vtrial(1,1)
334 #endif
335 
336  DBG_EXIT("COLL")
337 
338 end subroutine wvl_psitohpsi

ABINIT/wvl_tail_corrections [ Functions ]

[ Top ] [ 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

PARENTS

      afterscfloop

CHILDREN

      calculatetailcorrection,dcopy,wrtout,xmpi_allgatherv

SOURCE

469 subroutine wvl_tail_corrections(dtset, energies, etotal, mpi_enreg, psps, wvl, xcart)
470 
471 #if defined HAVE_BIGDFT
472   use BigDFT_API, only: CalculateTailCorrection
473 #endif
474 
475 !This section has been created automatically by the script Abilint (TD).
476 !Do not modify the following lines by hand.
477 #undef ABI_FUNC
478 #define ABI_FUNC 'wvl_tail_corrections'
479 !End of the abilint section
480 
481  implicit none
482 
483 !Arguments ------------------------------------
484 !scalars
485  real(dp),intent(out) :: etotal
486  type(MPI_type),intent(in) :: mpi_enreg
487  type(dataset_type),intent(in) :: dtset
488  type(energies_type),intent(inout) :: energies
489  type(pseudopotential_type),intent(in) :: psps
490  type(wvl_data),intent(inout) :: wvl
491 !arrays
492  real(dp),intent(in) :: xcart(3,dtset%natom)
493 
494 !Local variables-------------------------------
495 #if defined HAVE_BIGDFT
496 !scalars
497  integer :: ierr,me,nbuf,nproc,nsize,spaceComm
498  real(dp) :: ekin_sum,epot_sum,eproj_sum
499  logical :: parallel
500  character(len=500) :: message
501 !arrays
502  integer :: ntails(3)
503  real(dp) :: atails(3)
504 #endif
505 
506 ! *************************************************************************
507 
508 #if defined HAVE_BIGDFT
509 
510  spaceComm=mpi_enreg%comm_wvl
511  me=xmpi_comm_rank(spaceComm)
512  nproc=xmpi_comm_size(spaceComm)
513  parallel = (nproc > 1)
514 
515 !Write a message with the total energy before tail corrections.
516  etotal = energies%e_kinetic + energies%e_hartree + energies%e_xc + &
517 & energies%e_localpsp + energies%e_corepsp + energies%e_fock+&
518 & energies%e_entropy + energies%e_elecfield + energies%e_magfield+&
519 & energies%e_ewald + energies%e_chempot + energies%e_vdw_dftd
520  if (dtset%usepaw==0) etotal = etotal + energies%e_nonlocalpsp
521  if (dtset%usepaw/=0) etotal = etotal + energies%e_paw
522  write(message,'(a,2x,e19.12)') ' Total energy before tail correction', etotal
523  call wrtout(std_out, message, 'COLL')
524 
525 !Calculate kinetic energy correction due to boundary conditions
526  nbuf = nint(dtset%tl_radius / dtset%wvl_hgrid)
527  ntails = (/ wvl%descr%Glr%d%n1, wvl%descr%Glr%d%n2, wvl%descr%Glr%d%n3 /) + 2 * nbuf
528  atails = real(ntails, dp) * dtset%wvl_hgrid
529  write(message,'(a,a,i6,a,A,A,3F12.6,A,A,3I12,A)') ch10,&
530 & ' Tail requires ',nbuf,' additional grid points around cell.', ch10, &
531 & '  | new acell:', atails, ch10, &
532 & '  | new box size for wavelets:', ntails, ch10
533  call wrtout(std_out,message,'COLL')
534  call wrtout(ab_out,message,'COLL')
535 
536 
537 !Calculate energy correction due to finite size effects
538 !---reformat potential
539  nsize = wvl%descr%Glr%d%n1i * wvl%descr%Glr%d%n2i
540  ABI_ALLOCATE(wvl%den%denspot%pot_work, (nsize * wvl%descr%Glr%d%n3i * dtset%nsppol))
541 
542  if (parallel) then
543    call xmpi_allgatherv(wvl%den%denspot%rhov, &
544 &   nsize * wvl%den%denspot%dpbox%nscatterarr(me, 2), &
545 &   wvl%den%denspot%pot_work, nsize * wvl%den%denspot%dpbox%ngatherarr(:,2), &
546 &   nsize * wvl%den%denspot%dpbox%ngatherarr(:,3),spaceComm,ierr)
547  else
548    call dcopy(wvl%descr%Glr%d%n1i * wvl%descr%Glr%d%n2i * &
549 &   wvl%descr%Glr%d%n3i * dtset%nsppol,wvl%den%denspot%rhov,1,wvl%den%denspot%pot_work,1)
550  end if
551 
552  if(dtset%usepaw==1) then
553    call CalculateTailCorrection(me, nproc, wvl%descr%atoms, dtset%tl_radius, &
554 &   wvl%wfs%ks%orbs, wvl%wfs%ks%lzd%Glr, wvl%projectors%nlpsp, dtset%tl_nprccg, &
555 &   wvl%den%denspot%pot_work, dtset%wvl_hgrid, xcart, psps%gth_params%radii_cf, &
556 &   dtset%wvl_crmult, dtset%wvl_frmult, dtset%nsppol, &
557 &   wvl%wfs%ks%psi, .false., ekin_sum, epot_sum, eproj_sum,&
558 &   wvl%projectors%G,wvl%descr%paw)
559  else
560    call CalculateTailCorrection(me, nproc, wvl%descr%atoms, dtset%tl_radius, &
561 &   wvl%wfs%ks%orbs, wvl%wfs%ks%lzd%Glr, wvl%projectors%nlpsp, dtset%tl_nprccg, &
562 &   wvl%den%denspot%pot_work, dtset%wvl_hgrid, xcart, psps%gth_params%radii_cf, &
563 &   dtset%wvl_crmult, dtset%wvl_frmult, dtset%nsppol, &
564 &   wvl%wfs%ks%psi, .false., ekin_sum, epot_sum, eproj_sum)
565  end if
566 
567  ABI_DEALLOCATE(wvl%den%denspot%pot_work)
568 
569  energies%e_kinetic = ekin_sum
570  energies%e_localpsp = epot_sum - two * energies%e_hartree
571  energies%e_nonlocalpsp = eproj_sum
572  energies%e_corepsp = zero
573  energies%e_chempot = zero
574 #if defined HAVE_BIGDFT
575  energies%e_localpsp = energies%e_localpsp - wvl%e%energs%evxc
576 #endif
577 
578  write(message,'(a,3(1x,e18.11))') ' ekin_sum,epot_sum,eproj_sum',  &
579  ekin_sum,epot_sum,eproj_sum
580  call wrtout(std_out, message, 'COLL')
581  write(message,'(a,2(1x,e18.11))') ' ehart,eexcu', &
582 & energies%e_hartree,energies%e_xc
583  call wrtout(std_out, message, 'COLL')
584 
585  etotal = energies%e_kinetic + energies%e_hartree + energies%e_xc + &
586 & energies%e_localpsp + energies%e_corepsp + energies%e_fock+&
587 & energies%e_entropy + energies%e_elecfield + energies%e_magfield+&
588 & energies%e_ewald + energies%e_vdw_dftd
589  if (dtset%usepaw==0) etotal = etotal + energies%e_nonlocalpsp
590  if (dtset%usepaw/=0) etotal = etotal + energies%e_paw
591 
592  write(message,'(a,2x,e19.12)') ' Total energy with tail correction', etotal
593  call wrtout(std_out, message, 'COLL')
594 
595 !--- End if of tail calculation
596 
597 #else
598  BIGDFT_NOTENABLED_ERROR()
599  if (.false.) write(std_out,*) etotal,mpi_enreg%nproc,dtset%nstep,energies%e_ewald,psps%npsp,&
600 & wvl%wfs%ks,xcart(1,1)
601 #endif
602 
603 end subroutine wvl_tail_corrections