mirror of
https://github.com/cp2k/cp2k.git
synced 2026-07-28 14:15:19 -04:00
Use Newton's minimization method for MP and MV smearing (#4985)
This commit is contained in:
SY Wang
GitHub
parent
82afa516d5
commit
079df36cf9
2 changed files with 285 additions and 57 deletions
|
|
@ -91,7 +91,8 @@ MODULE bibliography
|
|||
Sander2015, Schreiber2008, vanSetten2015, Setyawan2010, Ahart2024, Knysh2024, &
|
||||
Schambeck2024, Mewes2018, Sertcan2024, Drautz2019, Lysogorskiy2021, Bochkarev2024, &
|
||||
VazdaCruz2021, Chen2025, Hernandez2025, Marek2025, Hehn2022, Hehn2024, Pasquier2025, &
|
||||
Hanasaki2025, Tan2025
|
||||
Hanasaki2025, Tan2025, FuHo1983, MethfesselPaxton1989, Marzari1999, dosSantos2023, &
|
||||
Mermin1965
|
||||
|
||||
CONTAINS
|
||||
|
||||
|
|
@ -2012,6 +2013,37 @@ CONTAINS
|
|||
source="J. Phys. Chem. A", volume="129", pages="7313-7344", &
|
||||
year=2025, doi="10.1021/acs.jpca.5c02969")
|
||||
|
||||
CALL add_reference(key=Mermin1965, &
|
||||
authors=s2a("N. D. Mermin"), &
|
||||
title="Thermal Properties of the Inhomogeneous Electron Gas", &
|
||||
source="Phys. Rev.", volume="137", pages="A1441-A1443", &
|
||||
year=1965, doi="10.1103/PhysRev.137.A1441")
|
||||
|
||||
CALL add_reference(key=FuHo1983, &
|
||||
authors=s2a("C.-L. Fu", "K.-M. Ho"), &
|
||||
title="First-principles calculation of the equilibrium ground-state "// &
|
||||
"properties of transition metals: Applications to Nb and Mo", &
|
||||
source="Phys. Rev. B", volume="28", pages="5480-5486", &
|
||||
year=1983, doi="10.1103/PhysRevB.28.5480")
|
||||
|
||||
CALL add_reference(key=MethfesselPaxton1989, &
|
||||
authors=s2a("M. Methfessel", "A. T. Paxton"), &
|
||||
title="High-precision sampling for Brillouin-zone integration in metals", &
|
||||
source="Phys. Rev. B", volume="40", pages="3616-3621", &
|
||||
year=1989, doi="10.1103/PhysRevB.40.3616")
|
||||
|
||||
CALL add_reference(key=Marzari1999, &
|
||||
authors=s2a("N. Marzari", "D. Vanderbilt", "A. De Vita", "M. C. Payne"), &
|
||||
title="Thermal Contraction and Disordering of the Al(110) Surface", &
|
||||
source="Phys. Rev. Lett.", volume="82", pages="3296-3299", &
|
||||
year=1999, doi="10.1103/PhysRevLett.82.3296")
|
||||
|
||||
CALL add_reference(key=dosSantos2023, &
|
||||
authors=s2a("F. J. dos Santos", "N. Marzari"), &
|
||||
title="Fermi energy determination for advanced smearing techniques", &
|
||||
source="Phys. Rev. B", volume="107", pages="195122", &
|
||||
year=2023, doi="10.1103/PhysRevB.107.195122")
|
||||
|
||||
END SUBROUTINE add_all_references
|
||||
|
||||
END MODULE bibliography
|
||||
|
|
|
|||
|
|
@ -23,6 +23,12 @@
|
|||
! **************************************************************************************************
|
||||
MODULE smearing_utils
|
||||
|
||||
USE bibliography, ONLY: FuHo1983,&
|
||||
Marzari1999,&
|
||||
Mermin1965,&
|
||||
MethfesselPaxton1989,&
|
||||
cite_reference,&
|
||||
dosSantos2023
|
||||
USE input_constants, ONLY: smear_fermi_dirac,&
|
||||
smear_gaussian,&
|
||||
smear_mp,&
|
||||
|
|
@ -41,11 +47,35 @@ MODULE smearing_utils
|
|||
! Unified interface (method as parameter)
|
||||
PUBLIC :: SmearOcc, SmearFixed, SmearFixedDeriv, SmearFixedDerivMV
|
||||
PUBLIC :: Smearkp, Smearkp2
|
||||
PRIVATE :: cite_smearing
|
||||
|
||||
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'smearing_utils'
|
||||
INTEGER, PARAMETER, PRIVATE :: BISECT_MAX_ITER = 400
|
||||
INTEGER, PARAMETER, PRIVATE :: NEWTON_MAX_ITER = 50
|
||||
|
||||
CONTAINS
|
||||
! **************************************************************************************************
|
||||
!> \brief Citation of Smearing methods
|
||||
!> \param method ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE cite_smearing(method)
|
||||
INTEGER, INTENT(IN) :: method
|
||||
|
||||
SELECT CASE (method)
|
||||
CASE (smear_fermi_dirac)
|
||||
CALL cite_reference(Mermin1965)
|
||||
CASE (smear_gaussian)
|
||||
CALL cite_reference(FuHo1983)
|
||||
CASE (smear_mp)
|
||||
CALL cite_reference(FuHo1983)
|
||||
CALL cite_reference(MethfesselPaxton1989)
|
||||
CALL cite_reference(dosSantos2023)
|
||||
CASE (smear_mv)
|
||||
CALL cite_reference(FuHo1983)
|
||||
CALL cite_reference(Marzari1999)
|
||||
CALL cite_reference(dosSantos2023)
|
||||
END SELECT
|
||||
END SUBROUTINE cite_smearing
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Returns occupations and smearing correction for a given set of
|
||||
|
|
@ -260,7 +290,12 @@ CONTAINS
|
|||
!> electron count equals N, for a given smearing method (Gamma point).
|
||||
!> Brackets mu by expanding outward from [min(e), max(e)] in steps
|
||||
!> of sigma, then bisects to machine precision.
|
||||
!> Could fail if mu lies far outside the eigenvalue range (to be fixed).
|
||||
!>
|
||||
!> For MP-1 and MV: the occupation function is non-monotonic, so it's
|
||||
!> possible that pure bisection find a spurious root.
|
||||
!> We first bisect with Gaussian smearing to get a reliable initial mu,
|
||||
!> then refine with Newton's method using the actual method's dN/dmu.
|
||||
!> (dos Santos & Marzari, PRB 2023)
|
||||
!>
|
||||
!> \param f occupations (output)
|
||||
!> \param mu chemical potential found by bisection (output)
|
||||
|
|
@ -281,9 +316,10 @@ CONTAINS
|
|||
INTEGER, INTENT(IN), OPTIONAL :: estate
|
||||
REAL(KIND=dp), INTENT(IN), OPTIONAL :: festate
|
||||
|
||||
INTEGER :: iter, my_estate
|
||||
REAL(KIND=dp) :: mu_max, mu_min, mu_now, my_festate, &
|
||||
N_now, N_tmp
|
||||
INTEGER :: iter, my_estate, Nstate
|
||||
REAL(KIND=dp) :: Gsum, mu_max, mu_min, mu_now, &
|
||||
my_festate, N_now, N_tmp
|
||||
REAL(KIND=dp), ALLOCATABLE :: gvec(:)
|
||||
|
||||
IF (PRESENT(estate) .AND. PRESENT(festate)) THEN
|
||||
my_estate = estate
|
||||
|
|
@ -293,53 +329,114 @@ CONTAINS
|
|||
my_festate = my_estate
|
||||
END IF
|
||||
|
||||
! Bracket mu from below
|
||||
mu_min = MINVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_min, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (N_tmp > N .OR. iter > 20) THEN
|
||||
mu_min = mu_min - sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
Nstate = SIZE(e)
|
||||
|
||||
! Bracket mu from above
|
||||
mu_max = MAXVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_max, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (N_tmp < N .OR. iter > 20) THEN
|
||||
mu_max = mu_max + sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
CALL cite_smearing(method)
|
||||
|
||||
! Bisection
|
||||
iter = 0
|
||||
DO WHILE (mu_max - mu_min > EPSILON(mu)*MAX(1.0_dp, ABS(mu_max), ABS(mu_min)))
|
||||
iter = iter + 1
|
||||
mu_now = (mu_max + mu_min)/2.0_dp
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu_now, sigma, maxocc, method, my_estate, my_festate)
|
||||
SELECT CASE (method)
|
||||
|
||||
IF (N_now <= N) THEN
|
||||
mu_min = mu_now
|
||||
ELSE
|
||||
mu_max = mu_now
|
||||
END IF
|
||||
! Non-monotonic methods: Gaussian bisection + Newton refinement
|
||||
CASE (smear_mp, smear_mv)
|
||||
! Step 1: Gaussian bisection for a reliable initial mu
|
||||
mu_min = MINVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_min, sigma, maxocc, smear_gaussian, my_estate, my_festate)
|
||||
IF (N_tmp > N .OR. iter > 20) THEN
|
||||
mu_min = mu_min - sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
IF (iter > BISECT_MAX_ITER) THEN
|
||||
CPWARN("SmearFixed: maximum bisection iterations reached")
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
mu_max = MAXVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_max, sigma, maxocc, smear_gaussian, my_estate, my_festate)
|
||||
IF (N_tmp < N .OR. iter > 20) THEN
|
||||
mu_max = mu_max + sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu, sigma, maxocc, method, my_estate, my_festate)
|
||||
iter = 0
|
||||
DO WHILE (mu_max - mu_min > EPSILON(mu)*MAX(1.0_dp, ABS(mu_max), ABS(mu_min)))
|
||||
iter = iter + 1
|
||||
mu_now = (mu_max + mu_min)/2.0_dp
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu_now, sigma, maxocc, smear_gaussian, my_estate, my_festate)
|
||||
IF (N_now <= N) THEN
|
||||
mu_min = mu_now
|
||||
ELSE
|
||||
mu_max = mu_now
|
||||
END IF
|
||||
IF (iter > BISECT_MAX_ITER) EXIT
|
||||
END DO
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
|
||||
! Step 2: Newton refinement with the actual method
|
||||
ALLOCATE (gvec(Nstate))
|
||||
DO iter = 1, NEWTON_MAX_ITER
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (ABS(N_now - N) < N*1.0e-12_dp) EXIT
|
||||
CALL compute_gvec(gvec, f, e, mu, sigma, maxocc, Nstate, method, my_estate, my_festate)
|
||||
Gsum = accurate_sum(gvec)
|
||||
IF (ABS(Gsum) < EPSILON(Gsum)) EXIT
|
||||
mu = mu + (N - N_now)/Gsum
|
||||
END DO
|
||||
DEALLOCATE (gvec)
|
||||
|
||||
! Final evaluation
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu, sigma, maxocc, method, my_estate, my_festate)
|
||||
|
||||
! Monotonic methods (FD, Gaussian): pure bisection
|
||||
CASE DEFAULT
|
||||
mu_min = MINVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_min, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (N_tmp > N .OR. iter > 20) THEN
|
||||
mu_min = mu_min - sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
mu_max = MAXVAL(e)
|
||||
iter = 0
|
||||
DO
|
||||
iter = iter + 1
|
||||
CALL SmearOcc(f, N_tmp, kTS, e, mu_max, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (N_tmp < N .OR. iter > 20) THEN
|
||||
mu_max = mu_max + sigma
|
||||
ELSE
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
iter = 0
|
||||
DO WHILE (mu_max - mu_min > EPSILON(mu)*MAX(1.0_dp, ABS(mu_max), ABS(mu_min)))
|
||||
iter = iter + 1
|
||||
mu_now = (mu_max + mu_min)/2.0_dp
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu_now, sigma, maxocc, method, my_estate, my_festate)
|
||||
IF (N_now <= N) THEN
|
||||
mu_min = mu_now
|
||||
ELSE
|
||||
mu_max = mu_now
|
||||
END IF
|
||||
IF (iter > BISECT_MAX_ITER) THEN
|
||||
CPWARN("SmearFixed: maximum bisection iterations reached")
|
||||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
CALL SmearOcc(f, N_now, kTS, e, mu, sigma, maxocc, method, my_estate, my_festate)
|
||||
|
||||
END SELECT
|
||||
|
||||
END SUBROUTINE SmearFixed
|
||||
|
||||
|
|
@ -350,6 +447,9 @@ CONTAINS
|
|||
!> or sigma*ln[(1-eps)/eps] for Fermi-Dirac, reflecting the different
|
||||
!> tail decay rates.
|
||||
!>
|
||||
!> For MP-1 and MV: Gaussian bisection + Newton refinement
|
||||
!> (dos Santos & Marzari, PRB 2023).
|
||||
!>
|
||||
!> \param f occupations (nmo x nkp, output)
|
||||
!> \param mu chemical potential (output)
|
||||
!> \param kTS smearing correction (output)
|
||||
|
|
@ -372,10 +472,25 @@ CONTAINS
|
|||
|
||||
REAL(KIND=dp), PARAMETER :: epsocc = 1.0e-12_dp
|
||||
|
||||
INTEGER :: iter
|
||||
REAL(KIND=dp) :: de, mu_max, mu_min, N_now
|
||||
INTEGER :: bisect_method, ik, is, iter, nkp, nmo
|
||||
REAL(KIND=dp) :: de, dNdmu, expu2, expx2, mu_max, mu_min, &
|
||||
N_now, u, x
|
||||
|
||||
nmo = SIZE(e, 1)
|
||||
nkp = SIZE(e, 2)
|
||||
|
||||
CALL cite_smearing(method)
|
||||
|
||||
! Choose bisection method: Gaussian for MP/MV, actual method for FD/Gaussian
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp, smear_mv)
|
||||
bisect_method = smear_gaussian
|
||||
CASE DEFAULT
|
||||
bisect_method = method
|
||||
END SELECT
|
||||
|
||||
! Initial bracket
|
||||
SELECT CASE (bisect_method)
|
||||
CASE (smear_fermi_dirac)
|
||||
de = sigma*LOG((1.0_dp - epsocc)/epsocc)
|
||||
CASE DEFAULT
|
||||
|
|
@ -383,13 +498,14 @@ CONTAINS
|
|||
END SELECT
|
||||
de = MAX(de, 0.5_dp)
|
||||
|
||||
! Bisection with bisect_method
|
||||
mu_min = MINVAL(e) - de
|
||||
mu_max = MAXVAL(e) + de
|
||||
iter = 0
|
||||
DO WHILE (mu_max - mu_min > EPSILON(mu)*MAX(1.0_dp, ABS(mu_max), ABS(mu_min)))
|
||||
iter = iter + 1
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
CALL Smear2(f, N_now, kTS, e, mu, wk, sigma, maxocc, method)
|
||||
CALL Smear2(f, N_now, kTS, e, mu, wk, sigma, maxocc, bisect_method)
|
||||
|
||||
IF (ABS(N_now - nel) < nel*epsocc) EXIT
|
||||
|
||||
|
|
@ -404,8 +520,38 @@ CONTAINS
|
|||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
|
||||
! Newton refinement for non-monotonic methods
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp, smear_mv)
|
||||
DO iter = 1, NEWTON_MAX_ITER
|
||||
CALL Smear2(f, N_now, kTS, e, mu, wk, sigma, maxocc, method)
|
||||
IF (ABS(N_now - nel) < nel*epsocc) EXIT
|
||||
|
||||
! Compute dN/dmu = sum_{ik} wk * g_i(k) inline
|
||||
dNdmu = 0.0_dp
|
||||
DO ik = 1, nkp
|
||||
DO is = 1, nmo
|
||||
x = (e(is, ik) - mu)/sigma
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp)
|
||||
expx2 = EXP(-x*x)
|
||||
dNdmu = dNdmu + maxocc*(3.0_dp - 2.0_dp*x*x)/(2.0_dp*sigma*rootpi)*expx2*wk(ik)
|
||||
CASE (smear_mv)
|
||||
u = x + sqrthalf
|
||||
expu2 = EXP(-u*u)
|
||||
dNdmu = dNdmu + maxocc*(2.0_dp + sqrt2*x)/(sigma*rootpi)*expu2*wk(ik)
|
||||
END SELECT
|
||||
END DO
|
||||
END DO
|
||||
|
||||
IF (ABS(dNdmu) < EPSILON(dNdmu)) EXIT
|
||||
mu = mu + (nel - N_now)/dNdmu
|
||||
END DO
|
||||
END SELECT
|
||||
|
||||
! Final evaluation with the actual method
|
||||
CALL Smear2(f, N_now, kTS, e, mu, wk, sigma, maxocc, method)
|
||||
|
||||
END SUBROUTINE Smearkp
|
||||
|
|
@ -415,6 +561,8 @@ CONTAINS
|
|||
!> chemical potential across both spin channels).
|
||||
!> Asserts that the third dimension of f and e is exactly 2.
|
||||
!>
|
||||
!> For MP-1 and MV: Gaussian bisection + Newton refinement.
|
||||
!>
|
||||
!> \param f occupations (nmo x nkp x 2, output)
|
||||
!> \param mu chemical potential (output)
|
||||
!> \param kTS smearing correction (output)
|
||||
|
|
@ -436,13 +584,26 @@ CONTAINS
|
|||
|
||||
REAL(KIND=dp), PARAMETER :: epsocc = 1.0e-12_dp
|
||||
|
||||
INTEGER :: iter
|
||||
REAL(KIND=dp) :: de, kTSa, kTSb, mu_max, mu_min, N_now, &
|
||||
na, nb
|
||||
INTEGER :: bisect_method, ik, is, ispin, iter, nkp, &
|
||||
nmo
|
||||
REAL(KIND=dp) :: de, dNdmu, expu2, expx2, kTSa, kTSb, &
|
||||
mu_max, mu_min, N_now, na, nb, u, x
|
||||
|
||||
CPASSERT(SIZE(f, 3) == 2 .AND. SIZE(e, 3) == 2)
|
||||
|
||||
nmo = SIZE(e, 1)
|
||||
nkp = SIZE(e, 2)
|
||||
|
||||
CALL cite_smearing(method)
|
||||
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp, smear_mv)
|
||||
bisect_method = smear_gaussian
|
||||
CASE DEFAULT
|
||||
bisect_method = method
|
||||
END SELECT
|
||||
|
||||
SELECT CASE (bisect_method)
|
||||
CASE (smear_fermi_dirac)
|
||||
de = sigma*LOG((1.0_dp - epsocc)/epsocc)
|
||||
CASE DEFAULT
|
||||
|
|
@ -450,14 +611,15 @@ CONTAINS
|
|||
END SELECT
|
||||
de = MAX(de, 0.5_dp)
|
||||
|
||||
! Bisection with bisect_method
|
||||
mu_min = MINVAL(e) - de
|
||||
mu_max = MAXVAL(e) + de
|
||||
iter = 0
|
||||
DO WHILE (mu_max - mu_min > EPSILON(mu)*MAX(1.0_dp, ABS(mu_max), ABS(mu_min)))
|
||||
iter = iter + 1
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
CALL Smear2(f(:, :, 1), na, kTSa, e(:, :, 1), mu, wk, sigma, 1.0_dp, method)
|
||||
CALL Smear2(f(:, :, 2), nb, kTSb, e(:, :, 2), mu, wk, sigma, 1.0_dp, method)
|
||||
CALL Smear2(f(:, :, 1), na, kTSa, e(:, :, 1), mu, wk, sigma, 1.0_dp, bisect_method)
|
||||
CALL Smear2(f(:, :, 2), nb, kTSb, e(:, :, 2), mu, wk, sigma, 1.0_dp, bisect_method)
|
||||
N_now = na + nb
|
||||
|
||||
IF (ABS(N_now - nel) < nel*epsocc) EXIT
|
||||
|
|
@ -473,8 +635,42 @@ CONTAINS
|
|||
EXIT
|
||||
END IF
|
||||
END DO
|
||||
|
||||
mu = (mu_max + mu_min)/2.0_dp
|
||||
|
||||
! Newton refinement for non-monotonic methods
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp, smear_mv)
|
||||
DO iter = 1, NEWTON_MAX_ITER
|
||||
CALL Smear2(f(:, :, 1), na, kTSa, e(:, :, 1), mu, wk, sigma, 1.0_dp, method)
|
||||
CALL Smear2(f(:, :, 2), nb, kTSb, e(:, :, 2), mu, wk, sigma, 1.0_dp, method)
|
||||
N_now = na + nb
|
||||
IF (ABS(N_now - nel) < nel*epsocc) EXIT
|
||||
|
||||
! dN/dmu across both spin channels (maxocc=1 per spin)
|
||||
dNdmu = 0.0_dp
|
||||
DO ispin = 1, 2
|
||||
DO ik = 1, nkp
|
||||
DO is = 1, nmo
|
||||
x = (e(is, ik, ispin) - mu)/sigma
|
||||
SELECT CASE (method)
|
||||
CASE (smear_mp)
|
||||
expx2 = EXP(-x*x)
|
||||
dNdmu = dNdmu + (3.0_dp - 2.0_dp*x*x)/(2.0_dp*sigma*rootpi)*expx2*wk(ik)
|
||||
CASE (smear_mv)
|
||||
u = x + sqrthalf
|
||||
expu2 = EXP(-u*u)
|
||||
dNdmu = dNdmu + (2.0_dp + sqrt2*x)/(sigma*rootpi)*expu2*wk(ik)
|
||||
END SELECT
|
||||
END DO
|
||||
END DO
|
||||
END DO
|
||||
|
||||
IF (ABS(dNdmu) < EPSILON(dNdmu)) EXIT
|
||||
mu = mu + (nel - N_now)/dNdmu
|
||||
END DO
|
||||
END SELECT
|
||||
|
||||
! Final evaluation with the actual method
|
||||
CALL Smear2(f(:, :, 1), na, kTSa, e(:, :, 1), mu, wk, sigma, 1.0_dp, method)
|
||||
CALL Smear2(f(:, :, 2), nb, kTSb, e(:, :, 2), mu, wk, sigma, 1.0_dp, method)
|
||||
kTS = kTSa + kTSb
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue