From 293f8179f027b5e1420c5dbeb8a20b14c8526661 Mon Sep 17 00:00:00 2001 From: HGAQ <893018484@qq.com> Date: Fri, 23 Oct 2026 12:50:53 -0400 Subject: [PATCH 1/4] feat(ri): add configurable smooth Coulomb cutoff for RPA exports Add a log-erfc radial cutoff and an independent BvK inscribed-radius option. Preserve hard/legacy defaults and apply the smooth Fourier factor through the post-SCF cut-Coulomb path, including RI fits using that metric. Register and validate the three INPUT parameters, log the selected radius, and document consistent downstream RI/Vcut/W use. Regenerate parameters.yaml and input-main.md from the production executable; add a usage guide and focused kernel/INPUT tests. Validation: - Fresh C++11 build of the kernel GoogleTests; CTest 1/1 passed (2 test cases). - Production executable --version and -h rpa_cut_coulomb passed. - Fresh --check-input: hard/legacy, smooth/legacy and smooth/inscribed accepted; invalid mode, radius and width rejected (6/6). - All 11 changed source/build/test files byte-match the compiled production source copy. Full application build and Si W/BSE runs were performed in that production copy, which also has pre-existing local interface changes; a full application build of this isolated checkout was not repeated. - git diff --cached --check passed; staged governance has no errors. Its sole header warning is required: std::array appears in the public shortest-vector function signature, so must be included. - New kernel/header/test files pass the repository quality scorer (3/3). Based on production baseline c2fa3ac19539bd10817d66e99ac5ff1cfd2156dd. Unrelated local BSE/interface changes and numerical output are excluded. --- docs/advanced/input_files/input-main.md | 21 +++ docs/advanced/interface/index.rst | 1 + .../interface/librpa-coulomb-cutoff.md | 94 +++++++++++++ docs/parameters.yaml | 24 ++++ source/source_hamilt/module_xc/exx_info.cpp | 19 +++ .../module_parameter/input_parameter.h | 3 + .../module_parameter/read_inp_exx_dftu.cpp | 52 +++++++ .../test_serial/read_input_item_test.cpp | 24 ++++ source/source_lcao/module_ri/CMakeLists.txt | 1 + .../module_ri/conv_coulomb_pot_k.cpp | 26 +++- .../source_lcao/module_ri/exx_lri_detail.cpp | 24 ++++ .../module_ri/smooth_cut_coulomb.cpp | 129 ++++++++++++++++++ .../module_ri/smooth_cut_coulomb.h | 15 ++ .../source_lcao/module_ri/test/CMakeLists.txt | 4 + .../test/smooth_cut_coulomb_test.cpp | 34 +++++ 15 files changed, 466 insertions(+), 5 deletions(-) create mode 100644 docs/advanced/interface/librpa-coulomb-cutoff.md create mode 100644 source/source_lcao/module_ri/smooth_cut_coulomb.cpp create mode 100644 source/source_lcao/module_ri/smooth_cut_coulomb.h create mode 100644 source/source_lcao/module_ri/test/smooth_cut_coulomb_test.cpp diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md index 8952ec2d060..3f74a8c86fa 100644 --- a/docs/advanced/input_files/input-main.md +++ b/docs/advanced/input_files/input-main.md @@ -338,6 +338,9 @@ - [exx\_real\_number](#exx_real_number) - [exx\_singularity\_correction](#exx_singularity_correction) - [rpa\_ccp\_rmesh\_times](#rpa_ccp_rmesh_times) + - [rpa\_cut\_coulomb](#rpa_cut_coulomb) + - [rpa\_cut\_radius](#rpa_cut_radius) + - [rpa\_cut\_width](#rpa_cut_width) - [exx\_symmetry\_realspace](#exx_symmetry_realspace) - [out\_ri\_cv](#out_ri_cv) - [Exact Exchange (PW)](#exact-exchange-pw) @@ -3378,6 +3381,24 @@ - **Description**: How many times larger the radial mesh required is to that of atomic orbitals in the postprocess calculation of the bare Coulomb matrix for RPA, GW, etc. - **Default**: 10 +### rpa_cut_coulomb + +- **Type**: String +- **Description**: Post-SCF cut-Coulomb kernel: hard preserves Spencer truncation; log_erfc uses erfc(log(r/Rc)/log(rpa_cut_width))/2r. Applies to compressed and uncompressed cut exports; full Coulomb is unchanged. Forces and stresses for log_erfc are not validated. +- **Default**: hard + +### rpa_cut_radius + +- **Type**: String +- **Description**: legacy preserves the existing cutoff radius. inscribed uses half the shortest nonzero BvK superlattice vector, including the lattice scale in Bohr. Independent of rpa_cut_coulomb; intended for controlled radius-versus-smoothing comparisons. +- **Default**: legacy + +### rpa_cut_width + +- **Type**: Real +- **Description**: Width factor for log_erfc. Must be finite and greater than 1 and at most 1.25. The default matches the logged Si FHI-aims width; it is not a universal FHI-aims constant. The finite radial integration range must resolve the smooth tail. +- **Default**: 1.091868692359 + ### exx_symmetry_realspace - **Type**: Boolean diff --git a/docs/advanced/interface/index.rst b/docs/advanced/interface/index.rst index 6622e821b7e..7f0d91b1596 100644 --- a/docs/advanced/interface/index.rst +++ b/docs/advanced/interface/index.rst @@ -19,3 +19,4 @@ Interfaces to Other Softwares ShengBTE candela TB2J + librpa-coulomb-cutoff diff --git a/docs/advanced/interface/librpa-coulomb-cutoff.md b/docs/advanced/interface/librpa-coulomb-cutoff.md new file mode 100644 index 00000000000..3e2f6e77d59 --- /dev/null +++ b/docs/advanced/interface/librpa-coulomb-cutoff.md @@ -0,0 +1,94 @@ +# Smooth Coulomb cutoff for LibRPA exports + +ABACUS can export a smooth truncated Coulomb interaction for post-SCF RPA/GW +workflows. The default remains the existing hard spherical (Spencer) cutoff +with the legacy radius convention. + +## Usage + +In an existing LCAO RPA export calculation, use: + +```text +rpa 1 +rpa_out_vel 1 +rpa_outdir ./OUT.librpa/ +out_unshrinked_v 1 + +rpa_cut_coulomb log_erfc +rpa_cut_radius legacy +rpa_cut_width 1.091868692359 +rpa_ccp_rmesh_times 10 +``` + +Keep the structure, pseudopotentials, orbital/auxiliary bases and k mesh from +the intended calculation. The width and radial-range setting above were used +for the silicon comparison; they are not universal convergence settings. +The width must be finite and satisfy `1 < rpa_cut_width <= 1.25`. + +Check standard output for `RPA cut Coulomb: mode=log_erfc radius=legacy` and +the actual `Rcut_bohr` and `width`. Then regenerate the downstream LibRPA +response and screened interaction. Consumers must use matching RI coefficients, +cut Coulomb matrices and screened interaction from the same cutoff workflow. +Changing only a downstream BSE input does not convert previously exported data. + +The existing LibRPA-interface options in the example request velocity and +uncompressed Coulomb output for downstream optical and matrix diagnostics. +They do not enable smoothing by themselves. + +To recover the hard-cutoff control, set `rpa_cut_coulomb hard` and +`rpa_cut_radius legacy`. Use separate output directories for each variant. + +## Kernel and radius + +The smooth radial kernel is + +$$ +v_{\mathrm{cut}}(r)=\frac{1}{2r} +\operatorname{erfc}\!\left(\frac{\ln(r/R_c)}{\ln w}\right). +$$ + +The implementation evaluates its Fourier factor as a positive Gaussian +average of hard-sphere factors. With $\sigma=\ln w$, + +$$ +\frac{k^2 v_{\mathrm{cut}}(k)}{4\pi} +=\frac{1}{\sqrt{\pi}}\int_{-\infty}^{\infty}e^{-t^2} +\left[1-\cos\!\left(kR_c e^{\sigma t}\right)\right]dt. +$$ + +This tends to the hard-sphere factor as $w\to1^+$. Its finite zero-wavevector +kernel limit is $2\pi R_c^2\exp(\sigma^2)$. The stored multiplicative factor +is zero at exactly $k=0$ because it represents $k^2v(k)/(4\pi)$. + +`legacy` retains the existing radius rule. The optional `inscribed` convention +uses half the shortest nonzero Born-von Karman superlattice vector: + +$$ +R_c=\frac12\min_{\mathbf n\in\mathbb Z^3\setminus\{0\}} +\left|\sum_i n_i N_i\mathbf a_i\right|. +$$ + +Lattice scale is included and lengths are in Bohr. The shortest-vector search +includes integer combinations, not just the three supplied lattice vectors. +Radius selection and kernel smoothing are independent controls. The smooth +kernel has a tail beyond $R_c$; the inscribed radius is not a strict support +boundary. `rpa_ccp_rmesh_times` controls a numerical radial range, not this +physical radius. Check radial-range convergence for the smooth tail. + +## Scope + +The option is passed explicitly into the post-SCF truncated-Coulomb path, +including compressed and uncompressed cut exports. The full-Coulomb kernel +definition is unchanged. RI fits using the cut metric also change, so this +is not generally a fixed-response replacement of the outer Coulomb factors +in the screened interaction. Downstream Fourier folding and interpolation +rules are unchanged. + +The default width was taken from a silicon FHI-aims log; this implementation +does not reproduce the complete FHI-aims short-/long-range split. Forces and +stresses for the smooth kernel have not been validated. + +Focused kernel tests are registered as `MODULE_RI_smooth_cut_coulomb` in an +MPI-enabled LCAO/LibRI test build. They cover the small-wavevector moment, +hard-cutoff limit, positivity, invalid arguments and cubic/fcc/skew-cell +radius selection. INPUT validation tests also cover the new switches. diff --git a/docs/parameters.yaml b/docs/parameters.yaml index 4c097ada35e..c339436657c 100644 --- a/docs/parameters.yaml +++ b/docs/parameters.yaml @@ -4641,6 +4641,30 @@ parameters: default_value: "10" unit: "" availability: "" + - name: rpa_cut_coulomb + category: Exact Exchange (LCAO) + type: String + description: | + Post-SCF cut-Coulomb kernel: hard preserves Spencer truncation; log_erfc uses erfc(log(r/Rc)/log(rpa_cut_width))/2r. Applies to compressed and uncompressed cut exports; full Coulomb is unchanged. Forces and stresses for log_erfc are not validated. + default_value: hard + unit: "" + availability: "" + - name: rpa_cut_radius + category: Exact Exchange (LCAO) + type: String + description: | + legacy preserves the existing cutoff radius. inscribed uses half the shortest nonzero BvK superlattice vector, including the lattice scale in Bohr. Independent of rpa_cut_coulomb; intended for controlled radius-versus-smoothing comparisons. + default_value: legacy + unit: "" + availability: "" + - name: rpa_cut_width + category: Exact Exchange (LCAO) + type: Real + description: | + Width factor for log_erfc. Must be finite and greater than 1 and at most 1.25. The default matches the logged Si FHI-aims width; it is not a universal FHI-aims constant. The finite radial integration range must resolve the smooth tail. + default_value: "1.091868692359" + unit: "" + availability: "" - name: exx_symmetry_realspace category: Exact Exchange (LCAO) type: Boolean diff --git a/source/source_hamilt/module_xc/exx_info.cpp b/source/source_hamilt/module_xc/exx_info.cpp index 923fef63d44..b6e0d142e71 100644 --- a/source/source_hamilt/module_xc/exx_info.cpp +++ b/source/source_hamilt/module_xc/exx_info.cpp @@ -7,6 +7,7 @@ #include #include #include +#include #include //---------------------------------------------------------- @@ -228,6 +229,24 @@ bool init_exx_info(Exx_Info& exx_info, const Input_para& inp) // nspin=4/SOC via restore_dm + restore_HR_nspin4), so symmetry=1 is honored here. exx_info.sync_from_global(); + // Explicit INPUT -> RI configuration; only init_spencer consumes these + // tags, so the full-Coulomb/ground-state exchange kernel is unchanged. + if (inp.rpa) + { + for (auto& cut_family : exx_info.info_ri.coulomb_param) + { + if (cut_family.first != Conv_Coulomb_Pot_K::Coulomb_Type::Fock) { continue; } + for (auto& cut_parameter : cut_family.second) + { + cut_parameter["rpa_cut_coulomb"] = inp.rpa_cut_coulomb; + cut_parameter["rpa_cut_radius"] = inp.rpa_cut_radius; + std::ostringstream cut_width_text; + cut_width_text.precision(17); + cut_width_text << inp.rpa_cut_width; + cut_parameter["rpa_cut_width"] = cut_width_text.str(); + } + } + } } return generate_opt_orb; diff --git a/source/source_io/module_parameter/input_parameter.h b/source/source_io/module_parameter/input_parameter.h index 01f28ff0040..197009af8d7 100644 --- a/source/source_io/module_parameter/input_parameter.h +++ b/source/source_io/module_parameter/input_parameter.h @@ -621,6 +621,9 @@ struct Input_para bool exx_symmetry_realspace = true; ///< whether to reduce the real-space sector in when using symmetry=1 in EXX calculation double rpa_ccp_rmesh_times = 10.0; ///< how many times larger the radial mesh required for ///< calculating Columb potential is to that of atomic orbitals + std::string rpa_cut_coulomb = "hard"; ///< post-SCF cut kernel: hard or log_erfc + std::string rpa_cut_radius = "legacy"; ///< legacy radius or BvK inscribed radius + double rpa_cut_width = 1.091868692359; ///< multiplicative log-erfc width double exx_cs_inv_thr = -1; ///< threshold to inverse Vq in abfs for generating Cs double shrink_abfs_pca_thr = -1; ///< threshold to shrink auxiliary basis for GW/RPA double shrink_LU_inv_thr = 1e-6; ///< threshold to get inverse of overlap matrix by LU decomposition in auxiliary basis representation diff --git a/source/source_io/module_parameter/read_inp_exx_dftu.cpp b/source/source_io/module_parameter/read_inp_exx_dftu.cpp index 13515d403ee..bc5a54874e4 100644 --- a/source/source_io/module_parameter/read_inp_exx_dftu.cpp +++ b/source/source_io/module_parameter/read_inp_exx_dftu.cpp @@ -2,6 +2,7 @@ #include "source_base/tool_quit.h" #include "read_input.h" #include "read_input_tool.h" +#include namespace ModuleIO { void ReadInput::item_exx() @@ -519,6 +520,57 @@ void ReadInput::item_exx() }; this->add_item(item); } + { + Input_Item item("rpa_cut_coulomb"); + item.annotation = "post-SCF truncated Coulomb kernel for RPA/GW export"; + item.category = "Exact Exchange (LCAO)"; + item.type = "String"; + item.description = "Post-SCF cut-Coulomb kernel: hard preserves Spencer truncation; log_erfc uses erfc(log(r/Rc)/log(rpa_cut_width))/2r. Applies to compressed and uncompressed cut exports; full Coulomb is unchanged. Forces and stresses for log_erfc are not validated."; + item.default_value = "hard"; + item.unit = ""; + read_sync_string(input.rpa_cut_coulomb); + item.check_value = [](const Input_Item& item, const Parameter& para) { + if (para.input.rpa_cut_coulomb != "hard" && para.input.rpa_cut_coulomb != "log_erfc") + { + ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_coulomb must be hard or log_erfc"); + } + }; + this->add_item(item); + } + { + Input_Item item("rpa_cut_radius"); + item.annotation = "radius convention for post-SCF cut Coulomb"; + item.category = "Exact Exchange (LCAO)"; + item.type = "String"; + item.description = "legacy preserves the existing cutoff radius. inscribed uses half the shortest nonzero BvK superlattice vector, including the lattice scale in Bohr. Independent of rpa_cut_coulomb; intended for controlled radius-versus-smoothing comparisons."; + item.default_value = "legacy"; + item.unit = ""; + read_sync_string(input.rpa_cut_radius); + item.check_value = [](const Input_Item& item, const Parameter& para) { + if (para.input.rpa_cut_radius != "legacy" && para.input.rpa_cut_radius != "inscribed") + { + ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_radius must be legacy or inscribed"); + } + }; + this->add_item(item); + } + { + Input_Item item("rpa_cut_width"); + item.annotation = "multiplicative width of the log-erfc cutoff"; + item.category = "Exact Exchange (LCAO)"; + item.type = "Real"; + item.description = "Width factor for log_erfc. Must be finite and greater than 1 and at most 1.25. The default matches the logged Si FHI-aims width; it is not a universal FHI-aims constant. The finite radial integration range must resolve the smooth tail."; + item.default_value = "1.091868692359"; + item.unit = ""; + read_sync_double(input.rpa_cut_width); + item.check_value = [](const Input_Item& item, const Parameter& para) { + if (!std::isfinite(para.input.rpa_cut_width) || para.input.rpa_cut_width <= 1.0 || para.input.rpa_cut_width > 1.25) + { + ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_width must be finite and 1 < width <= 1.25"); + } + }; + this->add_item(item); + } { Input_Item item("exx_symmetry_realspace"); item.annotation = "whether to reduce real-space sector in Hexx calculation"; diff --git a/source/source_io/test_serial/read_input_item_test.cpp b/source/source_io/test_serial/read_input_item_test.cpp index c5815f858ae..4552e93b05f 100644 --- a/source/source_io/test_serial/read_input_item_test.cpp +++ b/source/source_io/test_serial/read_input_item_test.cpp @@ -1674,6 +1674,30 @@ TEST_F(InputTest, Item_test2) output = testing::internal::GetCapturedStdout(); EXPECT_THAT(output, testing::HasSubstr("NOTICE")); } + { // independent smoothing and radius switches + const auto cut_mode_item = find_label("rpa_cut_coulomb", readinput.input_lists); + const auto cut_radius_item = find_label("rpa_cut_radius", readinput.input_lists); + const auto cut_width_item = find_label("rpa_cut_width", readinput.input_lists); + ASSERT_NE(cut_mode_item, readinput.input_lists.end()); + ASSERT_NE(cut_radius_item, readinput.input_lists.end()); + ASSERT_NE(cut_width_item, readinput.input_lists.end()); + EXPECT_EQ(param.input.rpa_cut_coulomb, "hard"); + EXPECT_EQ(param.input.rpa_cut_radius, "legacy"); + param.input.rpa_cut_coulomb = "log_erfc"; + param.input.rpa_cut_radius = "inscribed"; + cut_mode_item->second.check_value(cut_mode_item->second, param); + cut_radius_item->second.check_value(cut_radius_item->second, param); + cut_width_item->second.check_value(cut_width_item->second, param); + param.input.rpa_cut_coulomb = "invalid"; + EXPECT_EXIT(cut_mode_item->second.check_value(cut_mode_item->second, param), ::testing::ExitedWithCode(1), ""); + param.input.rpa_cut_radius = "invalid"; + EXPECT_EXIT(cut_radius_item->second.check_value(cut_radius_item->second, param), ::testing::ExitedWithCode(1), ""); + param.input.rpa_cut_width = 1.0; + EXPECT_EXIT(cut_width_item->second.check_value(cut_width_item->second, param), ::testing::ExitedWithCode(1), ""); + param.input.rpa_cut_coulomb = "hard"; + param.input.rpa_cut_radius = "legacy"; + param.input.rpa_cut_width = 1.091868692359; + } { // berry_phase auto it = find_label("berry_phase", readinput.input_lists); param.input.berry_phase = true; diff --git a/source/source_lcao/module_ri/CMakeLists.txt b/source/source_lcao/module_ri/CMakeLists.txt index c45ca9751a4..854e792d567 100644 --- a/source/source_lcao/module_ri/CMakeLists.txt +++ b/source/source_lcao/module_ri/CMakeLists.txt @@ -13,6 +13,7 @@ if (ENABLE_LIBRI) if(ENABLE_LCAO) list(APPEND objects conv_coulomb_pot_k.cpp + smooth_cut_coulomb.cpp exx_abfs_ctor_orbs.cpp exx_abfs_io.cpp exx_abfs_jle.cpp diff --git a/source/source_lcao/module_ri/conv_coulomb_pot_k.cpp b/source/source_lcao/module_ri/conv_coulomb_pot_k.cpp index 3970b98f41d..e59d5bc896e 100644 --- a/source/source_lcao/module_ri/conv_coulomb_pot_k.cpp +++ b/source/source_lcao/module_ri/conv_coulomb_pot_k.cpp @@ -1,4 +1,5 @@ #include "conv_coulomb_pot_k.h" +#include "smooth_cut_coulomb.h" #include "source_base/constants.h" #include "source_io/module_parameter/parameter.h" #include "source_basis/module_ao/orb_atomic_lm.h" @@ -162,11 +163,26 @@ namespace Conv_Coulomb_Pot_K case Conv_Coulomb_Pot_K::Coulomb_Type::Fock: { for (const auto& param: param_list.second) { - psik2_ccp = psik2_ccp - + std::stod(param.at("alpha")) - * cal_psi_fock_spencer(orbs.get_psif(), - orbs.get_k_radial(), - std::stod(param.at("Rcut"))); + const auto cut_mode = param.find("rpa_cut_coulomb"); + if (cut_mode != param.end() && cut_mode->second == "log_erfc") + { + const double cut_radius = std::stod(param.at("Rcut")); + const double cut_width = std::stod(param.at("rpa_cut_width")); + const double cut_alpha = std::stod(param.at("alpha")); + for (size_t cut_ik = 0; cut_ik < psik2_ccp.size(); ++cut_ik) + { + psik2_ccp[cut_ik] += cut_alpha * ModuleBase::FOUR_PI * orbs.get_psif()[cut_ik] + * log_erfc_cut_factor(orbs.get_k_radial()[cut_ik], cut_radius, cut_width); + } + } + else + { + psik2_ccp = psik2_ccp + + std::stod(param.at("alpha")) + * cal_psi_fock_spencer(orbs.get_psif(), + orbs.get_k_radial(), + std::stod(param.at("Rcut"))); + } } break; } diff --git a/source/source_lcao/module_ri/exx_lri_detail.cpp b/source/source_lcao/module_ri/exx_lri_detail.cpp index ff8046159b6..8df99e35500 100644 --- a/source/source_lcao/module_ri/exx_lri_detail.cpp +++ b/source/source_lcao/module_ri/exx_lri_detail.cpp @@ -5,6 +5,9 @@ #include "exx_lri_detail.h" +#include "smooth_cut_coulomb.h" +#include +#include #include "ri_util.h" #include "source_cell/klist.h" @@ -36,12 +39,33 @@ CoulombParam build_center2_cut_coulomb_param(const CoulombParam& coulomb_param, } for (auto& param: param_list.second) { + const auto cut_radius_policy = param.find("rpa_cut_radius"); + if (cut_radius_policy != param.end() && cut_radius_policy->second == "inscribed") + { + const std::array, 3> cut_bvk_lattice = {{ + {{ucell.lat0 * kv.nmp[0] * ucell.a1.x, ucell.lat0 * kv.nmp[0] * ucell.a1.y, ucell.lat0 * kv.nmp[0] * ucell.a1.z}}, + {{ucell.lat0 * kv.nmp[1] * ucell.a2.x, ucell.lat0 * kv.nmp[1] * ucell.a2.y, ucell.lat0 * kv.nmp[1] * ucell.a2.z}}, + {{ucell.lat0 * kv.nmp[2] * ucell.a3.x, ucell.lat0 * kv.nmp[2] * ucell.a3.y, ucell.lat0 * kv.nmp[2] * ucell.a3.z}} + }}; + std::ostringstream cut_radius_text; + cut_radius_text.precision(17); + cut_radius_text << Conv_Coulomb_Pot_K::bvk_inscribed_radius(cut_bvk_lattice); + param["Rcut"] = cut_radius_text.str(); + } auto rcut_it = param.find("Rcut"); if (rcut_it == param.end() || rcut_it->second.empty()) { param["Rcut"] = ModuleBase::GlobalFunc::TO_STRING(fallback_rcut); used_fallback_rcut = true; } + const auto cut_mode = param.find("rpa_cut_coulomb"); + if (cut_mode != param.end()) + { + std::cout << "RPA cut Coulomb: mode=" << cut_mode->second + << " radius=" << param.at("rpa_cut_radius") + << " Rcut_bohr=" << param.at("Rcut") + << " width=" << param.at("rpa_cut_width") << std::endl; + } } } diff --git a/source/source_lcao/module_ri/smooth_cut_coulomb.cpp b/source/source_lcao/module_ri/smooth_cut_coulomb.cpp new file mode 100644 index 00000000000..ea414aecb4c --- /dev/null +++ b/source/source_lcao/module_ri/smooth_cut_coulomb.cpp @@ -0,0 +1,129 @@ +#include "smooth_cut_coulomb.h" + +#include +#include +#include +#include + +namespace Conv_Coulomb_Pot_K +{ +double log_erfc_cut_factor(const double wave_number, + const double cutoff_radius, + const double width_factor) +{ + if (!std::isfinite(wave_number) || wave_number < 0.0 + || !std::isfinite(cutoff_radius) || cutoff_radius <= 0.0 + || !std::isfinite(width_factor) || width_factor <= 1.0 || width_factor > 1.25) + { + throw std::invalid_argument("log-erfc cutoff requires k>=0, Rc>0 and 1 1000000.0) + { + throw std::runtime_error("log-erfc cutoff quadrature exceeds one million panels"); + } + const int cut_panel_count = std::max(128, static_cast(cut_panel_estimate)); + const double cut_step = 16.0 / cut_panel_count; + double cut_integral = 0.0; + for (int cut_panel = 0; cut_panel < cut_panel_count; ++cut_panel) + { + const double cut_coordinate = -8.0 + (cut_panel + 0.5) * cut_step; + const double cut_phase_half = 0.5 * cut_argument * std::exp(cut_log_width * cut_coordinate); + const double cut_sine = std::sin(cut_phase_half); + cut_integral += std::exp(-cut_coordinate * cut_coordinate) * 2.0 * cut_sine * cut_sine; + } + return cut_integral * cut_step / std::sqrt(cut_pi); +} + +double bvk_inscribed_radius(const std::array, 3>& bvk_lattice) +{ + std::array, 3> bvk_cofactor{}; + double bvk_shortest_squared = std::numeric_limits::max(); + for (int bvk_axis = 0; bvk_axis < 3; ++bvk_axis) + { + double bvk_length_squared = 0.0; + for (int cart_axis = 0; cart_axis < 3; ++cart_axis) + { + if (!std::isfinite(bvk_lattice[bvk_axis][cart_axis])) + { + throw std::invalid_argument("non-finite BvK lattice"); + } + bvk_length_squared += bvk_lattice[bvk_axis][cart_axis] * bvk_lattice[bvk_axis][cart_axis]; + const int next_row = (bvk_axis + 1) % 3; + const int last_row = (bvk_axis + 2) % 3; + const int next_col = (cart_axis + 1) % 3; + const int last_col = (cart_axis + 2) % 3; + bvk_cofactor[bvk_axis][cart_axis] + = bvk_lattice[next_row][next_col] * bvk_lattice[last_row][last_col] + - bvk_lattice[next_row][last_col] * bvk_lattice[last_row][next_col]; + } + bvk_shortest_squared = std::min(bvk_shortest_squared, bvk_length_squared); + } + double bvk_determinant = 0.0; + for (int cart_axis = 0; cart_axis < 3; ++cart_axis) + { + bvk_determinant += bvk_lattice[0][cart_axis] * bvk_cofactor[0][cart_axis]; + } + if (bvk_determinant == 0.0 || !std::isfinite(bvk_determinant) || bvk_shortest_squared <= 0.0) + { + throw std::invalid_argument("singular BvK lattice"); + } + // For d=n*B, |n_i| <= |d|*||(B^-1)[:,i]||. These bounds include + // every vector shorter than the shortest supplied basis vector. + std::array bvk_bounds{}; + double bvk_candidate_count = 1.0; + for (int bvk_axis = 0; bvk_axis < 3; ++bvk_axis) + { + double bvk_inverse_norm_squared = 0.0; + for (int cart_axis = 0; cart_axis < 3; ++cart_axis) + { + const double bvk_inverse_entry = bvk_cofactor[bvk_axis][cart_axis] / bvk_determinant; + bvk_inverse_norm_squared += bvk_inverse_entry * bvk_inverse_entry; + } + const double bvk_bound = std::floor(std::sqrt(bvk_shortest_squared * bvk_inverse_norm_squared) + 1.e-10); + if (!std::isfinite(bvk_bound) || bvk_bound > 1000000.0) + { + throw std::runtime_error("BvK shortest-vector search needs a reduced lattice"); + } + bvk_bounds[bvk_axis] = static_cast(bvk_bound); + bvk_candidate_count *= 2.0 * bvk_bounds[bvk_axis] + 1.0; + } + if (bvk_candidate_count > 10000000.0) + { + throw std::runtime_error("BvK shortest-vector search exceeds ten million candidates"); + } + for (int first_cell = -bvk_bounds[0]; first_cell <= bvk_bounds[0]; ++first_cell) + { + for (int second_cell = -bvk_bounds[1]; second_cell <= bvk_bounds[1]; ++second_cell) + { + for (int third_cell = -bvk_bounds[2]; third_cell <= bvk_bounds[2]; ++third_cell) + { + if (first_cell == 0 && second_cell == 0 && third_cell == 0) { continue; } + double bvk_distance_squared = 0.0; + for (int cart_axis = 0; cart_axis < 3; ++cart_axis) + { + const double bvk_displacement = first_cell * bvk_lattice[0][cart_axis] + + second_cell * bvk_lattice[1][cart_axis] + third_cell * bvk_lattice[2][cart_axis]; + bvk_distance_squared += bvk_displacement * bvk_displacement; + } + bvk_shortest_squared = std::min(bvk_shortest_squared, bvk_distance_squared); + } + } + } + return 0.5 * std::sqrt(bvk_shortest_squared); +} +} diff --git a/source/source_lcao/module_ri/smooth_cut_coulomb.h b/source/source_lcao/module_ri/smooth_cut_coulomb.h new file mode 100644 index 00000000000..4b070872b23 --- /dev/null +++ b/source/source_lcao/module_ri/smooth_cut_coulomb.h @@ -0,0 +1,15 @@ +#ifndef ABACUS_SMOOTH_CUT_COULOMB_H +#define ABACUS_SMOOTH_CUT_COULOMB_H + +#include + +namespace Conv_Coulomb_Pot_K +{ +// k^2 v(k)/(4*pi) for f(r)/r, f=erfc(log(r/Rc)/log(width))/2. +double log_erfc_cut_factor(double wave_number, double cutoff_radius, double width_factor); + +// Half the shortest nonzero superlattice vector; row vectors are in Bohr. +double bvk_inscribed_radius(const std::array, 3>& bvk_lattice); +} + +#endif diff --git a/source/source_lcao/module_ri/test/CMakeLists.txt b/source/source_lcao/module_ri/test/CMakeLists.txt index ff75006e5b6..ebed5c8b185 100644 --- a/source/source_lcao/module_ri/test/CMakeLists.txt +++ b/source/source_lcao/module_ri/test/CMakeLists.txt @@ -1,6 +1,10 @@ abacus_disable_feature_definitions(__MLALGO) abacus_disable_feature_definitions(__CUDA) abacus_disable_feature_definitions(__ROCM) +AddTest( + TARGET MODULE_RI_smooth_cut_coulomb + SOURCES smooth_cut_coulomb_test.cpp ../smooth_cut_coulomb.cpp +) AddTest( TARGET MODULE_RI_dm_mixing_test LIBS parameter base device diff --git a/source/source_lcao/module_ri/test/smooth_cut_coulomb_test.cpp b/source/source_lcao/module_ri/test/smooth_cut_coulomb_test.cpp new file mode 100644 index 00000000000..f14a803d747 --- /dev/null +++ b/source/source_lcao/module_ri/test/smooth_cut_coulomb_test.cpp @@ -0,0 +1,34 @@ +#include "../smooth_cut_coulomb.h" +#include +#include +#include + +TEST(SmoothCutCoulomb, GammaMomentAndSharpLimit) +{ + const double cut_width = 1.091868692359; + const double cut_sigma = std::log(cut_width); + const double cut_z = 1.e-7; + EXPECT_NEAR(Conv_Coulomb_Pot_K::log_erfc_cut_factor(cut_z, 1.0, cut_width) / (cut_z * cut_z), + 0.5 * std::exp(cut_sigma * cut_sigma), 1.e-13); + EXPECT_DOUBLE_EQ(Conv_Coulomb_Pot_K::log_erfc_cut_factor(0.0, 32.0, cut_width), 0.0); + EXPECT_NEAR(Conv_Coulomb_Pot_K::log_erfc_cut_factor(0.2, 32.0, 1.00000001), + 1.0 - std::cos(6.4), 1.e-13); + for (int cut_index = 1; cut_index <= 100; ++cut_index) + { + const double cut_factor = Conv_Coulomb_Pot_K::log_erfc_cut_factor(cut_index, 20.0, cut_width); + EXPECT_GE(cut_factor, 0.0); + EXPECT_LE(cut_factor, 2.0 + 1.e-14); + } +} + +TEST(SmoothCutCoulomb, InscribedRadiusUsesShortestCombination) +{ + const std::array, 3> cut_cubic = {{{{8., 0., 0.}}, {{0., 10., 0.}}, {{0., 0., 12.}}}}; + const std::array, 3> cut_fcc = {{{{0., 4., 4.}}, {{4., 0., 4.}}, {{4., 4., 0.}}}}; + const std::array, 3> cut_skew = {{{{10., 0., 0.}}, {{9.9, 0.1, 0.}}, {{0., 0., 5.}}}}; + EXPECT_DOUBLE_EQ(Conv_Coulomb_Pot_K::bvk_inscribed_radius(cut_cubic), 4.0); + EXPECT_NEAR(Conv_Coulomb_Pot_K::bvk_inscribed_radius(cut_fcc), std::sqrt(8.0), 1.e-14); + EXPECT_NEAR(Conv_Coulomb_Pot_K::bvk_inscribed_radius(cut_skew), std::sqrt(0.02) / 2., 1.e-14); + EXPECT_THROW(Conv_Coulomb_Pot_K::log_erfc_cut_factor(1., 1., 1.), std::invalid_argument); + EXPECT_THROW(Conv_Coulomb_Pot_K::log_erfc_cut_factor(-1., 1., 1.1), std::invalid_argument); +} From 0133ed49093c0b044ac6c7bc8cb153a10e962e35 Mon Sep 17 00:00:00 2001 From: HGAQ <893018484@qq.com> Date: Fri, 23 Oct 2026 13:01:16 -0400 Subject: [PATCH 2/4] docs: remove smooth cutoff tutorial from published interface docs --- docs/advanced/interface/index.rst | 1 - .../interface/librpa-coulomb-cutoff.md | 94 ------------------- 2 files changed, 95 deletions(-) delete mode 100644 docs/advanced/interface/librpa-coulomb-cutoff.md diff --git a/docs/advanced/interface/index.rst b/docs/advanced/interface/index.rst index 7f0d91b1596..6622e821b7e 100644 --- a/docs/advanced/interface/index.rst +++ b/docs/advanced/interface/index.rst @@ -19,4 +19,3 @@ Interfaces to Other Softwares ShengBTE candela TB2J - librpa-coulomb-cutoff diff --git a/docs/advanced/interface/librpa-coulomb-cutoff.md b/docs/advanced/interface/librpa-coulomb-cutoff.md deleted file mode 100644 index 3e2f6e77d59..00000000000 --- a/docs/advanced/interface/librpa-coulomb-cutoff.md +++ /dev/null @@ -1,94 +0,0 @@ -# Smooth Coulomb cutoff for LibRPA exports - -ABACUS can export a smooth truncated Coulomb interaction for post-SCF RPA/GW -workflows. The default remains the existing hard spherical (Spencer) cutoff -with the legacy radius convention. - -## Usage - -In an existing LCAO RPA export calculation, use: - -```text -rpa 1 -rpa_out_vel 1 -rpa_outdir ./OUT.librpa/ -out_unshrinked_v 1 - -rpa_cut_coulomb log_erfc -rpa_cut_radius legacy -rpa_cut_width 1.091868692359 -rpa_ccp_rmesh_times 10 -``` - -Keep the structure, pseudopotentials, orbital/auxiliary bases and k mesh from -the intended calculation. The width and radial-range setting above were used -for the silicon comparison; they are not universal convergence settings. -The width must be finite and satisfy `1 < rpa_cut_width <= 1.25`. - -Check standard output for `RPA cut Coulomb: mode=log_erfc radius=legacy` and -the actual `Rcut_bohr` and `width`. Then regenerate the downstream LibRPA -response and screened interaction. Consumers must use matching RI coefficients, -cut Coulomb matrices and screened interaction from the same cutoff workflow. -Changing only a downstream BSE input does not convert previously exported data. - -The existing LibRPA-interface options in the example request velocity and -uncompressed Coulomb output for downstream optical and matrix diagnostics. -They do not enable smoothing by themselves. - -To recover the hard-cutoff control, set `rpa_cut_coulomb hard` and -`rpa_cut_radius legacy`. Use separate output directories for each variant. - -## Kernel and radius - -The smooth radial kernel is - -$$ -v_{\mathrm{cut}}(r)=\frac{1}{2r} -\operatorname{erfc}\!\left(\frac{\ln(r/R_c)}{\ln w}\right). -$$ - -The implementation evaluates its Fourier factor as a positive Gaussian -average of hard-sphere factors. With $\sigma=\ln w$, - -$$ -\frac{k^2 v_{\mathrm{cut}}(k)}{4\pi} -=\frac{1}{\sqrt{\pi}}\int_{-\infty}^{\infty}e^{-t^2} -\left[1-\cos\!\left(kR_c e^{\sigma t}\right)\right]dt. -$$ - -This tends to the hard-sphere factor as $w\to1^+$. Its finite zero-wavevector -kernel limit is $2\pi R_c^2\exp(\sigma^2)$. The stored multiplicative factor -is zero at exactly $k=0$ because it represents $k^2v(k)/(4\pi)$. - -`legacy` retains the existing radius rule. The optional `inscribed` convention -uses half the shortest nonzero Born-von Karman superlattice vector: - -$$ -R_c=\frac12\min_{\mathbf n\in\mathbb Z^3\setminus\{0\}} -\left|\sum_i n_i N_i\mathbf a_i\right|. -$$ - -Lattice scale is included and lengths are in Bohr. The shortest-vector search -includes integer combinations, not just the three supplied lattice vectors. -Radius selection and kernel smoothing are independent controls. The smooth -kernel has a tail beyond $R_c$; the inscribed radius is not a strict support -boundary. `rpa_ccp_rmesh_times` controls a numerical radial range, not this -physical radius. Check radial-range convergence for the smooth tail. - -## Scope - -The option is passed explicitly into the post-SCF truncated-Coulomb path, -including compressed and uncompressed cut exports. The full-Coulomb kernel -definition is unchanged. RI fits using the cut metric also change, so this -is not generally a fixed-response replacement of the outer Coulomb factors -in the screened interaction. Downstream Fourier folding and interpolation -rules are unchanged. - -The default width was taken from a silicon FHI-aims log; this implementation -does not reproduce the complete FHI-aims short-/long-range split. Forces and -stresses for the smooth kernel have not been validated. - -Focused kernel tests are registered as `MODULE_RI_smooth_cut_coulomb` in an -MPI-enabled LCAO/LibRI test build. They cover the small-wavevector moment, -hard-cutoff limit, positivity, invalid arguments and cubic/fcc/skew-cell -radius selection. INPUT validation tests also cover the new switches. From 73d698b649bb0bd91f6638fd0ccbee6490ae3e02 Mon Sep 17 00:00:00 2001 From: HGAQ <893018484@qq.com> Date: Fri, 23 Oct 2026 16:22:02 -0400 Subject: [PATCH 3/4] fix(ri): validate cutoff inputs before parsing and reject unsupported kernels Declare cutoff mode/radius choices and width bounds through the shared INPUT domain validator. Missing values, extra tokens, invalid choices and non-finite/out-of-range widths now fail before reader callbacks/broadcasts. Reject non-default cutoff requests for RPA configurations without a Fock component, including default HSE. Preserve existing Erfc kernel behavior and document the Fock-only scope in the INPUT metadata. Add parser-level regressions for the original malformed token sequences and configuration coverage for default HSE, rejected HSE cutoff requests, and accepted PBE cutoff settings. Numerical kernels and caching are unchanged. Validation: MODULE_IO_input_domain_reader, MODULE_IO_read_item_serial and MODULE_HAMILT_RpaCutoffConfig all passed (3/3 CTest targets). Staged diff whitespace and governance checks reported no errors. Generated YAML and Markdown were checked locally for enum choices, width bounds and scope; they remain untracked build artifacts under current develop policy. --- source/source_hamilt/module_xc/exx_info.cpp | 9 ++++ .../module_xc/test/CMakeLists.txt | 6 +++ .../module_xc/test/test_exx_info.cpp | 41 ++++++++++++++++++ .../module_parameter/read_inp_exx_dftu.cpp | 26 +++--------- .../source_io/test_serial/test_read_input.cpp | 42 +++++++++++++++++++ 5 files changed, 103 insertions(+), 21 deletions(-) create mode 100644 source/source_hamilt/module_xc/test/test_exx_info.cpp diff --git a/source/source_hamilt/module_xc/exx_info.cpp b/source/source_hamilt/module_xc/exx_info.cpp index 31b5daf7d5c..56da05e41be 100644 --- a/source/source_hamilt/module_xc/exx_info.cpp +++ b/source/source_hamilt/module_xc/exx_info.cpp @@ -9,6 +9,7 @@ #include #include #include +#include #include //---------------------------------------------------------- @@ -249,6 +250,14 @@ bool init_exx_info(Exx_Info& exx_info, const Input_para& inp) // tags, so the full-Coulomb/ground-state exchange kernel is unchanged. if (inp.rpa) { + const bool cut_nondefault = inp.rpa_cut_coulomb != "hard" || inp.rpa_cut_radius != "legacy"; + if (cut_nondefault + && exx_info.info_ri.coulomb_param.count(Conv_Coulomb_Pot_K::Coulomb_Type::Fock) == 0) + { + throw std::invalid_argument( + "Non-default rpa_cut_coulomb/rpa_cut_radius requires a Fock Coulomb component; " + "Erfc-only configurations such as default HSE are unsupported."); + } for (auto& cut_family : exx_info.info_ri.coulomb_param) { if (cut_family.first != Conv_Coulomb_Pot_K::Coulomb_Type::Fock) { continue; } diff --git a/source/source_hamilt/module_xc/test/CMakeLists.txt b/source/source_hamilt/module_xc/test/CMakeLists.txt index 82e8188923c..4a51c38880d 100644 --- a/source/source_hamilt/module_xc/test/CMakeLists.txt +++ b/source/source_hamilt/module_xc/test/CMakeLists.txt @@ -1,3 +1,9 @@ +AddTest( + TARGET MODULE_HAMILT_RpaCutoffConfig + LIBS parameter base device + SOURCES test_exx_info.cpp ../exx_info.cpp +) + AddTest( TARGET MODULE_HAMILT_XCTest_PBE LIBS parameter MPI::MPI_CXX Libxc::xc diff --git a/source/source_hamilt/module_xc/test/test_exx_info.cpp b/source/source_hamilt/module_xc/test/test_exx_info.cpp new file mode 100644 index 00000000000..393f836a503 --- /dev/null +++ b/source/source_hamilt/module_xc/test/test_exx_info.cpp @@ -0,0 +1,41 @@ +#include "../exx_info.h" +#include "source_io/module_parameter/input_parameter.h" +#include +#include + +TEST(RpaCutoffConfig, HseDefaultsRemainSupportedButNondefaultCutoffsAreRejected) +{ + Input_para cut_input; + cut_input.basis_type = "lcao"; + cut_input.rpa = true; + cut_input.dft_functional = "hse"; + // Prepared HSE defaults, after INPUT reset hooks have resolved "default". + cut_input.exx_fock_alpha.clear(); + cut_input.exx_fock_lambda.clear(); + cut_input.exx_erfc_alpha = {"0.25"}; + cut_input.exx_erfc_omega = {"0.11"}; + cut_input.exx_singularity_correction = "limits"; + cut_input.exx_real_number = "1"; + cut_input.exx_ccp_rmesh_times = "2"; + Exx_Info cut_default_info; + EXPECT_NO_THROW(init_exx_info(cut_default_info, cut_input)); + EXPECT_EQ(cut_default_info.info_ri.coulomb_param.count(Conv_Coulomb_Pot_K::Coulomb_Type::Fock), 0); + EXPECT_EQ(cut_default_info.info_ri.coulomb_param.count(Conv_Coulomb_Pot_K::Coulomb_Type::Erfc), 1); + + cut_input.rpa_cut_coulomb = "log_erfc"; + Exx_Info cut_smooth_info; + EXPECT_THROW(init_exx_info(cut_smooth_info, cut_input), std::invalid_argument); + cut_input.rpa_cut_coulomb = "hard"; + cut_input.rpa_cut_radius = "inscribed"; + Exx_Info cut_inscribed_info; + EXPECT_THROW(init_exx_info(cut_inscribed_info, cut_input), std::invalid_argument); + + cut_input.dft_functional = "pbe"; + cut_input.rpa_cut_coulomb = "log_erfc"; + Exx_Info cut_pbe_info; + ASSERT_NO_THROW(init_exx_info(cut_pbe_info, cut_input)); + const auto& cut_fock = cut_pbe_info.info_ri.coulomb_param.at(Conv_Coulomb_Pot_K::Coulomb_Type::Fock); + ASSERT_EQ(cut_fock.size(), 1); + EXPECT_EQ(cut_fock.front().at("rpa_cut_coulomb"), "log_erfc"); + EXPECT_EQ(cut_fock.front().at("rpa_cut_radius"), "inscribed"); +} diff --git a/source/source_io/module_parameter/read_inp_exx_dftu.cpp b/source/source_io/module_parameter/read_inp_exx_dftu.cpp index 913f8087917..f4fe601c5dc 100644 --- a/source/source_io/module_parameter/read_inp_exx_dftu.cpp +++ b/source/source_io/module_parameter/read_inp_exx_dftu.cpp @@ -2,7 +2,6 @@ #include "source_base/tool_quit.h" #include "read_input.h" #include "read_input_tool.h" -#include namespace ModuleIO { void ReadInput::item_exx() @@ -575,16 +574,11 @@ Use the same scheme for consistent SCF and NSCF bands. A source configuration mi item.annotation = "post-SCF truncated Coulomb kernel for RPA/GW export"; item.category = "Exact Exchange (LCAO)"; item.type = "String"; - item.description = "Post-SCF cut-Coulomb kernel: hard preserves Spencer truncation; log_erfc uses erfc(log(r/Rc)/log(rpa_cut_width))/2r. Applies to compressed and uncompressed cut exports; full Coulomb is unchanged. Forces and stresses for log_erfc are not validated."; + item.description = "Post-SCF cut-Coulomb kernel: hard preserves Spencer truncation; log_erfc uses erfc(log(r/Rc)/log(rpa_cut_width))/2r. Applies to Fock components in compressed and uncompressed cut exports; Erfc components retain their existing kernel. Non-default kernel/radius requests require a Fock component. Full Coulomb is unchanged. Forces and stresses for log_erfc are not validated."; item.default_value = "hard"; item.unit = ""; read_sync_string(input.rpa_cut_coulomb); - item.check_value = [](const Input_Item& item, const Parameter& para) { - if (para.input.rpa_cut_coulomb != "hard" && para.input.rpa_cut_coulomb != "log_erfc") - { - ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_coulomb must be hard or log_erfc"); - } - }; + item.domain.set_choices({"hard", "log_erfc"}); this->add_item(item); } { @@ -592,16 +586,11 @@ Use the same scheme for consistent SCF and NSCF bands. A source configuration mi item.annotation = "radius convention for post-SCF cut Coulomb"; item.category = "Exact Exchange (LCAO)"; item.type = "String"; - item.description = "legacy preserves the existing cutoff radius. inscribed uses half the shortest nonzero BvK superlattice vector, including the lattice scale in Bohr. Independent of rpa_cut_coulomb; intended for controlled radius-versus-smoothing comparisons."; + item.description = "legacy preserves the existing cutoff radius. inscribed uses half the shortest nonzero BvK superlattice vector, including the lattice scale in Bohr. Applies only to Fock components; inscribed is rejected if the RPA Coulomb configuration has no Fock component. Independent of rpa_cut_coulomb; intended for controlled radius-versus-smoothing comparisons."; item.default_value = "legacy"; item.unit = ""; read_sync_string(input.rpa_cut_radius); - item.check_value = [](const Input_Item& item, const Parameter& para) { - if (para.input.rpa_cut_radius != "legacy" && para.input.rpa_cut_radius != "inscribed") - { - ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_radius must be legacy or inscribed"); - } - }; + item.domain.set_choices({"legacy", "inscribed"}); this->add_item(item); } { @@ -613,12 +602,7 @@ Use the same scheme for consistent SCF and NSCF bands. A source configuration mi item.default_value = "1.091868692359"; item.unit = ""; read_sync_double(input.rpa_cut_width); - item.check_value = [](const Input_Item& item, const Parameter& para) { - if (!std::isfinite(para.input.rpa_cut_width) || para.input.rpa_cut_width <= 1.0 || para.input.rpa_cut_width > 1.25) - { - ModuleBase::WARNING_QUIT("ReadInput", "rpa_cut_width must be finite and 1 < width <= 1.25"); - } - }; + item.domain.set_rule("value>1 and value<=1.25"); this->add_item(item); } { diff --git a/source/source_io/test_serial/test_read_input.cpp b/source/source_io/test_serial/test_read_input.cpp index be1c52f1b28..b9eec3310ab 100644 --- a/source/source_io/test_serial/test_read_input.cpp +++ b/source/source_io/test_serial/test_read_input.cpp @@ -120,6 +120,48 @@ TEST_F(InputDomainReaderTest, PreservesLegacyEmptyValues) EXPECT_TRUE(input.inp.out_pchg.empty()); } +TEST_F(InputDomainReaderTest, RejectsMalformedRpaCutoffTokens) +{ + for (const std::string cut_name : {"rpa_cut_coulomb", "rpa_cut_radius", "rpa_cut_width"}) + { + for (const std::string cut_tokens : {"\nscf_nmax 20\n", "", " invalid\n", " 1 2\n"}) + { + const std::string cut_input = cut_name + cut_tokens; + const std::string cut_output = reject_input(cut_input); + const std::string cut_error = "Invalid INPUT parameter '" + cut_name + "'"; + EXPECT_NE(cut_output.find(cut_error), std::string::npos) << cut_input; + } + } + for (const std::string cut_input : {"rpa_cut_coulomb log_erfc extra\n", + "rpa_cut_radius inscribed extra\n", + "rpa_cut_width 1\n", "rpa_cut_width 1.2501\n", + "rpa_cut_width nan\n", "rpa_cut_width inf\n"}) + { + const std::string cut_output = reject_input(cut_input); + EXPECT_NE(cut_output.find("Invalid INPUT parameter 'rpa_cut_"), std::string::npos) << cut_input; + } +} + +TEST_F(InputDomainReaderTest, AcceptsRpaCutoffChoicesAndWidthUpperBound) +{ + for (const std::string cut_mode : {"hard", "log_erfc"}) + { + for (const std::string cut_radius : {"legacy", "inscribed"}) + { + const std::string cut_input = "rpa_cut_coulomb " + cut_mode + "\nrpa_cut_radius " + + cut_radius + "\nrpa_cut_width 1.25\n"; + write_input(cut_input); + ModuleIO::ReadInput cut_reader(0); + cut_reader.check_ntype_flag = false; + Parameter cut_parameters; + cut_reader.read_parameters(cut_parameters, filename); + EXPECT_EQ(cut_parameters.inp.rpa_cut_coulomb, cut_mode); + EXPECT_EQ(cut_parameters.inp.rpa_cut_radius, cut_radius); + EXPECT_DOUBLE_EQ(cut_parameters.inp.rpa_cut_width, 1.25); + } + } +} + TEST_F(InputDomainReaderTest, KeepsTrailingCommentsAndAbsolutePaths) { write_input("scf_thr 1e-7 // about iteration\n" From c4ba0f50c0c294b1809c67b515836640c3570653 Mon Sep 17 00:00:00 2001 From: HGAQ <893018484@qq.com> Date: Fri, 30 Oct 2026 16:51:45 -0400 Subject: [PATCH 4/4] ci: retrigger Read the Docs builds