fixed vvhv bug by reimplementing proto_hv construction#367
Open
v-agamshayit wants to merge 2 commits intomainfrom
Open
fixed vvhv bug by reimplementing proto_hv construction#367v-agamshayit wants to merge 2 commits intomainfrom
v-agamshayit wants to merge 2 commits intomainfrom
Conversation
v-agamshayit
requested review from
ConradJohnston,
nabbelbabbel,
sobolevnrm and
wavefunction91
as code owners
|
@v-agamshayit please read the following Contributor License Agreement(CLA). If you agree with the CLA, please reply with the following information.
Contributor License AgreementContribution License AgreementThis Contribution License Agreement (“Agreement”) is agreed to by the party signing below (“You”),
|
Contributor
There was a problem hiding this comment.
Pull request overview
This PR addresses non-determinism and correctness issues in the VVHV hard-virtual localization workflow, particularly around proto hard-virtual (proto-HV) construction under overlap degeneracies, and updates/expands regression tests to lock in the new deterministic behavior.
Changes:
- Reimplemented proto-HV construction logic and added an option to use Eigen’s eigensolver for deterministic orthonormalization.
- Fixed “weighted orthogonalization” to weight by inverse orbital spreads instead of spreads.
- Updated VVHV localization tests to use tighter reference metrics and added regression tests using reordered (“scrambled”) basis shells.
Reviewed changes
Copilot reviewed 2 out of 2 changed files in this pull request and generated 6 comments.
| File | Description |
|---|---|
cpp/src/qdk/chemistry/algorithms/microsoft/localization/vvhv.cpp |
Adds deterministic eigensolver option and rewrites proto-HV construction; adjusts weighted orthogonalization behavior. |
cpp/tests/test_localization.cpp |
Updates VVHV metric assertions and adds tests to validate stability under basis-shell reordering. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
Comment on lines
536
to
541
| this->orthonormalization( | ||
| num_atomic_orbitals_al_ori, num_atomic_orbitals_al_ori, | ||
| overlap_ori_al.data(), temp.data(), C_psi.data(), 1e-6, 0, | ||
| "initial orthonormalization of original basis on atom " + | ||
| std::to_string(atom_index) + " angular momentum " + | ||
| std::to_string(l)); |
There was a problem hiding this comment.
In proto_hv(), the initial orthonormalization of C_psi still uses the default LAPACK eigensolver (use_Eigen defaults to false). Since this basis feeds into the rest of the proto-HV construction, this undermines the intended determinism fix—pass use_Eigen through to this orthonormalization call (and use the same deterministic path used later for muBar).
| Eigen::MatrixXd pBar = Eigen::MatrixXd::Zero(num_atomic_orbitals_al_ori, num_atomic_orbitals_al_min); | ||
| this->orthonormalization( | ||
| num_atomic_orbitals_al_ori, num_atomic_orbitals_al_min, | ||
| overlap_ori_al.data(), pTilde.data(), pBar.data(), 1e-6, 0, ""); |
There was a problem hiding this comment.
proto_hv(): the orthonormalization used to build pBar currently uses the default LAPACK eigensolver (no use_Eigen argument provided). If determinism is required for handling degeneracies, this step should also use the deterministic eigensolver (pass use_Eigen) so pBar is reproducible across runs/platforms.
Suggested change
| overlap_ori_al.data(), pTilde.data(), pBar.data(), 1e-6, 0, ""); | |
| overlap_ori_al.data(), pTilde.data(), pBar.data(), 1e-6, 0, "", /*use_Eigen=*/true); |
Comment on lines
+863
to
+867
| Eigen::VectorXd spreads_hv(nhv); | ||
| this->calculate_orbital_spreads(C_hard_virtuals, spreads_hv); | ||
| // Weight each orbital by spread | ||
| for (int orb = 0; orb < nhv; ++orb) | ||
| C_hard_virtuals.col(orb) *= spreads_hv(orb); | ||
| C_hard_virtuals.col(orb) /= spreads_hv(orb); |
There was a problem hiding this comment.
Weighted orthogonalization now divides by spreads_hv(orb), but calculate_orbital_spreads() clamps spreads to >= 0 and can return exactly 0.0. This can introduce division-by-zero/infs and destabilize the subsequent orthonormalization; clamp spreads to a small epsilon (or skip weighting for near-zero spreads) before dividing.
Comment on lines
1012
to
1016
| if (use_Eigen) { | ||
| Eigen::SelfAdjointEigenSolver<Eigen::MatrixXd> eig(S); | ||
| eigenvalues = eig.eigenvalues(); | ||
| S = eig.eigenvectors(); // S now contains eigenvectors U | ||
| } |
There was a problem hiding this comment.
When use_Eigen is enabled, SelfAdjointEigenSolver is used without checking eig.info(). If the decomposition fails, eigenvalues/eigenvectors will be invalid and downstream math will silently corrupt results. Check eig.info() and throw (ideally including error_label) on failure.
cpp/tests/test_localization.cpp
Outdated
| // sort shells by angular momentum and exponent | ||
| std::stable_sort(shells.begin(), shells.end(), shell_comparator_ascending); | ||
|
|
||
| return std::make_shared<BasisSet>(basis->get_name() + "_scrambled", shells, basis->get_structure()); |
There was a problem hiding this comment.
scramble_basis_shells() constructs a new BasisSet without preserving the original basis' AOType (spherical vs cartesian). The constructor defaults to AOType::Spherical, so this helper can silently change the basis representation and invalidate the SCF/localization behavior being tested. Pass basis->get_atomic_orbital_type() into the BasisSet constructor.
Suggested change
| return std::make_shared<BasisSet>(basis->get_name() + "_scrambled", shells, basis->get_structure()); | |
| return std::make_shared<BasisSet>(basis->get_name() + "_scrambled", | |
| shells, | |
| basis->get_structure(), | |
| basis->get_atomic_orbital_type()); |
Comment on lines
+705
to
+708
| std::vector<size_t> indices(shell.get_num_primitives()); | ||
| std::iota(indices.begin(), indices.end(), 0); | ||
| std::sort(indices.begin(), indices.end(), [&shell](size_t a, size_t b) { | ||
| return shell.exponents(a) > shell.exponents(b); |
There was a problem hiding this comment.
scramble_basis_shells() uses std::iota/std::sort/std::stable_sort but this translation unit doesn't include the standard headers that guarantee those declarations (e.g., for std::iota and for sort/stable_sort). Add the required includes to avoid relying on transitive headers.
📊 Coverage Summary
Detailed Coverage ReportsC++ Coverage DetailsPython Coverage DetailsPybind11 Coverage Details |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Changed the formation of the proto HV's to address degeneracies. When the overlap matrix of the hard-virtual space has degenerate eigenvectors, proto HV's can be arbitrarily rotated into each other, resulting in slightly different hard virtual localizations. This is addressed by using Eigen's deterministic eigensolver when orthonormalizing the proto HV's.
Fixed the weighted orthogonalization of the hard virtuals to use the inverse orbital spreads as weights (rather than the spreads)