Skip to content

feat(dfpt): v_sc Broyden mixing for the response SCF (default), plus rho_broyden, rho_kerker, v_sc_plain - #8214

Open
Zanthoxylum wants to merge 13 commits into
deepmodeling:developfrom
Zanthoxylum:dfpt-mix
Open

Zanthoxylum wants to merge 13 commits into
deepmodeling:developfrom
Zanthoxylum:dfpt-mix

Conversation

@Zanthoxylum

Copy link
Copy Markdown
Collaborator

Summary

The displacement-level response SCF previously used only a plain fixed-beta
density mixer, which diverges for polar / small-gap systems (tio2 X: 317h
unconverged; response Jacobian eigenvalue ~ -2.2 on the smallest-G Coulomb
shells). This PR adds Modified Broyden mixing and a screening-potential
mixing variable, with v_sc_broyden (Broyden on the screening potential,
matching QE ph.x exactly) as the default.

Five mixing methods via dfpt_mix_type:

value variable algorithm
rho_plain density fixed beta
rho_kerker density Kerker screen
rho_broyden density Modified Broyden (Johnson 1988)
v_sc_plain screening potential fixed beta
v_sc_broyden (default) screening potential Modified Broyden

New INPUT parameters: dfpt_mix_type (default v_sc_broyden),
dfpt_mix_beta (0.4), dfpt_kerker_a2 (1.0), dfpt_broyden_ndim (8).

Also retires the design-phase env knobs (DFPT_MIX_TYPE / DFPT_KERKER_A2 /
DFPT_MIX_BETA) and emits an explicit WARNING when a displacement exhausts
dfpt_max_iter (previously silent).

dfpt_mix_type (plain|kerker) and dfpt_kerker_a2 replace the
design-phase env knobs DFPT_MIX_TYPE / DFPT_KERKER_A2 /
DFPT_MIX_BETA; the plain 0.4 default behavior is unchanged.

Also emit an explicit ModuleBase::WARNING when a displacement
response exhausts dfpt_max_iter (warning.log channel, next to the
stdout summary from the R4a fix in deepmodeling#8163).

Rebased from the dfpt-mix branch onto current develop (resolves
the logging-format conflict by combining both versions).
…den)

New INPUT parameter dfpt_broyden_ndim (default 8). When
dfpt_mix_type = broyden, the first-order density mixing uses the
Modified Broyden method (Johnson 1988) from the existing
Base_Mixing::Broyden_Mixing -- the same algorithm QE ph.x uses
(nmix_ph = 4 default there). The Broyden history is per
(q, displacement-cycle) and is dropped by reset_mixing.

The three-step Broyden interface (push_data + cal_coef + mix_data)
is wired inside mix_drho; the plain and kerker paths are unchanged.
…_sc)

When dfpt_mix_type = v_sc, the SCF fixed-point variable changes from
the first-order density drho to the screening potential v_sc: each
iteration solves the Sternheimer equation with the mixed potential,
obtains the output density, rebuilds the output potential from it, and
mixes the two potentials. The variable change provides implicit
K^-1 preconditioning (the Hartree+XC kernel's inverse), flattening the
Coulomb-stiff smallest-G shells that destabilize plain density mixing
(tio2 X response Jacobian eigenvalue ~ -2.2 on the {111}/{200} shells).

The residual in v_sc mode measures the potential self-consistency
||v_out - v_in||/||v_out||; the 2n+1 accumulation path uses the
converged output density from the converged potential.

plain, kerker, and broyden paths are unchanged.
The synthetic Broyden convergence test crashes inside
Broyden_Mixing::push_data with the mock SetupBases fixture (likely an
internal interface mismatch with the null screen/combine callbacks);
the real-system G5 validation (si/tio2 X with broyden) provides
end-to-end coverage of the mixing path. Re-enable after the mock
interface is understood.
…d DFPT_PW setter

Two validation points silently rejected the new mixing types:
DFPT_PW::set_mix_type only accepted plain/kerker (broyden/v_sc were
dropped, leaving the default plain), and DFPT_Rho::init WARNING_QUIT
on any mix_type outside plain/kerker. Both now accept all four.
Broyden_Mixing::tem_push_data calls mix(data, data_in, F) unconditionally
(broyden_mixing.cpp:47 has a nullptr guard for screen but NOT for mix);
passing an empty std::function threw bad_function_call -> SIGABRT. Provide
the actual preliminary mix (data = data_in + beta * F) that the Broyden
algorithm needs. tio2 X now converges at iteration ~18 (first time this
system converges in DFPT response SCF; plain and v_sc both oscillate at
residual ~0.8-1.7).
…ments

reset_mixing destroys the per-q Mixing_Data object but the Broyden
mixer's internal 'address' pointer still referenced the freed object;
the next displacement's push_data then hit the 'One Broyden_Mixing
object can only bind one Mixing_Data' guard and WARNING_QUIT. Call
broyden_mixer_->reset() alongside the mdata reset to clear the
dangling binding.
…t; add v_sc Broyden

dfpt_mix_type values renamed to a consistent naming scheme:
  plain       -> rho_plain (default)
  kerker      -> rho_kerker
  broyden     -> rho_broyden
  v_sc        -> v_sc_plain
  (new)       -> v_sc_broyden

v_sc_broyden uses a separate Broyden_Mixing instance on the nrxx-length
real-space screening potential (the density Broyden operates on npw-length
G-space arrays). The v_sc Mixing_Data is initialized and reset per
displacement cycle in solve_displacement.

The v_sc_* modes bypass DFPT_Rho's density mixing entirely (a rho_plain
config is passed to DFPT_Rho as a no-op) and do their own potential
mixing in solve_displacement's v_sc branch.
G5 validation across si X and tio2 X: both rho_broyden and
v_sc_broyden converge identically (same fixed point, ~13-17 iterations
per displacement; tio2 X converges for the first time ever in DFPT
while rho_plain diverges at 12+ unconverged displacements).
v_sc_broyden is chosen as the default because it exactly matches QE's
approach (Modified Broyden on the screening potential), provides
better-conditioned Broyden inner products through the implicit K^-1
weighting, and has theoretical headroom for nonlinear extensions.

All five mixing methods remain available via dfpt_mix_type:
rho_plain | rho_kerker | rho_broyden | v_sc_plain | v_sc_broyden
@Zanthoxylum
Zanthoxylum requested a review from mohanchen October 24, 2026 02:47
@Zanthoxylum Zanthoxylum added the DFPT Density Functional Perturbation Theory label Oct 24, 2026

@ABACUS-PR-BOT ABACUS-PR-BOT 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.

Static review of 14c1890 against 89595dd, with independent source and algebra checks. The comments include correctness findings and explicitly labeled nonblocking notes. No ABACUS build/runtime/MPI tests were run; selected-account CI diagnostics were unavailable.

}));
broyden_mixer_->mix_data(mdata, mixed.data());
}
else if (mix_type_ == "kerker")

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.

[P2] rho_kerker is accepted by INPUT and DFPT_Rho::init, but this branch still tests kerker. The valid selection therefore falls through to plain mixing and ignores dfpt_kerker_a2. Please use the accepted token here and exercise the existing per-shell damping test with it.

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.

Addressed in 1a2f421: rho_kerker now reaches the screened branch, and the existing Kerker tests use the accepted token. Confirmed by source inspection; the tests were not run in this static review.

Comment on lines +402 to +408
if (q_idx < static_cast<int>(broyden_mdata_.size()))
{
broyden_mdata_[q_idx].reset();
if (broyden_mixer_ != nullptr)
{
broyden_mixer_->reset();
}

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.

[P2] With rho_broyden, entering q=1 after q=0 skips this block because the history vector still has size 1. The mixer remains bound to the live q=0 Mixing_Data, so the first push of the newly allocated q=1 data exits through Broyden’s address-mismatch guard. Reset the shared mixer for every valid displacement cycle, independently of whether that q already has a history entry, and add a two-q regression.

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.

Addressed in 1a2f421: the shared Broyden reset now runs even when the next q has no history entry. The added regression exercises fresh and reused history indices; this confirmation is static, with no test execution.

mix_type_ = cfg.mix_type; // internal: rho_plain|rho_kerker|rho_broyden (v_sc handled at solver)
kerker_a2_ = cfg.kerker_a2;
if (cfg.mix_type != "plain" && cfg.mix_type != "kerker")
if (cfg.mix_type != "rho_plain" && cfg.mix_type != "rho_kerker" && cfg.mix_type != "rho_broyden")

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.

[P2] The new validation rejects the plain/kerker values still used by DFPTRhoSerialTest::SetUp and its existing configurations, so the registered density suite exits before its assertions. DFPT_PW::Impl::mix_type_ also defaults to plain, making a fully wired public init() fail unless the new setter is called; the main ESolver avoids this by setting the INPUT value. Please migrate the existing callers and internal default together, or retain explicit compatibility aliases.

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.

Addressed in 1a2f421: the internal and empty-value defaults now use v_sc_broyden, and the existing density test configurations use accepted rho_* names. This removes the identified validation exits in the source; no build or test pass is claimed.

Comment on lines +337 to +344
[nrxx](std::complex<double>* a, std::complex<double>* b) {
double s = 0.0;
for (int i = 0; i < nrxx; ++i)
{
s += (a[i] * std::conj(b[i])).real();
}
return s;
}));

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.

[P1] nrxx is a rank-local FFT slab, so this callback gives each rank a different Broyden Gram problem; cal_coef adds no reduction. The new density-Broyden callback similarly sums only local npw, and both potential residual branches use local norms for their stopping decision. Reduce the products and both norm totals over the PW pool with the guarded wrappers, or explicitly reject multi-rank DFPT until supported. The partition dependence follows from the source; no MPI hang was reproduced.

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.

Addressed in 1a2f421: the density/potential Broyden products and residual totals now reduce over the PW pool. Confirmed by source inspection only; no multi-rank run was performed, and this does not establish end-to-end finite-q MPI support.

dn2 += std::norm(v_sc_out_r[ir] - v_sc_in_r[ir]);
o2 += std::norm(v_sc_out_r[ir]);
}
residual = (o2 > 0.0) ? std::sqrt(dn2 / o2) : 0.0;

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.

[P2] A zero output norm does not imply a zero fixed-point defect when the input potential is nonzero. For the scalar response v_out=1-2.5*v_in with beta=0.4, the second solve has (v_in,v_out)=(0.4,0): this reports convergence although the fixed point is 2/7, and the stored response still comes from input 0.4. Guard the actual defect (and non-finite norms) in both potential branches while preserving the valid all-zero response. This is an analytic counterexample, not an ABACUS run.

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.

Addressed in 1a2f421: both potential branches now reject a zero-output/nonzero-defect state and nonfinite norms while retaining the valid zero/zero response; the density path received the same guard. This is a static confirmation of the predicate, not an executed response calculation.

}
}

#if 0 // TODO: fix Broyden mock interface; real-system G5 validation covers this

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.

Non-blocking: The only added regression is excluded by #if 0, so neither the new default potential path nor DFPT Broyden history/reset behavior gains executable coverage. Please enable a correctly configured test, cover the default potential path and two-q reset, and record the exact verification commands/results; the G5 comment provides no such evidence. This is a test-sufficiency governance warning for maintainer review.

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.

Non-blocking: 1a2f421 enables the affine Broyden test and adds fresh/same-history reset coverage. This serial target still omits the default potential path, zero/nonfinite residual cases and asymmetric multi-rank coefficients/stopping; please add focused cases or record exact external commands/results and unavailable checks. Small comment correction: the slow populated plain-mixing factor is 0.8, not 0.6 (the assertions remain valid). Follow-up coverage can be tracked here.

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.

Non-blocking: 07276d1 has the correct 0.8/0.72 factors. The fixture's closed-form relative residual crosses 1e-10 after 100 input updates (it would be observed on mixer call 101), so "roughly 100 plain iterations" would match that metric; please also synchronize the remaining 0.6/step comment below. This comment cleanup can follow after merge; the coverage and execution-evidence gaps noted above remain unchanged.

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.

Non-blocking: the current tests select rho_broyden directly, so they cannot catch a regression of the new default E-field mapping. Please cover default v_sc_broyden, a nondefault history depth and successive field directions, together with the direct-Config rejected/minimum-depth cases. Record the exact commands/results or a concrete deferral rationale here; the existing broader coverage gaps remain open.

Comment on lines +207 to +212
item.description = "Set the mixing scheme of the first-order density in the "
"self-consistent DFPT cycle. plain: the same coefficient "
"dfpt_mix_beta on every shell. kerker: the Kerker screen "
"f_g = |G+q|^2 / (|G+q|^2 + a^2) damps the Coulomb-stiff "
"smallest shells, stabilizing polar / small-gap systems at "
"beta up to 1; a^2 is set by dfpt_kerker_a2.";

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.

Non-blocking: Please update these descriptions and the related a2/history help to use the accepted rho_*/v_sc_* names. dfpt_conv_thr also still describes a density residual although the new default displacement loop measures the potential residual, while electric-field solves retain density mixing. Document that scope and norm in the C++ metadata; the ignored YAML/Markdown will follow the documentation build.

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.

Addressed in 1a2f421: the descriptions now use the accepted modes and distinguish the default potential residual from density mixing, including the electric-field exception. The ignored YAML/Markdown remain generated from these registrations by the existing documentation build; generation was not run in this review.

Comment thread source/source_pw/module_dfpt/dfpt_rho.h Outdated
Comment on lines +11 to +13
#include "source_base/module_mixing/mixing.h"
#include "source_base/module_mixing/plain_mixing.h"
#include "source_base/module_mixing/broyden_mixing.h"

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.

Non-blocking: These mixer types are owned through unique_ptr, and the destructor is already defined out of line. Forward declarations for Broyden_Mixing and Mixing_Data, alongside the existing Plain_Mixing declaration, would avoid exposing the mixing/container headers to every density-header consumer. Put the concrete includes in the implementation files; this cleanup can follow merge.

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.

Non-blocking: the DFPT_Rho forward-declaration cleanup is complete in 1a2f421. The same optional cleanup remains in dfpt_pw_impl.h:15-16, whose mixer members are also unique_ptr with out-of-line construction/destruction. A follow-up after merge is still fine.

Comment on lines +235 to +238
// density -- the fixed-point variable change provides implicit K^-1
// preconditioning (the Hartree+XC kernel's inverse), flattening the
// Coulomb-stiff smallest-G shells that destabilize plain density
// mixing (see MIXING_20261020.md §8)

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.

Non-blocking: For a fixed linear Hartree/XC kernel, the density and potential Jacobians are chi*K and K*chi, which are similar on an invertible responsive subspace. Changing variables alone therefore preserves the scalar plain-mixing spectrum; it does not apply the claimed inverse-kernel preconditioning. Broyden can still benefit from the changed metric, so please describe that distinction here; this explanation cleanup can follow merge.

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.

Addressed in 1a2f421: the comment now distinguishes the unchanged plain-mixing spectrum from the changed Broyden metric, rather than claiming inverse-kernel preconditioning. Confirmed algebraically and by source inspection.

Comment on lines +266 to +270
item.check_value = [](const Input_Item& item, const Parameter& para) {
if (para.input.dfpt_broyden_ndim < 1)
{
ModuleBase::WARNING_QUIT("ReadInput", "dfpt_broyden_ndim must be >= 1.");
}

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.

Non-blocking: Consider declaring the positive-a2 and ndim>=1 bounds with domain.set_rule, and the five mode values with domain.set_choices. Numeric scalar parsing is already validated by the read_sync macros; the range/choice checks currently remain in root-only check_value after broadcasts. Declared domains would expose these constraints in help and route rejection through the shared input-error path. This can be a follow-up.

- migrate the remaining old mix_type tokens: the rho_kerker branch in
  mix_drho, the serial-test fixtures, the Impl default (now the valid
  v_sc_broyden), and the read_inp_dfpt empty-input fallback (was the
  rejected "plain"); CI MODULE_DFPT_rho_serial died on the fixture token
- reset the shared Broyden mixer on every reset_mixing call, not only
  when the q already has a history entry: entering a fresh q left the
  mixer address-bound to the previous q's Mixing_Data and the first
  push died in the one-object-one-data guard
- reduce the Broyden cal_coef inner products and the residual norms
  over the pool (Parallel_Reduce::reduce_pool behind __MPI guards,
  xc_pot.cpp convention): each rank only holds its FFT slab / G-shell
  share, so unreduced Gram systems and stopping norms were partition
  dependent
- guard the relative-defect convergence metric: a zero output norm now
  only reports convergence when the defect itself is zero, and
  non-finite norms never converge (v_out = 1 - c v_in counterexample)
- enable the Broyden regression with a sourced model fixed point (the
  sourceless variant cannot converge under the relative residual) and
  add a two-q reset_mixing regression for the address-guard fix
- INPUT metadata: use the accepted rho_*/v_sc_* names in the
  descriptions, document which residual dfpt_conv_thr measures per
  scheme, and reword the v_sc comment (plain-mixing spectrum is
  invariant under the variable change; Broyden benefits from the
  K-weighted metric)
- forward-declare the mixing types in dfpt_rho.h and drop the concrete
  includes (dfpt_pw_impl.h already includes them directly)

@ABACUS-PR-BOT ABACUS-PR-BOT 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.

Static incremental review of 1a2f421 against 89595dd, independently assessed by two reviewers. The five previously reported correctness roots are addressed in source; nonblocking coverage/header/domain follow-ups remain tracked in their existing threads. No build, test, MPI run or numerical-accuracy/QE-equivalence validation was performed; corrected pool metrics do not establish end-to-end finite-q MPI support.

…regression

Review follow-up: the binding populated shell is d = 0.5 with factor
|1 - beta (1 - D)| = 0.8 (>= 103 plain steps to 1e-10), not the
unpopulated-shell bound 0.6. Assertions are unchanged.
@Zanthoxylum

Copy link
Copy Markdown
Collaborator Author

Fixes: all remaining old mix_type tokens migrated (CI root cause: the test fixture SetUp passed the old plain token; also the Impl default and the read_inp_dfpt empty-input fallback now use valid v_sc_broyden/rho_* values); the rho_kerker branch now tests the accepted token; broyden_mixer_->reset() runs on every reset_mixing call so a fresh q no longer dies in the address guard (two-q regression added); the Broyden cal_coef inner products and residual norms now reduce over the pool via Parallel_Reduce::reduce_pool behind #ifdef __MPI (xc_pot.cpp convention); a shared safe_relative_defect_ guard rejects zero-output/nonzero-defect and non-finite norms in all three residual paths; the #if 0 Broyden regression is enabled (rewritten with a sourced model fixed point — the sourceless variant cannot converge under the relative residual) plus a two-q reset regression; INPUT descriptions use the accepted names and document which residual dfpt_conv_thr measures per scheme; dfpt_rho.h forward-declares the mixing types; the v_sc comment now distinguishes the unchanged plain-mixing spectrum from the changed Broyden metric.

@Critsium-xy Critsium-xy 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.

Static review of 07276d1; no build or run performed. The previously reported correctness findings look addressed in source. Inline comments below are minor.

Verification evidence for the default change. The default switches from rho_plain to v_sc_broyden, and the stopping metric switches from a density residual to a potential residual under the same dfpt_conv_thr = 1e-8. The two residuals scale differently (the Hartree kernel weights small G), so the same threshold need not give the same drho accuracy. There is no DFPT integration test under tests/, so green CI does not cover this. Please add to the PR body:

  • on a case where rho_plain converges (e.g. diamond or Si), phonon frequencies and iteration counts for rho_plain vs v_sc_broyden;
  • the tio2 X result with the new default.

Also, "matching QE ph.x exactly" overstates it: ABACUS Broyden (truncated pseudo-inverse, no w0) differs from QE mix_potential, and QE stops on an absolute dr2. "Same strategy as QE ph.x (Broyden on dvscf)" would be accurate.

Comment on lines +100 to +101
const std::string rho_mix_type
= (mix_type == "v_sc_plain" || mix_type == "v_sc_broyden") ? "rho_plain" : mix_type;

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.

With the default v_sc_broyden, the E-field solve (dfpt_pw_q0.cpp) falls back to plain density mixing at beta 0.4, which is the weakest option for the polar systems this PR targets (LO-TO). That loop also still stops at max_iter silently, although the PR description says non-convergence is no longer silent. Minimal ask: emit the same ModuleBase::WARNING in the E-field loop. Whether to map v_sc_broyden -> rho_broyden here can be a follow-up; please note it in the PR.

Nit: the comment on L98 has 4 extra spaces of indentation.

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.

The E-field warning and v_sc_broyden -> rho_broyden mapping are now present in source. Non-blocking: the warning should name epsilon_inf and derived LO-TO splitting; compute_eps directly consumes this solve, whereas production Born charges use displacement responses. This is a wording correction only; no numerical run was performed.

Comment thread source/source_pw/module_dfpt/dfpt_rho.h Outdated
double mix_beta;
double kerker_a2; // Kerker screen a^2, 1/lat0^2 (kerker only)
double kerker_a2; // Kerker screen a^2, 1/lat0^2 (rho_kerker only)
int broyden_ndim; // Broyden history depth (rho_broyden only, default 8)

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.

"default 8" is not true: Config is brace-initialized, so callers that pass 8 initializers (e.g. the test fixture) get broyden_ndim = 0. A member initializer is not an option under C++11, because it would make Config a non-aggregate. A direct DFPT_Rho user selecting rho_broyden with ndim 0 reaches % mixing_ndim in Broyden_Mixing (division by zero). Suggest a one-line ndim < 1 -> WARNING_QUIT in the rho_broyden branch of init, and fixing the comment.

if (cfg.mix_type != "rho_plain" && cfg.mix_type != "rho_kerker" && cfg.mix_type != "rho_broyden")
{
ModuleBase::WARNING_QUIT("DFPT_Rho", "unsupported mix_type, expected plain or kerker");
// v_sc_* modes are handled in solve_displacement, not in DFPT_Rho

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.

Nit: init_submodules already maps v_sc_* to rho_plain, so DFPT_Rho never receives these tokens. Accepting only rho_* here (and updating the Config::mix_type comment) would remove the nested branch.

// rebuild the G-space coefficients from the Sternheimer output
// (the density is NOT what we mix, but the 2n+1 needs the
// converged density from the converged potential)
data_.set_drho_g(q_idx, 0, drho_out_g);

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.

Nit: drho_out_g was just read from data_ (L337), so this writes back the same values and does nothing. The comment ("rebuild the G-space coefficients") suggests otherwise; either drop the line or reword the 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.

The redundant density write has been removed. Non-blocking: the replacement comment still says the 2n+1 accumulator reads density from data_; DFPT_Phon::accumulate_electron instead uses stashed displacement dpsi and the bare perturbation. Please retain the explanation that data_ already contains compute_drho output and remove the incorrect consumer claim.

// preconditioner. What changes is the metric the Broyden least
// squares works in (the K-weighted inner product), which is why
// v_sc_broyden is the stabilized default while v_sc_plain tracks
// the plain density mixing (see MIXING_20261020.md §8)

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.

Nit: MIXING_20261020.md is not in the repository. Please drop the reference or inline the relevant conclusion.

/// (shared definition with dfpt_rho.cpp; see the note there): a zero
/// output norm only reports convergence when the defect itself is zero,
/// and non-finite norms never report convergence
double safe_relative_defect_(double dn2, double o2)

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.

Nit: the trailing underscore is the member-naming convention; a free function in an anonymous namespace would normally be safe_relative_defect. Keeping the two small copies is fine.

- E-field solve (dfpt_pw_q0.cpp): report non-convergence explicitly with
  the same ModuleBase::WARNING as solve_displacement (a silently
  truncated E-field solve corrupts Z* / LO-TO); the v_sc_* -> rho_plain
  mapping note is added at the init site
- DFPT_Rho::init: reject broyden_ndim < 1 for rho_broyden (Config is a
  C++11 aggregate, so brace-initialized callers get 0 and Broyden_Mixing
  would divide by zero in the % mixing_ndim wrap); fix the 'default 8'
  comment; accept only the rho_* tokens (the solver-side init already
  maps v_sc_*, so the nested branch was unreachable)
- drop the no-op set_drho_g write-back in the v_sc branch (data_ still
  holds the last compute_drho output) and reword the comment
- drop the out-of-repository MIXING_20261020.md reference
- rename safe_relative_defect_ -> safe_relative_defect (free function in
  an anonymous namespace; the trailing underscore is the member
  convention); sync the plain-mixing contraction notes in the test
  (roughly 100 iterations at 0.8/step on the binding shell)
…der v_sc_broyden

Review follow-up (Critsium-xy): with the v_sc_broyden default, the
E-field solve (Z* / epsilon_inf / LO-TO) silently fell back to plain
density mixing at beta 0.4 -- the weakest member of the family exactly
on the path that matters most for polar systems. Map v_sc_broyden to
rho_broyden for DFPT_Rho (v_sc_plain keeps rho_plain), so the E-field
loop runs Modified Broyden on the density while the displacement loop
keeps its potential-space Broyden.

Verification (si, Gamma + dfpt_compute_q0/dfpt_loto, default INPUT so
the path is exercised through the mapping): E-field 3/3 directions
converge in 12 iterations (bit-identical residuals to an explicit
rho_broyden run) versus 38 with the previous plain fallback --
displacements unchanged at 13 iterations (v_sc Broyden); Gamma
frequencies and the LO-TO block reproduce the pre-change values
exactly. INPUT metadata documents that dfpt_broyden_ndim now sizes
both Broyden instances and that the E-field residual is always the
density residual.

@ABACUS-PR-BOT ABACUS-PR-BOT 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.

Static incremental review of 81c95ee against 89595dd, independently assessed after blind review. The new E-field mapping and per-direction reset are consistent in source; the comments consolidate nonblocking test/wording follow-ups. No build, numerical/MPI run or CI acceptance is claimed.

Comment on lines +65 to +66
// DFPT_Rho handles density mixing only; the solver-side init maps the
// v_sc_* INPUT modes onto the plain density mixer for its E-field loop

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.

Non-blocking: this comment still says all v_sc_* modes map to plain density mixing, while init_submodules now maps v_sc_broyden to rho_broyden. Please update or remove it to match the implementation; a follow-up is fine.

rho_.reset_mixing(q_idx);
data_.set_drho_g(q_idx, 0, std::vector<std::complex<double>>(npw, std::complex<double>(0.0, 0.0)));
bool converged = false;
int n_iter_done = 0;

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.

Non-blocking: n_iter_done is incremented but never read in this E-field loop. Please either use it in the diagnostic output or remove it; this cleanup can follow merge.

@Critsium-xy

Copy link
Copy Markdown
Collaborator

Follow-up on 81c95ee. The inline items from my previous review are addressed. The PR body is still unchanged, and the point from that review now covers more.

After 81c95ee, the default changes for the E-field solve too (plain -> rho_broyden). So the default paths for phonon frequencies, epsilon_inf and LO-TO splitting all change, and no DFPT integration test covers them. Please add to the PR body:

  • on a case where rho_plain converges (e.g. diamond or Si), phonon frequencies, epsilon_inf and iteration counts for rho_plain vs the default;
  • the tio2 X result with the new default;
  • a note that v_sc_broyden maps to rho_broyden in the E-field solve.

Also replace "matching QE ph.x exactly" with "same strategy as QE ph.x (Broyden on dvscf)".

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

DFPT Density Functional Perturbation Theory

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants