Repository navigation
Conversation
e5a8e55 to
cd8346b
Compare
|
Thanks a lot for your hard work on this PR — the overall direction looks great! Below are a few suggestions, ordered from the easiest to the more involved ones. Most should be quick to address, and I'm happy to discuss any of them.
|
|
@mohanchen For the last point (7.), no new ESolver is added: esolver_lr_grad.cpp is just a separate .cpp file from esolver_lr_lcao_tddft.cpp, corresponding to esolver_lr_lcao_tddft.h. |
Port the LR-Grad series (44 commits, 7c6398d..6aea1fc, all by maki49) onto current develop as a single squashed change, adapted to develop's large-scale refactor. Path/layout migration: module_lr -> source_lcao/module_lr module_ri -> source_lcao/module_ri (Exx_LRI.h -> exx_lri.h) esolver_lrtd_lcao -> source_esolver/esolver_lr_lcao_tddft ESolver_LR now lives in namespace ModuleESolver, not LR. Interface adaptation: - Grid integration: the Gint_Gamma/Gint_k objects and TGint<T> selector are gone. All call sites now use the free functions in ModuleGint (cal_gint_rho / cal_gint_vl / cal_gint_fvl); the gint pointer has been removed from LR_Force, LR_Density, OperatorLRHxc, HamiltLR/HamiltULR, the Z-equation Hamiltonians and the multiplier helpers. - PulayForceStress::cal_pulay_fs (grid overload) keeps the explicit `nspin` argument the gradient code needs; develop's own callers pass PARAM.inp.nspin. - Parallel_Orbitals::get_row_size(iat)/get_col_size(iat) -> get_nrow_atom(iat)/get_ncol_atom(iat). - elecstate::Potential::init_pot(istep, chg) -> init_pot(chg); get_effective_v -> get_eff_v; Structure_Factor::setup_structure_factor -> setup; Charge::allocate now takes kin_den. - ModuleIO::output_single_R takes a SparseWriteOptions struct. - ModuleIO::read_wfc_nao uses develop's signature (ik2iktot/nkstot/nspin). - ModuleBase::timer::tick -> start/end pairs. - LR_Util::gather_2d_to_full keeps develop's explicit-dimension form; scatter_full_to_2d is added for the Z-equation solver. - get_DMR_real_imag_part / set_HR_real_imag_part are now the templated, nat-free versions; module_bse call sites updated, the now-dead lr_util_hcontainer.cpp is removed. - DensityMatrix::cal_DMR is const (_DMR is mutable). Known gap: LR_Force::cal_force_exx_gs_dm_relaxed_diff still uses RI::LR from LibRI, which no longer derives from RI::Exx in the current LibRI master (commit 5c6c262 refactored it into a k-space CVCX/BSE helper). That single function does not compile; everything else does.
…X force kernel LibRI commit 5c6c262 repurposed RI::LR: it no longer derives from RI::Exx and lost set_Ds/cal_Hs/cal_force/force. Reintroduce just the piece the LR-TDDFT gradient needs as LR::ExxForceTwoDM, a thin RI::Exx subclass whose cal_force overload overwrites post_2D.saves["Ds_<suffix>"] with the left density matrix before delegating to Exx::cal_force. This reproduces the old RI::LR semantics (D_IJ from Ds_left, D_KL from set_Ds) without pinning LibRI. Also drop dm_trans/dmr_complex.cpp from the build: develop moved the DensityMatrix<complex,complex>::cal_DMR specialization into density_matrix.cpp, so keeping both gave a duplicate definition at link time.
…ose HSolverLR timer
Two crashes found running the H2/TDHF/SZ gradient case:
1. PulayForceStress::cal_pulay_fs (the LR overload in pulay_force_hcontainer.h)
passed the density matrix's spin count as `nspin` to ModuleGint::cal_gint_fvl
while supplying a single potential channel. Gint_fvl::cal_fvl_svl_ indexes
vr_eff_[is] for is<nspin, so with a two-channel DM (the test_force path feeds
in the ground-state DM) it read past the vector and segfaulted inside the
OpenMP region. The LR kernel produces one potential channel, so the grid
integrals now run with nspin=1 on the leading DM channel -- exactly what the
old `Gint_inout(is=0, ...)` calls did.
2. HSolverLR::solve started a timer that was never ended, so the second
(triplet) call threw "timer::start HSolverLR::solve". Added the matching
ModuleBase::timer::end on the common exit path.
Verified on H2 / TDHF / SZ (nspin=2, gamma-only):
- SCF reproduces the pre-port reference bit-for-bit (E_KohnSham -2.3110142036 Ry,
TOTAL-FORCE -14.7806682608 eV/Ang).
- test_force cross-check: GS Hxc force via potential -28.4711400896 vs via the
LR kernel -28.4711414920 eV/Ang (7 significant digits).
- The reproduced EXX force (31.3817781891) is exactly 2x the SCF EXX force
(15.6908890935), confirming ExxForceTwoDM replaces RI::LR correctly.
- Excitations 22.484 eV (singlet) / 13.813 eV (triplet); excited-state
gradients -0.6313 / -2.01497 eV/Ang, antisymmetric with ~1e-14 transverse
components.
solve_Z_lapack built a replicated Hessian via HamiltLR::matrix() but handed LAPACK the pX-distributed right-hand side `R` directly. lapack_linear_solver then copies n_global*nstates entries out of a buffer that only holds ld*nstates = nk*pX[0].get_local_size()*nstates -- an out-of-bounds read as soon as the local size is smaller than the global one, and a segfault on any rank whose local size is 0. Gather R into a global buffer first, mirroring the scatter of Z_full that already follows the solve. Serial results are bit-for-bit unchanged (verified on H2/TDHF/SZ: excited-state gradients -0.631306 / -2.01497 eV/Ang, all per-term forces identical to the pre-fix run).
`PotHxcLR::cal_v_eff` dominates an LR run: 74% of the total for 07_LiF/pbe (949 s of 1289 s over 541 calls). Roughly half of that was FFT; the rest was serial element-wise work and allocator traffic. Five changes, no algorithmic change: 1. OpenMP on the element-wise loops. Everything these loops call (`grad_rho`, `grad_dot`, the Hartree kernel) was already parallel; the integrands here were the only part still running on one thread. 2. Raw pointers instead of `.at()` (43 sites). The bounds check is a branch per access, ~10 per grid point, and it blocks vectorization. The outer-vector lookups (`drho_gs.at(0)`) are hoisted out of the loops, where they keep their bounds check for free. 3. A scratch pool instead of per-call temporaries. `drho`, `gdot_terms` and `vxc_tmp` were allocated and value-initialized on every call and then overwritten before being read -- ~450 MB of pointless memset per call at a 200^3 grid, plus the page faults. `grad_dot` assigns rather than accumulates, so even `vxc_tmp`'s zero-fill was dead. The pool is shared by all instances (a run holds three `PotHxcLR`) so it does not multiply. 4. `add_v_hartree` replaces `H_Hartree_pw::v_hartree`, which (a) re-did the forward FFT of rho^X that the GGA branch needs anyway -- 1 of the 10 FFTs per call was pure duplication, (b) reduced a Hartree "energy" of the transition density through `Parallel_Reduce::reduce_pool` every call, a collective nobody reads that also clobbers the global `H_Hartree_pw::hartree_energy`, and (c) returned a `matrix` by value, to which `v_eff += 2 * (...)` added a second full-size temporary. 5. The nspin=2 singlet/triplet combinations `v2rho2_uu -+ v2rho2_ud` and `2*vsigma_uu -+ vsigma_ud` are pre-contracted once per potential instead of being rebuilt at every grid point of every call, which also replaces two strided reads with one contiguous one. Costs 8-16 B/point, and only for closed-shell nspin=2. `PotGradXCLR::cal_v_eff` (the g^xc branch feeding the Z-vector RHS) gets 1-3 of the same treatment, for -20% to -28%. The shared pool matters more there: a `PotGradXCLR` is constructed inside the loop over excited states, so per-object buffers would never be reused at all. Measured (single node, 16 threads), analytic gradients unchanged: 01_Si/lda 93 s -> 51 s (-45%) cal_v_eff 66 -> 45 ms/call 01_Si/pbe 229 s -> 125 s (-46%) cal_v_eff 356 -> 268 ms/call 02_Li2/lda 165 s -> 86 s (-48%) cal_v_eff 559 -> 278 ms/call 02_Li2/pbe 530 s -> 253 s (-52%) cal_v_eff 2651 -> 1553 ms/call (The total also benefits from solving the Z-vector equation once instead of three times, a separate fix; the per-call figures above isolate this commit.) Excitation energies are identical to every printed digit; the nspin=1 cases agree to 1e-14, and of 60 force components per nspin=2 case, 2-3 differ by exactly one unit in the last printed digit (the Hartree prefactor is now associated differently). Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01QMV7xFZ69hhc3HnVBHbD96
Add two-step four-rank BeH2 average and CH4 JT fixtures with gauge-independent subspace checks. Restore force checks for supported LDA/PBE Gamma cases using the CI image's LibXC version rebuilt with KXC and cached. Time the complete follow_target_state_ function, including its early return. Validation: cmake --build build -j2; python3 tests/integrate/tools/test_root_check.py (8 passed); four-rank Autotest on five LR cases (15 property checks passed) with KXC-enabled LibXC 5.1.7 runtime; workflow YAML and shell syntax; staged governance and diff checks. Timer instrumentation and test-only changes do not change INPUT behavior; no parameter documentation update is needed.
Resolve RI build-list conflicts by retaining both exx_lr_ws and develop's extracted ri_util in CMake and Makefile.Objects. Preserve the existing reviewed fix SHAs and the develop CI updates alongside KXC support. Validation: cmake --build build -j2 (GNU full executable passed); serial lr and three root-tracking CTest targets passed; GNU Makefile ri_util.o and exx_lr_ws.o passed using local MKL FFTW headers; eight root-checker tests passed; five four-rank LR integration cases passed all 15 property checks with OMP_NUM_THREADS=1 and KXC-enabled LibXC 5.1.7 runtime. Workflow YAML/shell syntax, governance and diff checks passed without blocking findings. No new INPUT behavior is introduced by conflict resolution.
Adapt LR relaxation validation to develop's validation_tools directory and modular collectors. Retain excited-state force extraction in props_basic, move the root checker and its tests, and update the CI test path. Preserve all previously reviewed fix SHAs. Validation: eight root-checker unit tests passed; recollecting existing logs for five LR cases matches their previous result.out exactly excluding timing; missing ionic step produces CompareRootTracking_pass 1; workflow YAML and collector shell syntax passed; no unresolved paths; conflict-resolution files pass whitespace checks and staged governance has no blocking findings. Physical calculations and C++ builds were not rerun for this script conflict resolution. Existing PR and incoming develop whitespace warnings were left unchanged.
Report the resolved full KS band count and the selected LR window count instead of reading eig_ks.nc before its allocation. This only changes diagnostic output; INPUT and numerical behavior are unchanged, so parameter documentation needs no update. Validation: cmake --build build --target abacus_std_para -j4 passed; four-rank MPI controls with OMP_NUM_THREADS=1 and lr_solver=lapack passed for both ks-lr and standalone lr, asserting KS=6 and LR window=5 and nonempty absorption.dat. The initial small-window Davidson controls hit its dimension assertion; LAPACK was used for the log checks. Staged governance and whitespace checks passed, with test/documentation review warnings addressed by this evidence and rationale.
Replace the unconditional KXC build reminder with a capability check on initialized functionals only when third derivatives are needed. Check independent mixed components recursively and name unsupported functionals with the loaded LibXC version. Release initialized functionals before reporting the normal ABACUS fatal error. Validation: cmake --build build -j2 passed; MODULE_LR_xc_kxc passed five tests with KXC-enabled LibXC and the system library; three four-rank LDA/PBE integration cases passed all 13 property checks with OMP_NUM_THREADS=1; a four-rank control loading system LibXC without KXC failed early with the expected functional/version/build guidance; C++11 syntax checks with and without LibXC passed. Staged governance has no blocking findings. No INPUT parameter behavior changes; this replaces an unconditional reminder and improves diagnostics for an already unsupported configuration, so parameter documentation updates are not required.
The updated write_vdata_palgrid interface requires a running log stream and data description. Pass the stream explicitly through LR density output and supply an electron-hole density description for both closed-shell and unrestricted callers, fixing the shared build failure across GNU, Intel, CUDA, serial, no-ELPA and Makefile jobs. Validation: cmake --build build -j4 passed; affected ESolver source passed GNU syntax compilation with serial and no-ELPA configuration flags; Intel Makefile esolver_lr_lcao_tddft.o passed; four-rank closed-shell LDA and unrestricted HF out_chg=1 controls generated 4 and 6 complete finite cube files with matching log descriptions, and both force totals passed existing references within 1e-6. OMP_NUM_THREADS=1. Conflict-resolution whitespace and staged governance checks have no blocking findings; inherited develop whitespace warnings are unchanged. No new INPUT behavior, so parameter docs need no update. Remote CI results remain pending.
OK, thanks. I will look into your details soon. |
Keep the existing dimensionless 1e-6 diagnostic threshold and explain its meaning in the warning. Rename state tracking sources and checker without changing numerical behavior. Test basenames use grp/chk abbreviations to retain the requested length limit; no INPUT or documentation behavior changes.
Separate orbital distributions and ground-state input, KS initialization and refresh, and dimension/kernel configuration into esolver_lr_orb.cpp, esolver_lr_init.cpp, and esolver_lr_cfg.cpp. Keep member definitions under source_esolver; update CMake and Make. Preserve numerical behavior and global dependency counts. No INPUT changes. Verified: full GNU MPI builds after each extraction; Intel Make object compilation; C++11 checks without MPI and ELPA; four-rank BeH2 average, CH4 JT, and unrestricted HF regressions.
|
@mohanchen The files named |
Apply the flattened UnitCell symmetry-analysis interface to esolver_lr_init.cpp while preserving the ESolver split and original feature SHAs. Verified full GNU MPI build, C++11 initialization checks without MPI and without ELPA, and four-rank two-step BeH2 state tracking. No INPUT behavior changes.
Use vector initialization for iat2it and iat2ia after the develop UnitCell refactor. Fix both the overlap block fixture and the cross-geometry AO overlap fixture. Verified both MPI-linked test targets, the four-rank state-overlap suite, and the serial state-overlap suite. No production behavior or reference changes.
Redistribute occupied and virtual MO-overlap subblocks to descriptors starting at (1,1) on the amplitude BLACS grid before PBLAS contractions. Preserve the existing block size and local distribution without gathering matrices. Remove the unit-block restriction while retaining dimension and layout checks. Regression tests cover block sizes 2, 4 and 32, misaligned occupied windows, multiple blocks, complex adjoints and empty ranks against independent contractions. Verified full GNU MPI build; 19 overlap tests in serial, MPI single-rank and four-rank runs; projection tests on a six-rank non-square grid; two-step BeH2 tracking with nb2d=32 on one rank and nb2d=2 on four ranks. No INPUT defaults or reference changes.
Redistribute selected old/new MO coefficient columns before AO-to-MO contraction, constructing only oo and vv overlaps once per spin. Reuse these aligned blocks for every tracked reference; projection now performs only the two matrix multiplications. Reuse the coefficient work buffer and avoid storing unused occupied-virtual overlaps. Preserve block sizes, numerical formulas and INPUT behavior. Verified full GNU MPI build; 21 overlap tests in serial and MPI single/four-rank runs; distributed tests on six ranks; two-step BeH2 tracking with nb2d=32 on one rank and nb2d=2 on four ranks. Added independent occupied/virtual construction checks including complex and misaligned large blocks. No reference changes or user-facing documentation changes required.
Merge cd2c7ce, including the upstream consistent AO radial spline values/derivatives correction (deepmodeling#8119). Resolve the unrestricted HF reference conflict by retaining the LR force check and adopting the verified new spectrum. Refresh force sums in the five Gamma LR cases; multi-k spectra and BeH2/CH4 state-following references remain unchanged. Keep existing thresholds and timing references. Verification: cmake --build build -j8 passed. Run all eight tests/08_RI/lr_* integration cases with OMP_NUM_THREADS=1, MKL_NUM_THREADS=1 and mpirun -np 4, using LibXC 5.1.7 with KXC as in CI; collect with catch_properties.sh. All 28 reference checks pass. Repeat unrestricted HF on four ranks (same force sum) and one rank (3.12e-7 deviation, below 1e-6). ABACUS v3.11.0-beta10.
|
If this PR is ready for review again, let me know. @maki49 |
Yes, ready. |
There was a problem hiding this comment.
Let's keep this file as simple as possible, we don't have exx, ofdft, sdft or rt-tddft information here, so I think the codes related to LR-TDDFT can be moved to somewhere else?
There was a problem hiding this comment.
how about 'exited_state_check' instead of state_check?
There was a problem hiding this comment.
this file name still contains 'root'?
There was a problem hiding this comment.
this file name contains 'root' ?
|
Did you introduce KXC and LXC in this PR? If so, please make sure to mention it. I noticed that compiling KXC increases the overall CI/CD build time. |
Yes, kxc is necessary for lda/pbe LR force and relax tests. lxc is not introduced. I've mentioned it in the description. |
# Conflicts: # source/source_lcao/module_lr/utils/test/CMakeLists.txt # source/source_lcao/module_ri/test/CMakeLists.txt # source/source_lcao/module_ri/unittests/test_exx_lr_ws.cpp
Features
Refactors
Other Behavior changes
1e-6to1e-14, and the gradient cutoff from1e-10to1e-20. This can change excitation energies even when forces are disabled. A four-rank A/B run reproduces the LDA gamma reference changes. The periodic Si multi-k case has a minimum spin density of approximately1.29e-3, above the old cutoff; its four excitation energies agree at 17-digit output precision between the two settings.nspin=1singlet response uses the factor-two Hartree/XC contribution, while ground-stateS1_gsretains weight one. Existing single-spin excitation results can therefore change. BSE grid initialization usesS1_gsbecauseHamiltBSEalready applies the singlet factor; the existing RI BSE path is unaffected by this grid correction.1e-6when evaluating XC derivatives. This can change HSE06 response energies and gradients near vanishing density gradients; it does not apply this clamp to ordinary PBE.kxc, rather than allowing the last functional to overwrite earlier contributions. This affects the XC contribution to analytic excited-state forces.lr_nvirtcontinues to select the response basis. The gradient can depend on the available ground-state band window even when the response basis is held fixed.test_forceruns.