Skip to content

Commit c5d5eab

Browse files
linpeizePeizeLin
andauthored
update Conv_Coulomb_Pot_K::Ccp_Type (#5702)
Co-authored-by: linpz <[email protected]>
1 parent fd3a6d8 commit c5d5eab

File tree

6 files changed

+39
-42
lines changed

6 files changed

+39
-42
lines changed

source/module_io/input_conv.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -377,7 +377,7 @@ void Input_Conv::Convert()
377377
} else if (dft_functional_lower == "hse") {
378378
GlobalC::exx_info.info_global.cal_exx = true;
379379
GlobalC::exx_info.info_global.ccp_type
380-
= Conv_Coulomb_Pot_K::Ccp_Type::Hse;
380+
= Conv_Coulomb_Pot_K::Ccp_Type::Erfc;
381381
} else if (dft_functional_lower == "opt_orb") {
382382
GlobalC::exx_info.info_global.cal_exx = false;
383383
Exx_Abfs::Jle::generate_matrix = true;
@@ -391,12 +391,12 @@ void Input_Conv::Convert()
391391
else if ( dft_functional_lower == "wp22" )
392392
{
393393
GlobalC::exx_info.info_global.cal_exx = true;
394-
GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::erf; // use the error function erf(w|r-r'|), exx just has the long-range part
394+
GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erf; // use the error function erf(w|r-r'|), exx just has the long-range part
395395
}
396396
else if ( dft_functional_lower == "cwp22" )
397397
{
398398
GlobalC::exx_info.info_global.cal_exx = true;
399-
GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hse; // use the erfc(w|r-r'|), exx just has the short-range part
399+
GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erfc; // use the erfc(w|r-r'|), exx just has the short-range part
400400
}
401401
else {
402402
GlobalC::exx_info.info_global.cal_exx = false;

source/module_lr/esolver_lrtd_lcao.cpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -226,7 +226,7 @@ LR::ESolver_LR<T, TR>::ESolver_LR(ModuleESolver::ESolver_KS_LCAO<T, TR>&& ks_sol
226226
{
227227
// set ccp_type according to the xc_kernel
228228
if (xc_kernel == "hf") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hf; }
229-
else if (xc_kernel == "hse") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hse; }
229+
else if (xc_kernel == "hse") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erfc; }
230230
this->exx_lri = std::make_shared<Exx_LRI<T>>(exx_info.info_ri);
231231
this->exx_lri->init(MPI_COMM_WORLD, this->kv, ks_sol.orb_);
232232
this->exx_lri->cal_exx_ions(input.out_ri_cv);
@@ -396,7 +396,7 @@ LR::ESolver_LR<T, TR>::ESolver_LR(const Input_para& inp, UnitCell& ucell) : inpu
396396
{
397397
// set ccp_type according to the xc_kernel
398398
if (xc_kernel == "hf") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hf; }
399-
else if (xc_kernel == "hse") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hse; }
399+
else if (xc_kernel == "hse") { exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erfc; }
400400
this->exx_lri = std::make_shared<Exx_LRI<T>>(exx_info.info_ri);
401401
this->exx_lri->init(MPI_COMM_WORLD, this->kv, orb);
402402
this->exx_lri->cal_exx_ions(input.out_ri_cv);

source/module_ri/Exx_LRI.hpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -78,7 +78,7 @@ void Exx_LRI<Tdata>::init(const MPI_Comm &mpi_comm_in, const K_Vectors &kv_in, c
7878
const double hf_Rcut = std::pow(0.75 * this->p_kv->get_nkstot_full()/nspin0 * GlobalC::ucell.omega / (ModuleBase::PI), 1.0/3.0);
7979
return {{"hf_Rcut", hf_Rcut}};
8080
}
81-
case Conv_Coulomb_Pot_K::Ccp_Type::Hse:
81+
case Conv_Coulomb_Pot_K::Ccp_Type::Erfc:
8282
return {{"hse_omega", this->info.hse_omega}};
8383
default:
8484
throw std::domain_error(std::string(__FILE__)+" line "+std::to_string(__LINE__)); break;

source/module_ri/conv_coulomb_pot_k.cpp

Lines changed: 26 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -11,9 +11,8 @@ namespace Conv_Coulomb_Pot_K
1111
const std::vector<double> & psif)
1212
{
1313
std::vector<double> psik2_ccp(psif.size());
14-
for( size_t ik=0; ik<psif.size(); ++ik ) {
15-
psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik];
16-
}
14+
for( size_t ik=0; ik<psif.size(); ++ik )
15+
{ psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik]; }
1716
return psik2_ccp;
1817
}
1918

@@ -25,22 +24,20 @@ namespace Conv_Coulomb_Pot_K
2524
const double hf_Rcut)
2625
{
2726
std::vector<double> psik2_ccp(psif.size());
28-
for (size_t ik = 0; ik < psif.size(); ++ik) {
29-
psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * (1 - std::cos(k_radial[ik] * hf_Rcut));
30-
}
27+
for (size_t ik = 0; ik < psif.size(); ++ik)
28+
{ psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * (1 - std::cos(k_radial[ik] * hf_Rcut)); }
3129
return psik2_ccp;
3230
}
3331

3432

35-
std::vector<double> cal_psi_hse(
33+
std::vector<double> cal_psi_erfc(
3634
const std::vector<double> & psif,
3735
const std::vector<double> & k_radial,
3836
const double hse_omega)
3937
{
4038
std::vector<double> psik2_ccp(psif.size());
41-
for( size_t ik=0; ik<psif.size(); ++ik ) {
42-
psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * (1-std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega)));
43-
}
39+
for( size_t ik=0; ik<psif.size(); ++ik )
40+
{ psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * (1-std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega))); }
4441
return psik2_ccp;
4542
}
4643

@@ -53,10 +50,9 @@ namespace Conv_Coulomb_Pot_K
5350
const double hf_Rcut)
5451
{
5552
std::vector<double> psik2_ccp(psif.size());
56-
for( size_t ik=0; ik<psif.size(); ++ik ) {
57-
psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * ( std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega)) - std::cos(k_radial[ik] * hf_Rcut) );
58-
}
59-
// psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * ( std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega)) );
53+
for( size_t ik=0; ik<psif.size(); ++ik )
54+
{ psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * ( std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega)) - std::cos(k_radial[ik] * hf_Rcut) ); }
55+
// { psik2_ccp[ik] = ModuleBase::FOUR_PI * psif[ik] * ( std::exp(-(k_radial[ik]*k_radial[ik])/(4*hse_omega*hse_omega)) ); }
6056
return psik2_ccp;
6157
}
6258

@@ -75,30 +71,28 @@ namespace Conv_Coulomb_Pot_K
7571
psik2_ccp = cal_psi_ccp( orbs.get_psif() ); break;
7672
case Ccp_Type::Hf:
7773
psik2_ccp = cal_psi_hf( orbs.get_psif(), orbs.get_k_radial(), parameter.at("hf_Rcut")); break;
78-
case Ccp_Type::Hse:
79-
psik2_ccp = cal_psi_hse( orbs.get_psif(), orbs.get_k_radial(), parameter.at("hse_omega") ); break;
80-
case Ccp_Type::erf:
74+
case Ccp_Type::Erfc:
75+
psik2_ccp = cal_psi_erfc( orbs.get_psif(), orbs.get_k_radial(), parameter.at("hse_omega") ); break;
76+
case Ccp_Type::Erf:
8177
psik2_ccp = cal_psi_erf( orbs.get_psif(), orbs.get_k_radial(), parameter.at("hse_omega"), parameter.at("hf_Rcut") ); break;
8278
default:
83-
throw( ModuleBase::GlobalFunc::TO_STRING(__FILE__)+" line "+ModuleBase::GlobalFunc::TO_STRING(__LINE__) ); break;
79+
throw( std::string(__FILE__) + " line " + std::to_string(__LINE__) ); break;
8480
}
8581

8682
const double dr = orbs.get_rab().back();
8783
const int Nr = (static_cast<int>(orbs.getNr()*rmesh_times)) | 1;
84+
8885
std::vector<double> rab(Nr);
89-
for( size_t ir=0; ir<std::min(orbs.getNr(),Nr); ++ir ) {
90-
rab[ir] = orbs.getRab(ir);
91-
}
92-
for( size_t ir=orbs.getNr(); ir<Nr; ++ir ) {
93-
rab[ir] = dr;
94-
}
86+
for( size_t ir=0; ir<std::min(orbs.getNr(),Nr); ++ir )
87+
{ rab[ir] = orbs.getRab(ir); }
88+
for( size_t ir=orbs.getNr(); ir<Nr; ++ir )
89+
{ rab[ir] = dr; }
90+
9591
std::vector<double> r_radial(Nr);
96-
for( size_t ir=0; ir<std::min(orbs.getNr(),Nr); ++ir ) {
97-
r_radial[ir] = orbs.getRadial(ir);
98-
}
99-
for( size_t ir=orbs.getNr(); ir<Nr; ++ir ) {
100-
r_radial[ir] = orbs.get_r_radial().back() + (ir - orbs.getNr() + 1) * dr;
101-
}
92+
for( size_t ir=0; ir<std::min(orbs.getNr(),Nr); ++ir )
93+
{ r_radial[ir] = orbs.getRadial(ir); }
94+
for( size_t ir=orbs.getNr(); ir<Nr; ++ir )
95+
{ r_radial[ir] = orbs.get_r_radial().back() + (ir - orbs.getNr() + 1) * dr; }
10296

10397
Numerical_Orbital_Lm orbs_ccp;
10498
orbs_ccp.set_orbital_info(
@@ -126,9 +120,8 @@ namespace Conv_Coulomb_Pot_K
126120
{
127121
for(int ir=orbs.getNr()-1; ir>=0; --ir)
128122
{
129-
if(std::abs(orbs.getPsi(ir))>=psi_threshold) {
130-
return static_cast<double>(ir)/orbs.getNr();
131-
}
123+
if(std::abs(orbs.getPsi(ir))>=psi_threshold)
124+
{ return static_cast<double>(ir)/orbs.getNr(); }
132125
}
133126
return 0.0;
134127
}

source/module_ri/conv_coulomb_pot_k.h

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -10,8 +10,8 @@ namespace Conv_Coulomb_Pot_K
1010
enum class Ccp_Type{ // parameter:
1111
Ccp, //
1212
Hf, // "hf_Rcut"
13-
Hse, // "hse_omega"
14-
erf}; // "hse_omega"
13+
Erfc, // "hse_omega"
14+
Erf}; // "hse_omega", "hf_Rcut"
1515

1616
template<typename T> T cal_orbs_ccp(
1717
const T &orbs,

source/module_ri/exx_lip.hpp

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -334,13 +334,17 @@ void Exx_Lip<T, Device>::qkg2_exp(const int ik, const int iq)
334334
this->sum2_factor += this->recip_qkg2[ig] * std::exp(-info.lambda * qkg2);
335335
this->recip_qkg2[ig] = sqrt(this->recip_qkg2[ig]);
336336
}
337-
else if (Conv_Coulomb_Pot_K::Ccp_Type::Hse == info.ccp_type)
337+
else if (Conv_Coulomb_Pot_K::Ccp_Type::Erfc == info.ccp_type)
338338
{
339339
if (std::abs(qkg2) < 1e-10)
340340
{ this->recip_qkg2[ig] = 1.0 / (2 * info.hse_omega); }
341341
else
342342
{ this->recip_qkg2[ig] = sqrt((1 - std::exp(-qkg2 / (4 * info.hse_omega * info.hse_omega))) / qkg2); }
343343
}
344+
else
345+
{
346+
throw( std::string(__FILE__) + " line " + std::to_string(__LINE__) );
347+
}
344348
}
345349
ModuleBase::timer::tick("Exx_Lip", "qkg2_exp");
346350
}

0 commit comments

Comments
 (0)