1! 2! Copyright (C) 2004 PWSCF group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8! 9!-------------------------------------------------------------------------- 10subroutine descreening 11 !-------------------------------------------------------------------------- 12 ! 13 ! This routine descreens the local potential and the ddd 14 ! coefficients (the latter only in the US case) 15 ! The charge density is computed with the test configuration, 16 ! not the one used to generate the pseudopotential 17 ! 18 use kinds, only: dp 19 use io_global, only : stdout, ionode, ionode_id 20 use mp, only : mp_bcast 21 use mp_world, only : world_comm 22 use radial_grids, only: ndmx 23 use ld1_parameters, only: nwfsx 24 use ld1inc, only: grid, nlcc, vxt, lsd, vpstot, vpsloc, file_screen, & 25 vh, enne, rhoc, latt, rhos, enl, & 26 nbeta, bmat, qvan, qvanl, jjs, lls, ikk, pseudotype, & 27 nwfts, enlts, octs, llts, jjts, phits, nstoaets, lpaw, & 28 which_augfun 29 implicit none 30 31 integer :: & 32 ns, & ! counter on pseudo functions 33 ns1, & ! counter on pseudo functions 34 ib,jb, & ! counter on beta functions 35 lam ! the angular momentum 36 37 real(DP) :: & 38 vaux(ndmx,2) ! work space 39 40 real(DP), external :: int_0_inf_dr ! the integral function 41 42 real(DP), parameter :: & 43 thresh= 1.e-12_dp ! threshold for selfconsistency 44 45 integer :: & 46 n, nst, iwork(nwfsx), ios, nerr 47 ! 48 ! descreening the local potential: NB: this descreening is done with 49 ! the occupation of the test configuration. This is required 50 ! for pseudopotentials with semicore states. In the other cases 51 ! a test configuration equal to the one used for pseudopotential 52 ! generation is strongly suggested 53 ! 54 do n=1,nwfts 55 enlts(n)=enl(nstoaets(n)) 56 enddo 57 ! 58 ! compute the pseudowavefunctions in the test configuration 59 ! 60 call ascheqps_drv(vpsloc, 1, thresh, .false., nerr) 61 ! 62 ! descreening the D coefficients 63 ! 64 if (pseudotype.eq.3) then 65 do ib=1,nbeta 66 do jb=1,ib 67 if (lls(ib).eq.lls(jb).and.abs(jjs(ib)-jjs(jb)).lt.1.e-7_dp) then 68 lam=lls(ib) 69 nst=(lam+1)*2 70 IF (which_augfun=='PSQ') then 71 do n=1,ikk(ib) 72 vaux(n,1)=qvanl(n,ib,jb,0)*vpsloc(n) 73 enddo 74 ELSE 75 do n=1,ikk(ib) 76 vaux(n,1)=qvan(n,ib,jb)*vpsloc(n) 77 enddo 78 ENDIF 79 bmat(ib,jb)= bmat(ib,jb) & 80 - int_0_inf_dr(vaux(1,1),grid,ikk(ib),nst) 81 endif 82 bmat(jb,ib)=bmat(ib,jb) 83 enddo 84 enddo 85 write(stdout,'(/5x,'' The ddd matrix'')') 86 do ns1=1,nbeta 87 write(stdout,'(6f12.5)') (bmat(ns1,ns),ns=1,nbeta) 88 enddo 89 endif 90 ! 91 ! descreening the local pseudopotential 92 ! 93 iwork=1 94 call chargeps(rhos,phits,nwfts,llts,jjts,octs,iwork) 95 96 call new_potential(ndmx,grid%mesh,grid,0.0_dp,vxt,lsd,nlcc,latt,enne,& 97 rhoc,rhos,vh,vaux,1) 98 99 do n=1,grid%mesh 100 vpstot(n,1)=vpsloc(n) 101 vpsloc(n)=vpsloc(n)-vaux(n,1) 102 enddo 103 104 if (file_screen .ne.' ') then 105 if (ionode) & 106 open(unit=20,file=file_screen, status='unknown', iostat=ios, err=100 ) 107100 call mp_bcast(ios, ionode_id, world_comm) 108 call errore('descreening','opening file'//file_screen,abs(ios)) 109 if (ionode) then 110 do n=1,grid%mesh 111 write(20,'(i5,7e12.4)') n,grid%r(n), vpsloc(n)+vaux(n,1), & 112 vpsloc(n), vaux(n,1), rhos(n,1) 113 enddo 114 close(20) 115 endif 116 endif 117 118 return 119end subroutine descreening 120