Skip to content

DFT+U refactor, step 8 - #7912

Open
mohanchen wants to merge 33 commits into
deepmodeling:developfrom
mohanchen:2026-09-06-line2
Open

DFT+U refactor, step 8#7912
mohanchen wants to merge 33 commits into
deepmodeling:developfrom
mohanchen:2026-09-06-line2

Conversation

@mohanchen

@mohanchen mohanchen commented Sep 6, 2026

Copy link
Copy Markdown
Collaborator

Summary
This PR refactors the DFT+U module to reduce cross-layer coupling, clarify the PW/LCAO split, and clean up legacy interfaces — while preserving all existing physics and INPUT behavior.

Architectural refactoring: decouple DFT+U from concrete types
Hold Plus_U_Base* (instead of the concrete DFTU type) in OperatorDFTU, ESolver_KS, HamiltLCAO, force_stress, and spar_u interfaces, demoting ~10 public APIs (pot_onsite, pot_uterm_, cal_energy_correction, cal_occ_mat_, getForceStress, output_mat_sparse, ctrl_scf_lcao, LCAO_domain::set_pot, etc.) from concrete-DFTU parameters to base-class pointers.
Make Plus_U_Base destructor virtual to enable safe polymorphic use.
Extract cal_occ_pw from Plus_U_Base into a free function; remove the Plus_U_Base::nspin member variable (mutable workflow state).
PW-side reorganization
Split dftu_pw into dftu_pw_tools (pure functions, formerly dftu_base_tools) and dftu_pw (orchestration), deleting the obsolete dftu_base_occ.cpp.
Move the corresponding unit tests to source_pw/module_pwdft/test/ and remove the #define private public hack from test_dftu_pw_tools.cpp.
Naming and semantics cleanup
Rename cal_type to a UForm enum with explicit formalism semantics.
Rename occupancy-matrix state accessors to ready/stale semantics.
Build system
Compile module_dftu sources only under ENABLE_LCAO, and drop now-redundant #ifdef _LCAO guards in module_dftu headers/sources.
Guard MODULE_IO_numerical_basis_test behind ENABLE_LCAO.
Fix a bug in Makefile.Objects.
Bug fixes
Restore the INPUT parameter name orbital_corr (the internal C++ rename to l_channel had leaked into the user-facing INPUT), and fix an LCAO dftu
null-pointer dereference.
Sync docs/parameters.yaml with the C++ source for l_channel.
Other cleanups
Pass Input_para explicitly to ctrl_output_fp (reducing global-state reads).
Update occ_matrix/occ_mixer and their tests to match the new accessors.
Verification: unit tests updated for occ_mixer, dftu_base, dftu_pw_tools, and read_input; no change to numerical results is intended

abacus_fixer added 2 commits September 6, 2026 11:04
- Rename member variable `orbital_corr` to `l_channel` in Plus_U_Base
- Rename member variable `orbital_corr_` to `l_channel_` in YukawaScreening and OccMatMixer
- Rename accessor `get_orbital_corr(int)` to `get_l_channel(int)`
- Rename accessor `get_orbital_corr_vec()` to `get_l_channel_vec()`
- Rename predicate `has_correlated_orbital(int)` to `has_l_channel(int)`
- Rename local variable `orbital_corr_tmp` to `l_channel_tmp` in onsite_proj_tools_force.cpp
- Update all call sites across PW, LCAO, IO, and test modules

The new name `l_channel` directly indicates the angular momentum quantum
number `l` of the channel receiving the +U correction, avoiding confusion
with "correlated orbital" (many-body physics concept) and the ambiguous
abbreviation `corr`.
@mohanchen mohanchen added DFT+U Issues related to DFT plus U function Refactor Refactor ABACUS codes The Absolute Zero Reduce the "entropy" of the code to 0 labels Sep 6, 2026
@mohanchen
mohanchen requested a review from lanshuyue September 6, 2026 06:48
abacus_fixer added 24 commits September 6, 2026 14:58
Move Plus_U_Base::cal_occ_pw to DFTU_BASE::cal_occ_pw in
dftu_base_occ.cpp, with all dependencies passed explicitly as
parameters. This decouples the PW occupation-matrix calculation
from the Plus_U_Base class, making it directly unit-testable.

Changes:
- dftu_base_tools.h: declare DFTU_BASE::cal_occ_pw with 15 params
- dftu_base_occ.cpp: implement as free function, remove member version
- dftu_base.h: remove cal_occ_pw declaration, add public accessors
  (get_nspin, get_device, get_u_current_vec, get_pot_uterm_pw_index,
  get_pot_uterm_pw, energy_ref) needed by callers
- setup_dftu_pw.cpp: call DFTU_BASE::cal_occ_pw via accessors
- dftu_base.cpp: update comment to reflect free-function status
Remove the mutable nspin member from Plus_U_Base to avoid hidden
workflow state. The value is now passed explicitly:

- init_base keeps its nspin parameter (needed for layout setup)
- get_pot_uterm_pw_spin / get_size_pot_uterm_pw_spin take nspin
  as a function argument
- Call sites pass PARAM.inp.nspin at the PW operator boundaries
  where signatures are fixed by the framework

Also add missing parameter.h include to onsite_proj_force_stress.cpp.
… semantics

Replace the opaque int cal_type (always 3, never written) in Plus_U_Base
with a type-safe enum class UForm that explicitly names the DFT+U
formalism and double-counting scheme:

- lich_fll = 1: Lichtenstein rotationally invariant + FLL DC
- lich_amf = 2: Lichtenstein rotationally invariant + AMF DC
- dud_fll  = 3: Dudarev simplified + FLL DC (default, implemented)

Changes:
- dftu_base.h: add UForm enum, replace member and getter
- dftu_nao_pots.cpp: switch dispatch to enum cases
- dftu_nao_pots.h: update doc comment
- dftu_core_test.cpp: update comments

The enum preserves the existing placeholder values for future
implementation of the Lichtenstein formalism paths.

Also rename the per-atom DFT+U effective-potential buffer from
pot_uterm_pw/pot_uterm_pw_index to uterm_mat/uterm_mat_index and the
related accessors (get_uterm_mat*, get_size_uterm_mat*). The old name
suggested a basis-expanded potential; it is actually the small (2l+1)^2
projector coefficient matrix, analogous to the nonlocal-pseudopotential
D_ij coefficients, and is basis-independent despite living in the
PW/LCAO-shared base class.
Update docs/parameters.yaml to use l_channel instead of orbital_corr,
matching the parameter name change in the C++ source (commit eabfabc).
Remove PARAM.inp dependency from ModuleIO::ctrl_output_fp by adding a
const Input_para& parameter, and update the single call site in
ESolver_FP::after_scf to forward *this->inp_.
Rename the DFT+U occupation-matrix state accessors to shorter, clearer names:
- is_occ_mat_initialized() -> is_occmat_ready()
- mark_occ_mat_initialized() -> set_occmat_ready()
- mark_occ_mat_dirty() -> set_occmat_stale()
- member occ_mat_initialized -> occmat_ready_

Also unify the base-class internal write paths to use set_occmat_ready()
instead of direct member assignment, and update all call sites and comments
across PW, LCAO, and DFPT modules.
… (orchestration)

- dftu_pw_tools.h/cpp: 7 pure functions in namespace pw, unit-testable
- dftu_pw.h/cpp: cal_occ_pw + accumulate_occ_one_k in namespace DFTU_BASE
- test_dftu_pw_tools.cpp: tests for pure functions, links only dftu_pw_tools.cpp
- Update CMakeLists.txt accordingly
- Replace PARAM.inp.nspin with local variable c.nspin in EnergyWeightsAllNspin
- Remove unnecessary parameter.h include
The three DFTU_LCAO free functions in dftu_nao_pots only use base-class
accessors (occmat, get_l_channel, get_form, yukawa, get_u_current), so
their dftu parameter can take Plus_U_Base& instead of the LCAO-derived
Plus_U&. Callers still pass Plus_U, which upconverts implicitly, so this
compiles standalone. This is step A1 of decoupling PW DFT+U from the
LCAO-derived type.
The four pot_uterm free functions only use the base-class
is_occmat_ready() plus the already-demoted pot_onsite_* helpers, so
their dftu parameter can take Plus_U_Base&. Callers pass Plus_U, which
upconverts implicitly. Step A2 of decoupling PW DFT+U from Plus_U.
cal_energy_correction only uses base-class accessors (is_occmat_ready,
occmat, get_l_channel, use_yukawa, yukawa, get_u_current, set_energy)
plus the already-demoted get_onsite_pot, so its dftu parameter can take
Plus_U_Base&. set_energy is mutable so the reference stays non-const.
Step A3 of decoupling PW DFT+U from Plus_U.
cal_occ_mat, cal_occ_mat_k and cal_occ_mat_gamma only use base-class
accessors (occmat, get_l_channel_vec, has_occ_mixer, occ_mixer,
is_occmat_ready, set_occmat_ready), so their dftu parameter can take
Plus_U_Base&. Step A4 of decoupling PW DFT+U from Plus_U.
force_stress only used Plus_U for get_orb_cutoff(); all other state goes
through base-class accessors (pot_onsite_*, get_l_channel_vec, occmat).
The per-type orbital cutoff is LCAO-specific (PW uses onsite_radius), so
pass it explicitly as a new orb_cutoff parameter instead of moving it to
the base class. The internal cal_force_k/cal_stress_k/cal_stress_gamma
already take orb_cutoff as a parameter. Caller passes orb.cutoffs().
Step A5 (B-beta) of decoupling PW DFT+U from Plus_U.
cal_HR_dftu and cal_HR_dftu_soc only forward dftu to the already-demoted
pot_uterm_HR_real/complex, so their parameter can take Plus_U_Base&.
Step B1 of decoupling PW DFT+U from Plus_U.
Both operator class templates now store and accept Plus_U_Base* instead
of the LCAO-derived Plus_U*. OperatorDFTU only uses base-class interfaces
(pot_uterm_*), so it demotes cleanly. DFTU still needs the LCAO-only
get_dmr() for the density matrix, so the five get_dmr() call sites use a
static_cast<Plus_U*> back to the derived type; the cast in cal_fs_nao_r
is hoisted to a single local at function entry.

This confines the base-to-derived cast to the operators' DM access, so
upstream callers (HamiltLCAO, getForceStress) can pass base-class
pointers. Step B (operators) of decoupling PW DFT+U from Plus_U.
The HamiltLCAO constructor only forwards p_dftu to the DFTU/OperatorDFTU
operator constructors, which now accept Plus_U_Base*. Demote the parameter
and switch the forward declaration from Plus_U to Plus_U_Base. Callers pass
&this->dftu (a Plus_U) which upcasts implicitly.
getForceStress only forwards dftu to the already-demoted force_stress and
passes &dftu to the DFTU operator constructor (now Plus_U_Base*), so its
parameter can take Plus_U_Base& with no cast needed.
The p_dftu parameter is unused inside output_mat_sparse (only forwarded
between the two overloads), so demote it to Plus_U_Base* across both
declarations, both definitions, and the four explicit instantiations.
Change the ctrl_scf_lcao template and its three explicit
specializations to take Plus_U_Base& instead of Plus_U&. The dftu
argument is only forwarded (to output_mat_sparse, already on
Plus_U_Base*), so no cast is needed. Switch the header include from
dftu_nao.h to dftu_base.h.

Verified: make -j 30 in build_max_para_test builds abacus_max_para.
set_pot only needs the LCAO-specific init(ucell, ..., &orb), which lives
on the derived Plus_U. Take Plus_U_Base& at the interface and cast to
Plus_U& at the single init call site, matching the existing set_dmr
pattern. The header now forward-declares Plus_U_Base instead of
including dftu_nao.h; the .cpp includes dftu_nao.h for the cast.

Verified: make -j 30 in build_max_para_test builds abacus_max_para.
Plus_U adds members (orb_cutoff_, dm_in_dftu_*) that require proper
destruction. The upcoming esolver change holds the object through a
Plus_U_Base* and deletes it in the base-class destructor, which is
undefined behavior without a virtual destructor. Make it virtual now as
a standalone, separately-compilable step.

Verified: make -j 30 in build_max_para_test builds abacus_max_para.
Replace the value member `Plus_U dftu` with a base-class pointer
`Plus_U_Base* dftu_` in ESolver_KS, so that the PW/LCAO layering no
longer requires the base class to know the LCAO-derived Plus_U type.

- esolver_ks.h: member -> Plus_U_Base* dftu_ = nullptr; include
  dftu_base.h instead of dftu_nao.h.
- esolver_ks.cpp: delete dftu_ in ~ESolver_KS (virtual dtor on
  Plus_U_Base makes this safe); update u_converged() call.
- esolver_ks_pw.cpp: allocate `new Plus_U_Base()` in ctor (PW basis).
- esolver_ks_lcao.cpp: allocate `new Plus_U()` in ctor (LCAO basis);
  include dftu_nao.h for the derived type.
- Adapt all call sites to the pointer/reference forms expected by the
  already-base-classified helpers (Hamilt*, setup_pot, chgmixing,
  ctrl_scf_lcao, init/finish_dftu_lcao void* API, force/stress/dfpt).

Derived classes (TDDFT, DoubleXC, DFPT_PW, lcao_others) inherit the
already-constructed dftu_ via the ESolver_KS_PW / ESolver_KS_LCAO
constructors, so each object is allocated exactly once.

Build: make -j 30 in build_max_para_test, abacus_max_para links OK.
All sources in source_lcao/module_dftu implement the LCAO-specific
Plus_U derived class and the DFTU_LCAO operators; the PW DFT+U
implementation lives in source_pw/module_pwdft/dftu_base.* and is
compiled unconditionally. Move the seven previously-unconditional
sources (dftu_nao, dftu_nao_fs_k, dftu_nao_folding, dftu_nao_pots,
dftu_nao_occ, dftu_nao_energy, dftu_hamilt) into the existing
if(ENABLE_LCAO) block, following the module_rt precedent.

Because the dftu OBJECT library now has no sources when LCAO is off,
wrap add_library/add_coverage/test in if(ENABLE_LCAO) and move the
`dftu` entry in the final link list from the unconditional block into
the if(ENABLE_LCAO) block in source/CMakeLists.txt, mirroring how
`tddft` is handled.

Verification:
- LCAO build (build_max_para_test): cmake reconfigure + make -j 30,
  abacus_max_para links OK.
- PW-only build: cmake -DENABLE_LCAO=OFF -DBUILD_TESTING=OFF,
  make -j 30 abacus_pw_para links OK with no dftu target.
abacus_fixer added 4 commits September 7, 2026 14:56
These seven .cpp files are now compiled only when ENABLE_LCAO is on
(see the previous CMakeLists change), so the #ifdef __LCAO guards
inside them are always true and can be removed. Only the preprocessor
lines are deleted; all function bodies are unchanged, and the
pre-existing #ifdef __MPI blocks (ScalapackConnector calls) are kept
intact.

The header files keep their #ifdef __LCAO guards: dftu_nao.h is still
included by PW-compiled code (esolver_ks.h) and must let Plus_U
degrade to an empty shell when LCAO is off.

Verification:
- LCAO build (build_max_para_test): make -j 30 abacus_max_para links OK.
- PW-only build (ENABLE_LCAO=OFF): abacus_pw_para links OK.

Note: code_quality_score flags dftu_nao_fs_k.cpp (6) and
dftu_nao_occ.cpp (49) for pre-existing debt (file length, parameter
count, cyclomatic complexity) untouched by this change.
These six headers are only included from LCAO-only translation units
(verified: no source_pw/ references), and the whole dftu module is now
compiled only under ENABLE_LCAO. The __LCAO guard inside them is dead
in both build flavors, so remove it; only preprocessor lines are dropped,
declarations are unchanged. dftu_nao.h keeps its guards because the
PW-shared esolver_ks.h still includes it.

Verified: abacus_max_para links (ENABLE_LCAO=ON) and abacus_pw_para
links (ENABLE_LCAO=OFF). The unrelated MODULE_IO_numerical_basis_test
fails to link under ENABLE_LCAO=OFF because its CMakeLists links
numerical_atomic_orbitals/orb outside if(ENABLE_LCAO) -- a pre-existing
issue independent of this change.
This test links the LCAO-only libraries numerical_atomic_orbitals and orb
and compiles source_lcao/center2orb.cpp, so it cannot build under
ENABLE_LCAO=OFF. Wrap its AddTest in if(ENABLE_LCAO), matching the
neighboring MODULE_IO_to_qo_test.

Verified: ENABLE_LCAO=OFF full build completes (abacus_pw_para links) and
the test still builds under ENABLE_LCAO=ON.
…u_ null pointer

The rename commit eabfabc changed the INPUT parameter from orbital_corr
to l_channel, but test cases (scf_u_spin2, scf_u_spin2_old, scf_u_yukawa)
still use the old name. Restore the INPUT parameter name to orbital_corr
while keeping the internal C++ member variable as l_channel.

Additionally, the refactor to base-class pointer (156a5d8) left
ESolver_KS_LCAO without initializing dftu_ (unlike PW path which news
Plus_U_Base), causing a segfault on DFT+U LCAO calculations. Add the
missing allocation in the ESolver_KS_LCAO constructor.

Verified: scf_u_spin2 passes with exit code 0, SCF converges,
final energy -11033.13 eV.

@lanshuyue lanshuyue left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

dftu_ is owning storage allocated with new in both ESolver_KS_PW and ESolver_KS_LCAO, but the current ESolver_KS destructor does not delete it. This also contradicts the member comment and commit 156a5d8, both of which state that dftu_ is deleted in ~ESolver_KS. Please fix the ownership before merging. Preferably make this a std::unique_ptr<Plus_U_Base> so the ownership is explicit.

abacus_fixer and others added 3 commits September 8, 2026 16:32
ESolver_KS::dftu_ was allocated with new in both ESolver_KS_PW and
ESolver_KS_LCAO constructors but never deleted in ~ESolver_KS,
contradicting the member comment and commit 156a5d8. Convert dftu_
to std::unique_ptr<Plus_U_Base> so ownership is explicit and the
object is released automatically on destruction.

Call sites that pass raw pointers to downstream functions use .get();
C++11 baseline requires reset(new ...) instead of std::make_unique.
Merge from develop reverted the raw pointer back to unique_ptr member
in init_dftu_lcao call site, causing compile error on void* parameter.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

DFT+U Issues related to DFT plus U function Refactor Refactor ABACUS codes The Absolute Zero Reduce the "entropy" of the code to 0

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants