Fix: split LibXC GGA threshold masks for vrho and vsigma following QE convention - #8017
Conversation
Add cal_sgn_vxc returning separate masks: exc/vrho are evaluated down to rho_threshold_lda (1e-10), while only the vsigma gradient term is suppressed below rho_threshold_gga (1e-6) / grho_threshold_gga (1e-10), matching Quantum ESPRESSO's libxc interface (XClib/xc_wrapper_gga.f90).
Cover the two-tier mask logic in libxc_tools.cpp: vrho kept down to rho_threshold_vrho while vsigma is suppressed below rho_threshold_vsigma / grho_threshold_vsigma, GGA vs LDA behavior, and joint spin-channel masking for nspin=2.
There was a problem hiding this comment.
Copilot review overview
🔵 Needs a closer look
The change affects low-density LibXC behavior across multiple calculation paths and warrants final human review.
Review effort: Lite
Findings: None
What changed in this PR
This PR aligns LibXC GGA low-density masking with QE’s two-tier convention, preserving vrho while selectively suppressing unstable vsigma terms.
Changes:
- Added separate
vrhoandvsigmamasks. - Updated LibXC potential calculations and thresholds.
- Added focused unit tests and CMake registration.
| File | Summary |
|---|---|
source/source_hamilt/module_xc/test/test_libxc_tools.cpp |
Tests density, gradient, LDA, and spin masking. |
source/source_hamilt/module_xc/test/CMakeLists.txt |
Registers the new test target. |
source/source_hamilt/module_xc/libxc_tools.cpp |
Implements split threshold masks. |
source/source_hamilt/module_xc/libxc_pot.cpp |
Applies the masks during LibXC evaluation. |
source/source_hamilt/module_xc/libxc_abacus.h |
Declares updated interfaces. |
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
|
End-to-end verification against issue #4695 (the H2O case whose unoccupied-state mismatch motivated #7996): with this PR, all five bands agree with QE-LibXC to within 0.0002 eV (the unoccupied state moves from 0.051 eV above QE to -0.9224 eV vs QE's -0.9226 eV). Detailed comparison table: #4695 (comment) |
The two-tier vrho/vsigma threshold masks in v_xc_libxc shift the HSE energies, forces and stresses of the H2O-based EXX cases, which have low-density regions. Regenerate the six affected result.ref files with the new code (values identical to the CI run of this PR); totaltimeref entries are kept unchanged.
…ask change Merging upstream PR deepmodeling#8017 applies LibXC vrho/exc down to rho=1e-10 (QE convention) instead of masking at 1e-6. The direct effect on the initial density is ~3e-10 eV, but this deliberately-minimal metallic EXX+SOC case never reaches a fixed point and stops mid-oscillation, so the shifted vxc seed moves its endpoint by ~0.11 eV; converged HSE cases shift only by ~1e-5 eV. Regenerate result.ref from the bit-reproducible np=1 run (etot -2926.9134396250 eV, stress sum 3851.624753). Symmetry and magnetic-group assertions (O_h / C_4h / nksibz=3) are unchanged. Slim the threshold file to a single energy override (1e-5 eV): measured np=4 spread is within 9e-7 eV, while force/stress/fatal fit the global defaults (1e-4 / 1e-3 / 1). Verified: np=4 three runs and np=1 one run, 8/8 checks each.
Merging upstream PR deepmodeling#8017 applies LibXC vrho/exc down to rho=1e-10 (QE convention) instead of masking at 1e-6. This NSCF case reads a fixed density matrix and is fully deterministic: the new value is bit-identical across np=1 and np=4 runs (spread 3e-12 eV), and matches the CI cal value to all printed digits (-429.67342983 eV). Regenerate result.ref from the np=1 run (etot -429.6734298284445 eV). No threshold override is needed; the default 1e-7 eV tolerance leaves ample margin for the ~1e-11 eV MPI spread. Verified: np=1 two runs and np=4 two runs, 2/2 checks each.
Linked Issue
Fix #7996
What's Changed
In
XC_Functional_Libxc::v_xc_libxc(source/source_hamilt/module_xc/libxc_pot.cpp),xc_func_set_dens_threshold(&func, 1E-6)made LibXC zero all outputs —exc,vrhoandvsigma— below rho = 1E-6. QE's LibXC interface (XClib/xc_wrapper_gga.f90) instead setsdens_threshold = small (1E-10)and applies two-tier masks:exc/vrhoare kept down to 1E-10, while onlyvsigmais suppressed below rho = 1E-6 or sqrt(|sigma|) = 1E-10. ABACUS-LibXC was the only one of the four PBE code paths (QE built-in / QE LibXC / ABACUS built-in / ABACUS LibXC) that truncatedvrhoat 1E-6, shifting diffuse unoccupied states of finite systems up by ~0.1 eV (see the issue for the full analysis).This PR mirrors QE's two-tier mask convention, exactly as proposed in the issue:
xc_func_set_dens_threshold(&func, 1E-10)instead of 1E-6;XC_Functional_Libxc::cal_sgn_vxc(libxc_tools.cpp) returns two masks:sgn_vrhozeroesexc/vrhoonly where rho <= 1E-10 (for nspin=2 both channels are masked jointly when either spin density fails, matching QE), andsgn_vsigmazeroes onlyvsigmawhere rho <= 1E-6 or sqrt(|sigma|) <= 1E-10;convert_vtxc_vandcal_dhnow take the two masks separately:sgn_vrhomultiplies thevrhocontribution,sgn_vsigmamultiplies the gradient (vsigma) contribution.Scope / side effects
v_xc_meta(mGGA), ingcxc_libxc/gcxc_spin_libxc(stress path), and in the othercal_sgncallers (write_libxc_r.cpp,module_lr/potentials/xc_kernel.cpp). These are left unchanged here and can be addressed separately, as noted in the issue.docs/parameters.yamlordocs/advanced/input_files/input-main.mdis required.Tests
source/source_hamilt/module_xc/test/test_libxc_tools.cpp(targetMODULE_HAMILT_XCTest_LIBXC_TOOLS) coverscal_sgn_vxc: the two-tier density thresholds for nspin=1 GGA, the sqrt(|sigma|) trigger, the LDA case (vsigma mask stays 1), joint spin-channel masking for nspin=2, and that the up-down cross component of sigma does not trigger the vsigma mask.Verification
Executable:
./build/rel/abacus_std_para --version-> ABACUS version v3.11.0-beta10 (GCC, RelWithDebInfo, ENABLE_LIBXC=ON, ENABLE_MLALGO=OFF).OMP_NUM_THREADS=1 ctest --test-dir build/rel -R MODULE_HAMILT_XCTest --output-on-failure: 9/9 passed, including the newMODULE_HAMILT_XCTest_LIBXC_TOOLSand the existingMODULE_HAMILT_XCTest_VXC(which exercisesv_xc_libxcend-to-end).OMP_NUM_THREADS=1, 2 MPI ranks,tests/integrate/Autotest.sh -r <regex>:tests/01_PW: 14 matched cases, 0 failed — including202_PW_ONCV_Libxcand204_PW_SY.tests/07_OFDFT: 10 matched cases, including22_OF_LibxcPBEand26_OF_MD_LibxcPBEpassed; no reference updates were needed anywhere.ENABLE_MLALGO=OFF(CI builds with-DENABLE_MLALGO=ON,.github/workflows/test.yml):104_PW_Gene_Descriptors: the SCF part converged and reproduces the reference etot to 1.7E-9 eV (threshold 1E-7 eV); only the MLKEDF descriptor.npyoutputs are not generated without MLALGO.06_OF_KE_MPN: aborts at input check with "ML KEDF requires ENABLE_MLALGO option" before any calculation. The system is periodic fcc Al (no low-density region), so this PR cannot change its result; CI coverage is needed.python3 tools/03_code_analysis/agent_governance_check.py --base upstream/develop --head HEAD --format text: only the documentation-sync warning remains, addressed by the no-docs-needed statement above.python3 tools/03_code_analysis/code_quality_score.pyon the changed files:libxc_pot.cppscores 32, up from 31 for the unmodifiedupstream/developversion — the sub-60 score is pre-existing historical debt (file length,v_xc_metacomplexity), not introduced by this PR; the other changed files score 68/71/95.