[CP2K-user] [CP2K:22445] Re: DDEC6 Charges, Valence Density Mismatch
tkuehne
tkuehne at gmail.com
Wed Oct 14 06:28:55 UTC 2026
Dear Nick,
the mismatch comes from interpreting the GAPW soft density as a
valence-only density. Your log already accounts for approximately
102 electrons in the complete electronic density, so it does not
indicate an extra electron or an electron-counting error in the SCF.
In your all-electron calculation,
rho_soft + rho1_hard - rho1_soft
is the complete density, including core and valence electrons. It
should therefore integrate to 102 electrons, not 76. The GAPW
hard/soft decomposition is not a chemical core/valence partition.
NGAUSS is used by XRD_INTERFACE and does not remove core electrons
from the SOFT_DENSITY cube.
For DDEC6, I recommend exporting the Gaussian-orbital wavefunction
in WFX format rather than supplying the soft-density cube. Chargemol
can then evaluate and partition the complete electronic density
using its own density treatment. For an all-electron WFX file, do
not declare the 26 core electrons as missing electrons, since they
are already included. Rescaling the soft cube to 76 electrons would
not restore the missing density shape.
CP2K already has a CHARGEMOL/WFX exporter, but it currently rejects
OT wavefunctions. I have opened a PR to enable OT export, including
UKS and ADMM, using private copies of the occupied orbitals without
changing the SCF state:
https://github.com/cp2k/cp2k/pull/6276
With a build containing this patch, you can add the following under
FORCE_EVAL/DFT/PRINT to export every tenth MD step:
&CHARGEMOL ON
FILENAME chargemol
&EACH
MD 10
&END EACH
&END CHARGEMOL
The PR is not merged yet. Without the patch, diagonalization-based
SCF remains an alternative for WFX export.
I tested the export with all-electron OT calculations, including
UKS and ADMM. For a simplified PBE/6-31G* calculation using your
24-atom molecular geometry, Chargemol accounted for all 102 electrons
and completed DDEC6 normally. This validates the export and electron
accounting, not the numerical charges of your original hybrid
QM/MM trajectory, which I have not reproduced.
For the cutoff question, WFX avoids the cube representation of
sharp core-density peaks, but does not replace convergence tests
of the underlying calculation. Check energies, forces and the
charges of interest against CUTOFF, REL_CUTOFF and the GAPW atomic
grids.
Finally, since the gold slab is classical MM, this analysis covers
the QM ion pair only. It cannot assign DDEC6 charges to the gold
atoms or describe electron transfer involving metallic electronic
states that are absent from the model.
Best regards,
Thomas
On Wednesday, May 7, 2025 at 7:48:17 PM UTC+2 Nicholas Laws wrote:
> Hi CP2K Community,
>
> I am trying to obtain time-dependent DDEC6 charges for a mixed QM/MM
> simulation in which an EMI-BF4 ion pair (76 valence e, 26 core e)
> collides with a gold slab at various kinetic energies. The QM region is
> treated with ωB97X-V/aug-cc-pVTZ and all-electron potentials in GAPW.
> During production (input *ESC**ollision.inp**,* attached) I print a cube
> with:
>
> &PRINT
>
> &E_DENSITY_CUBE
>
> DENSITY_INCLUDE SOFT_DENSITY
>
> NGAUSS 0
>
> STRIDE 1 1 1
>
> FILENAME valence_density
>
> &END
>
> &END
>
> and a 5-level PW grid of:
>
> &MGRID
>
> NGRIDS 5
>
> CUTOFF 800
>
> REL_CUTOFF 80
>
> COMMENSURATE .TRUE.
>
> &END
>
> The log shows:
>
> Electronic density on regular grids: -76.7658266487
> <(765)%20826-6487> 25.2341733513
>
> Core density on regular grids: 102.0000000000
> -0.0000000000
>
> Hard and soft densities (Lebedev): -694.5370264067 -669.
> 3028566617 <(302)%20856-6617>
>
> Total Rho_soft + Rho1_hard - Rho1_soft -101.9999963937
>
> Total charge density (r-space): 0.0000036063
>
> Total Rho_soft + Rho0_soft (g-space): 0.0000033449
>
> Overlap energy of the core charge distribution:
> 0.00002215241243
>
> Self energy of the core charge distribution:
> -720.29674329623936
>
> Core Hamiltonian energy:
> -75.86850888340894
>
> Hartree energy:
> 209.46847518535864
>
> Exchange-correlation energy:
> -26.02788481604846
>
> Hartree-Fock Exchange energy:
> -28.87249823090302
>
> Dispersion energy:
> 0.25782644744755
>
> so the cube integrates to ≈77 e instead of the expected 76 e. Chargemol
> therefore exits with
>
> ncore = 26.0000
>
> nvalence = 76.0000
>
> pixelvolume = 1.3079E-03
>
> numerically integrated valence density = 7.6767E+01
>
> sum_valence_occupancy_correction = 0.0000E+00
>
> checkme = 7.6659E-01
>
> The electrons are not properly accounted for. Either the grid in your
> electron density input file is too coarse, you have specified the incorrect
> net charge in the chargemol_job.m file, or the number of core electrons for
> each atom has not been set up correctly.
>
> Program will terminate.
>
> My questions are:
>
> 1.
>
> Is there any keyword combination in GAPW that prints only the valence
> density (ρ_soft + ρ₁,hard − ρ₁,soft) so the integral matches 76 e, or must
> I run a GPW single-point for every MD frame?
> 2.
>
> Is there a published benchmark for the PW cutoff required to converge
> ωB97X-V + aug-cc-pVTZ?
>
> Any hints, keywords, code pointers, or literature would be greatly
> appreciated!
>
> Many thanks in advance,
>
> Nick
>
--
You received this message because you are subscribed to the Google Groups "cp2k" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cp2k+unsubscribe at googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/cp2k/afbc0155-eefe-4125-800e-d4053666036an%40googlegroups.com.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20261013/7fd35675/attachment-0001.htm>
More information about the CP2K-user
mailing list