Skip to content

Commit ab49ad9

Browse files
1041176461jiyuangQianruipkuESROAMERmaki49
authored
Feature: long-range corrected hybrid functionals implementation with Spencer-Alavi truncation Coulomb potential (#6351)
* Use template to reconstruct parse_expression * Feature: output R matrix at each MD step * Modify'matrix_HS' to 'matrix' for R matrix output * Merge branches 'develop' and 'develop' of https://github.com/1041176461/abacus-develop into develop * Fix: modify index in parse_expression * Fix: add regfree for parse_expression * Doc: update phonopy doc * Doc: update phonopy doc * fix tdos plot for nspin=2 * optimize dosplot for nspin=2 * fix legend for dosplot * Add files via upload * Update cal_edm_tddft.cpp * Refactor: modify exx loop for relax/md * Update result.ref * Fix wrong parameters in integrate test * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * add exx when istep>0 * Update dos.py * Update esolver_sdft_pw.cpp * Update lcao_before_scf.cpp * Update Exx_LRI_interface.h * Update Exx_LRI_interface.hpp * Fix: compile error * Fix: compile error * Fix: change HSE relax/md result.ref for new framework * Fix: compile error * compatible with exx_iter_finish * Add files via upload * Update esolver_ks_lcao_tddft.h * Update esolver_ks_lcao_tddft.cpp * Fix: support negative value in parse_expression * [pre-commit.ci lite] apply automatic fixes * Fix: parse_expression for scientific notation * Update input_conv.h * add complex erf function * refactor ccp framework * fix compile * fix compile * fix compile * fix compile * fix compile * fix compile * fix compile * fix empty error * add revised spencer * add cal_psi_erfc_spencer * add exx_singularity_correction * add exx_singularity_correction test * add lc-corrected hybrid functionals * fix compile * fix compile * fix compile error * fix libxc error * fix libxc error * fix libxc error * add klist.h in RI_Util.h * fix include * fix ptest for hse * fix lc_pbe name for vdw autoset * fix test error * fix cam_pbeh test * fix fock settings * support rsh functionals * fix libxc param * add input * fix input --------- Co-authored-by: jiyuang <jiyuyang@mail.ustc.com> Co-authored-by: Qianrui <76200646+Qianruipku@users.noreply.github.com> Co-authored-by: HTZhao <104255052+ESROAMER@users.noreply.github.com> Co-authored-by: maki49 <1579492865@qq.com> Co-authored-by: pre-commit-ci-lite[bot] <117423508+pre-commit-ci-lite[bot]@users.noreply.github.com>
1 parent 4b15663 commit ab49ad9

File tree

24 files changed

+487
-77
lines changed

24 files changed

+487
-77
lines changed

docs/advanced/input_files/input-main.md

Lines changed: 31 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -2769,37 +2769,41 @@ Support for hybrid functionals in the *pw [basis_type](#basis_type)* is under ac
27692769

27702770
The following parameters apply to *[basis_type](#basis_type)==lcao/lcao_in_pw/pw*. For basis specific parameters, see the sections *[Exact Exchange (LCAO/LCAO in PW)](#exact-exchange-lcaolcao-in-pw)* and *[Exact Exchange (PW)](#exact-exchange-pw)*.
27712771

2772-
**Availablity**: *[dft_functional](#dft_functional)==hse/hf/pbe0/scan0/opt_orb* or *[rpa](#rpa)==True*.
2772+
#### Hybrid Functional Parameters {#hybrid_func_params}
2773+
|[dft_functional](#dft_functional)|[exx_fock_alpha](#exx_fock_alpha)|[exx_erfc_alpha](#exx_erfc_alpha)|[exx_erfc_omega](#exx_erfc_omega)|
2774+
|:-:|:-:|:-:|:-:|
2775+
|hf|1|0|0|
2776+
|lc_pbe|1|-1|0.33|
2777+
|lc_wpbe|1|-1|0.4|
2778+
|lrc_wpbe|1|-1|0.3|
2779+
|lrc_wpbeh|1|-0.8|0.2|
2780+
|pbe0|0.25|0|0|
2781+
|hse|0.0|0.25|0.11|
2782+
|scan0|0.25|0|0|
2783+
|cam_pbeh|0.2|0.8|0.7|
2784+
|b3lyp|0.2|0|0|
2785+
|muller|1|0|0|
2786+
|power|1|0|0|
2787+
|wp22|1|-1|0.11|
2788+
|cwp22|0|1|0.11|
27732789

27742790
### exx_fock_alpha
27752791

27762792
- **Type**: Real \[Real...\](optional)
2777-
- **Description**: Fraction of Fock exchange 1/r in hybrid functionals, so that $E_{X} = \alpha E_{X} + (1-\alpha)E_{X,\text{LDA/GGA}}$
2778-
- **Default**:
2779-
- 1: if *[dft_functional](#dft_functional)==hf*
2780-
- 0.25: if *[dft_functional](#dft_functional)==pbe0*
2781-
- 0.2: if *[dft_functional](#dft_functional)==b3lyp*
2782-
- 0.25: if *[dft_functional](#dft_functional)==scan0*
2783-
- 1: if *[dft_functional](#dft_functional)==muller*
2784-
- 1: if *[dft_functional](#dft_functional)==power*
2785-
- 1: if *[dft_functional](#dft_functional)==wp22*
2786-
- 0: else
2793+
- **Description**: Fraction of full-ranged Fock exchange 1/r ($\alpha$) in range-separated hybrid funtionals, so that $E_{X} = \alpha E_{X}^\text{HF-LR}+(\alpha+\beta) E_{X}^\text{HF-SR}+(1-\alpha)E_{X}^\text{KS-LR}+[1-(\alpha+\beta)]E_{X}^\text{KS-SR}$.
2794+
- **Default**: see [hybrid_func_params](#hybrid_func_params)
27872795

27882796
### exx_erfc_alpha
27892797

27902798
- **Type**: Real \[Real...\](optional)
2791-
- **Description**: Fraction of exchange erfc(wr)/r in hybrid functionals, so that $E_{X} = \alpha E_{X}^{\text{SR}} + (1-\alpha)E_{X,\text{LDA/GGA}}^{\text{SR}} + E_{X,\text{LDA/GGA}}^{\text{LR}}$
2792-
- **Default**:
2793-
- 0.25: if *[dft_functional](#dft_functional)==hse*
2794-
- 1: if *[dft_functional](#dft_functional)==cwp22*
2795-
- -1: if *[dft_functional](#dft_functional)==wp22*
2796-
- 0: else
2799+
- **Description**: Fraction of short-ranged Fock exchange erfc(wr)/r ($\beta$) in range-separated hybrid funtionals, so that $E_{X} = \alpha E_{X}^\text{HF-LR}+(\alpha+\beta) E_{X}^\text{HF-SR}+(1-\alpha)E_{X}^\text{KS-LR}+[1-(\alpha+\beta)]E_{X}^\text{KS-SR}$.
2800+
- **Default**: see [hybrid_func_params](#hybrid_func_params)
27972801

27982802
### exx_erfc_omega
27992803

28002804
- **Type**: Real \[Real...\](optional)
28012805
- **Description**: Range-separation parameter in exchange, such that $1/r=\text{erfc}(\omega r)/r+\text{erf}(\omega r)/r$
2802-
- **Default**: 0.11
2806+
- **Default**: see [hybrid_func_params](#hybrid_func_params)
28032807

28042808
### exx_separate_loop
28052809

@@ -2934,6 +2938,15 @@ These variables are relevant when using hybrid functionals with *[basis_type](#b
29342938
- True: if gamma_only
29352939
- False: else
29362940

2941+
### exx_singularity_correction
2942+
2943+
- **Type**: String
2944+
- **Description**:
2945+
- spencer: see Phys. Rev. B 77, 193110 (2008).
2946+
- revised_spencer: see Phys. Rev. Mater. 5, 013807 (2021).
2947+
Set the scheme of Coulomb singularity correction.
2948+
- **Default**: default
2949+
29372950
### rpa_ccp_rmesh_times
29382951

29392952
- **Type**: Real

source/module_parameter/input_parameter.h

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -523,6 +523,7 @@ struct Input_para
523523
std::vector<std::string> exx_erfc_omega = {"default"}; ///< range-separation parameter in HSE functional
524524
bool exx_separate_loop = true; ///< if 1, a two-step method is employed, else it will start
525525
///< with a GGA-Loop, and then Hybrid-Loop
526+
std::string exx_singularity_correction = "default"; ///< set the scheme of Coulomb singularity correction
526527
int exx_hybrid_step = 100; ///< the maximal electronic iteration number in
527528
///< the evaluation of Fock exchange
528529
double exx_mixing_beta = 1.0; ///< mixing_beta for outer-loop when exx_separate_loop=1

source/source_hamilt/module_vdw/vdwd3_autoset_xcname.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -229,7 +229,7 @@ const std::map<std::string, std::string> xcname_libxc_xc_ = {
229229
{"XC_HYB_GGA_XC_WHPBE0", "whpbe0"},
230230
{"XC_HYB_GGA_XC_LC_BLYP_EA", "lc_blyp_ea"},
231231
{"XC_HYB_GGA_XC_LC_BOP", "lc_bop"},
232-
{"XC_HYB_GGA_XC_LC_PBEOP", "lc_pbeop"},
232+
{"XC_HYB_GGA_XC_LC_PBEOP", "lc_pbe"},
233233
{"XC_HYB_GGA_XC_LC_BLYPR", "lc_blypr"},
234234
{"XC_HYB_GGA_XC_MCAM_B3LYP", "mcam_b3lyp"},
235235
{"XC_MGGA_XC_VCML_RVV10", "vcml_rvv10"},

source/source_hamilt/module_xc/exx_info.h

Lines changed: 11 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -18,14 +18,14 @@ struct Exx_Info
1818

1919
// Fock:
2020
// "alpha": "0"
21-
// "Rcut_type": "limits" / "spencer"
21+
// "singularity_correction": "limits" / "spencer" / "revised_spencer"
2222
// "lambda": "0.3"
23-
// //"Rcut"
23+
// "Rcut"
2424
// Erfc:
2525
// "alpha": "0"
2626
// "omega": "0.11"
27-
// "Rcut_type": "limits"
28-
// //"Rcut"
27+
// "singularity_correction": "limits" / "spencer" / "revised_spencer"
28+
// "Rcut"
2929

3030
Conv_Coulomb_Pot_K::Ccp_Type ccp_type;
3131
double hybrid_alpha = 0.25;
@@ -51,10 +51,7 @@ struct Exx_Info
5151

5252
struct Exx_Info_RI
5353
{
54-
std::map<Conv_Coulomb_Pot_K::Coulomb_Method,
55-
std::pair<bool,
56-
std::map<Conv_Coulomb_Pot_K::Coulomb_Type,
57-
std::vector<std::map<std::string,std::string>>>>> coulomb_settings;
54+
const std::map<Conv_Coulomb_Pot_K::Coulomb_Type, std::vector<std::map<std::string,std::string>>> &coulomb_param;
5855

5956
bool real_number = false;
6057

@@ -72,10 +69,15 @@ struct Exx_Info
7269
double kmesh_times = 4;
7370

7471
int abfs_Lmax = 0; // tmp
72+
73+
Exx_Info_RI(const Exx_Info::Exx_Info_Global& info_global)
74+
: coulomb_param(info_global.coulomb_param)
75+
{
76+
}
7577
};
7678
Exx_Info_RI info_ri;
7779

78-
Exx_Info() : info_lip(this->info_global)
80+
Exx_Info() : info_lip(this->info_global), info_ri(this->info_global)
7981
{
8082
}
8183
};

source/source_hamilt/module_xc/xc_functional.cpp

Lines changed: 38 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -29,9 +29,10 @@ void XC_Functional::set_xc_first_loop(const UnitCell& ucell)
2929
the first scf iteration only calculate the functional without exact
3030
exchange. but in "nscf" calculation, there is no need of "two-level"
3131
method. */
32-
if (ucell.atoms[0].ncpp.xc_func == "HF"
33-
|| ucell.atoms[0].ncpp.xc_func == "PBE0"
34-
|| ucell.atoms[0].ncpp.xc_func == "HSE") {
32+
if (ucell.atoms[0].ncpp.xc_func == "HF" || ucell.atoms[0].ncpp.xc_func == "HSE"
33+
|| ucell.atoms[0].ncpp.xc_func == "PBE0"|| ucell.atoms[0].ncpp.xc_func == "LC_PBE"
34+
|| ucell.atoms[0].ncpp.xc_func == "LC_WPBE" || ucell.atoms[0].ncpp.xc_func == "LRC_WPBEH"
35+
|| ucell.atoms[0].ncpp.xc_func == "CAM_PBEH") {
3536
XC_Functional::set_xc_type("pbe");
3637
}
3738
else if (ucell.atoms[0].ncpp.xc_func == "SCAN0") {
@@ -149,6 +150,36 @@ void XC_Functional::set_xc_type(const std::string xc_func_in)
149150
func_type = 5;
150151
use_libxc = true;
151152
}
153+
else if( xc_func == "LC_PBE")
154+
{
155+
func_id.push_back(XC_HYB_GGA_XC_LC_PBEOP);
156+
func_type = 4;
157+
use_libxc = true;
158+
}
159+
else if( xc_func == "LC_WPBE")
160+
{
161+
func_id.push_back(XC_HYB_GGA_XC_LC_WPBE);
162+
func_type = 4;
163+
use_libxc = true;
164+
}
165+
else if( xc_func == "LRC_WPBE")
166+
{
167+
func_id.push_back(XC_HYB_GGA_XC_LRC_WPBE);
168+
func_type = 4;
169+
use_libxc = true;
170+
}
171+
else if( xc_func == "LRC_WPBEH")
172+
{
173+
func_id.push_back(XC_HYB_GGA_XC_LRC_WPBEH);
174+
func_type = 4;
175+
use_libxc = true;
176+
}
177+
else if( xc_func == "CAM_PBEH")
178+
{
179+
func_id.push_back(XC_HYB_GGA_XC_CAM_PBEH);
180+
func_type = 4;
181+
use_libxc = true;
182+
}
152183
#endif
153184
else if( xc_func == "HF")
154185
{
@@ -271,9 +302,11 @@ void XC_Functional::set_xc_type(const std::string xc_func_in)
271302

272303
#ifndef USE_LIBXC
273304
if(xc_func == "SCAN" || xc_func == "HSE" || xc_func == "SCAN0"
274-
|| xc_func == "MULLER" || xc_func == "POWER" || xc_func == "WP22" || xc_func == "CWP22")
305+
|| xc_func == "MULLER" || xc_func == "POWER" || xc_func == "WP22" || xc_func == "CWP22" ||
306+
xc_func == "LC_PBE" || xc_func == "LC_WPBE" || xc_func == "LRC_WPBE" ||
307+
xc_func == "LRC_PBEH" || xc_func == "CAM_PBEH")
275308
{
276-
ModuleBase::WARNING_QUIT("set_xc_type","to use SCAN, SCAN0, or HSE, LIBXC is required");
309+
ModuleBase::WARNING_QUIT("set_xc_type","to use SCAN, SCAN0, HSE, long-range corrected (LC_PBE, LC_WPBE...) or CAM_PBEH LIBXC is required");
277310
}
278311
use_libxc = false;
279312
#endif

source/source_hamilt/module_xc/xc_functional_libxc.cpp

Lines changed: 44 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -175,6 +175,50 @@ const std::vector<double> in_built_xc_func_ext_params(const int id)
175175
case XC_GGA_C_LYPR:
176176
return {0.04918, 0.132, 0.2533, 0.349,
177177
0.35/2.29, 2.0/2.29, GlobalC::exx_info.info_global.hse_omega};
178+
// Long-range corrected functionals:
179+
case XC_HYB_GGA_XC_LC_PBEOP: // LC version of PBE
180+
{
181+
// This is a range-separated hybrid functional with range-separation constant 0.330,
182+
// and 0.0% short-range and 100.0% long-range exact exchange,
183+
// using the error function kernel.
184+
return { GlobalC::exx_info.info_global.hse_omega }; //Range separation constant: 0.33
185+
}
186+
case XC_HYB_GGA_XC_LC_WPBE: // Long-range corrected PBE (LC-wPBE) by Vydrov and Scuseria
187+
{
188+
// This is a range-separated hybrid functional with range-separation constant 0.400,
189+
// and 0.0% short-range and 100.0% long-range exact exchange,
190+
// using the error function kernel.
191+
return { std::stod(PARAM.inp.exx_fock_alpha[0]), //Fraction of Hartree-Fock exchange: 1.0
192+
std::stod(PARAM.inp.exx_erfc_alpha[0]), //Fraction of short-range exact exchange: -1.0
193+
GlobalC::exx_info.info_global.hse_omega }; //Range separation constant: 0.4
194+
}
195+
case XC_HYB_GGA_XC_LRC_WPBE: // Long-range corrected PBE (LRC-wPBE) by by Rohrdanz, Martins and Herbert
196+
{
197+
// This is a range-separated hybrid functional with range-separation constant 0.300,
198+
// and 0.0% short-range and 100.0% long-range exact exchange,
199+
// using the error function kernel.
200+
return { std::stod(PARAM.inp.exx_fock_alpha[0]), //Fraction of Hartree-Fock exchange: 1.0
201+
std::stod(PARAM.inp.exx_erfc_alpha[0]), //Fraction of short-range exact exchange: -1.0
202+
GlobalC::exx_info.info_global.hse_omega }; //Range separation constant: 0.3
203+
}
204+
case XC_HYB_GGA_XC_LRC_WPBEH: // Long-range corrected short-range hybrid PBE (LRC-wPBEh) by Rohrdanz, Martins and Herbert
205+
{
206+
// This is a range-separated hybrid functional with range-separation constant 0.200,
207+
// and 20.0% short-range and 100.0% long-range exact exchange,
208+
// using the error function kernel.
209+
return { std::stod(PARAM.inp.exx_fock_alpha[0]), //Fraction of Hartree-Fock exchange: 1.0
210+
std::stod(PARAM.inp.exx_erfc_alpha[0]), //Fraction of short-range exact exchange: -0.8
211+
GlobalC::exx_info.info_global.hse_omega }; //Range separation constant: 0.2
212+
}
213+
case XC_HYB_GGA_XC_CAM_PBEH: // CAM hybrid screened exchange PBE version
214+
{
215+
// This is a range-separated hybrid functional with range-separation constant 0.700,
216+
// and 100.0% short-range and 20.0% long-range exact exchange,
217+
// using the error function kernel.
218+
return { std::stod(PARAM.inp.exx_fock_alpha[0]), //Fraction of Hartree-Fock exchange: 0.2
219+
std::stod(PARAM.inp.exx_erfc_alpha[0]), //Fraction of short-range exact exchange: 0.8
220+
GlobalC::exx_info.info_global.hse_omega }; //Range separation constant: 0.7
221+
}
178222
#endif
179223
default:
180224
return std::vector<double>{};

source/source_io/input_conv.cpp

Lines changed: 21 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -297,7 +297,12 @@ void Input_Conv::Convert()
297297
if (dft_functional_lower == "hf" || dft_functional_lower == "pbe0"
298298
|| dft_functional_lower == "hse"
299299
|| dft_functional_lower == "opt_orb"
300-
|| dft_functional_lower == "scan0") {
300+
|| dft_functional_lower == "scan0"
301+
|| dft_functional_lower == "lc_pbe"
302+
|| dft_functional_lower == "lc_wpbe"
303+
|| dft_functional_lower == "lrc_wpbe"
304+
|| dft_functional_lower == "lrc_wpbeh"
305+
|| dft_functional_lower == "cam_pbeh") {
301306
GlobalC::restart.info_load.load_charge = true;
302307
GlobalC::restart.info_load.load_H = true;
303308
}
@@ -322,10 +327,15 @@ void Input_Conv::Convert()
322327
dft_functional_lower.begin(),
323328
tolower);
324329
if (dft_functional_lower == "hf"
325-
|| dft_functional_lower == "pbe0" || dft_functional_lower == "b3lyp" || dft_functional_lower == "hse"
326-
|| dft_functional_lower == "scan0"
327-
|| dft_functional_lower == "muller" || dft_functional_lower == "power"
328-
|| dft_functional_lower == "cwp22" || dft_functional_lower == "wp22")
330+
|| dft_functional_lower == "pbe0" || dft_functional_lower == "b3lyp" || dft_functional_lower == "hse"
331+
|| dft_functional_lower == "scan0"
332+
|| dft_functional_lower == "muller" || dft_functional_lower == "power"
333+
|| dft_functional_lower == "cwp22" || dft_functional_lower == "wp22"
334+
|| dft_functional_lower == "lc_pbe"
335+
|| dft_functional_lower == "lc_wpbe"
336+
|| dft_functional_lower == "lrc_wpbe"
337+
|| dft_functional_lower == "lrc_wpbeh"
338+
|| dft_functional_lower == "cam_pbeh")
329339
{
330340
GlobalC::exx_info.info_global.cal_exx = true;
331341

@@ -352,15 +362,13 @@ void Input_Conv::Convert()
352362
{
353363
if(PARAM.inp.basis_type == "lcao")
354364
{
355-
std::map<Conv_Coulomb_Pot_K::Coulomb_Type, std::vector<std::map<std::string,std::string>>> coulomb_param;
356-
coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock].resize(fock_alpha.size());
365+
GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock].resize(fock_alpha.size());
357366
for(std::size_t i=0; i<fock_alpha.size(); ++i)
358367
{
359-
coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock][i] = {{
368+
GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock] = {{
360369
{"alpha", ModuleBase::GlobalFunc::TO_STRING(fock_alpha[i])},
361-
{"Rcut_type", "spencer"} }};
370+
{"singularity_correction", PARAM.inp.exx_singularity_correction} }};
362371
}
363-
GlobalC::exx_info.info_ri.coulomb_settings[Conv_Coulomb_Pot_K::Coulomb_Method::Center2] = std::make_pair(true, coulomb_param);
364372
}
365373
else if(PARAM.inp.basis_type == "lcao_in_pw")
366374
{
@@ -392,16 +400,14 @@ void Input_Conv::Convert()
392400
assert(erfc_alpha.size() == PARAM.inp.exx_erfc_omega.size());
393401
if(PARAM.inp.basis_type == "lcao")
394402
{
395-
std::map<Conv_Coulomb_Pot_K::Coulomb_Type, std::vector<std::map<std::string,std::string>>> coulomb_param;
396-
coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc].resize(erfc_alpha.size());
403+
GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc].resize(erfc_alpha.size());
397404
for(std::size_t i=0; i<erfc_alpha.size(); ++i)
398405
{
399-
coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{
406+
GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{
400407
{"alpha", ModuleBase::GlobalFunc::TO_STRING(erfc_alpha[i])},
401408
{"omega", ModuleBase::GlobalFunc::TO_STRING(PARAM.inp.exx_erfc_omega[i])},
402-
{"Rcut_type", "limits"} }};
409+
{"singularity_correction", PARAM.inp.exx_singularity_correction} }};
403410
}
404-
GlobalC::exx_info.info_ri.coulomb_settings[Conv_Coulomb_Pot_K::Coulomb_Method::Center2] = std::make_pair(true, coulomb_param);
405411
}
406412
}
407413
}

0 commit comments

Comments
 (0)