Conversation
Add configurable linear, Pulay, and Anderson polarization mixing with adaptive controls and diagnostics. Use a periodic mass-weighted ionic center consistently for PCC2D moments, energies, potentials, and forces, and extend relax, INPUT, and focused unit-test coverage.
Use existing Parallel_Global lifecycle wrappers in new SCCS tests, name intermediate potential/adapter arguments, and forward-declare ChargeReduction in PCC operator headers. No documentation update required: INPUT behavior and numerical formulas are unchanged. Verified an isolated GNU build and all 23 surchem CTest targets with OMP_NUM_THREADS=1.
…s-develop into sccs-pcc-port-7972
pcc_2d fixed the open direction to Cartesian +y: the geometry demanded b along +y with a and c in the x-z plane, and the moments, potential, force and ionic center read y directly. Describe the slab instead by the open lattice vector (axis 0, 1 or 2) and its unit normal, which must be perpendicular to the two periodic vectors, and measure every coordinate along that normal. The moments and parameters lose their _y suffixes and the debug labels name the axis. SurchemParameters::pcc_2d_axis stays 1 here, so the geometry is the old one: all eight SCCS/PCC integration cases give identical energies and forces to the last printed digit.
Critsium-xy
reviewed
Sep 30, 2026
|
|
||
| - **Type**: Real | ||
| - **Availability**: *[`imp_sol`](#imp_sol)==2* | ||
| - **Description**: Low-pass filter of the SCCS switching-function derivatives, as Environ deriv_lowpass_p1 with deriv_method fft: when sccs_lowpass_p1 and sccs_lowpass_p2 are both positive, every Fourier derivative of the switching function is multiplied by 0.5 erfc(p1 G^2/Gcut^2 - p2), Gcut^2 being the ecutrho sphere, and the electronic potential becomes the exact derivative of the discrete SCCS energy, so forces agree with energy differences. Only with assume_isolated pcc_0d or pcc_2d. The default -1 turns it off and reproduces Environ deriv_method fft (continuum cavity potential). 10 with sccs_lowpass_p2 5 was validated at ecutrho 300-500 Ry; the filter changes the model energy (about 10 meV for H3O+). |
Collaborator
There was a problem hiding this comment.
Non-blocking documentation suggestion: could you explicitly state that, with lowpass disabled (the default), analytical forces may differ from finite differences of the self-consistent energy? For geometry optimization with PCC, please recommend considering lowpass and checking force accuracy, while retaining the existing note that it changes the energy. Please update the Input_Item description and regenerate the docs. No algorithm or default change is requested.
Author
There was a problem hiding this comment.
- The description now says that with lowpass off (the default), analytic forces may differ from finite differences of the self-consistent energy.
- It recommends considering lowpass for geometry optimization with PCC and checking the force accuracy, and it keeps the note that the filter changes the energy.
Add the INPUT key pcc_2d_axis (0, 1 or 2; default 2, like efield_dir) for the lattice vector along which assume_isolated=pcc_2d is open. The Gamma-only k-point check follows it. The existing pcc_2d integration cases and examples, which are open along y, now set pcc_2d_axis 1 and reproduce their previous energies and forces exactly. Tests constructing SurchemParameters for a y-open slab state the axis explicitly instead of relying on a struct default.
The rho and force symmetrization assume the Hamiltonian is invariant under every analyzed operation, but the PCC correction is not invariant under all operations of the structure. For pcc_2d the potential depends on the coordinate along the open axis, so an operation that mixes that lattice vector with a periodic one, or a fractional translation along it, would symmetrize the density with a symmetry the correction lacks. pcc_0d is likewise not invariant under fractional translations of a non-primitive cell. After the symmetry analysis, stop with a message that names the violation and suggests a primitive cell or symmetry 0/-1. Operations act on direct coordinates, so the pcc_2d test reads row and column pcc_2d_axis of each rotation directly. The eight SCCS/PCC integration cases are compatible and still pass; a lone atom in a cubic pcc_2d cell is now rejected.
A unit test builds the same standalone pcc_2d slab open along y and, after the cyclic relabeling (x, y, z) -> (z, x, y), open along z. The PCC energy agrees to 1e-11 relative and the forces agree component by component after the permutation. The integration case 222_PW_SCCS_PCC2D_AXIS_Z is case 214 relabeled the same way with pcc_2d_axis 2. Its energies and permuted forces agree with case 214 to within the SCF convergence; at scf_thr 1e-11 the energy agrees to 1e-8 eV and the forces to 3e-5 eV/Angstrom.
The assume_isolated description now explains why any lattice vector perpendicular to the periodic plane can be the pcc_2d open one, that a tilted open vector needs the equivalent perpendicular cell, that sampling along it must be Gamma only, and which symmetry operations stop the run for pcc_2d and pcc_0d. pcc_2d_axis refers to these requirements.
Critsium-xy
reviewed
Sep 30, 2026
| // chgmixing_ks broadcasts drho only later; decide the delayed SCCS start | ||
| // on the root value so every rank enters the SCCS reductions together. | ||
| double sccs_start_drho = this->drho; | ||
| Parallel_Common::bcast_double(sccs_start_drho); |
Collaborator
There was a problem hiding this comment.
Could you guard this broadcast so it runs only while SCCS is enabled and awaiting activation? It currently adds an MPI_COMM_WORLD broadcast to every iteration of ordinary non-SCCS calculations, and continues after SCCS has activated. The guard must be consistent across ranks; please retain synchronized activation, the mixing reset, and the convergence veto for that iteration.
As an optional readability improvement, an explicit convergence-eligibility flag would express the later scf_thr = -1.0 workaround more clearly. I am not requesting a broader SCF refactor or claiming a demonstrated numerical/performance regression.
Author
There was a problem hiding this comment.
The broadcast now runs only when solvent.uses_sccs() && !solvent.sccs_is_active(). Both terms come from the configuration and from the activation decision, which is itself taken on the broadcast value, so every rank evaluates the guard identically.
…heck v_correction_pcc checked that the open lattice vector is perpendicular to the periodic plane to a relative 1e-10, which rejects a rotated cell whose vectors were typed with six digits. Use 1e-6: a tilt of that size is physically irrelevant, and axis-aligned cells are unaffected.
Critsium-xy
reviewed
Sep 30, 2026
Critsium-xy
left a comment
Collaborator
There was a problem hiding this comment.
The inline comment addresses the new static PCC energy transfer.
Non-blocking follow-ups, which can be handled after this PR is merged: common/pw_grid.h can forward-declare Matrix3, and sccs/sccs_charge.h can forward-declare ChargeReduction, with the complete includes in their implementation files. Broader isolation of SCCS details from surchem.h is also optional: the current value members require complete types, so simply removing sccs_driver.h would not be sufficient. These header cleanups should not delay merging.
| static double Ael; | ||
| // PCC open-boundary correction energy (Ry), reported separately from | ||
| // the solvation terms Ael and Acav; zero without assume_isolated pcc_*. | ||
| static double Epcc; |
Collaborator
There was a problem hiding this comment.
Please avoid introducing a new mutable process-wide energy through static Epcc. The instance already stores pcc_energy_rydberg_, but ElecState::get_pcc_energy() reads this shared copy without identifying the owning solvent instance or checking its validity. Interleaved solver instances could therefore overwrite each other's energy. Could you expose the instance energy through a const accessor/result and pass it explicitly to the corresponding electronic state, covering every potential-update path used before energy evaluation? The existing Ael/Acav statics are historical debt and need not be refactored in this PR.
The pcc_2d monopole constant was ENVIRON's -pi*q/(3*L). The correction is the open planar kernel -2*pi*|u|/A minus the zero-mean periodic kernel; with the open kernel zero on the plane of the charge, the constant is -pi*q*L/(3*A), and the two agree only for A = L^2. The ENVIRON value adds 0.5*q^2*(pi*L/(3*A) - pi/(3*L)) to the energy of a charged slab, a term linear in the cell length. With the open planar constant, a charged slab's energy converges with the vacuum size. A new unit test shows that the periodic energy plus the PCC of a charged Gaussian layer equals its analytic open-boundary energy at L = 12, 20 and 32 bohr. For an H3O+ layer at L = 16..60 bohr, the energy now spreads by 0.23 meV in vacuum and 0.32 meV with SCCS. Before, it grew by 0.151 and 0.0019 eV/bohr, which is the extra term above. Neutral slabs and all forces are unaffected. Drop the warnings that charged pcc_2d energies at different cell lengths are not comparable. Document that charged-slab and solvation energies are referenced to zero potential on the plane of the charge and differ from ENVIRON unless A = L^2. Case 217_PW_SCCS_PCC2D_H3O is regenerated. Its E_pcc changes by -1.2665 eV, the predicted -0.046542 Ha for A = 144 and L = 20 bohr, and its total energy changes by -0.016 eV, that amount screened by the solvent (1/78.3). A finite-difference test uses a larger, still exact, step because the larger constant raised its roundoff.
Critsium-xy
reviewed
Sep 30, 2026
| { | ||
| if (!this->uses_pcc() || !this->pcc_result_valid_) | ||
| { | ||
| throw std::logic_error("correction summary requires a current SCCS or PCC result"); |
Collaborator
There was a problem hiding this comment.
Non-blocking robustness suggestion: write_iteration() calls this only on the output rank, with no surrounding exception handler. If the invariant fails, that rank terminates while peers may continue into MPI collectives, leaving job cleanup to the launcher. Could validation happen in the computation path with consistent failure handling across ranks, so printing does not introduce a rank-local fatal path? Output-rank-only WARNING_QUIT() also just exits that process. I have not identified a normal execution path that triggers this; this is not a reproduced deadlock.
Critsium-xy
reviewed
Sep 30, 2026
| 212_PW_USPP_BPCG | ||
| 213_PW_USPP_pchg_wfc | ||
| 214_PW_SCAN_NLCC | ||
| 214_PW_SCCS_PCC2D_CF |
Collaborator
There was a problem hiding this comment.
Non-blocking naming suggestion: 214 is already used by 214_PW_SCAN_NLCC immediately above. Could you rename this new case to an unused number (e.g. 220_PW_SCCS_PCC2D_CF) and update this list and any references? The runner uses full directory names, so this does not cause a test collision, but unique numbers would avoid ambiguous filtering and references.
Author
There was a problem hiding this comment.
Renamed to consecutive unused numbers 220–226
Critsium-xy
reviewed
Sep 30, 2026
Comment on lines
+35
to
+37
| positions[ir] = ModuleBase::Vector3<double>(fx * a1.x + fy * a2.x + fz * a3.x, | ||
| fx * a1.y + fy * a2.y + fz * a3.y, | ||
| fx * a1.z + fy * a2.z + fz * a3.z); |
Collaborator
There was a problem hiding this comment.
Please compute these coordinates in named local variables before constructing Vector3, to follow AGENTS.md rule 14. The same pattern appears in the new PCC helpers, e.g. scaled(add(add(a, b), c), 0.5) and the inline gradient components in pcc_0d.cpp. Please apply this convention to the new call sites, preserving the arithmetic order.
Author
There was a problem hiding this comment.
Done for all new call sites, with the arithmetic order preserved
Critsium-xy
reviewed
Sep 30, 2026
| #ifndef SURCHEM_PW_GRID_H | ||
| #define SURCHEM_PW_GRID_H | ||
|
|
||
| #include "source_base/matrix3.h" |
Collaborator
There was a problem hiding this comment.
Non-blocking dependency cleanup (fine as a follow-up after merge): Matrix3 is only used by reference here, so it can be forward-declared; the same applies to ChargeReduction in sccs_charge.h. Keep the required definitions in the implementation files. Conversely, please include <algorithm> directly in read_inp_sccs.cpp for std::find, and source_base/constants.h in sccs_response.cpp for the constants it uses, rather than relying on transitive includes.
Broadcast the density residual only while SCCS awaits activation. Preserve synchronized activation, mixing-history reset, and the convergence veto for the activation iteration. Document that forces with lowpass disabled may differ from finite differences, recommend checking force accuracy for PCC geometry optimization, and regenerate the parameter documentation. Algorithms and defaults are unchanged. Verification: incremental abacus_std_para build passed; OMP_NUM_THREADS=1 python3 /tmp/test_sccs_broadcast_guard.py passed four MPI cases with mpirun -np 4 (ordinary, immediate, delayed, activation veto); ctest --test-dir /home/lyt/DFT/abacus_sccs/build_clean_20260929 --output-on-failure -R '^MODULE_HAMILT_surchem_h_corr_sccs$' passed 1/1. CLI help and documentation generation passed. git diff --check passed. Staged governance check passed with a test-evidence warning; focused runtime evidence is recorded here.
…ollectively Remove static Epcc and expose the valid energy of each solvent instance through its owning Potential. Pass that energy explicitly to ElecState::cal_energies at all production call sites, including energy-threshold and converged-potential updates. Preserve historical Ael and Acav statics. Add regressions for interleaved solvent and electronic-state instances and invalidated PCC results. Validate PCC summary readiness collectively in the computation path when SCCS/PCC debug output is enabled. Every rank takes the same failure branch; printing skips unavailable results and avoids throwing result accessors. Add a four-rank regression with an invalid result on rank zero only. This is defensive hardening, not a reproduced normal-run deadlock or performance regression. Verification (OMP_NUM_THREADS=1; toolchain/install/setup sourced; MPI outside sandbox): - cmake --build /home/lyt/DFT/abacus_sccs/build_clean_20260929 --target abacus_std_para MODULE_HAMILT_surchem_h_corr_sccs -j 12: passed after final edits. - Built MODULE_ESTATE_elecstate_energy, MODULE_ESTATE_potentials_new and MODULE_HAMILT_surchem_sol_force; focused CTest selection passed 4/4 before the printing follow-up. After the follow-up, ctest --test-dir /home/lyt/DFT/abacus_sccs/build_clean_20260929 --output-on-failure -R '^MODULE_HAMILT_surchem_h_corr_sccs$' passed 1/1. - mpirun -np 4 /home/lyt/DFT/abacus_sccs/build_clean_20260929/source/source_hamilt/module_surchem/test/MODULE_HAMILT_surchem_h_corr_sccs --gtest_filter=HCorrSccs.InvalidIterationResultIsReportedOnEveryRank: passed on every rank. - python3 /tmp/test_pcc_instance_energy.py: five four-rank SCF cases passed after final edits (ordinary, immediate SCCS, delayed SCCS, energy threshold, activation convergence veto). Immediate and delayed energies agreed within 1e-5 eV. - git diff --check: passed. - python3 tools/03_code_analysis/agent_governance_check.py --staged: no blocking findings; two warnings explained below. No INPUT behavior, default, or algorithm changes; no parameter documentation update required. The RDMFT PARAM usage is an unchanged call argument whose line was extended for the explicit PCC energy; added=1, removed=1, net change=0.
Assign the seven PW SCCS/PCC cases unique consecutive IDs 220-226. Update CASES_CPU.txt and README cross-references while preserving existing upstream cases. Verification: python3 /tmp/check_sccs_case_renumbering.py passed (42 renamed files match HEAD apart from intended README reference updates; unique case IDs; list/directory consistency; no old names in tracked files). python3 tools/03_code_analysis/agent_governance_check.py --staged reported no findings. git diff --cached --check passed. Calculations were not rerun because only directory names, list entries and README references changed; calculation inputs and references are byte-identical.
Compute grid coordinates, lattice-vector components, PCC gradients and forces, and nested geometry/moment helper results in named local variables before passing them to calls, following AGENTS.md rule 14. Preserve arithmetic grouping, including adding a and b, adding c, then scaling by 0.5 for the cell origin. Verification with OMP_NUM_THREADS=1 and toolchain/install/setup sourced: - cmake --build /home/lyt/DFT/abacus_sccs/build_clean_20260929 --target abacus_std_para MODULE_HAMILT_surchem_pw_grid MODULE_HAMILT_surchem_pcc_0d MODULE_HAMILT_surchem_pcc_2d MODULE_HAMILT_surchem_h_corr_sccs MODULE_HAMILT_surchem_sol_force -j 12: passed. Rebuilt affected production and correction/force targets after the final ionic-center helper edit; passed. - ctest --test-dir /home/lyt/DFT/abacus_sccs/build_clean_20260929 --output-on-failure -R '^(MODULE_HAMILT_surchem_pw_grid|MODULE_HAMILT_surchem_pcc_0d|MODULE_HAMILT_surchem_pcc_2d|MODULE_HAMILT_surchem_h_corr_sccs|MODULE_HAMILT_surchem_sol_force)$': 5/5 passed after final edits, outside the sandbox. - git diff --check and git diff --cached --check: passed. - python3 tools/03_code_analysis/agent_governance_check.py --staged: no blocking findings; warnings request test evidence and documentation rationale, provided here. This is a mechanical readability cleanup; existing numerical tests were rerun rather than adding tests that duplicate implementation. No INPUT behavior, algorithm or default change; no documentation update required. Runtime SCF calculations were not rerun.
Forward-declare Matrix3 in pw_grid.h and ChargeReduction in sccs_charge.h. Keep complete definitions in the implementation files and tests that require them. Include algorithm directly for std::find in read_inp_sccs.cpp and source_base/constants.h for the constants used by sccs_response.cpp. Verification: after sourcing toolchain/install/setup, OMP_NUM_THREADS=1 cmake --build /home/lyt/DFT/abacus_sccs/build_clean_20260929 -j 4 passed for all targets. The first build exposed an indirect ChargeReduction include in sccs_pw_nonel.cpp; that and dependent users were given explicit includes before the successful build. git diff --cached --check passed. The staged governance checker found no blockers; no INPUT behavior, algorithm or default changes, so no docs update is required. The complete CTest suite is scheduled to run with -j 4 after this commit, as requested by the user; test results are not claimed in this commit message.
Apply AGENTS.md rule 14 to the remaining new SCCS call sites that the PCC cleanup did not cover: cavity switching function, ionic normalization and far-field Gauss-law checks, Gaussian-ion force, surface-gradient norm, sqrt-CG residual norms, the deferred-SCCS and PCC-force buffers, and two exception messages. Each computed argument is assigned to a named local first; the arithmetic grouping is unchanged. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced, MPI outside the sandbox): - cmake --build /home/lyt/DFT/abacus_sccs/build_clean_20260929 -j 12: passed. - ctest -R "MODULE_HAMILT_surchem|MODULE_IO_read_input_serial| MODULE_IO_input_test_para|MODULE_IO_input_help_test|MODULE_IO_read_item| MODULE_ESTATE_elecstate_energy|MODULE_ESTATE_elecstate_print| MODULE_ESTATE_potentials_new|MODULE_ESTATE_charge_mixing": 30/30 passed. - Autotest.sh -n 4 -o 1 for 01_PW '^22[0-6]_PW_SCCS' (41 checks), 02_NAO_Gamma and 03_NAO_multik scf_sccs_pcc2d (6 checks each): passed. Per-iteration E_KohnSham/E_Harris, E_sol_el, E_sol_cav, E_pcc, final energies and forces of all nine cases are bit-identical to the four-rank run at f48579d. - git diff --cached --check and agent_governance_check.py --staged: no blockers. No INPUT behavior, default or algorithm change; no documentation update is required, and the unchanged numerical results are the test evidence.
… PCC symmetry check pcc_symmetry_violation also checks the primitive-cell translations that rhog_symmetry applies when pricell_loop is set, and the antiunitary operations of magnetic groups. Add a regression for both: an in-plane supercell translation is accepted, half the open vector is rejected only with pricell_loop, and an antiunitary operation that exchanges the open vector with a periodic one is rejected while a mirror normal to it is accepted. The static pricell_loop flag is restored before the assertions. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced): - cmake --build /home/lyt/DFT/abacus_sccs/build_clean_20260929 -j 12: passed. - MODULE_HAMILT_surchem_h_corr_sccs --gtest_filter='SurchemInput.*': 13/13 passed; the surchem/INPUT/elecstate ctest selection: 30/30 passed. - git diff --cached --check and agent_governance_check.py --staged: no blockers. Test-only change; no production code, INPUT or documentation change.
… check Finish AGENTS.md rule 14 in the two remaining PCC files: the minimum-image displacement and wrapped center of pcc_2d_system_center, and the translation and identity errors of the symmetry check, are assigned to named locals before they are passed on. The arithmetic is unchanged. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced, MPI outside the sandbox): cmake --build of all targets passed; the surchem/INPUT/elecstate ctest selection passed 30/30; Autotest.sh -n 4 -o 1 for the seven PW SCCS/PCC cases and both NAO scf_sccs_pcc2d cases passed, with per-iteration energies, energy terms, final energies and forces bit-identical to the run at f48579d. No INPUT, default or documentation change.
… PCC code - Pcc2dGeometry defaulted to the y axis from before pcc_2d_axis existed; default it to axis 2 like the INPUT keyword, and state that the multipoles are taken about the origin with the coordinate cut half a cell away from it. Production always builds the geometry with pcc_2d_geometry(); the unit-test helper that relied on the old default now selects y explicitly. - The zero-dimensional cubic-cell error named SCCS although standalone PCC reaches it too. - evaluate_pw_sccs validated the PCC geometry right before constructing the Coulomb operators, whose constructors validate it again. - ADD_SCCS_REAL_ITEM used the full parameter description as the INPUT.info annotation, so every sccs_* line of INPUT.info carried hundreds of characters; give each item a one-line annotation like the other keywords. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced, MPI outside the sandbox): cmake --build of all targets passed; the surchem/INPUT/elecstate ctest selection passed 30/30; abacus_std_para --generate-parameters-yaml is identical to docs/parameters.yaml (annotations are not part of it), so no documentation update is needed. No INPUT behavior or default change.
- Suites follow the file under test: SccsPcc2d -> Pcc2d, SccsPcc -> Pcc0d, SccsPwCharge -> PwGrid and SccsPeriodic -> SccsResponse. PCC and the grid helper are not SCCS code, and pw_grid replaced the old pw_charge helper. - Test names and comments from the y-only implementation now refer to the open axis or the cell boundary (AccumulatesOnlyNormalMoments, ReducesNormalMoments, L instead of L_y); the moment-reduction fragment selects its y open axis explicitly instead of relying on a default. - Drop an unused reversed-normal lattice: a reversed open vector has been accepted since the open axis became configurable. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced): cmake --build of all targets passed; ctest -R MODULE_HAMILT_surchem passed 22/22, and --gtest_list_tests shows the new suite names. Test-only change; no production code, INPUT or documentation change.
- assume_isolated: name the point-counter-charge correction, state that the pcc_2d multipoles are taken about the mass-weighted ionic center and that the plane half a cell away must lie in vacuum, state the symmetry rules in one place (a mirror normal to the open vector may carry a fractional translation; primitive-cell translations are checked as well), move the text shared by pcc_0d and pcc_2d out of the pcc_2d bullet, list the calculation types the INPUT check allows, and cite Andreussi and Marzari (2014). - imp_sol: a charged solute needs PCC for any water preset, not only when sccs_epsilon exceeds 1; list the nspin, device, field, Makov-Payne and DFT-1/2 restrictions that check_solvation enforces; cite the SCCS paper. - pcc_2d_axis: drop the efield_dir comparison, efield is incompatible with PCC. - sccs_maxiter: drop the obsolete remark about a discrete adjoint. - sccs_start_drho: "at the start of the run" instead of "on a cold start", and note that the SCF does not stop in the activation iteration. - sccs_lowpass_p1: the continuum cavity potential uses the FFT gradient of the PCC-corrected potential, which oscillates around the potential step at the cell boundary (test SccsPcc2dSqrtCg.LayeredCavityMatchesOpenOneDimensionalField measures it); advise keeping the dielectric transition away from it. - docs/advanced/scf/advanced.md described only imp_sol 1; add a section on SCCS and the PCC open-boundary corrections. docs/parameters.yaml was regenerated with abacus_std_para --generate-parameters-yaml and input-main.md with docs/generate_input_main.py. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced): cmake --build of all targets passed; ctest -R "MODULE_IO_read_input_serial| MODULE_IO_input_test_para|MODULE_IO_input_help_test|MODULE_IO_read_item| MODULE_HAMILT_surchem_h_corr_sccs" passed 6/6; abacus_std_para -h imp_sol and -h sccs_maxiter print the new text. Descriptions only; no INPUT behavior or default change.
…e READMEs - The two SCCS examples used a custom sccs_epsilon of 1.1, a numerical test setting rather than a solvent; use sccs_preset water-neutral and drop the keywords that repeat defaults or are overridden by the preset. - examples/27_imp_sol/README described only imp_sol 1: add the SCCS and PCC keywords, how to take the solvation energy with the same assume_isolated in both runs, and the SCCS and PCC references; rewrite the two SCCS example entries without the pre-pcc_2d_axis "(+y)" wording. - Shorten the SCCS/PCC integration-case READMEs to what each case checks. The 226 README still pointed to case 214, which is 220 after the renumbering. Verification (OMP_NUM_THREADS=1, toolchain/install/setup sourced, 8 MPI ranks): both new example INPUTs, and the same inputs with imp_sol 0, converged (#SCF IS CONVERGED#). Solvation energies: PW -0.339 eV, LCAO -0.400 eV for one water molecule per 10 x 10 bohr in-plane cell at ecutwfc 60 Ry. git grep finds no remaining old case numbers. Integration-case inputs and references are unchanged.
After the merge with develop, the cell object library contains mdcell_reader.cpp, which calls DomainDecomposition. The three surchem tests that link cell objects (sccs_gaussian_ion, h_corr_sccs, sol_force) failed to link in CI with undefined references to DomainDecomposition.
Author
|
@Critsium-xy Thanks. The y restriction was a scope limit of the first implementation. It was not required by the numerical formulation and was not inherited from ENVIRON, which selects the axis with |
…c boundaries A charged solute in a dielectric solvent with assume_isolated none used to stop the run. The periodic Poisson solver drops the G = 0 component of the net charge, so such a calculation carries a cell-size error in an inhomogeneous dielectric, but it is well defined: the screening reduces the leading Madelung error by about 1/epsilon, and ENVIRON and VASPsol also allow it. It now runs with a warning that recommends pcc_0d or pcc_2d for converged charged energies. The imp_sol description and the SCCS section of advanced.md are updated, and docs/parameters.yaml and input-main.md are regenerated.
This branch has not been deployed
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
3 participants
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.
Summary
This PR adds the self-consistent continuum solvation model (SCCS, Andreussi, Dabo & Marzari 2012) and point-counter-charge (PCC) open-boundary corrections for 0D and 2D systems to ABACUS. The implementation follows QE-Environclosely, so results can be checked against it directly.
Technical approach
sccs_rho_min,sccs_rho_max,sccs_epsilon). The water presets (water-neutral,water-cation,water-anion) are those of Environ.generalized_sqrt). The preconditioner's Coulomb operator is periodic, or periodic plus the PCC term, so a charged solute keeps its physical gauge and its screening charge. The solute charge is the electron density plus Gaussian ions (spread 0.5 bohr, as Environ).fullmode.assume_isolated pcc_0d: cubic cell, Madelung constant, origin at the mass-weighted ionic center.assume_isolated pcc_2d: slab periodic.imp_sol 0or2).sccs_solvent_mode full: Environsolvent_mode='full'.sccs_lowpass_p1/p2: Environderiv_lowpass, with the exact derivative of the discrete energy for the cavity potential.sccs_start_drho/sccs_start_nmax: delayed start of the solvent.sccs_debug: diagnostic output.imp_sol 2selects SCCS;imp_sol 1remains the legacy solvent.1. PCC: convergence with cell size
10.1103/PhysRevB.77.115139:
The test uses pyridazine (neutral) and pyridazinium (+1), LCAO APNS DZP, Gamma point, no solvent, with cubic cells of 12–30 bohr.
makov_payne.cpptakes the electronic dipole and quadrupole integrals at voxel centers, (i+0.5)/N.2. Comparison with QE-Environ's SCCS
fullmode, PW): 0.001 meV. ABACUS and QE both give −72.894 kcal/mol, against an experimental −73.4.electronicmode fails for S/Cl anions in both codes. The QE crash was "too many bands are not converged"; the reason is explained in section 3.Force check. O is displaced by ±0.01 Å. The table shows the analytic force minus the finite difference, for the solvation part only (solvent minus vacuum + PCC), in eV/Å.
full, corespread 1.53. Solvation energies: pseudopotentials and basis sets
Literature:
The test set is 34 ions (15 cations, 19 anions) from 10.1063/1.4832475. All runs use LCAO DZP, PBE, PCC0D, and the water-cation or water-anion preset.
electronicfullWhen to use
sccs_solvent_mode fullsccs_rho_max. For example, APNS gives 1.2×10⁻³ for S and 4.5×10⁻³ for Cl. The presets userho_max= 3.5×10⁻³ (cation) and 1.55×10⁻² (anion). Inelectronicmode, dielectric then appears inside the atom and the SCF diverges.electronic, 6 of the 34 ions do not converge and 1 gives an unreasonable energy.fullmode adds Environ's core-electron Gaussians to the cavity density only. With it, all 34 ions converge in both codes.fullwhen the pseudo-valence density at a nucleus falls belowsccs_rho_max.fullchanges the energy by 3×10⁻⁸ eV.Basis set: DZP to TZP
Small anions are sensitive to the LCAO basis, because the extra electron is diffuse in vacuum. For ethanethiolate, the PBE HOMO in vacuum lies at +0.2 eV.
full)The DZP error comes mostly from the vacuum reference. Its vacuum total energy is 1.9 eV above PW, against 1.6 eV above PW in solvent. TZP reduces the error from 6.9 to 1.4 kcal/mol. For anions, use TZP or PW.
When to use
sccs_lowpass_p1/p2Below is an example that applied ccqn from atst-tools to search for TS (with sccs_lowpass_p1/p2=10/5) and then obtained sp energy with tzp level basis (without sccs_lowpass and with sccs_solvent_mode=full). And we got similar results compared with the literature.
Tests and verification
source/source_hamilt/module_surchem/test/, plus INPUT and elecstate tests. The 28 relevant ctest entries pass.tests/01_PW/214–219(SCCS with PCC2D, periodic, PCC0D, charged H3O⁺, lowpass);tests/02_NAO_Gamma/scf_sccs_pcc2d;tests/03_NAO_multik/scf_sccs_pcc2d.