From e6dadc430712cbe8103ac5615002dd165efcb5cd Mon Sep 17 00:00:00 2001 From: linpz Date: Fri, 20 Jun 2025 00:50:28 +0800 Subject: [PATCH 1/3] change INPUT exx_ for linear combination of coulomb_param --- docs/advanced/input_files/input-main.md | 56 +++++--- .../module_xc/xc_functional_libxc.cpp | 4 +- source/module_io/input_conv.cpp | 133 +++++++++++++----- source/module_io/read_input.cpp | 1 + source/module_io/read_input_item_exx_dftu.cpp | 129 ++++++++++++----- source/module_parameter/input_parameter.h | 11 +- source/module_rdmft/rdmft.h | 1 - source/module_rdmft/rdmft_tools.h | 1 - source/source_main/driver.cpp | 1 + 9 files changed, 242 insertions(+), 95 deletions(-) diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index bfcd443bea..3ca3b37fd4 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -262,13 +262,15 @@ - [block\_up](#block_up) - [block\_height](#block_height) - [Exact Exchange (Common)](#exact-exchange-common) - - [exx\_hybrid\_alpha](#exx_hybrid_alpha) - - [exx\_hse\_omega](#exx_hse_omega) + - [exx\_fock\_alpha](#exx_fock_alpha) + - [exx\_erfc\_alpha](#exx_erfc_alpha) + - [exx\_erfc\_omega](#exx_erfc_omega) - [exx\_separate\_loop](#exx_separate_loop) - - [Exact Exchange (LCAO/LCAO in PW)](#exact-exchange-lcaolcao-in-pw) - [exx\_hybrid\_step](#exx_hybrid_step) - [exx\_mixing\_beta](#exx_mixing_beta) - - [exx\_lambda](#exx_lambda) + - [Exact Exchange (LCAO in PW)](#exact-exchange-lcao-in-pw) + - [exx\_erfc\_lambda](#exx_erfc_lambda) + - [Exact Exchange (LCAO)](#exact-exchange-lcao) - [exx\_pca\_threshold](#exx_pca_threshold) - [exx\_c\_threshold](#exx_c_threshold) - [exx\_v\_threshold](#exx_v_threshold) @@ -2754,18 +2756,34 @@ The following parameters apply to *[basis_type](#basis_type)==lcao/lcao_in_pw/pw **Availablity**: *[dft_functional](#dft_functional)==hse/hf/pbe0/scan0/opt_orb* or *[rpa](#rpa)==True*. -### exx_hybrid_alpha +### exx_fock_alpha -- **Type**: Real -- **Description**: fraction of Fock exchange in hybrid functionals, so that $E_{X}=\alpha E_{X}+(1-\alpha)E_{X,\text{LDA/GGA}}$ +- **Type**: Real \[Real...\](optional) +- **Description**: fraction of Fock exchange 1/r in hybrid functionals, so that $E_{X} = \alpha E_{X} + (1-\alpha)E_{X,\text{LDA/GGA}}$ - **Default**: - 1: if *[dft_functional](#dft_functional)==hf* - - 0.25: else + - 0.25: if *[dft_functional](#dft_functional)==pbe0* + - 0.2: if *[dft_functional](#dft_functional)==b3lyp* + - 0.25: if *[dft_functional](#dft_functional)==scan0* + - 1: if *[dft_functional](#dft_functional)==muller* + - 1: if *[dft_functional](#dft_functional)==power* + - 1: if *[dft_functional](#dft_functional)==wp22* + - 0: else -### exx_hse_omega +### exx_erfc_alpha -- **Type**: Real -- **Description**: range-separation parameter in HSE functional, such that $1/r=\text{erfc}(\omega r)/r+\text{erf}(\omega r)/r$ +- **Type**: Real \[Real...\](optional) +- **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}}$ +- **Default**: + - 0.25: if *[dft_functional](#dft_functional)==hse* + - 1: if *[dft_functional](#dft_functional)==cwp22* + - -1: if *[dft_functional](#dft_functional)==wp22* + - 0: else + +### exx_erfc_omega + +- **Type**: Real \[Real...\](optional) +- **Description**: range-separation parameter in exchange, such that $1/r=\text{erfc}(\omega r)/r+\text{erf}(\omega r)/r$ - **Default**: 0.11 ### exx_separate_loop @@ -2776,10 +2794,6 @@ The following parameters apply to *[basis_type](#basis_type)==lcao/lcao_in_pw/pw - True: A two-step method is employed, i.e. in the inner iterations, density matrix is updated, while in the outer iterations, $H_{exx}$ is calculated based on density matrix that converges in the inner iteration. - **Default**: True -## Exact Exchange (LCAO/LCAO in PW) - -These variables are relevant when using hybrid functionals with *[basis_type](#basis_type)==lcao/lcao_in_pw*. - ### exx_hybrid_step - **Type**: Integer @@ -2794,13 +2808,21 @@ These variables are relevant when using hybrid functionals with *[basis_type](#b - **Description**: mixing_beta for densty matrix in each iteration of the outer-loop - **Default**: 1.0 -### exx_lambda +## Exact Exchange (LCAO in PW) -- **Type**: Real +These variables are relevant when using hybrid functionals with *[basis_type](#basis_type)==lcao/lcao_in_pw*. + +### exx_fock_lambda + +- **Type**: Real \[Real...\](optional) - **Availability**: *[basis_type](#basis_type)==lcao_in_pw* - **Description**: It is used to compensate for divergence points at G=0 in the evaluation of Fock exchange using *lcao_in_pw* method. - **Default**: 0.3 +## Exact Exchange (LCAO) + +These variables are relevant when using hybrid functionals with *[basis_type](#basis_type)==lcao/lcao_in_pw*. + ### exx_pca_threshold - **Type**: Real diff --git a/source/module_hamilt_general/module_xc/xc_functional_libxc.cpp b/source/module_hamilt_general/module_xc/xc_functional_libxc.cpp index f3ea0fbbd0..e5dc28f100 100644 --- a/source/module_hamilt_general/module_xc/xc_functional_libxc.cpp +++ b/source/module_hamilt_general/module_xc/xc_functional_libxc.cpp @@ -170,11 +170,11 @@ const std::vector in_built_xc_func_ext_params(const int id) GlobalC::exx_info.info_global.hse_omega}; // short-range of B88_X case XC_GGA_X_ITYH: - return {PARAM.inp.exx_hse_omega}; + return {GlobalC::exx_info.info_global.hse_omega}; // short-range of LYP_C case XC_GGA_C_LYPR: return {0.04918, 0.132, 0.2533, 0.349, - 0.35/2.29, 2.0/2.29, PARAM.inp.exx_hse_omega}; + 0.35/2.29, 2.0/2.29, GlobalC::exx_info.info_global.hse_omega}; #endif default: return std::vector{}; diff --git a/source/module_io/input_conv.cpp b/source/module_io/input_conv.cpp index be89b612a8..06d618e9a9 100644 --- a/source/module_io/input_conv.cpp +++ b/source/module_io/input_conv.cpp @@ -321,53 +321,113 @@ void Input_Conv::Convert() PARAM.inp.dft_functional.end(), dft_functional_lower.begin(), tolower); + if (dft_functional_lower == "hf" + || dft_functional_lower == "pbe0" || dft_functional_lower == "b3lyp" || dft_functional_lower == "hse" + || dft_functional_lower == "scan0" + || dft_functional_lower == "muller" || dft_functional_lower == "power" + || dft_functional_lower == "cwp22" || dft_functional_lower == "wp22") + { + GlobalC::exx_info.info_global.cal_exx = true; + + GlobalC::exx_info.info_global.hybrid_alpha = 0; + std::vector fock_alpha(PARAM.inp.exx_fock_alpha.size()); + for(std::size_t i=0; i erfc_alpha(PARAM.inp.exx_erfc_alpha.size()); + for(std::size_t i=0; i0); + for(std::size_t i=0; irank == 0) { diff --git a/source/module_io/read_input_item_exx_dftu.cpp b/source/module_io/read_input_item_exx_dftu.cpp index ea9a9c3927..9a37d52fe0 100644 --- a/source/module_io/read_input_item_exx_dftu.cpp +++ b/source/module_io/read_input_item_exx_dftu.cpp @@ -8,56 +8,127 @@ void ReadInput::item_exx() { // EXX { - Input_Item item("exx_hybrid_alpha"); - item.annotation = "fraction of Fock exchange in hybrid functionals"; - read_sync_string(input.exx_hybrid_alpha); + Input_Item item("exx_fock_alpha"); + item.annotation = "fraction of Fock exchange 1/r in hybrid functionals"; + item.read_value = [](const Input_Item& item, Parameter& para) + { + para.input.exx_fock_alpha = item.str_values; + }; item.reset_value = [](const Input_Item& item, Parameter& para) { - if (para.input.exx_hybrid_alpha == "default") + if (para.input.exx_fock_alpha.size()==1 && para.input.exx_fock_alpha[0]=="default") { std::string& dft_functional = para.input.dft_functional; std::string dft_functional_lower = dft_functional; std::transform(dft_functional.begin(), dft_functional.end(), dft_functional_lower.begin(), tolower); if (dft_functional_lower == "hf") { - para.input.exx_hybrid_alpha = "1"; + para.input.exx_fock_alpha = {"1"}; } - else if (dft_functional_lower == "pbe0" || dft_functional_lower == "hse" - || dft_functional_lower == "scan0") + else if (dft_functional_lower == "pbe0" || dft_functional_lower == "scan0") { - para.input.exx_hybrid_alpha = "0.25"; + para.input.exx_fock_alpha = {"0.25"}; } - // added by jghan 2024-07-06 - else if (dft_functional_lower == "muller" || dft_functional_lower == "power" - || dft_functional_lower == "wp22" || dft_functional_lower == "cwp22") + else if (dft_functional_lower == "b3lyp") { - para.input.exx_hybrid_alpha = "1"; + para.input.exx_fock_alpha = {"0.2"}; } - else if (dft_functional_lower == "b3lyp") + else if (dft_functional_lower == "muller" || dft_functional_lower == "power" ) + { + para.input.exx_fock_alpha = {"1"}; + } + else if (dft_functional_lower == "wp22") { - para.input.exx_hybrid_alpha = "0.2"; + para.input.exx_fock_alpha = {"1"}; } else - { // no exx in scf, but will change to non-zero in - // postprocess like rpa - para.input.exx_hybrid_alpha = "0"; + { // no exx in scf, but will change to non-zero in postprocess like rpa + para.input.exx_fock_alpha = {}; } } }; - item.check_value = [](const Input_Item& item, const Parameter& para) { - const double exx_hybrid_alpha_value = std::stod(para.input.exx_hybrid_alpha); - if (exx_hybrid_alpha_value < 0 || exx_hybrid_alpha_value > 1) + sync_stringvec(input.exx_fock_alpha, para.input.exx_fock_alpha.size(), ""); + this->add_item(item); + } + { + Input_Item item("exx_erfc_alpha"); + item.annotation = "fraction of exchange erfc(wr)/r in hybrid functionals"; + item.read_value = [](const Input_Item& item, Parameter& para) + { + para.input.exx_erfc_alpha = item.str_values; + }; + item.reset_value = [](const Input_Item& item, Parameter& para) + { + if (para.input.exx_erfc_alpha.size()==1 && para.input.exx_erfc_alpha[0]=="default") { - ModuleBase::WARNING_QUIT("ReadInput", - "The Hartree-Fock fraction (exx_hybrid_alpha) can only be in range [0, 1]"); + std::string& dft_functional = para.input.dft_functional; + std::string dft_functional_lower = dft_functional; + std::transform(dft_functional.begin(), dft_functional.end(), dft_functional_lower.begin(), tolower); + if (dft_functional_lower == "hse") + { + para.input.exx_erfc_alpha = {"0.25"}; + } + else if (dft_functional_lower == "cwp22") + { + para.input.exx_erfc_alpha = {"1"}; + } + else if (dft_functional_lower == "wp22") + { + para.input.exx_erfc_alpha = {"-1"}; + } + else + { // no exx in scf, but will change to non-zero in postprocess like rpa + para.input.exx_erfc_alpha = {}; + } } }; + sync_stringvec(input.exx_erfc_alpha, para.input.exx_erfc_alpha.size(), ""); this->add_item(item); } { - Input_Item item("exx_hse_omega"); - item.annotation = "range-separation parameter in HSE functional"; - read_sync_double(input.exx_hse_omega); + Input_Item item("exx_erfc_omega"); + item.annotation = "range-separation parameter erfc(wr)/r in hybrid functionals"; + item.read_value = [](const Input_Item& item, Parameter& para) + { + para.input.exx_erfc_omega = item.str_values; + }; + item.reset_value = [](const Input_Item& item, Parameter& para) + { + if (para.input.exx_erfc_omega.size()==1 && para.input.exx_erfc_omega[0]=="default") + { + std::string& dft_functional = para.input.dft_functional; + std::string dft_functional_lower = dft_functional; + std::transform(dft_functional.begin(), dft_functional.end(), dft_functional_lower.begin(), tolower); + if (dft_functional_lower == "hse" || dft_functional_lower == "cwp22" || dft_functional_lower == "wp22") + { + para.input.exx_erfc_omega = {"0.11"}; + } + else + { + para.input.exx_erfc_omega = {}; + } + } + }; + sync_stringvec(input.exx_erfc_omega, para.input.exx_erfc_omega.size(), ""); + this->add_item(item); + } + { + Input_Item item("exx_fock_lambda"); + item.annotation = "used to compensate for divergence points at G=0 in " + "the evaluation of Fock exchange using " + "lcao_in_pw method"; + item.read_value = [](const Input_Item& item, Parameter& para) + { + para.input.exx_fock_lambda = item.str_values; + }; + item.reset_value = [](const Input_Item& item, Parameter& para) + { + if (para.input.exx_fock_lambda.size()==1 && para.input.exx_fock_lambda[0]=="default") + { + para.input.exx_fock_lambda = std::vector(para.input.exx_fock_alpha.size(), "0.3"); + } + }; + sync_stringvec(input.exx_fock_lambda, para.input.exx_fock_lambda.size(), ""); this->add_item(item); } { @@ -87,14 +158,6 @@ void ReadInput::item_exx() read_sync_double(input.exx_mixing_beta); this->add_item(item); } - { - Input_Item item("exx_lambda"); - item.annotation = "used to compensate for divergence points at G=0 in " - "the evaluation of Fock exchange using " - "lcao_in_pw method"; - read_sync_double(input.exx_lambda); - this->add_item(item); - } { Input_Item item("exx_real_number"); item.annotation = "exx calculated in real or complex"; diff --git a/source/module_parameter/input_parameter.h b/source/module_parameter/input_parameter.h index fd3e29aade..d60839319a 100644 --- a/source/module_parameter/input_parameter.h +++ b/source/module_parameter/input_parameter.h @@ -515,15 +515,16 @@ struct Input_para // exx // Peize Lin add 2018-06-20 // ========================================================== - std::string exx_hybrid_alpha = "default"; ///< fraction of Fock exchange in hybrid functionals - double exx_hse_omega = 0.11; ///< range-separation parameter in HSE functional + std::vector exx_fock_alpha = {"default"}; ///< fraction of Fock exchange 1/r in hybrid functionals + std::vector exx_fock_lambda = {"default"}; ///< used to compensate for divergence points at G=0 in the + ///< evaluation of Fock exchange using lcao_in_pw method + std::vector exx_erfc_alpha = {"default"}; ///< fraction of exchange erfc(wr)/r in hybrid functionals + std::vector exx_erfc_omega = {"default"}; ///< range-separation parameter in HSE functional bool exx_separate_loop = true; ///< if 1, a two-step method is employed, else it will start ///< with a GGA-Loop, and then Hybrid-Loop int exx_hybrid_step = 100; ///< the maximal electronic iteration number in ///< the evaluation of Fock exchange double exx_mixing_beta = 1.0; ///< mixing_beta for outer-loop when exx_separate_loop=1 - double exx_lambda = 0.3; ///< used to compensate for divergence points at G=0 in the - ///< evaluation of Fock exchange using lcao_in_pw method std::string exx_real_number = "default"; ///< exx calculated in real or complex double exx_pca_threshold = 0.0001; ///< threshold to screen on-site ABFs in exx double exx_c_threshold = 0.0001; ///< threshold to screen C matrix in exx @@ -647,7 +648,7 @@ struct Input_para // RDMFT jghan added on 2024-07-06 bool rdmft = false; // rdmft, reduced density matrix funcional theory double rdmft_power_alpha = 0.656; // the alpha parameter of power-functional, g(occ_number) = occ_number^alpha - // double rdmft_wp22_omega; // the omega parameter of wp22-functional = exx_hse_omega + // double rdmft_wp22_omega; // the omega parameter of wp22-functional = exx_erfc_omega // ============== #Parameters (22.EXX PW) ===================== // EXX for planewave basis, rhx0820 2025-03-10 diff --git a/source/module_rdmft/rdmft.h b/source/module_rdmft/rdmft.h index bc02480637..6e0dd5b046 100644 --- a/source/module_rdmft/rdmft.h +++ b/source/module_rdmft/rdmft.h @@ -26,7 +26,6 @@ #include "module_ri/Exx_LRI.h" #include "module_ri/module_exx_symmetry/symmetry_rotation.h" // there are some operator reload to print data in different formats -#include "module_ri/test_code/test_function.h" #endif #include diff --git a/source/module_rdmft/rdmft_tools.h b/source/module_rdmft/rdmft_tools.h index a18e436c47..8075264c1f 100644 --- a/source/module_rdmft/rdmft_tools.h +++ b/source/module_rdmft/rdmft_tools.h @@ -30,7 +30,6 @@ #include "module_ri/RI_2D_Comm.h" #include "module_ri/Exx_LRI.h" // there are some operator reload to print data in different formats -#include "module_ri/test_code/test_function.h" #endif #include diff --git a/source/source_main/driver.cpp b/source/source_main/driver.cpp index 37f99d482a..6ee43a68f1 100644 --- a/source/source_main/driver.cpp +++ b/source/source_main/driver.cpp @@ -109,6 +109,7 @@ void Driver::print_start_info() void Driver::reading() { + ModuleBase::TITLE("Driver", "reading"); ModuleBase::timer::tick("Driver", "reading"); // temperarily GlobalV::MY_RANK = PARAM.globalv.myrank; From c0bdbfbbda017ca61fd3339d8085268811808b18 Mon Sep 17 00:00:00 2001 From: linpz Date: Fri, 20 Jun 2025 03:04:36 +0800 Subject: [PATCH 2/3] git revert a1452937d02b0fb18 --- .../module_xc/exx_info.h | 12 +- source/module_io/input_conv.cpp | 14 +- source/module_ri/Exx_LRI.h | 16 +- source/module_ri/Exx_LRI.hpp | 99 +-- source/module_ri/LRI_CV.h | 6 +- source/module_ri/LRI_CV.hpp | 29 +- source/module_ri/LRI_CV_Tools.h | 348 +++------- source/module_ri/LRI_CV_Tools.hpp | 614 ++++++------------ source/module_ri/RPA_LRI.h | 1 + source/module_ri/RPA_LRI.hpp | 13 +- source/module_ri/conv_coulomb_pot_k.h | 1 - 11 files changed, 376 insertions(+), 777 deletions(-) diff --git a/source/module_hamilt_general/module_xc/exx_info.h b/source/module_hamilt_general/module_xc/exx_info.h index 9aca0a4088..24acbdb9b2 100644 --- a/source/module_hamilt_general/module_xc/exx_info.h +++ b/source/module_hamilt_general/module_xc/exx_info.h @@ -15,10 +15,7 @@ struct Exx_Info { bool cal_exx = false; - std::unordered_map>>>> coulomb_settings; + std::unordered_map>> coulomb_param; // Fock: // "alpha": "0" // "Rcut_type": "limits" / "spencer" @@ -55,10 +52,7 @@ struct Exx_Info struct Exx_Info_RI { - const std::unordered_map>>>> &coulomb_settings; + const std::unordered_map>> &coulomb_param; bool real_number = false; @@ -80,7 +74,7 @@ struct Exx_Info int abfs_Lmax = 0; // tmp Exx_Info_RI(const Exx_Info::Exx_Info_Global& info_global) - : coulomb_settings(info_global.coulomb_settings) + : coulomb_param(info_global.coulomb_param) { } }; diff --git a/source/module_io/input_conv.cpp b/source/module_io/input_conv.cpp index d24e8dec47..e9c8582449 100644 --- a/source/module_io/input_conv.cpp +++ b/source/module_io/input_conv.cpp @@ -328,12 +328,10 @@ void Input_Conv::Convert() { GlobalC::exx_info.info_global.cal_exx = true; GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hf; - std::unordered_map>> coulomb_param; - coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock] = {{ + GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock] = {{ {"alpha", "1"}, {"Rcut_type", "spencer"}, {"lambda", std::to_string(PARAM.inp.exx_lambda)} }}; - GlobalC::exx_info.info_global.coulomb_settings[Conv_Coulomb_Pot_K::Coulomb_Method::Center2] = std::make_pair(true, coulomb_param); } // use the error function erf(w|r-r'|), exx just has the short-range part else if (dft_functional_lower == "hse" @@ -341,28 +339,24 @@ void Input_Conv::Convert() { GlobalC::exx_info.info_global.cal_exx = true; GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erfc; - std::unordered_map>> coulomb_param; - coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{ + GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{ {"alpha", "1"}, {"omega", std::to_string(PARAM.inp.exx_hse_omega)}, {"Rcut_type", "limits"} }}; - GlobalC::exx_info.info_global.coulomb_settings[Conv_Coulomb_Pot_K::Coulomb_Method::Center2] = std::make_pair(true, coulomb_param); } // use the error function erf(w|r-r'|), exx just has the long-range part else if ( dft_functional_lower == "wp22" ) { GlobalC::exx_info.info_global.cal_exx = true; GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Erf; - std::unordered_map>> coulomb_param; - coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock] = {{ + GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock] = {{ {"alpha", "1"}, {"Rcut_type", "spencer"}, {"lambda", std::to_string(PARAM.inp.exx_lambda)} }}; - coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{ + GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc] = {{ {"alpha", "-1"}, {"omega", std::to_string(PARAM.inp.exx_hse_omega)}, {"Rcut_type", "limits"} }}; - GlobalC::exx_info.info_global.coulomb_settings[Conv_Coulomb_Pot_K::Coulomb_Method::Center2] = std::make_pair(true, coulomb_param); } #ifdef __EXX else if (dft_functional_lower == "opt_orb") diff --git a/source/module_ri/Exx_LRI.h b/source/module_ri/Exx_LRI.h index 1b24261b12..63f90150e9 100644 --- a/source/module_ri/Exx_LRI.h +++ b/source/module_ri/Exx_LRI.h @@ -15,7 +15,6 @@ #include #include #include -#include #include #include @@ -38,15 +37,6 @@ class OperatorLREXX; } -template -class Exx_Obj -{ - // match with Conv_Coulomb_Pot_K::Coulomb_Method - public: - LRI_CV cv; - std::vector>> abfs_ccp; -}; - template class Exx_LRI { @@ -95,9 +85,9 @@ class Exx_LRI std::vector>> lcaos; std::vector>> abfs; - //std::vector>> abfs_ccp; - std::unordered_map> exx_objs; - //LRI_CV cv; + std::vector>> abfs_ccp; + + LRI_CV cv; RI::Exx exx_lri; void post_process_Hexx( std::map>> &Hexxs_io ) const; diff --git a/source/module_ri/Exx_LRI.hpp b/source/module_ri/Exx_LRI.hpp index 7d0b4d4c57..54d401c02e 100644 --- a/source/module_ri/Exx_LRI.hpp +++ b/source/module_ri/Exx_LRI.hpp @@ -48,18 +48,19 @@ void Exx_LRI::init(const MPI_Comm &mpi_comm_in, { this->abfs = Exx_Abfs::IO::construct_abfs( abfs_same_atom, orb, this->info.files_abfs, this->info.kmesh_times ); } Exx_Abfs::Construct_Orbs::print_orbs_size(ucell, this->abfs, GlobalV::ofs_running); + const std::unordered_map>> + coulomb_param_updated = RI_Util::update_coulomb_param(this->info.coulomb_param, ucell.omega, this->p_kv->get_nkstot_full()); + this->abfs_ccp = Conv_Coulomb_Pot_K::cal_orbs_ccp(this->abfs, coulomb_param_updated, this->info.ccp_rmesh_times); + for( size_t T=0; T!=this->abfs.size(); ++T ) { GlobalC::exx_info.info_ri.abfs_Lmax = std::max( GlobalC::exx_info.info_ri.abfs_Lmax, static_cast(this->abfs[T].size())-1 ); } - for(const auto &settings_list : this->info.coulomb_settings) - { - const std::unordered_map>> - coulomb_param_updated = RI_Util::update_coulomb_param(settings_list.second.second, ucell.omega, this->p_kv->get_nkstot_full()); - this->exx_objs[settings_list.first].abfs_ccp = Conv_Coulomb_Pot_K::cal_orbs_ccp(this->abfs, coulomb_param_updated, this->info.ccp_rmesh_times); - this->exx_objs[settings_list.first].cv.set_orbitals(ucell, orb, - this->lcaos, this->abfs, this->exx_objs[settings_list.first].abfs_ccp, - this->info.kmesh_times, this->info.ccp_rmesh_times ); - } + this->cv.set_orbitals( + ucell, + orb, + this->lcaos, this->abfs, this->abfs_ccp, + this->info.kmesh_times, this->info.ccp_rmesh_times ); + ModuleBase::timer::tick("Exx_LRI", "init"); } @@ -94,40 +95,26 @@ void Exx_LRI::cal_exx_ions(const UnitCell& ucell, const std::pair, std::vector>>>> list_As_Vs = RI::Distribute_Equally::distribute_atoms_periods(this->mpi_comm, atoms, period_Vs, 2, false); - std::map>> Vs; - std::map, Ndim>>> dVs; - for(const auto &settings_list : this->info.coulomb_settings) - { - std::map>> - Vs_temp = this->exx_objs[settings_list.first].cv.cal_Vs(ucell, - list_As_Vs.first, list_As_Vs.second[0], - {{"writable_Vws",true}}); - this->exx_objs[settings_list.first].cv.Vws = LRI_CV_Tools::get_CVws(ucell,Vs_temp); - Vs = Vs.empty() ? Vs_temp : LRI_CV_Tools::add(Vs, Vs_temp); - - if(PARAM.inp.cal_force || PARAM.inp.cal_stress) - { - std::map, Ndim>>> - dVs_temp = this->exx_objs[settings_list.first].cv.cal_dVs(ucell, - list_As_Vs.first, list_As_Vs.second[0], - {{"writable_dVws",true}}); - this->exx_objs[settings_list.first].cv.dVws = LRI_CV_Tools::get_dCVws(ucell,dVs_temp); - dVs = dVs.empty() ? dVs_temp : LRI_CV_Tools::add(dVs, dVs_temp); - } - - } + std::map>> + Vs = this->cv.cal_Vs(ucell, + list_As_Vs.first, list_As_Vs.second[0], + {{"writable_Vws",true}}); + this->cv.Vws = LRI_CV_Tools::get_CVws(ucell,Vs); if (write_cv && GlobalV::MY_RANK == 0) { LRI_CV_Tools::write_Vs_abf(Vs, PARAM.globalv.global_out_dir + "Vs"); } this->exx_lri.set_Vs(std::move(Vs), this->info.V_threshold); if(PARAM.inp.cal_force || PARAM.inp.cal_stress) { - std::array>>, Ndim> - dVs_order = LRI_CV_Tools::change_order(std::move(dVs)); - this->exx_lri.set_dVs(std::move(dVs_order), this->info.V_grad_threshold); + std::array>>,3> + dVs = this->cv.cal_dVs(ucell, + list_As_Vs.first, list_As_Vs.second[0], + {{"writable_dVws",true}}); + this->cv.dVws = LRI_CV_Tools::get_dCVws(ucell,dVs); + this->exx_lri.set_dVs(std::move(dVs), this->info.V_grad_threshold); if(PARAM.inp.cal_stress) { - std::array>>,3>,3> dVRs = LRI_CV_Tools::cal_dMRs(ucell,dVs_order); + std::array>>,3>,3> dVRs = LRI_CV_Tools::cal_dMRs(ucell,dVs); this->exx_lri.set_dVRs(std::move(dVRs), this->info.V_grad_R_threshold); } } @@ -136,47 +123,29 @@ void Exx_LRI::cal_exx_ions(const UnitCell& ucell, const std::pair, std::vector>>>> list_As_Cs = RI::Distribute_Equally::distribute_atoms_periods(this->mpi_comm, atoms, period_Cs, 2, false); - std::map>> Cs; - std::map, 3>>> dCs; - for(const auto &settings_list : this->info.coulomb_settings) - { - if(settings_list.second.first) - { - std::pair>>, - std::map, 3>>>> - Cs_dCs = this->exx_objs[settings_list.first].cv.cal_Cs_dCs( - ucell, - list_As_Cs.first, list_As_Cs.second[0], - {{"cal_dC",PARAM.inp.cal_force||PARAM.inp.cal_stress}, - {"writable_Cws",true}, {"writable_dCws",true}, {"writable_Vws",false}, {"writable_dVws",false}}); - std::map>> &Cs_temp = std::get<0>(Cs_dCs); - this->exx_objs[settings_list.first].cv.Cws = LRI_CV_Tools::get_CVws(ucell,Cs_temp); - Cs = Cs.empty() ? Cs_temp : LRI_CV_Tools::add(Cs, Cs_temp); - - if(PARAM.inp.cal_force || PARAM.inp.cal_stress) - { - std::map, 3>>> &dCs_temp = std::get<1>(Cs_dCs); - this->exx_objs[settings_list.first].cv.dCws = LRI_CV_Tools::get_dCVws(ucell,dCs_temp); - dCs = dCs.empty() ? dCs_temp : LRI_CV_Tools::add(dCs, dCs_temp); - } - } - } + std::pair>>, std::array>>,3>> + Cs_dCs = this->cv.cal_Cs_dCs( + ucell, + list_As_Cs.first, list_As_Cs.second[0], + {{"cal_dC",PARAM.inp.cal_force||PARAM.inp.cal_stress}, + {"writable_Cws",true}, {"writable_dCws",true}, {"writable_Vws",false}, {"writable_dVws",false}}); + std::map>> &Cs = std::get<0>(Cs_dCs); + this->cv.Cws = LRI_CV_Tools::get_CVws(ucell,Cs); if (write_cv && GlobalV::MY_RANK == 0) { LRI_CV_Tools::write_Cs_ao(Cs, PARAM.globalv.global_out_dir + "Cs"); } this->exx_lri.set_Cs(std::move(Cs), this->info.C_threshold); if(PARAM.inp.cal_force || PARAM.inp.cal_stress) { - std::array>>, Ndim> - dCs_order = LRI_CV_Tools::change_order(std::move(dCs)); - this->exx_lri.set_dCs(std::move(dCs_order), this->info.C_grad_threshold); + std::array>>,3> &dCs = std::get<1>(Cs_dCs); + this->cv.dCws = LRI_CV_Tools::get_dCVws(ucell,dCs); + this->exx_lri.set_dCs(std::move(dCs), this->info.C_grad_threshold); if(PARAM.inp.cal_stress) { - std::array>>,3>,3> dCRs = LRI_CV_Tools::cal_dMRs(ucell,dCs_order); + std::array>>,3>,3> dCRs = LRI_CV_Tools::cal_dMRs(ucell,dCs); this->exx_lri.set_dCRs(std::move(dCRs), this->info.C_grad_R_threshold); } } - ModuleBase::timer::tick("Exx_LRI", "cal_exx_ions"); } diff --git a/source/module_ri/LRI_CV.h b/source/module_ri/LRI_CV.h index 8e455c2bf7..bfca46208c 100644 --- a/source/module_ri/LRI_CV.h +++ b/source/module_ri/LRI_CV.h @@ -47,14 +47,14 @@ class LRI_CV const std::vector &list_A0, const std::vector &list_A1, const std::map &flags); // "writable_Vws" - inline std::map, 3>>> + inline std::array>>,3> cal_dVs( const UnitCell &ucell, const std::vector &list_A0, const std::vector &list_A1, const std::map &flags); // "writable_dVws" - std::pair>>, - std::map, 3>>>> + std::pair>>, + std::array>>,3>> cal_Cs_dCs( const UnitCell &ucell, const std::vector &list_A0, diff --git a/source/module_ri/LRI_CV.hpp b/source/module_ri/LRI_CV.hpp index 7b9d5aba1e..06e2b51685 100644 --- a/source/module_ri/LRI_CV.hpp +++ b/source/module_ri/LRI_CV.hpp @@ -146,14 +146,15 @@ auto LRI_CV::cal_dVs( const std::vector &list_A0, const std::vector &list_A1, const std::map &flags) // + "writable_dVws" --> std::map, 3>>> +-> std::array>>,3> { ModuleBase::TITLE("LRI_CV","cal_dVs"); const T_func_DPcal_data,3>> func_DPcal_dV = std::bind( &LRI_CV::DPcal_dV, this, std::placeholders::_1, std::placeholders::_2, std::placeholders::_3, std::placeholders::_4); - return this->cal_datas(ucell,list_A0, list_A1, flags, this->ccp_rmesh_times, func_DPcal_dV); + return LRI_CV_Tools::change_order( + this->cal_datas(ucell,list_A0, list_A1, flags, this->ccp_rmesh_times, func_DPcal_dV)); } template @@ -162,9 +163,7 @@ auto LRI_CV::cal_Cs_dCs( const std::vector &list_A0, const std::vector &list_A1, const std::map &flags) // "cal_dC" + "writable_Cws", "writable_dCws", "writable_Vws", "writable_dVws" --> std::pair< - std::map>>, - std::map, 3>>>> +-> std::pair>>, std::array>>,3>> { ModuleBase::TITLE("LRI_CV","cal_Cs_dCs"); const T_func_DPcal_data, std::array,3>>> @@ -174,16 +173,16 @@ auto LRI_CV::cal_Cs_dCs( std::map, std::array,3>>>> Cs_dCs_tmp = this->cal_datas(ucell,list_A0, list_A1, flags, std::min(1.0,this->ccp_rmesh_times), func_DPcal_C_dC); - std::map>> Cs; - std::map, 3>>> dCs; - for (auto& Cs_dCs_A: Cs_dCs_tmp) - for (auto& Cs_dCs_B: Cs_dCs_A.second) { - Cs[Cs_dCs_A.first][Cs_dCs_B.first] - = std::move(std::get<0>(Cs_dCs_B.second)); - if (flags.at("cal_dC")) - dCs[Cs_dCs_A.first][Cs_dCs_B.first] - = std::move(std::get<1>(Cs_dCs_B.second)); - } + std::map>> Cs; + std::array>>,3> dCs; + for(auto &Cs_dCs_A : Cs_dCs_tmp) + for(auto &Cs_dCs_B : Cs_dCs_A.second) + { + Cs[Cs_dCs_A.first][Cs_dCs_B.first] = std::move(std::get<0>(Cs_dCs_B.second)); + if(flags.at("cal_dC")) + for(int ix=0; ix<3; ++ix) + dCs[ix][Cs_dCs_A.first][Cs_dCs_B.first] = std::move(std::get<1>(Cs_dCs_B.second)[ix]); + } return std::make_pair(Cs, dCs); } diff --git a/source/module_ri/LRI_CV_Tools.h b/source/module_ri/LRI_CV_Tools.h index 28a7d6263d..43eea41ada 100644 --- a/source/module_ri/LRI_CV_Tools.h +++ b/source/module_ri/LRI_CV_Tools.h @@ -6,265 +6,113 @@ #ifndef LRI_CV_TOOLS_H #define LRI_CV_TOOLS_H - -#include "abfs.h" #include "source_base/abfs-vector3_order.h" #include "module_ri/abfs.h" -#include #include -#include + #include -#include +#include #include +#include namespace LRI_CV_Tools { -template -extern RI::Tensor cal_I(const RI::Tensor& m); -template -extern std::vector>> cal_I(const std::vector>>& ms); - -template -inline RI::Tensor transform_Rm(const RI::Tensor& V); -template -inline std::array, 3> transform_Rm(const std::array, 3>& dV); - -// template inline bool exist(const T &V); - -// template -// extern Treturn mul1(const T1 &t1, const T2 &t2); -// template -// extern Treturn mul2(const T1 &mat, const T2 &vec); - -template -inline bool exist(const RI::Tensor& V); -template -inline bool exist(const std::array& dV); - -template -extern RI::Tensor mul1(const RI::Tensor& t1, const RI::Tensor& t2); -template -extern std::array mul1(const std::array& t1, const T& t2); - -template -extern std::vector> mul2(const std::vector>>& mat, - const std::vector>& vec); -template -extern std::array mul2(const T1& t1, const std::array& t2); -template -extern RI::Tensor mul2(const T& t1, const RI::Tensor& t2); -template -extern std::map> mul2(const T& t1, const std::map>& t2); - -// template -// std::array operator-(const std::array &v1, const std::array -// &v2); template std::vector operator-(const std::vector &v1, -// const std::vector &v2); -template -extern std::vector> minus(const std::vector>& v1, - const std::vector>& v2); -template -extern std::array>, N> minus( - std::array>, N>& v1, - std::array>, N>& v2); -template -inline std::map>> minus( - std::map>>& v1, - std::map>>& v2); -template -extern std::map> minus(std::map>& v1, - std::map>& v2); - -template -extern std::vector> add(const std::vector>& v1, - const std::vector>& v2); -template -extern std::array>, N> add( - std::array>, N>& v1, - std::array>, N>& v2); -template -inline std::map>> add( - std::map>>& v1, - std::map>>& v2); -template -extern std::map> add(std::map>& v1, - std::map>& v2); - -template -extern std::array negative(const std::array& v_in); - -// template T transpose12(const T &c_in); -template -RI::Tensor transpose12(const RI::Tensor& c_in); -template -std::array transpose12(const std::array& c_in); - -template -extern std::array, N> change_order(std::vector>&& ds_in); -template -std::vector> change_order(std::array, N>&& ds_in); -template -extern std::array>, N> change_order(std::vector>>&& ds_in); -template -extern std::array>, N> change_order( - std::map>>&& ds_in); -template -extern std::map>> change_order( - std::array>, N>&& ds_in); - -template -extern std::array cal_latvec_range(const double& rcut_times, - const UnitCell& ucell, - const std::vector& orb_cutoff); - -template -extern std::map, RI::Tensor>>> get_CVws( - const UnitCell& ucell, - const std::map>, RI::Tensor>>& CVs); -template -extern std::map, std::array, 3>>>> get_dCVws( - const UnitCell& ucell, - const std::map>, std::array, 3>>>& dCVs); -template -extern std::array, RI::Tensor>>, 3>, 3> cal_dMRs( - const UnitCell& ucell, - const std::array, RI::Tensor>>, 3>& dMs); - -using TC = std::array; -using TAC = std::pair; -template -using TLRI = std::map>>; -template -TLRI read_Cs_ao(const std::string& file_path, const double& threshold = 1e-10); -template -void write_Cs_ao(const TLRI& Vs, const std::string& file_path); -template -TLRI read_Vs_abf(const std::string& file_path, const double& threshold = 1e-10); -template -void write_Vs_abf(const TLRI& Vs, const std::string& file_path); - -template -struct is_std_array : std::false_type -{ -}; -template -struct is_std_array> : std::true_type -{ -}; -template -struct is_tensor : std::false_type -{ -}; -template -struct is_tensor> : std::true_type -{ -}; - -template -struct TinType; - -template -struct TinType> -{ - using type = T; -}; - -template -struct TinType, N>> -{ - using type = T; -}; - -template ::value>> -inline void init_elem(Tdata& data, const size_t ndim0, const size_t ndim1) -{ - data = Tdata({ndim0, ndim1}); -}; -template -extern void init_elem(std::array, N>& data, const size_t ndim0, const size_t ndim1); - -template ::value && !is_tensor::value>> -inline void add_elem(Tdata& data, const Tdata& val, const Tdata& frac) -{ - data += frac * val; -}; -template -extern void add_elem(std::array& data, const T& val, const T& frac); -template ::value>> -inline void add_elem(const Tdata& data, - const int lmp, - const int lmq, - const typename TinType::type& val, - const typename TinType::type& frac) -{ - data(lmp, lmq) += frac * val; -}; -template -extern void add_elem(std::array, N>& data, - const int lmp, - const int lmq, - const std::array& val, - const T& frac); -template ::value>> -inline void add_elem(Tdata& data, - const int lmp0, - const int lmq0, - const Tdata& val, - const int lmp1, - const int lmq1, - const typename TinType::type& frac) -{ - data(lmp0, lmq0) += frac * val(lmp1, lmq1); -}; -template -extern void add_elem(std::array, N>& data, - const int lmp0, - const int lmq0, - const std::array, N>& val, - const int lmp1, - const int lmq1, - const T& frac); - -template -inline RI::Tensor convert(RI::Tensor&& data); -template -extern std::array, N> convert(std::array, N>&& data); - -// template -// typename std::enable_if::value, T>::type -// inline check_zero(T value) { -// return (std::abs(value) < 1e-8) ? static_cast(0) : value; -// } - -// template -// typename std::enable_if::value, T>::type -// inline check_zero(const T& value) { -// using RealType = typename T::value_type; -// RealType real_part = std::real(value); -// RealType imag_part = std::imag(value); - -// real_part = (std::abs(real_part) < 1e-8) ? 0 : real_part; -// imag_part = (std::abs(imag_part) < 1e-8) ? 0 : imag_part; - -// return std::complex(real_part, imag_part); -// } - -// template -// extern RI::Tensor check_zero(RI::Tensor&& data); -// template -// extern std::array, N> check_zero(std::array, N>&& data); - -template -struct plus -{ - T operator()(const T& lhs, const T& rhs) const - { - using namespace RI::Array_Operator; - return lhs + rhs; - } -}; -} // namespace LRI_CV_Tools + template extern RI::Tensor cal_I( const RI::Tensor &m ); + template extern std::vector>> cal_I( const std::vector>> &ms ); + + template inline RI::Tensor transform_Rm(const RI::Tensor &V ); + template inline std::array,3> transform_Rm(const std::array,3> &dV); + + //template inline bool exist(const T &V); + + //template + //extern Treturn mul1(const T1 &t1, const T2 &t2); + //template + //extern Treturn mul2(const T1 &mat, const T2 &vec); + + template inline bool exist(const RI::Tensor &V); + template inline bool exist(const std::array &dV); + + template + extern RI::Tensor mul1(const RI::Tensor &t1, const RI::Tensor &t2); + template + extern std::array mul1(const std::array &t1, const T &t2); + + template + extern std::vector> mul2(const std::vector>> &mat, + const std::vector> &vec); + template + extern std::array mul2(const T1 &t1, const std::array &t2); + + //template + //std::array operator-(const std::array &v1, const std::array &v2); + //template + //std::vector operator-(const std::vector &v1, const std::vector &v2); + template + extern std::vector> minus( + const std::vector> &v1, + const std::vector> &v2); + + template + extern std::array negative(const std::array &v_in); + + //template T transpose12(const T &c_in); + template RI::Tensor transpose12(const RI::Tensor &c_in); + template std::array transpose12(const std::array &c_in); + + template + extern std::array,N> + change_order( + std::vector> &&ds_in); + template + std::vector> + change_order( + std::array,N> &&ds_in); + template + extern std::array>,N> + change_order( + std::vector>> &&ds_in); + template + extern std::array>,N> + change_order( + std::map>> && ds_in); + + template + extern std::array cal_latvec_range(const double &rcut_times, + const UnitCell &ucell, + const std::vector& orb_cutoff); + + template + extern std::map,RI::Tensor>>> + get_CVws( + const UnitCell &ucell, + const std::map>,RI::Tensor>> &CVs); + template + extern std::map,std::array,3>>>> + get_dCVws( + const UnitCell &ucell, + const std::array>,RI::Tensor>>,3> &dCVs); + + template + extern std::array,RI::Tensor>>,3>,3> + cal_dMRs( + const UnitCell &ucell, + const std::array,RI::Tensor>>,3> &dMs); + + using TC = std::array; + using TAC = std::pair; + template + using TLRI = std::map>>; + template + TLRI read_Cs_ao(const std::string& file_path, const double& threshold = 1e-10); + template + void write_Cs_ao(const TLRI& Vs, const std::string& file_path); + template + TLRI read_Vs_abf(const std::string& file_path, const double& threshold = 1e-10); + template + void write_Vs_abf(const TLRI& Vs, const std::string& file_path); +} #include "LRI_CV_Tools.hpp" #include "write_ri_cv.hpp" diff --git a/source/module_ri/LRI_CV_Tools.hpp b/source/module_ri/LRI_CV_Tools.hpp index e38973f3d2..6e11d620de 100644 --- a/source/module_ri/LRI_CV_Tools.hpp +++ b/source/module_ri/LRI_CV_Tools.hpp @@ -6,367 +6,251 @@ #ifndef LRI_CV_TOOLS_HPP #define LRI_CV_TOOLS_HPP - -#include "../source_base/mathzone.h" -#include "Inverse_Matrix.h" #include "LRI_CV_Tools.h" -#include "RI_Util.h" - -#include -#include +#include "Inverse_Matrix.h" +#include "../source_base/mathzone.h" #include "../module_hamilt_pw/hamilt_pwdft/global.h" -template -RI::Tensor LRI_CV_Tools::cal_I(const RI::Tensor& m) { - Inverse_Matrix I; - I.input(m); - I.cal_inverse(Inverse_Matrix::Method::potrf); - return I.output(); +template +RI::Tensor +LRI_CV_Tools::cal_I( const RI::Tensor &m ) +{ + Inverse_Matrix I; + I.input(m); + I.cal_inverse( Inverse_Matrix::Method::potrf ); + return I.output(); } -template +template std::vector>> - LRI_CV_Tools::cal_I(const std::vector>>& ms) { - Inverse_Matrix I; - I.input(ms); - I.cal_inverse(Inverse_Matrix::Method::potrf); - return I.output({ms[0][0].shape[0], ms[1][0].shape[0]}, - {ms[0][0].shape[1], ms[0][1].shape[1]}); +LRI_CV_Tools::cal_I( const std::vector>> &ms ) +{ + Inverse_Matrix I; + I.input(ms); + I.cal_inverse( Inverse_Matrix::Method::potrf ); + return I.output({ms[0][0].shape[0], ms[1][0].shape[0]}, {ms[0][0].shape[1], ms[0][1].shape[1]}); } -template -RI::Tensor LRI_CV_Tools::transform_Rm(const RI::Tensor& V) { - return V.transpose(); + + +template +RI::Tensor LRI_CV_Tools::transform_Rm(const RI::Tensor &V) +{ + return V.transpose(); } -template -std::array, 3> - LRI_CV_Tools::transform_Rm(const std::array, 3>& dV) { - return std::array, 3>{-dV[0].transpose(), - -dV[1].transpose(), - -dV[2].transpose()}; +template +std::array,3> LRI_CV_Tools::transform_Rm(const std::array,3> &dV) +{ + return std::array,3>{-dV[0].transpose(), -dV[1].transpose(), -dV[2].transpose()}; } -template -bool LRI_CV_Tools::exist(const RI::Tensor& V) { - return !V.empty(); +template +bool LRI_CV_Tools::exist(const RI::Tensor &V) +{ + return !V.empty(); } -template -bool LRI_CV_Tools::exist(const std::array& dV) { - for (size_t i = 0; i < 3; ++i) - if (!dV[i].empty()) - return true; - return false; +template +bool LRI_CV_Tools::exist(const std::array &dV) +{ + for(size_t i=0; i<3; ++i) + { + if(!dV[i].empty()) + return true; + } + return false; } -template -RI::Tensor LRI_CV_Tools::mul1(const RI::Tensor& t1, - const RI::Tensor& t2) { - const size_t sa0 = t1.shape[0], sa1 = t2.shape[0], sl0 = t2.shape[1], - sl1 = t2.shape[2]; - return (t1 * t2.reshape({sa1, sl0 * sl1})).reshape({sa0, sl0, sl1}); + +template +RI::Tensor LRI_CV_Tools::mul1( + const RI::Tensor &t1, + const RI::Tensor &t2) +{ + const size_t sa0=t1.shape[0], sa1=t2.shape[0], sl0=t2.shape[1], sl1=t2.shape[2]; + return (t1 * t2.reshape({sa1,sl0*sl1})).reshape({sa0,sl0,sl1}); } -template -std::array LRI_CV_Tools::mul1(const std::array& t1, const T& t2) { - return std::array{mul1(t1[0], t2), mul1(t1[1], t2), mul1(t1[2], t2)}; +template +std::array LRI_CV_Tools::mul1( + const std::array &t1, + const T &t2) +{ + return std::array{ + mul1(t1[0],t2), mul1(t1[1],t2), mul1(t1[2],t2) }; } /* template std::array LRI_CV_Tools::mul1( - const T &t1, - const std::array &t2) + const T &t1, + const std::array &t2) { - return std::array{ - mul1(t1,t2[0]), mul1(t1,t2[1]), mul1(t1,t2[2]) }; + return std::array{ + mul1(t1,t2[0]), mul1(t1,t2[1]), mul1(t1,t2[2]) }; } */ -template -std::vector> - LRI_CV_Tools::mul2(const std::vector>>& mat, - const std::vector>& vec) { - const size_t sa0 = vec[0].shape[0], sa1 = vec[1].shape[0], - sl0 = vec[0].shape[1], sl1 = vec[0].shape[2]; - const RI::Tensor vec0 = vec[0].reshape({sa0, sl0 * sl1}), - vec1 = vec[1].reshape({sa1, sl0 * sl1}); - return std::vector>{ - (mat[0][0] * vec0 + mat[0][1] * vec1).reshape({sa0, sl0, sl1}), - (mat[1][0] * vec0 + mat[1][1] * vec1).reshape({sa1, sl0, sl1})}; +template +std::vector> LRI_CV_Tools::mul2( + const std::vector>> &mat, + const std::vector> &vec) +{ + const size_t sa0=vec[0].shape[0], sa1=vec[1].shape[0], sl0=vec[0].shape[1], sl1=vec[0].shape[2]; + const RI::Tensor vec0=vec[0].reshape({sa0,sl0*sl1}), vec1=vec[1].reshape({sa1,sl0*sl1}); + return std::vector> + {( mat[0][0]*vec0 + mat[0][1]*vec1 ).reshape({sa0,sl0,sl1}), + ( mat[1][0]*vec0 + mat[1][1]*vec1 ).reshape({sa1,sl0,sl1})}; } /* template std::array LRI_CV_Tools::mul2( - const std::array &t1, - const T2 &t2) + const std::array &t1, + const T2 &t2) { - return std::array{ - mul2(t1[0],t2), mul2(t1[1],t2), mul2(t1[2],t2) }; + return std::array{ + mul2(t1[0],t2), mul2(t1[1],t2), mul2(t1[2],t2) }; } */ -template -std::array LRI_CV_Tools::mul2(const T1& t1, - const std::array& t2) { - return std::array{mul2(t1, t2[0]), mul2(t1, t2[1]), mul2(t1, t2[2])}; -} - -template -RI::Tensor LRI_CV_Tools::mul2(const T& t1, const RI::Tensor& t2) { - return t1 * t2; -} - -template -std::map> - LRI_CV_Tools::mul2(const T& t1, - const std::map>& t2) { - std::map> res; - for (const auto& outerPair: t2) { - const TkeyA keyA = outerPair.first; - const std::map& innerMap = outerPair.second; - std::map newInnerMap; - - for (const auto& innerPair: innerMap) { - const TkeyB keyB = innerPair.first; - const Tvalue value = innerPair.second; - newInnerMap[keyB] = mul2(t1, value); - } - - res[keyA] = newInnerMap; - } - - return res; +template +std::array LRI_CV_Tools::mul2( + const T1 &t1, + const std::array &t2) +{ + return std::array{ + mul2(t1,t2[0]), mul2(t1,t2[1]), mul2(t1,t2[2]) }; } /* template -std::array LRI_CV_Tools::operator-(const std::array &v1, const -std::array &v2) +std::array LRI_CV_Tools::operator-(const std::array &v1, const std::array &v2) { - std::array v; - for(std::size_t i=0; i v; + for(std::size_t i=0; i -std::vector LRI_CV_Tools::operator-(const std::vector &v1, const -std::vector &v2) +std::vector LRI_CV_Tools::operator-(const std::vector &v1, const std::vector &v2) { - assert(v1.size()==v2.size()); - std::vector v(v1.size()); - for(std::size_t i=0; i v(v1.size()); + for(std::size_t i=0; i -std::vector> - LRI_CV_Tools::minus(const std::vector>& v1, - const std::vector>& v2) { - assert(v1.size() == v2.size()); - std::vector> v(v1.size()); - for (std::size_t i = 0; i < v.size(); ++i) - for (std::size_t j = 0; j < N; ++j) - v[i][j] = v1[i][j] - v2[i][j]; - return v; -} - -template -std::map>> LRI_CV_Tools::minus( - std::map>>& v1, - std::map>>& v2) { - std::array>, N> v1_order - = change_order(std::move(v1)); - std::array>, N> v2_order - = change_order(std::move(v2)); - auto dv = minus(v1_order, v2_order); - return change_order(std::move(dv)); -} - -template -std::array>, N> LRI_CV_Tools::minus( - std::array>, N>& v1, - std::array>, N>& v2) { - std::array>, N> dv; - for (size_t i = 0; i != N; ++i) - dv[i] = minus(v1[i], v2[i]); - return dv; -} - -template -std::map> - LRI_CV_Tools::minus(std::map>& v1, - std::map>& v2) { - assert(v1.size() == v2.size()); - using namespace RI::Map_Operator; - using namespace RI::Array_Operator; - - std::map> dv; - auto it1 = v1.begin(); - auto it2 = v2.begin(); - while (it1 != v1.end() && it2 != v2.end()) { - assert(it1->first == it2->first); - const TkeyA& keyA = it1->first; - const std::map& map1 = it1->second; - const std::map& map2 = it2->second; - dv[keyA] = map1 - map2; - ++it1; - ++it2; - } - return dv; -} - -template -std::vector> - LRI_CV_Tools::add(const std::vector>& v1, - const std::vector>& v2) { - assert(v1.size() == v2.size()); - std::vector> v(v1.size()); - for (std::size_t i = 0; i < v.size(); ++i) - for (std::size_t j = 0; j < N; ++j) - v[i][j] = v1[i][j] + v2[i][j]; - return v; +template +std::vector> LRI_CV_Tools::minus( + const std::vector> &v1, + const std::vector> &v2) +{ + assert(v1.size()==v2.size()); + std::vector> v(v1.size()); + for(std::size_t i=0; i -std::map>> LRI_CV_Tools::add( - std::map>>& v1, - std::map>>& v2) { - std::array>, N> v1_order - = change_order(std::move(v1)); - std::array>, N> v2_order - = change_order(std::move(v2)); - auto dv = add(v1_order, v2_order); - return change_order(std::move(dv)); -} -template -std::array>, N> LRI_CV_Tools::add( - std::array>, N>& v1, - std::array>, N>& v2) { - std::array>, N> dv; - for (size_t i = 0; i != N; ++i) - dv[i] = add(v1[i], v2[i]); - return dv; +template +std::array LRI_CV_Tools::negative(const std::array &v_in) +{ + std::array v_out; + for(std::size_t i=0; i -std::map> - LRI_CV_Tools::add(std::map>& v1, - std::map>& v2) { - assert(v1.size() == v2.size()); - using namespace RI::Map_Operator; - using namespace RI::Array_Operator; - - std::map> dv; - auto it1 = v1.begin(); - auto it2 = v2.begin(); - while (it1 != v1.end() && it2 != v2.end()) { - assert(it1->first == it2->first); - const TkeyA& keyA = it1->first; - const std::map& map1 = it1->second; - const std::map& map2 = it2->second; - dv[keyA] = map1 + map2; - ++it1; - ++it2; - } - return dv; -} -template -std::array LRI_CV_Tools::negative(const std::array& v_in) { - std::array v_out; - for (std::size_t i = 0; i < N; ++i) - v_out[i] = -v_in[i]; - return v_out; +template +RI::Tensor LRI_CV_Tools::transpose12(const RI::Tensor &c_in) +{ + RI::Tensor c_out({c_in.shape[0], c_in.shape[2], c_in.shape[1]}); + for(size_t i0=0; i0 -RI::Tensor LRI_CV_Tools::transpose12(const RI::Tensor& c_in) { - RI::Tensor c_out({c_in.shape[0], c_in.shape[2], c_in.shape[1]}); - for (size_t i0 = 0; i0 < c_in.shape[0]; ++i0) - for (size_t i1 = 0; i1 < c_in.shape[1]; ++i1) - for (size_t i2 = 0; i2 < c_in.shape[2]; ++i2) - c_out(i0, i2, i1) = c_in(i0, i1, i2); - return c_out; +template +std::array LRI_CV_Tools::transpose12(const std::array &c_in) +{ + std::array c_out; + for(size_t i=0; i -std::array LRI_CV_Tools::transpose12(const std::array& c_in) { - std::array c_out; - for (size_t i = 0; i < N; ++i) - c_out[i] = transpose12(c_in[i]); - return c_out; -} -template -std::array, N> - LRI_CV_Tools::change_order(std::vector>&& ds_in) { - std::array, N> ds; - for (int ix = 0; ix < N; ++ix) { - ds[ix].resize(ds_in.size()); - for (int iv = 0; iv < ds_in.size(); ++iv) - ds[ix][iv] = std::move(ds_in[iv][ix]); - } - return ds; +template +std::array,N> +LRI_CV_Tools::change_order(std::vector> &&ds_in) +{ + std::array,N> ds; + for(int ix=0; ix -std::vector> - LRI_CV_Tools::change_order(std::array, N>&& ds_in) { - std::vector> ds(ds_in[0].size()); - for (int ix = 0; ix < N; ++ix) { - assert(ds.size() == ds_in[ix].size()); - for (int iv = 0; iv < ds.size(); ++iv) - ds[iv][ix] = std::move(ds_in[ix][iv]); - } - return ds; +template +std::vector> +LRI_CV_Tools::change_order(std::array,N> &&ds_in) +{ + std::vector> ds(ds_in[0].size()); + for(int ix=0; ix -std::array>, N> LRI_CV_Tools::change_order( - std::vector>>&& ds_in) { - std::array>, N> ds; - for (int ix = 0; ix < N; ++ix) { - ds[ix].resize(ds_in.size()); - for (int i0 = 0; i0 < ds_in.size(); ++i0) { - ds[ix][i0].resize(ds_in[i0].size()); - for (int i1 = 0; i1 < ds_in[i0].size(); ++i1) - ds[ix][i0][i1] = std::move(ds_in[i0][i1][ix]); - } - } - return ds; +template +std::array>,N> +LRI_CV_Tools::change_order(std::vector>> &&ds_in) +{ + std::array>,N> ds; + for(int ix=0; ix -std::array>, N> - LRI_CV_Tools::change_order( - std::map>>&& ds_in) { - std::array>, N> ds; - for (auto& ds_A: ds_in) - for (auto& ds_B: ds_A.second) - for (int ix = 0; ix < N; ++ix) - ds[ix][ds_A.first][ds_B.first] = std::move(ds_B.second[ix]); - return ds; +template +std::array>,N> +LRI_CV_Tools::change_order(std::map>> && ds_in) +{ + std::array>,N> ds; + for(auto &ds_A : ds_in) + for(auto &ds_B : ds_A.second) + for(int ix=0; ix -std::map>> - LRI_CV_Tools::change_order( - std::array>, N>&& ds_in) { - std::map>> ds; - for (int ix = 0; ix < N; ++ix) - for (auto& ds_A: ds_in[ix]) - for (auto& ds_B: ds_A.second) - ds[ds_A.first][ds_B.first][ix] = std::move(ds_B.second); - return ds; -} -template -std::array LRI_CV_Tools::cal_latvec_range(const double& rcut_times, +template +std::array +LRI_CV_Tools::cal_latvec_range(const double &rcut_times, const UnitCell &ucell, - const std::vector& orb_cutoff) { - double Rcut_max = 0; - for(int T=0; T& orb_cutoff) +{ + double Rcut_max = 0; + for(int T=0; T proj = ModuleBase::Mathzone::latvec_projection( std::array,3>{ucell.a1, ucell.a2, ucell.a3}); @@ -405,115 +289,36 @@ LRI_CV_Tools::get_CVws( return CVws; } -template -std::map, std::array, 3>>>> LRI_CV_Tools:: - get_dCVws(const UnitCell& ucell, - const std::map>, std::array, 3>>>& dCVs) +template +std::map,std::array,3>>>> +LRI_CV_Tools::get_dCVws( + const UnitCell &ucell, + const std::array>,RI::Tensor>>,3> &dCVs) { - std::map, std::array, 3>>>> dCVws; - for (const auto& dCVs_A: dCVs) - { - const TA iat0 = dCVs_A.first; - const int it0 = ucell.iat2it[iat0]; - const int ia0 = ucell.iat2ia[iat0]; - const ModuleBase::Vector3 tau0 = ucell.atoms[it0].tau[ia0]; - for (const auto& dCVs_B: dCVs_A.second) - { - const TA iat1 = dCVs_B.first.first; - const int it1 = ucell.iat2it[iat1]; - const int ia1 = ucell.iat2ia[iat1]; - const std::array& cell1 = dCVs_B.first.second; - const ModuleBase::Vector3 tau1 = ucell.atoms[it1].tau[ia1]; - const Abfs::Vector3_Order R_delta - = -tau0 + tau1 + (RI_Util::array3_to_Vector3(cell1) * ucell.latvec); - dCVws[it0][it1][R_delta] = dCVs_B.second; - } - } - return dCVws; -} - -template -void LRI_CV_Tools::init_elem(std::array, N>& data, - const size_t ndim0, - const size_t ndim1) { - for (size_t i = 0; i < N; ++i) { - data[i] = RI::Tensor({ndim0, ndim1}); - } -} - -template -void LRI_CV_Tools::add_elem(std::array& data, - const T& val, - const T& frac) { - for (size_t i = 0; i < N; ++i) - data[i] += frac * val; -} - -template -void LRI_CV_Tools::add_elem(std::array, N>& data, - const int lmp, - const int lmq, - const std::array& val, - const T& frac) { - for (size_t i = 0; i < N; ++i) { - data[i](lmp, lmq) += frac * val[i]; - } -} - -template -void LRI_CV_Tools::add_elem(std::array, N>& data, - const int lmp0, - const int lmq0, - const std::array, N>& val, - const int lmp1, - const int lmq1, - const T& frac) { - for (size_t i = 0; i < N; ++i) { - data[i](lmp0, lmq0) += frac * val[i](lmp1, lmq1); - } -} - -template -RI::Tensor LRI_CV_Tools::convert(RI::Tensor&& data) { - return RI::Global_Func::convert(data); -} - -template -std::array, N> - LRI_CV_Tools::convert(std::array, N>&& data) { - std::array, N> out; - for (size_t i = 0; i != N; ++i) - out[i] = RI::Global_Func::convert(data[i]); - return out; + std::map,std::array,3>>>> dCVws; + for(int ix=0; ix<3; ++ix) + { + for(const auto &dCVs_A : dCVs[ix]) + { + const TA iat0 = dCVs_A.first; + const int it0 = ucell.iat2it[iat0]; + const int ia0 = ucell.iat2ia[iat0]; + const ModuleBase::Vector3 tau0 = ucell.atoms[it0].tau[ia0]; + for(const auto &dCVs_B : dCVs_A.second) + { + const TA iat1 = dCVs_B.first.first; + const int it1 = ucell.iat2it[iat1]; + const int ia1 = ucell.iat2ia[iat1]; + const std::array &cell1 = dCVs_B.first.second; + const ModuleBase::Vector3 tau1 = ucell.atoms[it1].tau[ia1]; + const Abfs::Vector3_Order R_delta = -tau0+tau1+(RI_Util::array3_to_Vector3(cell1)*ucell.latvec); + dCVws[it0][it1][R_delta][ix] = dCVs_B.second; + } + } + } + return dCVws; } -// template -// RI::Tensor LRI_CV_Tools::check_zero(RI::Tensor&& data) { -// RI::Tensor result(data.shape); - -// const std::size_t rows = data.shape[0]; -// const std::size_t cols = data.shape[1]; - -// for (std::size_t i = 0; i < rows; ++i) { -// for (std::size_t j = 0; j < cols; ++j) { -// result(i, j) = LRI_CV_Tools::check_zero(data(i, j)); -// } -// } - -// return result; -// } - -// template -// std::array, N> -// LRI_CV_Tools::check_zero(std::array, N>&& data) { -// std::array, N> result; - -// for (size_t i = 0; i != N; ++i) -// result[i] = LRI_CV_Tools::check_zero(std::move(data[i])); - -// return result; -// } - // dMRs[ipos0][ipos1] = \nabla_{ipos0} M R_{ipos1} template @@ -556,5 +361,4 @@ LRI_CV_Tools::cal_dMRs( } return dMRs; } - #endif diff --git a/source/module_ri/RPA_LRI.h b/source/module_ri/RPA_LRI.h index 60f1645c86..2007715019 100644 --- a/source/module_ri/RPA_LRI.h +++ b/source/module_ri/RPA_LRI.h @@ -69,6 +69,7 @@ template class RPA_LRI std::vector>> lcaos; std::vector>> abfs; + std::vector>> abfs_ccp; // Exx_LRI exx_postSCF_double(info); // LRI_CV cv; diff --git a/source/module_ri/RPA_LRI.hpp b/source/module_ri/RPA_LRI.hpp index 567bd9d78e..79c851fa3a 100644 --- a/source/module_ri/RPA_LRI.hpp +++ b/source/module_ri/RPA_LRI.hpp @@ -23,6 +23,7 @@ void RPA_LRI::init(const MPI_Comm& mpi_comm_in, const K_Vectors& kv_in this->orb_cutoff_ = orb_cutoff; this->lcaos = exx_lri_rpa.lcaos; this->abfs = exx_lri_rpa.abfs; + this->abfs_ccp = exx_lri_rpa.abfs_ccp; this->p_kv = &kv_in; // this->cv = std::move(exx_lri_rpa.cv); @@ -43,7 +44,7 @@ void RPA_LRI::cal_rpa_cv(const UnitCell& ucell) const std::pair, std::vector>>>> list_As_Vs = RI::Distribute_Equally::distribute_atoms(this->mpi_comm, atoms, period_Vs, 2, false); - std::map>> Vs = exx_lri_rpa.exx_objs[Conv_Coulomb_Pot_K::Coulomb_Method::Center2].cv.cal_Vs(ucell, + std::map>> Vs = exx_lri_rpa.cv.cal_Vs(ucell, list_As_Vs.first, list_As_Vs.second[0], { @@ -56,8 +57,8 @@ void RPA_LRI::cal_rpa_cv(const UnitCell& ucell) = RI::Distribute_Equally::distribute_atoms_periods(this->mpi_comm, atoms, period_Cs, 2, false); std::pair>>, - std::map, 3>>>> - Cs_dCs = exx_lri_rpa.exx_objs[Conv_Coulomb_Pot_K::Coulomb_Method::Center2].cv.cal_Cs_dCs(ucell, + std::array>>, 3>> + Cs_dCs = exx_lri_rpa.cv.cal_Cs_dCs(ucell, list_As_Cs.first, list_As_Cs.second[0], { @@ -330,7 +331,7 @@ void RPA_LRI::out_coulomb_k(const UnitCell &ucell) for (int I = 0; I != ucell.nat; I++) { mu_shift[I] = all_mu; - all_mu += exx_lri_rpa.exx_objs[Conv_Coulomb_Pot_K::Coulomb_Method::Center2].cv.get_index_abfs_size(ucell.iat2it[I]); + all_mu += exx_lri_rpa.cv.get_index_abfs_size(ucell.iat2it[I]); } const int nks_tot = PARAM.inp.nspin == 2 ? (int)p_kv->get_nks() / 2 : p_kv->get_nks(); std::stringstream ss; @@ -343,7 +344,7 @@ void RPA_LRI::out_coulomb_k(const UnitCell &ucell) for (auto& Ip: this->Vs_period) { auto I = Ip.first; - size_t mu_num = exx_lri_rpa.exx_objs[Conv_Coulomb_Pot_K::Coulomb_Method::Center2].cv.get_index_abfs_size(ucell.iat2it[I]); + size_t mu_num = exx_lri_rpa.cv.get_index_abfs_size(ucell.iat2it[I]); for (int ik = 0; ik != nks_tot; ik++) { @@ -372,7 +373,7 @@ void RPA_LRI::out_coulomb_k(const UnitCell &ucell) { auto iJ = vq_Jp.first; auto& vq_J = vq_Jp.second; - size_t nu_num = exx_lri_rpa.exx_objs[Conv_Coulomb_Pot_K::Coulomb_Method::Center2].cv.get_index_abfs_size(ucell.iat2it[iJ]); + size_t nu_num = exx_lri_rpa.cv.get_index_abfs_size(ucell.iat2it[iJ]); ofs << all_mu << " " << mu_shift[I] + 1 << " " << mu_shift[I] + mu_num << " " << mu_shift[iJ] + 1 << " " << mu_shift[iJ] + nu_num << std::endl; ofs << ik + 1 << " " << p_kv->wk[ik] / 2.0 * PARAM.inp.nspin << std::endl; diff --git a/source/module_ri/conv_coulomb_pot_k.h b/source/module_ri/conv_coulomb_pot_k.h index c9870e71c4..7daca7e1d7 100644 --- a/source/module_ri/conv_coulomb_pot_k.h +++ b/source/module_ri/conv_coulomb_pot_k.h @@ -14,7 +14,6 @@ namespace Conv_Coulomb_Pot_K Hf, // "hf_Rcut" Erfc, // "hse_omega" Erf}; // "hse_omega", "hf_Rcut" - enum class Coulomb_Method{Center2, Ewald}; // Different methods for constructing the Coulomb matrix. template extern T cal_orbs_ccp( const T &orbs, From 7980251085e035c82729fd4d7a982d2b18acd334 Mon Sep 17 00:00:00 2001 From: linpz Date: Fri, 20 Jun 2025 04:40:41 +0800 Subject: [PATCH 3/3] update read_input_ptest.cpp, read_input_item_test.cpp, integrate_tests/INPUT --- source/module_io/test/read_input_ptest.cpp | 14 +++--- source/module_io/test/support/INPUT | 14 +++--- .../test_serial/read_input_item_test.cpp | 43 ++++++------------- tests/08_EXX/07_KP_CR_HSE/INPUT | 1 - tests/08_EXX/08_KP_HSE_symm/INPUT | 1 - tests/08_EXX/09_KP_HSE_comp/INPUT | 1 - tests/08_EXX/10_KP_HF/INPUT | 1 - tests/08_EXX/11_KP_PBE0/INPUT | 1 - tests/08_EXX/12_KP_OXC/INPUT | 1 - 9 files changed, 25 insertions(+), 52 deletions(-) diff --git a/source/module_io/test/read_input_ptest.cpp b/source/module_io/test/read_input_ptest.cpp index 87de823b3a..cee426fe38 100644 --- a/source/module_io/test/read_input_ptest.cpp +++ b/source/module_io/test/read_input_ptest.cpp @@ -278,26 +278,24 @@ TEST_F(InputParaTest, ParaRead) EXPECT_EQ(param.inp.vdw_cutoff_period[0], 3); EXPECT_EQ(param.inp.vdw_cutoff_period[1], 3); EXPECT_EQ(param.inp.vdw_cutoff_period[2], 3); - EXPECT_EQ(std::stod(param.inp.exx_hybrid_alpha), 0.25); + EXPECT_DOUBLE_EQ(std::stod(param.inp.exx_fock_alpha[0]), 1); + EXPECT_DOUBLE_EQ(std::stod(param.inp.exx_erfc_alpha[0]), 0.25); EXPECT_EQ(param.inp.exx_real_number, "1"); - EXPECT_DOUBLE_EQ(param.inp.exx_hse_omega, 0.11); + EXPECT_DOUBLE_EQ(std::stod(param.inp.exx_erfc_omega[0]), 0.11); EXPECT_TRUE(param.inp.exx_separate_loop); EXPECT_EQ(param.inp.exx_hybrid_step, 100); - EXPECT_DOUBLE_EQ(param.inp.exx_lambda, 0.3); + EXPECT_DOUBLE_EQ(std::stod(param.inp.exx_fock_lambda[0]), 0.3); EXPECT_DOUBLE_EQ(param.inp.exx_mixing_beta, 1.0); EXPECT_DOUBLE_EQ(param.inp.exx_pca_threshold, 0); EXPECT_DOUBLE_EQ(param.inp.exx_c_threshold, 0); EXPECT_DOUBLE_EQ(param.inp.exx_v_threshold, 0); EXPECT_DOUBLE_EQ(param.inp.exx_dm_threshold, 0); - EXPECT_DOUBLE_EQ(param.inp.exx_schwarz_threshold, 0); - EXPECT_DOUBLE_EQ(param.inp.exx_cauchy_threshold, 0); EXPECT_DOUBLE_EQ(param.inp.exx_c_grad_threshold, 0); EXPECT_DOUBLE_EQ(param.inp.exx_v_grad_threshold, 0); - EXPECT_DOUBLE_EQ(param.inp.exx_cauchy_force_threshold, 0); - EXPECT_DOUBLE_EQ(param.inp.exx_cauchy_stress_threshold, 0); + EXPECT_DOUBLE_EQ(param.inp.exx_c_grad_r_threshold, 0); + EXPECT_DOUBLE_EQ(param.inp.exx_v_grad_r_threshold, 0); EXPECT_EQ(param.inp.exx_ccp_rmesh_times, "1.5"); EXPECT_DOUBLE_EQ(param.inp.rpa_ccp_rmesh_times, 10.0); - EXPECT_EQ(param.inp.exx_distribute_type, "htime"); EXPECT_EQ(param.inp.exx_opt_orb_lmax, 0); EXPECT_DOUBLE_EQ(param.inp.exx_opt_orb_ecut, 0.0); EXPECT_DOUBLE_EQ(param.inp.exx_opt_orb_tolerence, 0.0); diff --git a/source/module_io/test/support/INPUT b/source/module_io/test/support/INPUT index d366d5b0db..db194646d5 100644 --- a/source/module_io/test/support/INPUT +++ b/source/module_io/test/support/INPUT @@ -274,25 +274,23 @@ vdw_cn_thr_unit Bohr #unit of cn_thr, Bohr or Angstrom vdw_cutoff_period 3 3 3 #periods of periodic structure #Parameters (14.exx) -exx_hybrid_alpha default # -exx_hse_omega 0.11 # +exx_fock_alpha 1 # +exx_erfc_alpha 0.25 # +exx_erfc_omega 0.11 # exx_separate_loop 1 #0 or 1 exx_hybrid_step 100 # exx_mixing_beta 1.0 # -exx_lambda 0.3 # +exx_fock_lambda 0.3 # exx_real_number default #0 or 1 exx_pca_threshold 0 # exx_c_threshold 0 # exx_v_threshold 0 # exx_dm_threshold 0 # -exx_schwarz_threshold 0 # -exx_cauchy_threshold 0 # exx_c_grad_threshold 0 # exx_v_grad_threshold 0 # -exx_cauchy_force_threshold 0 # -exx_cauchy_stress_threshold 0 # +exx_c_grad_r_threshold 0 # +exx_v_grad_r_threshold 0 # exx_ccp_rmesh_times default # -exx_distribute_type htime #htime or kmeans1 or kmeans2 exx_opt_orb_lmax 0 # exx_opt_orb_ecut 0 # exx_opt_orb_tolerence 0 # diff --git a/source/module_io/test_serial/read_input_item_test.cpp b/source/module_io/test_serial/read_input_item_test.cpp index edba8327c4..b54b1d9602 100644 --- a/source/module_io/test_serial/read_input_item_test.cpp +++ b/source/module_io/test_serial/read_input_item_test.cpp @@ -1263,38 +1263,29 @@ TEST_F(InputTest, Item_test2) output = testing::internal::GetCapturedStdout(); EXPECT_THAT(output, testing::HasSubstr("NOTICE")); } - { // exx_hybrid_alpha - auto it = find_label("exx_hybrid_alpha", readinput.input_lists); - param.input.exx_hybrid_alpha = "default"; + { // exx_fock_alpha + auto it = find_label("exx_fock_alpha", readinput.input_lists); + param.input.exx_fock_alpha[0] = "default"; param.input.dft_functional = "HF"; it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.exx_hybrid_alpha, "1"); + EXPECT_EQ(param.input.exx_fock_alpha[0], "1"); - param.input.exx_hybrid_alpha = "default"; + param.input.exx_fock_alpha[0] = "default"; param.input.dft_functional = "PBE0"; it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.exx_hybrid_alpha, "0.25"); + EXPECT_EQ(param.input.exx_fock_alpha[0], "0.25"); - param.input.exx_hybrid_alpha = "default"; + param.input.exx_fock_alpha[0] = "default"; param.input.dft_functional = "SCAN0"; it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.exx_hybrid_alpha, "0.25"); - - param.input.exx_hybrid_alpha = "default"; + EXPECT_EQ(param.input.exx_fock_alpha[0], "0.25"); + } + { // exx_erfc_alpha + auto it = find_label("exx_erfc_alpha", readinput.input_lists); + param.input.exx_erfc_alpha[0] = "default"; param.input.dft_functional = "HSE"; it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.exx_hybrid_alpha, "0.25"); - - param.input.exx_hybrid_alpha = "default"; - param.input.dft_functional = "none"; - it->second.reset_value(it->second, param); - EXPECT_EQ(param.input.exx_hybrid_alpha, "0"); - - param.input.exx_hybrid_alpha = "-1"; - testing::internal::CaptureStdout(); - EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("NOTICE")); + EXPECT_EQ(param.input.exx_erfc_alpha[0], "0.25"); } { // exx_hybrid_step auto it = find_label("exx_hybrid_step", readinput.input_lists); @@ -1349,14 +1340,6 @@ TEST_F(InputTest, Item_test2) output = testing::internal::GetCapturedStdout(); EXPECT_THAT(output, testing::HasSubstr("NOTICE")); } - { // exx_distribute_type - auto it = find_label("exx_distribute_type", readinput.input_lists); - param.input.exx_distribute_type = "none"; - testing::internal::CaptureStdout(); - EXPECT_EXIT(it->second.check_value(it->second, param), ::testing::ExitedWithCode(1), ""); - output = testing::internal::GetCapturedStdout(); - EXPECT_THAT(output, testing::HasSubstr("NOTICE")); - } { // exx_opt_orb_lmax auto it = find_label("exx_opt_orb_lmax", readinput.input_lists); param.input.exx_opt_orb_lmax = -1; diff --git a/tests/08_EXX/07_KP_CR_HSE/INPUT b/tests/08_EXX/07_KP_CR_HSE/INPUT index f36f16b619..27e3eb61d8 100644 --- a/tests/08_EXX/07_KP_CR_HSE/INPUT +++ b/tests/08_EXX/07_KP_CR_HSE/INPUT @@ -31,7 +31,6 @@ exx_pca_threshold 1E-3 exx_c_threshold 1E-3 exx_v_threshold 1 exx_dm_threshold 1E-3 -exx_cauchy_threshold 1E-4 exx_ccp_rmesh_times 1 exx_separate_loop 1 exx_hybrid_step 2 diff --git a/tests/08_EXX/08_KP_HSE_symm/INPUT b/tests/08_EXX/08_KP_HSE_symm/INPUT index 233f905222..aa98760d7a 100644 --- a/tests/08_EXX/08_KP_HSE_symm/INPUT +++ b/tests/08_EXX/08_KP_HSE_symm/INPUT @@ -32,7 +32,6 @@ exx_pca_threshold 1E-4 exx_c_threshold 1E-4 exx_v_threshold 1 exx_dm_threshold 1E-4 -exx_cauchy_threshold 1E-6 exx_ccp_rmesh_times 1 exx_separate_loop 1 exx_hybrid_step 3 diff --git a/tests/08_EXX/09_KP_HSE_comp/INPUT b/tests/08_EXX/09_KP_HSE_comp/INPUT index e33adad4e4..c1bcc5c5cc 100644 --- a/tests/08_EXX/09_KP_HSE_comp/INPUT +++ b/tests/08_EXX/09_KP_HSE_comp/INPUT @@ -30,7 +30,6 @@ exx_pca_threshold 1E-3 exx_c_threshold 1E-3 exx_v_threshold 1 exx_dm_threshold 1E-3 -exx_cauchy_threshold 1E-4 exx_ccp_rmesh_times 1 exx_separate_loop 1 exx_hybrid_step 2 diff --git a/tests/08_EXX/10_KP_HF/INPUT b/tests/08_EXX/10_KP_HF/INPUT index 35cf8e13c4..9e7990ac1c 100644 --- a/tests/08_EXX/10_KP_HF/INPUT +++ b/tests/08_EXX/10_KP_HF/INPUT @@ -31,7 +31,6 @@ exx_pca_threshold 1E-3 exx_c_threshold 1E-3 exx_v_threshold 1 exx_dm_threshold 1E-3 -exx_cauchy_threshold 1E-4 exx_ccp_rmesh_times 2 exx_separate_loop 1 exx_hybrid_step 2 diff --git a/tests/08_EXX/11_KP_PBE0/INPUT b/tests/08_EXX/11_KP_PBE0/INPUT index c3495e2d41..fdf3f3c4db 100644 --- a/tests/08_EXX/11_KP_PBE0/INPUT +++ b/tests/08_EXX/11_KP_PBE0/INPUT @@ -31,7 +31,6 @@ exx_pca_threshold 1E-3 exx_c_threshold 1E-3 exx_v_threshold 1 exx_dm_threshold 1E-3 -exx_cauchy_threshold 1E-4 exx_ccp_rmesh_times 2 exx_separate_loop 1 exx_hybrid_step 2 diff --git a/tests/08_EXX/12_KP_OXC/INPUT b/tests/08_EXX/12_KP_OXC/INPUT index bfdada3d0d..1b3f014ab2 100644 --- a/tests/08_EXX/12_KP_OXC/INPUT +++ b/tests/08_EXX/12_KP_OXC/INPUT @@ -35,7 +35,6 @@ exx_pca_threshold 1E-2 exx_c_threshold 1E-2 exx_v_threshold 1 exx_dm_threshold 1E-2 -exx_cauchy_threshold 1E-3 exx_ccp_rmesh_times 1 exx_separate_loop 1 exx_hybrid_step 1