1! 2! Copyright (C) 2004-201 Quantum ESPRESSO 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 lderiv 11 !--------------------------------------------------------------- 12 ! 13 ! numerical integration of the radial schroedinger equation 14 ! computing logarithmic derivatives for Coulomb potential 15 ! 16 ! 17 USE kinds, ONLY : dp 18 USE radial_grids, ONLY : ndmx 19 USE io_global, ONLY : stdout 20 USE mp, ONLY : mp_bcast 21 USE ld1_parameters, ONLY : nwfsx 22 USE ld1inc, ONLY : file_logder, grid, vpot, rel, nspin, nld, zed, & 23 npte, deld, eminld, emaxld, rlderiv 24 25 IMPLICIT NONE 26 27 INTEGER :: & 28 lam, & ! the angular momentum 29 ikrld, & ! index of matching radius 30 nc, & ! counter on logarithmic derivatives 31 nin, & ! integer variable for lschps 32 is, & ! counter on spin 33 nstop, & ! integer to monitor errors 34 ios, & ! used for I/O control 35 n,ie ! generic counter 36 37 real(DP) :: & 38 aux(ndmx), & ! the square of the wavefunction 39 aux_dir(ndmx,2), & ! the square of the wavefunction 40 ze2, & ! the nuclear charge in Ry units 41 e, & ! the eigenvalue 42 j, & ! total angular momentum for log_der 43 thrdum = 0.0_dp ! real variable (not used) for lschps 44 45 real(DP), EXTERNAL :: compute_log 46 47 real(DP), ALLOCATABLE :: & 48 ene(:), & ! the energy grid 49 dlchi(:, :) ! the logarithmic derivative 50 51 CHARACTER(len=256) :: flld 52 53 IF (nld == 0 .or. file_logder == ' ') RETURN 54 IF (nld > nwfsx) CALL errore('lderiv','nld is too large',1) 55 56 ze2=-zed*2.0_dp 57 58 DO n=1,grid%mesh 59 IF (grid%r(n) > rlderiv) GOTO 10 60 ENDDO 61 CALL errore('lderiv','wrong rlderiv?',1) 6210 ikrld = n-1 63 WRITE(stdout,'(5x,''Computing logarithmic derivative in'',f10.5)') & 64 (grid%r(ikrld)+grid%r(ikrld+1))*0.5_dp 65 66 npte= (emaxld-eminld)/deld + 1 67 ALLOCATE ( dlchi(npte, nld) ) 68 ALLOCATE ( ene(npte) ) 69 DO ie=1,npte 70 ene(ie)= eminld+deld*(ie-1) 71 ENDDO 72 73 DO is=1,nspin 74 DO nc=1,nld 75 IF (rel < 2) THEN 76 lam=nc-1 77 j=0.0_dp 78 ELSE 79 lam=nc/2 80 IF (mod(nc,2)==0) j=lam-0.5_dp 81 IF (mod(nc,2)==1) j=lam+0.5_dp 82 ENDIF 83 DO ie=1,npte 84 e=ene(ie) 85 ! 86 ! integrate outward up to ikrld+1 87 ! 88 IF (rel == 1) THEN 89 nin = ikrld+5 90 CALL lschps (3, zed, thrdum, grid, nin, 1, lam, e, & 91 vpot(1,is), aux, nstop ) 92 ELSEIF (rel == 2) THEN 93 CALL dir_outward(ndmx,ikrld+5,lam,j,e,grid%dx,& 94 aux_dir,grid%r,grid%rab,vpot(1,is)) 95 aux(:)=aux_dir(:,1) 96 ELSE 97 CALL intref(lam,e,ikrld+5,grid,vpot(1,is),ze2,aux) 98 ENDIF 99 ! 100 ! compute the logarithmic derivative and save in dlchi 101 ! 102 dlchi(ie, nc) = compute_log(aux(ikrld-3),grid%r(ikrld),grid%dx) 103 ENDDO 104 ENDDO 105 106 IF (nspin == 2 .and. is == 1) THEN 107 flld = trim(file_logder)//'up' 108 ELSEIF (nspin == 2 .and. is == 2) THEN 109 flld = trim(file_logder)//'dw' 110 ELSE 111 flld = trim(file_logder) 112 ENDIF 113 114 CALL write_efun(flld,dlchi,ene,npte,nld) 115 ! 116 ENDDO 117 118 DEALLOCATE (ene) 119 DEALLOCATE (dlchi) 120 RETURN 121END SUBROUTINE lderiv 122 123