From ec95336946014fd95ba1753a06b404c38679a213 Mon Sep 17 00:00:00 2001 From: fabrizio22 Date: Fri, 18 Sep 2020 17:52:46 +0200 Subject: [PATCH] XClib - some external calls removed (slater, gcxc..) --- PP/src/ppacf.f90 | 72 +++++++++++++++++++++++++------------ XClib/qe_external_calls.f90 | 67 ---------------------------------- XClib/xc_interfaces.f90 | 29 +-------------- 3 files changed, 50 insertions(+), 118 deletions(-) diff --git a/PP/src/ppacf.f90 b/PP/src/ppacf.f90 index 883e1b985..3e7a3e7bb 100644 --- a/PP/src/ppacf.f90 +++ b/PP/src/ppacf.f90 @@ -56,8 +56,8 @@ PROGRAM do_ppacf USE funct, ONLY : get_iexch, get_icorr, get_igcx, get_igcc USE funct, ONLY : set_exx_fraction, set_auxiliary_flags, & enforce_input_dft, is_libxc - USE xc_interfaces, ONLY : xc, gcxc, gcx_spin, gcc_spin, slater, & - slater_spin, xclib_set_threshold + USE xc_interfaces, ONLY : xc, xc_gcx, & ! gcxc, gcx_spin, gcc_spin, & + xclib_set_threshold !USE xc_lda_lsda, ONLY : xc USE wvfct, ONLY : npw, npwx USE environment, ONLY : environment_start, environment_end @@ -119,7 +119,8 @@ PROGRAM do_ppacf ! REAL(DP) :: grho2(2), sx(1), sc(1), scp, scm, & etxcgc, vtxcgc, segno, rh, grh2(1) - REAL(DP) :: v1x(1), v2x(1), v1c(1), v2c(1), v1cs(1,2), v1xs(1,2), v2xs(1,2) + REAL(DP) :: v1x(1,2), v2x(1,2), v1c(1,2), v2c(1,2) + ! REAL(DP) :: dq0_dq ! The derivative of the saturated REAL(DP) :: grid_cell_volume ! @@ -392,11 +393,11 @@ PROGRAM do_ppacf arhox(1,1) = ABS(rhox) IF (arhox(1,1) > vanishing_charge) THEN rs = pi34 /arhox(1,1)**third - IF (iexch == 1) THEN - CALL slater( rs, ex(1), vx(1,1) ) ! \epsilon_x,\lambda[n]=\epsilon_x[n] - ELSE - CALL xc( 1, nspin, nspin, arhox(:,1:1), ex, ec, vx(:,1:1), vc(:,1:1) ) - ENDIF + !IF (iexch == 1) THEN + ! CALL slater( rs, ex(1), vx(1,1) ) ! \epsilon_x,\lambda[n]=\epsilon_x[n] + !ELSE + CALL xc( 1, nspin, nspin, arhox(:,1:1), ex, ec, vx(:,1:1), vc(:,1:1) ) + !ENDIF etx = etx + e2*ex(1)*rhox etxlda = etxlda + e2*ex(1)*rhox grho2(1) = grho(1,ir,1)**2 + grho(2,ir,1)**2 + grho(3,ir,1)**2 @@ -434,11 +435,18 @@ PROGRAM do_ppacf IF (grho2(1)>epsg .AND. igcc/=0) THEN segno = SIGN( 1.D0, rhoout(ir,1) ) ! - CALL gcxc( 1, arhox(:,1)/ccp3, grho2/ccp8, sx, sc, v1x, v2x, v1c, v2c ) + CALL xc_gcx( 1, 1, arhox(:,1:1)/ccp3, grho(:,ir:ir,1:1)/ccp4, sx, sc, v1x(:,1:1),& + v2x(:,1:1), v1c(:,1:1), v2c(:,1:1) ) scp = sc(1) - CALL gcxc( 1, arhox(:,1)/ccm3, grho2/ccm8, sx, sc, v1x, v2x, v1c, v2c ) + CALL xc_gcx( 1, 1, arhox(:,1:1)/ccm3, grho(:,ir:ir,1:1)/ccm4, sx, sc, v1x(:,1:1),& + v2x(:,1:1), v1c(:,1:1), v2c(:,1:1) ) scm = sc(1) ! + !CALL gcxc( 1, arhox(:,1)/ccp3, grho2/ccp8, sx, sc, v1x, v2x, v1c, v2c ) + !scp = sc(1) + !CALL gcxc( 1, arhox(:,1)/ccm3, grho2/ccm8, sx, sc, v1x, v2x, v1c, v2c ) + !scm = sc(1) + ! ecgc_l = (ccp2*scp*ccp3-ccm2*scm*ccm3)/dcc*0.5_DP etcgclambda = etcgclambda + e2*ecgc_l*segno ENDIF @@ -457,7 +465,10 @@ PROGRAM do_ppacf IF ( grho2(1) > epsg ) THEN segno = SIGN( 1.D0, rhoout(ir,1) ) ! - CALL gcxc( 1, arhox(:,1), grho2, sx, sc, v1x, v2x, v1c, v2c ) + CALL xc_gcx( 1, 1, arhox(:,1:1), grho(:,ir:ir,1:1), sx, sc, v1x(:,1:1),& + v2x(:,1:1), v1c(:,1:1), v2c(:,1:1) ) + ! + !CALL gcxc( 1, arhox(:,1), grho2, sx, sc, v1x, v2x, v1c, v2c ) ! etx = etx + e2*sx(1)*segno etxgc = etxgc + e2*sx(1)*segno @@ -490,16 +501,18 @@ PROGRAM do_ppacf rhoupdw(1,1) = (rho%of_r(ir,1) + rho%of_r(ir,2) + rho_core(ir))*0.5_DP rhoupdw(1,2) = (rho%of_r(ir,1) - rho%of_r(ir,2) + rho_core(ir))*0.5_DP IF (ABS(zeta(1)) > 1.D0) zeta(1) = SIGN(1.D0, zeta(1)) - IF (iexch == 1) THEN - CALL slater_spin( arhox(1,1), zeta(1), ex(1), vx(1,1), vx(1,2) ) - ELSE - CALL xc( 1, nspin, nspin, rhoupdw, ex, ec, vx, vc ) - ENDIF + !IF (iexch == 1) THEN + ! CALL slater_spin( arhox(1,1), zeta(1), ex(1), vx(1,1), vx(1,2) ) + !ELSE + CALL xc( 1, nspin, nspin, rhoupdw, ex, ec, vx, vc ) + !ENDIF etx = etx + e2*ex(1)*rhox etxlda = etxlda+e2*ex(1)*rhox - grh2(1) = ( grho(1,ir,1) + grho(1,ir,2) )**2 + & - ( grho(2,ir,1) + grho(2,ir,2) )**2 + & - ( grho(3,ir,1) + grho(3,ir,2) )**2 + ! + ! + !grh2(1) = ( grho(1,ir,1) + grho(1,ir,2) )**2 + & + ! ( grho(2,ir,1) + grho(2,ir,2) )**2 + & + ! ( grho(3,ir,1) + grho(3,ir,2) )**2 IF (cc > 0._DP) THEN ccp = cc+dcc ccm = cc-dcc @@ -528,10 +541,17 @@ PROGRAM do_ppacf ! IF (igcc /= 0) THEN arhox(1,1) = rhox - CALL gcc_spin( 1, arhox(:,1)/ccp3, zeta, grh2/ccp8, sc, v1cs, v2c ) + ! + CALL xc_gcx( 1, 2, arhox/ccp3, grho(:,ir:ir,:)/ccp4, sx, sc, v1x, v2x, v1c, v2c ) scp = sc(1) - CALL gcc_spin( 1, arhox(:,1)/ccm3, zeta, grh2/ccm8, sc, v1cs, v2c ) + CALL xc_gcx( 1, 2, arhox/ccm3, grho(:,ir:ir,:)/ccm4, sx, sc, v1x, v2x, v1c, v2c ) scm = sc(1) + ! + !CALL gcc_spin( 1, arhox(:,1)/ccp3, zeta, grh2/ccp8, sc, v1cs, v2c ) + !scp = sc(1) + !CALL gcc_spin( 1, arhox(:,1)/ccm3, zeta, grh2/ccm8, sc, v1cs, v2c ) + !scm = sc(1) + ! ecgc_l = (ccp2*scp*ccp3-ccm2*scm*ccm3)/dcc*0.5_DP etcgclambda = etcgclambda+e2*ecgc_l arhox(1,1) = ABS(rhox) @@ -552,12 +572,18 @@ PROGRAM do_ppacf ! r_v(1,:) = rhoout(ir,:) s2_v(1,:) = grho2(:) - CALL gcx_spin( 1, r_v, s2_v, sx, v1xs, v2xs ) + ! + CALL xc_gcx( 1, 2, rhoout(ir:ir,:), grho(:,ir:ir,:), sx, sc, v1x, v2x, v1c, v2c ) + ! + !CALL gcx_spin( 1, r_v, s2_v, sx, v1xs, v2xs ) ! etx = etx + e2*sx(1) etxgc = etxgc + e2*sx(1) r_v(1,1) = rhox - CALL gcc_spin( 1, r_v(1:1,1), zeta, grh2, sc, v1cs, v2c ) + ! + + ! + !CALL gcc_spin( 1, r_v(1:1,1), zeta, grh2, sc, v1cs, v2c ) etcgc = etcgc + e2*sc(1) etc = etc + e2*sc(1) IF ( icc == ncc ) THEN diff --git a/XClib/qe_external_calls.f90 b/XClib/qe_external_calls.f90 index efa04be7f..174ffe766 100644 --- a/XClib/qe_external_calls.f90 +++ b/XClib/qe_external_calls.f90 @@ -1,70 +1,3 @@ -!----------------------------------------------------------------------- -SUBROUTINE slater_ext( rs, ex, vx ) ! - !--------------------------------------------------------------------- - !! Slater exchange with alpha=2/3 - ! - USE kind_l, ONLY: DP - ! - IMPLICIT NONE - !! - REAL(DP), INTENT(IN) :: rs - !! Wigner-Seitz radius - REAL(DP), INTENT(OUT) :: ex - !! Exchange energy (per unit volume) - REAL(DP), INTENT(OUT) :: vx - !! Exchange potential - ! - ! ... local variables - ! - REAL(DP), PARAMETER :: f = -0.687247939924714_DP, alpha = 2.0_DP/3.0_DP - ! f = -9/8*(3/2pi)^(2/3) - ex = f * alpha / rs - vx = 4._DP / 3._DP * f * alpha / rs - ! - RETURN - ! -END SUBROUTINE slater_ext -! -!----------------------------------------------------------------------- -SUBROUTINE slater_spin_ext( rho, zeta, ex, vx_up, vx_dw ) ! - !----------------------------------------------------------------------- - !! Slater exchange with alpha=2/3, spin-polarized case. - ! - USE kind_l, ONLY : DP - ! - IMPLICIT NONE - ! - REAL(DP), INTENT(IN) :: rho - !! total charge density - REAL(DP), INTENT(IN) :: zeta - !! zeta = (rho_up - rho_dw) / rho_tot - REAL(DP), INTENT(OUT) :: ex - !! exchange energy - REAL(DP), INTENT(OUT) :: vx_up, vx_dw - !! exchange potential (up, down) - ! - ! ... local variables - ! - REAL(DP), PARAMETER :: f = -1.10783814957303361d0, alpha = 2.0d0/3.0d0 - ! f = -9/8*(3/pi)^(1/3) - REAL(DP), PARAMETER :: third = 1.d0/3.d0, p43 = 4.d0/3.d0 - REAL(DP) :: exup, exdw, rho13 - ! - ! - rho13 = ( (1.d0 + zeta)*rho )**third - exup = f * alpha * rho13 - vx_up = p43 * f * alpha * rho13 - ! - rho13 = ( (1.d0 - zeta)*rho )**third - exdw = f * alpha * rho13 - vx_dw = p43 * f * alpha * rho13 - ! - ex = 0.5d0 * ( (1.d0 + zeta)*exup + (1.d0 - zeta)*exdw) - ! - RETURN - ! -END SUBROUTINE slater_spin_ext -! ! !----------------------------------------------------------------------- SUBROUTINE pw_ext( rs, iflag, ec, vc ) ! diff --git a/XClib/xc_interfaces.f90 b/XClib/xc_interfaces.f90 index 041df21b0..8ebe2110a 100644 --- a/XClib/xc_interfaces.f90 +++ b/XClib/xc_interfaces.f90 @@ -7,7 +7,7 @@ MODULE xc_interfaces ! ! LDA PUBLIC :: XC, DMXC - PUBLIC :: SLATER, SLATER_SPIN, PW, PW_SPIN + PUBLIC :: PW, PW_SPIN ! GGA PUBLIC :: XC_GCX, DGCXC PUBLIC :: GCXC, GCX_SPIN, GCC_SPIN @@ -233,33 +233,6 @@ MODULE xc_interfaces ! !---PROVISIONAL .. for cases when functional routines are called outside xc-wrappers--- ! - INTERFACE SLATER - ! - SUBROUTINE slater_ext( rs, ex, vx ) - ! - USE kind_l, ONLY: DP - IMPLICIT NONE - REAL(DP), INTENT(IN) :: rs - REAL(DP), INTENT(OUT) :: ex - REAL(DP), INTENT(OUT) :: vx - ! - END SUBROUTINE slater_ext - ! - END INTERFACE - ! - INTERFACE SLATER_SPIN - ! - SUBROUTINE slater_spin_ext( rho, zeta, ex, vx_up, vx_dw ) - ! - USE kind_l, ONLY: DP - IMPLICIT NONE - REAL(DP), INTENT(IN) :: rho, zeta - REAL(DP), INTENT(OUT) :: ex, vx_up, vx_dw - ! - END SUBROUTINE slater_spin_ext - ! - END INTERFACE - ! ! INTERFACE PW !