From d11cffa60e249254714fda4c7cc06cd445421b25 Mon Sep 17 00:00:00 2001 From: maki49 <1579492865@qq.com> Date: Thu, 24 Sep 2026 12:43:17 -0400 Subject: [PATCH 1/3] Fix spin-channel mismatch in LR complex transition dipole (length gauge) LR_Spectrum>::cal_transition_dipole_istate_length looped over all nspin_x spin channels of the transition density matrix but called get_DMR_real_imag_part() without selecting a channel, so the function copied every channel of DM_trans (size nspin_x) into the single-channel DM_trans_real_imag. For nspin_x==2 (open-shell/spin- polarized, multi-k) this hits an assertion failure (or out-of-bounds access in release builds), and even when it doesn't crash it double- counts the merged density on each iteration of the outer loop. Add an overload of get_DMR_real_imag_part() that copies a single named spin channel into the (always single-channel) DMR_real, and use it from the two call sites in cal_transition_dipole_istate_length so each outer loop iteration only processes its own channel, matching the double (gamma-only) specialization's behavior. Verified: reproduced the original assertion failure with gdb on the open-shell (nspin=2) multi-k length-gauge path, confirmed the fix resolves it, and confirmed the fixed multi-k oscillator strength / transition dipoles match the gamma-only result exactly. Co-Authored-By: Claude Sonnet 5 --- source/source_lcao/module_lr/lr_spectrum.cpp | 4 +-- .../module_lr/utils/lr_util_hcontainer.cpp | 30 +++++++++++++++++++ .../module_lr/utils/lr_util_hcontainer.h | 7 +++++ 3 files changed, 39 insertions(+), 2 deletions(-) diff --git a/source/source_lcao/module_lr/lr_spectrum.cpp b/source/source_lcao/module_lr/lr_spectrum.cpp index 9e1368d0bec..abb5bd81765 100644 --- a/source/source_lcao/module_lr/lr_spectrum.cpp +++ b/source/source_lcao/module_lr/lr_spectrum.cpp @@ -100,13 +100,13 @@ ModuleBase::Vector3> LR::LR_Spectrum>: LR_Util::initialize_DMR(DM_trans_real_imag, this->pmat, this->ucell, this->gd_, this->orb_cutoff_); // real part - LR_Util::get_DMR_real_imag_part(DM_trans, DM_trans_real_imag, ucell.nat, 'R'); + LR_Util::get_DMR_real_imag_part(DM_trans, DM_trans_real_imag, ucell.nat, is, 'R'); ModuleBase::GlobalFunc::ZEROS(rho_trans_real[0], this->rho_basis.nrxx); ModuleGint::cal_gint_rho(DM_trans_real_imag.get_dmr_vec(), 1, rho_trans_real, false); // LR_Util::print_grid_nonzero(rho_trans_real[0], this->rho_basis.nrxx, 10, "rho_trans"); // imag part - LR_Util::get_DMR_real_imag_part(DM_trans, DM_trans_real_imag, ucell.nat, 'I'); + LR_Util::get_DMR_real_imag_part(DM_trans, DM_trans_real_imag, ucell.nat, is, 'I'); ModuleBase::GlobalFunc::ZEROS(rho_trans_imag[0], this->rho_basis.nrxx); ModuleGint::cal_gint_rho(DM_trans_real_imag.get_dmr_vec(), 1, rho_trans_imag, false); // LR_Util::print_grid_nonzero(rho_trans_imag[0], this->rho_basis.nrxx, 10, "rho_trans"); diff --git a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp index 97c75d8243f..6165353273c 100644 --- a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp +++ b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp @@ -32,6 +32,36 @@ namespace LR_Util } } + void get_DMR_real_imag_part(const elecstate::DensityMatrix, std::complex>& DMR, + elecstate::DensityMatrix, double>& DMR_real, + const int& nat, + const int& is, + const char& type) + { + assert(is < DMR.get_DMR_vector().size()); + assert(DMR_real.get_DMR_vector().size() == 1); + bool get_imag = (type == 'I' || type == 'i'); + auto dr = DMR.get_DMR_vector()[is]; //get_DMR_pointer() has bug when is=0 + auto dr_real = DMR_real.get_DMR_vector()[0]; + assert(dr != nullptr); + assert(dr_real != nullptr); + for (int ia = 0;ia < nat;ia++) { + for (int ja = 0;ja < nat;ja++) + { + auto ap = dr->find_pair(ia, ja); + auto ap_real = dr_real->find_pair(ia, ja); + for (int iR = 0;iR < ap->get_R_size();++iR) + { + // R index may be different between the two HContainers, find by R value instead of R-index + auto dR = ap->get_R_index(iR); + auto ptr = ap->get_HR_values(iR).get_pointer(); + auto ptr_real = ap_real->get_HR_values(dR.x, dR.y, dR.z).get_pointer(); + for (int i = 0;i < ap->get_size();++i) { ptr_real[i] = (get_imag ? ptr[i].imag() : ptr[i].real()); } + } + } + } + } + void set_HR_real_imag_part(const hamilt::HContainer& HR_real, hamilt::HContainer>& HR, const int& nat, diff --git a/source/source_lcao/module_lr/utils/lr_util_hcontainer.h b/source/source_lcao/module_lr/utils/lr_util_hcontainer.h index cb89c988cad..9147aa45149 100644 --- a/source/source_lcao/module_lr/utils/lr_util_hcontainer.h +++ b/source/source_lcao/module_lr/utils/lr_util_hcontainer.h @@ -42,6 +42,13 @@ namespace LR_Util module_dm::DensityMatrix, double>& DMR_real, const int& nat, const char& type = 'R'); + /// overload: only copy the `is`-th spin channel of DMR (source) into the (single-channel) DMR_real, + /// to avoid mixing/overlapping spin channels when DMR has more than one spin channel + void get_DMR_real_imag_part(const elecstate::DensityMatrix, std::complex>& DMR, + elecstate::DensityMatrix, double>& DMR_real, + const int& nat, + const int& is, + const char& type = 'R'); void set_HR_real_imag_part(const hamilt::HContainer& HR_real, hamilt::HContainer>& HR, const int& nat, From 2f9902cfaef9a0b43d34e6b29d6b038fe05b6ea1 Mon Sep 17 00:00:00 2001 From: maki49 <1579492865@qq.com> Date: Thu, 24 Sep 2026 12:53:50 -0400 Subject: [PATCH 2/3] Fix segfault in LR real/imag HContainer helpers under MPI-parallel runs MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit get_DMR_real_imag_part() (both overloads) and set_HR_real_imag_part() dereferenced the result of HContainer::find_pair(ia, ja) without a null check. Under real multi-rank MPI parallelism, the HContainer is distributed 2D block-cyclic, so find_pair() legitimately returns nullptr for atom pairs not owned by the current rank — dereferencing that pointer segfaults. Reproduced with gdb on tests/08_EXX/54_GO_ULR_HF (gamma_only=0, KPT 1 1 1, mpirun -np 2): SIGSEGV in get_DMR_real_imag_part, called from OperatorLRHxc::grid_calculation during the Casida eigenvalue solve — a different call site/crash than the spin-channel-mismatch bug fixed in the previous commit, and only reproducible with more than one MPI rank (a single-rank run owns every atom pair locally, masking the bug). Skip atom pairs not present on the local rank, matching the existing pattern used elsewhere in the codebase for MPI-parallel HContainer access. Verified: the same 54_GO_ULR_HF case now runs cleanly under mpirun -np 2 (gdb shows a clean exit, no crash) and its excitation energies exactly match both the gamma_only=1 run and result.ref. Co-Authored-By: Claude Sonnet 5 --- .../source_lcao/module_lr/utils/lr_util_hcontainer.cpp | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp index 6165353273c..0a4d0272acc 100644 --- a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp +++ b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp @@ -19,6 +19,9 @@ namespace LR_Util { auto ap = dr->find_pair(ia, ja); auto ap_real = dr_real->find_pair(ia, ja); + // under MPI-parallel (2D block-cyclic) HContainer, an atom pair not owned by this rank + // is absent from find_pair() and returns nullptr here; skip it + if (!ap || !ap_real) { continue; } for (int iR = 0;iR < ap->get_R_size();++iR) { // R index may be different between the two HContainers, find by R value instead of R-index @@ -50,6 +53,9 @@ namespace LR_Util { auto ap = dr->find_pair(ia, ja); auto ap_real = dr_real->find_pair(ia, ja); + // under MPI-parallel (2D block-cyclic) HContainer, an atom pair not owned by this rank + // is absent from find_pair() and returns nullptr here; skip it + if (!ap || !ap_real) { continue; } for (int iR = 0;iR < ap->get_R_size();++iR) { // R index may be different between the two HContainers, find by R value instead of R-index @@ -73,6 +79,9 @@ namespace LR_Util { auto ap = HR.find_pair(ia, ja); auto ap_real = HR_real.find_pair(ia, ja); + // under MPI-parallel (2D block-cyclic) HContainer, an atom pair not owned by this rank + // is absent from find_pair() and returns nullptr here; skip it + if (!ap || !ap_real) { continue; } for (int iR = 0;iR < ap->get_R_size();++iR) { // R index may be different between the two HContainers, find by R value instead of R-index From b664a54afa728a46cf0699c57ce791487ac2de48 Mon Sep 17 00:00:00 2001 From: maki49 <1579492865@qq.com> Date: Thu, 24 Sep 2026 22:51:28 -0400 Subject: [PATCH 3/3] Adapt new get_DMR_real_imag_part overload to module_dm refactor Upstream #8000 renamed elecstate::DensityMatrix to module_dm::DensityMatrix and get_DMR_vector/get_DMR_pointer to get_dmr_vec/get_dmr_ptr; update the spin-index overload added in this branch accordingly. Co-Authored-By: Claude Sonnet 5 --- .../module_lr/utils/lr_util_hcontainer.cpp | 12 ++++++------ .../source_lcao/module_lr/utils/lr_util_hcontainer.h | 4 ++-- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp index 0a4d0272acc..b655f00c99a 100644 --- a/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp +++ b/source/source_lcao/module_lr/utils/lr_util_hcontainer.cpp @@ -35,17 +35,17 @@ namespace LR_Util } } - void get_DMR_real_imag_part(const elecstate::DensityMatrix, std::complex>& DMR, - elecstate::DensityMatrix, double>& DMR_real, + void get_DMR_real_imag_part(const module_dm::DensityMatrix, std::complex>& DMR, + module_dm::DensityMatrix, double>& DMR_real, const int& nat, const int& is, const char& type) { - assert(is < DMR.get_DMR_vector().size()); - assert(DMR_real.get_DMR_vector().size() == 1); + assert(is < DMR.get_dmr_vec().size()); + assert(DMR_real.get_dmr_vec().size() == 1); bool get_imag = (type == 'I' || type == 'i'); - auto dr = DMR.get_DMR_vector()[is]; //get_DMR_pointer() has bug when is=0 - auto dr_real = DMR_real.get_DMR_vector()[0]; + auto dr = DMR.get_dmr_vec()[is]; //get_dmr_ptr() has bug when is=0 + auto dr_real = DMR_real.get_dmr_vec()[0]; assert(dr != nullptr); assert(dr_real != nullptr); for (int ia = 0;ia < nat;ia++) { diff --git a/source/source_lcao/module_lr/utils/lr_util_hcontainer.h b/source/source_lcao/module_lr/utils/lr_util_hcontainer.h index 9147aa45149..aaa9f8b742d 100644 --- a/source/source_lcao/module_lr/utils/lr_util_hcontainer.h +++ b/source/source_lcao/module_lr/utils/lr_util_hcontainer.h @@ -44,8 +44,8 @@ namespace LR_Util const char& type = 'R'); /// overload: only copy the `is`-th spin channel of DMR (source) into the (single-channel) DMR_real, /// to avoid mixing/overlapping spin channels when DMR has more than one spin channel - void get_DMR_real_imag_part(const elecstate::DensityMatrix, std::complex>& DMR, - elecstate::DensityMatrix, double>& DMR_real, + void get_DMR_real_imag_part(const module_dm::DensityMatrix, std::complex>& DMR, + module_dm::DensityMatrix, double>& DMR_real, const int& nat, const int& is, const char& type = 'R');