Commit 6a3fd460 authored by Samuel Poncé's avatar Samuel Poncé

Update of superconductivity

Issue related to ef and ef0 being different. Courtesy of R. Margine.
Also change from eps2/ry2meV to eps8.
Update of superconductivity reference.
Addition of epf17(:,:,:,:) = czero initialization in the q-loop.
parent 51e3f0d0
......@@ -225,6 +225,7 @@
WRITE(iudecayH,*) wslen(ir) * celldm (1) * bohr2ang, tmp
!
ENDDO
!
! RMDB
DO ir = 1, nrr
DO jbnd = 1, nbndsub
......@@ -234,7 +235,8 @@
ENDDO
ENDDO
ENDDO
CLOSE(300)
!
CLOSE(iudecayH)
!
ENDIF
CALL mp_barrier(inter_pool_comm)
......@@ -471,11 +473,11 @@
!DO ir = 1, nrr
! DO jbnd = 1, nbndsub
! DO ibnd = 1, nbndsub
! WRITE(300,'(5I5,6F12.6)') irvec(:,ir), ibnd, jbnd, cdmew(:,ibnd,jbnd,ir)
! WRITE(iudecayP,'(5I5,6F12.6)') irvec(:,ir), ibnd, jbnd, cdmew(:,ibnd,jbnd,ir)
! ENDDO
! ENDDO
!ENDDO
!
!
CLOSE(iudecayP)
ENDIF
CALL mp_barrier(inter_pool_comm)
......@@ -1171,6 +1173,7 @@
WRITE(iuwane, *) wslen(ir) * celldm(1) * bohr2ang, tmp
!
ENDDO
!
! RMDB
!DO ir = 1, nrr
! DO jbnd = 1, nbndsub
......@@ -1179,6 +1182,7 @@
! ENDDO
! ENDDO
!ENDDO
!
CLOSE(iuwane)
ENDIF
!
......
......@@ -98,11 +98,10 @@
! nkfs : nr. of irreducible k-points within the Fermi shell on the fine mesh
! nbndfs : nr. of electronic bands within the Fermi shell
!
INTEGER, ALLOCATABLE :: equivk(:), ixkff(:), ixkf(:), ixkqf(:,:), ixqfs(:,:), nqfs(:)
INTEGER, ALLOCATABLE :: ixkff(:), ixkf(:), ixkqf(:,:), ixqfs(:,:), nqfs(:)
!
! nkf = nr of irreducible k-points on the fine grid, if mp_mesh_k = .true.
! nkf = total nr of k-points on the fine grid, otherwise
! equivk : index of equivalent k-points on the fine k-grid equivk(nkf)
! ixkff : index of k-point on the full k-grid ixkff(nkftot)
! ixkf : index of k-point on the irreducible k-grid within the Fermi shell ixkf(nkf)
! ixkqf : index k+q or k-q on the irreducilble k-grid within the Fermi shell ixkqf(nkfs,nqftot)
......
......@@ -42,8 +42,8 @@
specfun_pl, lindabs, mob_maxiter, use_ws, &
epmatkqread, selecqread
USE noncollin_module, ONLY : noncolin
USE constants_epw, ONLY : ryd2ev, ryd2mev, one, two, eps2, zero, czero, &
twopi, ci, kelvin2eV, eps6
USE constants_epw, ONLY : ryd2ev, ryd2mev, one, two, zero, czero, &
twopi, ci, kelvin2eV, eps6, eps8
USE io_files, ONLY : prefix, diropn, tmp_dir
USE io_global, ONLY : stdout, ionode
USE io_epw, ONLY : lambda_phself, linewidth_phself, iunepmatwe, &
......@@ -248,8 +248,6 @@
!! External function to calculate the fermi energy
REAL(kind=DP), EXTERNAL :: efermig_seq
!! Same but in sequential
REAL(kind=DP), PARAMETER :: eps = 0.01/ryd2mev
!! Tolerence
REAL(kind=DP), ALLOCATABLE :: etf_all(:,:)
!! Eigen-energies on the fine grid collected from all pools in parallel case
REAL(kind=DP), ALLOCATABLE :: w2 (:)
......@@ -933,12 +931,21 @@
CLOSE(iunrestart)
ENDIF
CALL mp_bcast(iq_restart, ionode_id, world_comm )
#if defined(__MPI)
CALL MPI_BCAST( ind_tot, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
CALL MPI_BCAST( ind_totcb, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
CALL MPI_BCAST( lrepmatw2, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
CALL MPI_BCAST( lrepmatw4, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
CALL MPI_BCAST( lrepmatw5, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
CALL MPI_BCAST( lrepmatw6, 1, MPI_OFFSET, ionode_id, world_comm, ierr)
#else
CALL mp_bcast( ind_tot, ionode_id, world_comm )
CALL mp_bcast( ind_totcb, ionode_id, world_comm )
CALL mp_bcast( lrepmatw2, ionode_id, world_comm )
CALL mp_bcast( lrepmatw4, ionode_id, world_comm )
CALL mp_bcast( lrepmatw5, ionode_id, world_comm )
CALL mp_bcast( lrepmatw6, ionode_id, world_comm )
#endif
IF( ierr /= 0 ) CALL errore( 'ephwann_shuffle', 'error in MPI_BCAST',1 )
!
! Now, the iq_restart point has been done, so we need to do the next one except if last
......@@ -953,7 +960,7 @@
!
DO iqq = iq_restart, totq
! This needs to be uncommented.
!epf17(:,:,:,:) = czero
epf17(:,:,:,:) = czero
!
iq = selecq(iqq)
!
......@@ -1134,7 +1141,7 @@
!
CALL compute_umn_f( nbndsub, cufkk, cufkq, bmatf )
!
IF ( (abs(xxq(1)) > eps) .or. (abs(xxq(2)) > eps) .or. (abs(xxq(3)) > eps) ) THEN
IF ( (abs(xxq(1)) > eps8) .or. (abs(xxq(2)) > eps8) .or. (abs(xxq(3)) > eps8) ) THEN
!
CALL cryst_to_cart (1, xxq, bg, 1)
CALL rgd_blk_epw_fine(nq1, nq2, nq3, xxq, uf, epmatf, &
......@@ -1173,12 +1180,12 @@
IF (specfun_ph ) CALL spectral_func_ph( iqq, iq, totq )
IF (specfun_pl .and. .not. vme ) CALL spectral_func_pl_q( iqq, iq, totq )
IF (ephwrite) THEN
IF ( iqq == 1 ) THEN
IF ( iq == 1 ) THEN
CALL kmesh_fine
CALL kqmap_fine
ENDIF
CALL write_ephmat( iqq, iq, totq )
CALL count_kpoints( iqq )
CALL write_ephmat( iq )
CALL count_kpoints( iq )
ENDIF
!
IF (.NOT. scatread) THEN
......@@ -1205,7 +1212,7 @@
etf (ibnd, ikq) = etf (ibnd, ikq) + scissor
ENDDO
ENDDO
IF ( iqq == 1 ) THEN
IF ( iq == 1 ) THEN
WRITE(stdout, '(5x,"Applying a scissor shift of ",f9.5," eV to the conduction states")' ) scissor * ryd2ev
ENDIF
ENDIF
......@@ -1894,9 +1901,9 @@
!! Finds the Fermi energy - Gaussian Broadening
!! (see Methfessel and Paxton, PRB 40, 3616 (1989 )
!!
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
USE constants, ONLY : rytoev
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
USE constants_epw, ONLY : ryd2ev, eps10
!
implicit none
!
......@@ -1924,7 +1931,6 @@
!
! Local variables
!
real(DP), parameter :: eps= 1.0d-10
integer, parameter :: maxiter = 300
real(DP) :: Ef, Eup, Elw, sumkup, sumklw, sumkmid
real(DP), external:: sumkg_seq
......@@ -1945,15 +1951,15 @@
!
sumkup = sumkg_seq (et, nbnd, nks, wk, Degauss, Ngauss, Eup, is, isk)
sumklw = sumkg_seq (et, nbnd, nks, wk, Degauss, Ngauss, Elw, is, isk)
if ( (sumkup - nelec) < -eps .or. (sumklw - nelec) > eps ) &
if ( (sumkup - nelec) < -eps10 .or. (sumklw - nelec) > eps10 ) &
call errore ('efermig_seq', 'internal error, cannot bracket Ef', 1)
DO i = 1, maxiter
Ef = (Eup + Elw) / 2.d0
sumkmid = sumkg_seq (et, nbnd, nks, wk, Degauss, Ngauss, Ef, is, isk)
if (abs (sumkmid-nelec) < eps) then
if (abs (sumkmid-nelec) < eps10) then
efermig_seq = Ef
return
elseif ( (sumkmid-nelec) < -eps) then
elseif ( (sumkmid-nelec) < -eps10) then
Elw = Ef
else
Eup = Ef
......@@ -1962,7 +1968,7 @@
IF (is /= 0) WRITE(stdout, '(5x,"Spin Component #",i3)') is
WRITE( stdout, '(5x,"Warning: too many iterations in bisection"/ &
& 5x,"Ef = ",f10.6," sumk = ",f10.6," electrons")' ) &
Ef * rytoev, sumkmid
Ef * ryd2ev, sumkmid
!
efermig_seq = Ef
RETURN
......@@ -2121,7 +2127,7 @@
USE kinds, ONLY : DP
USE io_global, ONLY : stdout
USE elph2, ONLY : etf, nkf, wkf
USE constants_epw, ONLY : ryd2ev, bohr2ang, ang2cm
USE constants_epw, ONLY : ryd2ev, bohr2ang, ang2cm, eps5
USE noncollin_module, ONLY : noncolin
USE pwcom, ONLY : nelec
USE epwcom, ONLY : int_mob, nbndsub, ncarrier, &
......@@ -2182,8 +2188,6 @@
!! Electron carrier density
REAL(DP), PARAMETER :: maxarg = 200.d0
!! Maximum value for the argument of the exponential
REAL(DP), PARAMETER :: eps= 1.0d-5
!! Tolerence to be converged [relative]
!
Ef = 0.0d0
!
......@@ -2316,11 +2320,11 @@
! WRITE(stdout,*),'electron_density ',electron_density * (1.0d0/omega) * (bohr2ang * ang2cm )**(-3)
rel_err = (hole_density-electron_density) / hole_density
!
IF (abs (rel_err) < eps) THEN
IF (abs (rel_err) < eps5) THEN
fermi_exp = Ef
fermicarrier = (- log (fermi_exp) * temp) + evbm
return
ELSEIF( (rel_err) > eps) THEN
ELSEIF( (rel_err) > eps5) THEN
Elw = Ef
ELSE
Eup = Ef
......@@ -2360,11 +2364,11 @@
!
rel_err = (electron_density-ncarrier) / electron_density
!
IF (abs (rel_err) < eps) THEN
IF (abs (rel_err) < eps5) THEN
fermi_exp = Ef
fermicarrier = ecbm - ( log (fermi_exp) * temp)
return
ELSEIF( (rel_err) > eps) THEN
ELSEIF( (rel_err) > eps5) THEN
Eup = Ef
ELSE
Elw = Ef
......@@ -2404,11 +2408,11 @@
! In this case ncarrier is a negative number
rel_err = (hole_density-ABS(ncarrier)) / hole_density
!
IF (abs (rel_err) < eps) THEN
IF (abs (rel_err) < eps5) THEN
fermi_exp = Ef
fermicarrier = (- log (fermi_exp) * temp) + evbm
return
ELSEIF( (rel_err) > eps) THEN
ELSEIF( (rel_err) > eps5) THEN
Elw = Ef
ELSE
Eup = Ef
......
......@@ -43,8 +43,8 @@
nqf2, nqf3, mp_mesh_k, restart, ncarrier, plselfen, &
specfun_pl, use_ws, epmatkqread, selecqread
USE noncollin_module, ONLY : noncolin
USE constants_epw, ONLY : ryd2ev, ryd2mev, one, two, eps2, zero, czero, cone, &
twopi, ci, kelvin2eV
USE constants_epw, ONLY : ryd2ev, ryd2mev, one, two, zero, czero, cone, &
twopi, ci, kelvin2eV, eps8
USE io_files, ONLY : prefix, diropn, tmp_dir
USE io_global, ONLY : stdout, ionode
USE io_epw, ONLY : lambda_phself, linewidth_phself, iunepmatwe, &
......@@ -249,8 +249,6 @@
!! External function to calculate the fermi energy
REAL(kind=DP), EXTERNAL :: efermig_seq
!! Same but in sequential
REAL(kind=DP), PARAMETER :: eps = 0.01/ryd2mev
!! Tolerence
REAL(kind=DP), ALLOCATABLE :: etf_all(:,:)
!! Eigen-energies on the fine grid collected from all pools in parallel case
REAL(kind=DP), ALLOCATABLE :: w2 (:)
......@@ -927,6 +925,8 @@
! -----------------------------------------------------------------------------
!
DO iqq = iq_restart, totq
! This needs to be uncommented.
epf17(:,:,:,:) = czero
!
iq = selecq(iqq)
!
......@@ -1114,7 +1114,7 @@
!
CALL compute_umn_f( nbndsub, cufkk, cufkq, bmatf )
!
IF ( (abs(xxq(1)) > eps) .or. (abs(xxq(2)) > eps) .or. (abs(xxq(3)) > eps) ) THEN
IF ( (abs(xxq(1)) > eps8) .or. (abs(xxq(2)) > eps8) .or. (abs(xxq(3)) > eps8) ) THEN
!
CALL cryst_to_cart (1, xxq, bg, 1)
CALL rgd_blk_epw_fine_mem(imode, nq1, nq2, nq3, xxq, uf, epmatlrT(:,:,imode,ik), &
......@@ -1173,12 +1173,12 @@
IF (specfun_ph ) CALL spectral_func_ph( iqq, iq, totq )
IF (specfun_pl .and. .not. vme ) CALL spectral_func_pl_q( iqq, iq, totq )
IF (ephwrite) THEN
IF ( iqq == 1 ) THEN
IF ( iq == 1 ) THEN
CALL kmesh_fine
CALL kqmap_fine
ENDIF
CALL write_ephmat( iqq, iq, totq )
CALL count_kpoints( iqq )
CALL write_ephmat( iq )
CALL count_kpoints( iq )
ENDIF
!
! Conductivity ---------------------------------------------------------
......
......@@ -136,8 +136,6 @@
REAL(kind=DP), EXTERNAL :: w0gauss
!! This function computes the derivative of the Fermi-Dirac function
!! It is therefore an approximation for a delta function
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
REAL(kind=DP), ALLOCATABLE :: xkf_all(:,:)
!! Collect k-point coordinate from all pools in parallel case
REAL(kind=DP), ALLOCATABLE :: etf_all(:,:)
......
This diff is collapsed.
......@@ -28,7 +28,7 @@
xqf, zi_allvb, zi_allcb, xkf, wkf, dmef, vmef, nqf
USE transportcom, ONLY : transp_temp, lower_bnd
USE constants_epw, ONLY : zero, one, two, pi, ryd2mev, kelvin2eV, ryd2ev, &
eps6, eps10, bohr2ang, ang2cm
eps6, eps10, bohr2ang, ang2cm, eps4, eps8
USE io_files, ONLY : prefix, diropn, tmp_dir
USE mp, ONLY : mp_barrier, mp_sum, mp_bcast
USE mp_global, ONLY : world_comm, my_pool_id, npool
......@@ -193,10 +193,6 @@
!! Compute the approximate theta function. Here computes Fermi-Dirac
REAL(KIND=DP), EXTERNAL :: w0gauss
!! The derivative of wgauss: an approximation to the delta function
REAL(kind=DP), PARAMETER :: eps = 1.d-4
!! Tolerence parameter for the velocity
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
REAL(kind=DP) :: carrier_density, fnk, inv_cell
!
inv_cell = 1.0d0/omega
......@@ -317,8 +313,8 @@
!
! In case of q=\Gamma, then the short-range = the normal g. We therefore
! need to treat it like the normal g with abs(g).
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp, KIND=DP )
......@@ -375,7 +371,7 @@
!
! In this case we are also computing the scattering rate for another Fermi level position
! This is used to compute both the electron and hole mobility at the same time.
IF ( ABS(efcb(itemp)) > eps ) THEN
IF ( ABS(efcb(itemp)) > eps4 ) THEN
!
DO ibnd = 1, ibndmax-ibndmin+1
!
......@@ -418,8 +414,8 @@
!
! In case of q=\Gamma, then the short-range = the normal g. We therefore
! need to treat it like the normal g with abs(g).
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp, KIND=DP)
......
......@@ -29,7 +29,7 @@
USE epwcom, ONLY : nbndsub
USE elph2, ONLY : etf, ibndmin, ibndmax, nkqf, xqf, &
nkf, epf17, xkf, nkqtotf, wf
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero, eps8
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : inter_pool_comm
USE mp_world, ONLY : mpime
......@@ -88,8 +88,6 @@
!! g vectex accross all pools
REAL(kind=DP), ALLOCATABLE :: epc_sym(:,:,:)
!! Temporary g-vertex for each pool
REAL(kind=DP), PARAMETER :: eps = 0.01/ryd2mev
!! Tolerence to be degenerate
!
! find the bounds of k-dependent arrays in the parallel case in each pool
CALL fkbounds( nkqtotf/2, lower_bnd, upper_bnd )
......@@ -135,7 +133,7 @@
n = 0
DO mu = 1, nmodes
w_2 = wf(mu,iq)
IF ( abs(w_2-w_1).lt.eps ) THEN
IF ( abs(w_2-w_1) < eps8 ) THEN
n = n + 1
g2 = g2 + epc(ibnd,jbnd,mu,ik+lower_bnd-1)*epc(ibnd,jbnd,mu,ik+lower_bnd-1)
ENDIF
......@@ -155,7 +153,7 @@
n = 0
DO pbnd = 1, ibndmax-ibndmin+1
w_2 = etf (ibndmin-1+pbnd, ikk)
IF ( abs(w_2-w_1).lt.eps ) THEN
IF ( abs(w_2-w_1) < eps8 ) THEN
n = n + 1
g2 = g2 + epc(pbnd,jbnd,nu,ik+lower_bnd-1)*epc(pbnd,jbnd,nu,ik+lower_bnd-1)
ENDIF
......@@ -176,7 +174,7 @@
n = 0
DO pbnd = 1, ibndmax-ibndmin+1
w_2 = etf(ibndmin-1+pbnd, ikq)
IF ( abs(w_2-w_1).lt.eps ) then
IF ( abs(w_2-w_1) < eps8 ) then
n = n + 1
g2 = g2 + epc(ibnd,pbnd,nu,ik+lower_bnd-1)*epc(ibnd,pbnd,nu,ik+lower_bnd-1)
ENDIF
......
......@@ -23,7 +23,7 @@
USE elph2, ONLY : epmatq, zstar, epsi, bmat
USE epwcom, ONLY : lpolar
USE modes, ONLY : nmodes
USE constants_epw, ONLY : cone, czero, one, ryd2mev
USE constants_epw, ONLY : cone, czero, one, ryd2mev, eps8
USE pwcom, ONLY : nbnd, nks
USE ions_base, ONLY : amass, ityp
USE phcom, ONLY : nq1, nq2, nq3
......@@ -65,7 +65,6 @@
!
REAL(kind=DP) :: massfac
!! square root of mass
REAL(kind=DP), parameter :: eps = 0.01/ryd2mev
!
COMPLEX(kind=DP) :: eptmp( nmodes)
!! temporary e-p matrix elements
......@@ -178,7 +177,7 @@
epmatq_opt(ibnd, jbnd, ik, :), 1, czero, eptmp, 1 )
!
IF (lpolar) THEN
IF ( (abs(xq(1)).gt.eps) .or. (abs(xq(2)).gt.eps) .or. (abs(xq(3)).gt.eps) ) THEN
IF ( (abs(xq(1)) > eps8) .or. (abs(xq(2)) > eps8) .or. (abs(xq(3)) > eps8) ) THEN
CALL rgd_blk_epw (nq1, nq2, nq3, xq, cz2t, eptmp, &
nmodes, epsi, zstar, bmat(ibnd,jbnd,ik,iq), -one)
ENDIF
......
......@@ -43,7 +43,7 @@
sigmar_all, sigmai_all, sigmai_mode, zi_all, efnew
USE transportcom, ONLY : lower_bnd
USE control_flags, ONLY : iverbosity
USE constants_epw, ONLY : ryd2mev, one, ryd2ev, two, zero, pi, ci, eps6
USE constants_epw, ONLY : ryd2mev, one, ryd2ev, two, zero, pi, ci, eps6, eps8
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : inter_pool_comm
USE mp_world, ONLY : mpime
......@@ -135,8 +135,6 @@
REAL(kind=DP), external :: w0gauss
!! This function computes the derivative of the Fermi-Dirac function
!! It is therefore an approximation for a delta function
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
REAL(kind=DP), ALLOCATABLE :: xkf_all(:,:)
!! Collect k-point coordinate from all pools in parallel case
REAL(kind=DP), ALLOCATABLE :: etf_all(:,:)
......@@ -262,8 +260,8 @@
! with hbar = 1 and M already contained in the eigenmodes
! g2 is Ry^2, wkf must already account for the spin factor
!
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp )
......
......@@ -37,7 +37,7 @@
use elph2, ONLY : epf17, ibndmax, ibndmin, etf, wkf, xqf, wqf, nkqf, &
nkf, wf, nkqtotf, xqf, lambda_all, lambda_v_all, &
dmef, vmef, gamma_all,gamma_v_all, efnew
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero, pi, eps4, eps6
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero, pi, eps4, eps6, eps8
use mp, ONLY : mp_barrier, mp_sum
use mp_global, ONLY : inter_pool_comm
!
......@@ -153,8 +153,6 @@
!! It is therefore an approximation for a delta function
REAL(kind=DP), external :: efermig
!! Return the fermi energy
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
!
IF ( iq == 1 ) THEN
WRITE(stdout,'(/5x,a)') repeat('=',67)
......@@ -310,8 +308,8 @@
! with hbar = 1 and M already contained in the eigenmodes
! g2 is Ry^2, wkf must already account for the spin factor
!
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp )
......
......@@ -34,7 +34,7 @@
USE elph2, ONLY : etf, ibndmin, ibndmax, nkqf, xqf, &
epf17, wkf, nkf, wf, wqf, xkf, nkqtotf,&
esigmar_all, esigmai_all, a_all
USE constants_epw, ONLY : ryd2mev, one, ryd2ev, two, zero, pi, ci
USE constants_epw, ONLY : ryd2mev, one, ryd2ev, two, zero, pi, ci, eps8
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : me_pool, inter_pool_comm
USE division, ONLY : fkbounds
......@@ -103,8 +103,6 @@
!! Current frequency
REAL(kind=DP) :: dw
!! Frequency intervals
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
real(kind=DP) :: specfun_sum, esigmar0
real(kind=DP) :: fermi(nw_specfun)
real(kind=DP), external :: efermig, dos_ef, wgauss
......@@ -223,8 +221,8 @@
! with hbar = 1 and M already contained in the eigenmodes
! g2 is Ry^2, wkf must already account for the spin factor
!
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp )
......
......@@ -36,7 +36,7 @@
wkf, xqf, nkqf, &
nkf, wf, a_all, &
gammai_all,gammar_all, efnew
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero, pi, cone, ci
USE constants_epw, ONLY : ryd2mev, ryd2ev, two, zero, pi, cone, ci, eps8
USE mp_world, ONLY : mpime
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : inter_pool_comm, ionode_id
......@@ -114,8 +114,6 @@
REAL(kind=DP), external :: w0gauss
!! This function computes the derivative of the Fermi-Dirac function
!! It is therefore an approximation for a delta function
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
!
dw = ( wmax_specfun - wmin_specfun ) / dble (nw_specfun-1)
!
......@@ -228,8 +226,8 @@
! with hbar = 1 and M already contained in the eigenmodes
! g2 is Ry^2, wkf must already account for the spin factor
!
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp ) !* epsTF
......
......@@ -201,7 +201,7 @@
xqf, zi_allvb, zi_allcb
USE transportcom, ONLY : transp_temp, lower_bnd
USE constants_epw, ONLY : zero, one, two, pi, ryd2mev, kelvin2eV, ryd2ev, &
eps6
eps6, eps8
USE mp, ONLY : mp_barrier, mp_sum
USE mp_global, ONLY : world_comm
USE io_scattering, ONLY : scattering_write, tau_write, merge_read
......@@ -303,8 +303,6 @@
!! The derivative of wgauss: an approximation to the delta function
REAL(kind=DP), PARAMETER :: eps = 1.d-4
!! Tolerence parameter for the velocity
REAL(kind=DP), PARAMETER :: eps2 = 0.01/ryd2mev
!! Tolerence
!
CALL start_clock ( 'SCAT' )
!
......@@ -429,8 +427,8 @@
!
! In case of q=\Gamma, then the short-range = the normal g. We therefore
! need to treat it like the normal g with abs(g).
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp, KIND=DP )
......@@ -510,8 +508,8 @@
!
! In case of q=\Gamma, then the short-range = the normal g. We therefore
! need to treat it like the normal g with abs(g).
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps2 .OR. abs(xqf (2, iq))> eps2 &
.OR. abs(xqf (3, iq))> eps2 )) THEN
IF ( shortrange .AND. ( abs(xqf (1, iq))> eps8 .OR. abs(xqf (2, iq))> eps8 &
.OR. abs(xqf (3, iq))> eps8 )) THEN
! SP: The abs has to be removed. Indeed the epf17 can be a pure imaginary
! number, in which case its square will be a negative number.
g2 = REAL( (epf17 (jbnd, ibnd, imode, ik)**two)*inv_wq*g2_tmp, KIND=DP)
......
......@@ -150,7 +150,9 @@ clean:
epw_*/*_elcond_* epw_*/fort.* epw_*/decay.* epw_*/scattering_rate* \
epw_*/*.fc.* epw_*/*_band.dat epw_*/*F_restart* epw_*/*F_restart_CB* \
epw_*/*.mmn epw_*/*.epmatkq* epw_*/*.epmatkqcb* epw_*/sparse* \
epw_*/*.Fin_restart* epw_*/*.Fin_restartcb*
epw_*/*.Fin_restart* epw_*/*.Fin_restartcb* epw_*/*.acon_iso_* \
epw_*/*.imag_aniso* epw_*/*.imag_iso* epw_*/*.pade_* epw_*/*qdos_*
# Special cases for EPW
@rm -rf epw_base/save epw_super/save
......
......@@ -33,7 +33,7 @@
eps_acustic = 5.0 ! Lowest boundary for the
ephwrite = .true. ! Writes .ephmat files used when wliasberg = .true.
fsthick = 0.4 ! eV
fsthick = 20.0 ! eV
eptemp = 300 ! K
degaussw = 0.10 ! eV
degaussq = 0.5 ! meV
......@@ -50,8 +50,8 @@
wscut = 0.5 ! eV Upper limit over frequency integration/summation in the Elisashberg eq
nstemp = 1
tempsmin = 25.00
tempsmax = 30.00
tempsmin = 15.00
tempsmax = 20.00
nsiter = 500
......
......@@ -33,7 +33,7 @@
eps_acustic = 5.0 ! Lowest boundary for the
ephwrite = .true. ! Writes .ephmat files used when wliasberg = .true.
fsthick = 0.4 ! eV
fsthick = 20.0 ! eV
eptemp = 300 ! K
degaussw = 0.10 ! eV
degaussq = 0.5 ! meV
......@@ -52,8 +52,8 @@
wscut = 0.5 ! eV Upper limit over frequency integration/summation in the Elisashberg eq
nstemp = 1
tempsmin = 25.00
tempsmax = 30.00
tempsmin = 15.00
tempsmax = 20.00
nsiter = 500
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment