Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add energy contributions #36

Open
wants to merge 8 commits into
base: develop
Choose a base branch
from
+177 −20
Open

Add energy contributions #36

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
13d82b4
Add different energies to energies and inner-scf blocks.
Dec 9, 2023
60efdb4
Adjust tests and change output format for levelparser.
Dec 9, 2023
49b9657
Update cp2k_output_tools/blocks/energies.py
hdsassnick Dec 11, 2023
ede7500
Update cp2k_output_tools/blocks/energies.py
hdsassnick Dec 11, 2023
a055b49
Update cp2k_output_tools/blocks/energies.py
hdsassnick Dec 11, 2023
e74b640
Update cp2k_output_tools/blocks/energies.py
hdsassnick Dec 11, 2023
48e3554
Name INNER_SCF_CONV_RE pattern back to INNER_SCF_END_RE to avoid back…
Dec 11, 2023
3591c12
Reverse nested list for spans.
Dec 11, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
32 changes: 28 additions & 4 deletions cp2k_output_tools/blocks/energies.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,23 @@

from .common import FLOAT, BlockMatch

MAIN_ENERGY_RE = re.compile(
rf"""
^\s*Overlap\ energy\ of\ the\ core\ charge\ distribution:\s*(?P<overlap_core>{FLOAT})\n
\s*Self\ energy\ of\ the\ core\ charge\ distribution:\s*(?P<self_core>{FLOAT})\n
\s*Core\ Hamiltonian\ energy:\s*(?P<core_hamiltonian>{FLOAT})\n
\s*Hartree\ energy:\s*(?P<hartree>{FLOAT})\n
\s*Exchange-correlation\ energy:\s*(?P<xc>{FLOAT})\n
(\s*Electronic\ entropic\ energy:\s*(?P<electronic_entropic>{FLOAT})\n)?
(\s*Fermi\ energy:\s*(?P<fermi>{FLOAT})\n)?
(\s*Dispersion\ energy:\s*(?P<dispersion>{FLOAT})\n)?
hdsassnick marked this conversation as resolved.
Show resolved Hide resolved
\n
\s*Total\ energy:\s*(?P<total>{FLOAT})\n
""",
re.VERBOSE | re.MULTILINE,
)


FORCE_EVAL_ENERGY_RE = re.compile(
rf"""
^\s*ENERGY\|\ Total\ FORCE_EVAL [^:]+:\s*(?P<value>{FLOAT})\n
Expand All @@ -13,9 +30,16 @@


def match_energies(content: str) -> Optional[BlockMatch]:
spans = []
energies = {}
match = MAIN_ENERGY_RE.search(content)
if match:
energies = match.groupdict()
spans = match.spans(0)
hdsassnick marked this conversation as resolved.
Show resolved Hide resolved
match = FORCE_EVAL_ENERGY_RE.search(content)

if not match:
if match:
energies["total force_eval"] = match["value"]
spans += match.spans(0)
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
spans += match.spans(0)
spans += [match.spans(0)]

Copy link
Author

Choose a reason for hiding this comment

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

I reversed this change since it was resulting in the cli tests to fail. Checking the output I believe that my original implementation was actually already giving the correct output?

if len(energies) == 0:
hdsassnick marked this conversation as resolved.
Show resolved Hide resolved
return None

return BlockMatch({"energies": {"total force_eval": float(match["value"])}}, match.spans(0))
return BlockMatch({"energies": {key: float(val) for key, val in energies.items() if val is not None}}, spans)
25 changes: 20 additions & 5 deletions cp2k_output_tools/blocks/scf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@

from . import UREG
from .common import FLOAT, Level
from .energies import FORCE_EVAL_ENERGY_RE
from .energies import FORCE_EVAL_ENERGY_RE, MAIN_ENERGY_RE
from .mulliken import MULLIKEN_POPULATION_ANALYSIS_RE
from .warnings import Message, match_messages

Expand Down Expand Up @@ -71,7 +71,7 @@


INNER_SCF_START_RE = re.compile(r"^\s+ SCF\ WAVEFUNCTION\ OPTIMIZATION", re.VERBOSE | re.MULTILINE)
INNER_SCF_END_RE = re.compile(
INNER_SCF_CONV_RE = re.compile(
Copy link
Contributor

Choose a reason for hiding this comment

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

Please note that the blocks themselves can be used by other tools, changing the names will break them.
But if you think there is good reason to do so, please continue. The only thing I would then ask you is whether you can add a new entry in CHANGELOG.md for the next minor version (major.minor.patch) to make sure we don't forget to bump it on the next release for the backwards incompatible change.
The reason I called it SCF_END_RE is because it's notoriously difficult to reliably determine whether the CP2K (SCF) finished properly. The idea was to assume that the presence of that line indicates a finished SCF, everything else not.

Copy link
Author

@hdsassnick hdsassnick Dec 11, 2023 •
edited
Loading

Choose a reason for hiding this comment

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

Sorry, I didn't take that into consideration. I just thought that as the energies are parsed further below in the output file and SCF_END_RE would refer to the end of the inner SCF cycle it would make sense to rename it.

However, having backwards compatibility is probably more important in this case, thus I reversed this change.

r"^\s+ (?P<convtxt>\*\*\*\ SCF\ run\ converged\ in|Leaving\ inner\ SCF\ loop\ after\ reaching) \s+ (?P<nsteps>\d+)\ steps",
re.VERBOSE | re.MULTILINE,
)
Expand Down Expand Up @@ -132,6 +132,15 @@ class OuterSCF(Level):
class InnerSCF(Level):
converged: bool
nsteps: int
overlap_core_energy: Decimal
self_core_energy: Decimal
core_hamiltonian_energy: Decimal
hartree_energy: Decimal
xc_energy: Decimal
total_energy: Decimal
electronic_entropic_energy: Optional[Decimal] = None
fermi_energy: Optional[Decimal] = None
dispersion_energy: Optional[Decimal] = None


@dataclass
Expand Down Expand Up @@ -194,11 +203,17 @@ def match_scf(content: str, start: int = 0, end: int = sys.maxsize) -> Optional[
match = INNER_SCF_START_RE.search(content, start, end)
if match:
start = match.span()[1]
ematch = INNER_SCF_END_RE.search(content, start, end)
if ematch:
energy_match = MAIN_ENERGY_RE.search(content, start, end)
conv_match = INNER_SCF_CONV_RE.search(content, start, end)
if conv_match and energy_match:
start = match.span()[1]
kwargs = {
key + "_energy": Decimal(val) * UREG.hartree for key, val in energy_match.groupdict().items() if val is not None
}
sublevels.append(
InnerSCF(converged="SCF run converged" in ematch["convtxt"], nsteps=int(ematch["nsteps"]), sublevels=[])
InnerSCF(
converged="SCF run converged" in conv_match["convtxt"], nsteps=int(conv_match["nsteps"]), **kwargs, sublevels=[]
)
)

force_eval_energy: Optional[Decimal] = None
Expand Down
14 changes: 13 additions & 1 deletion tests/test_energies.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,4 +9,16 @@ def test_energies():
result = next(parse_iter(fhandle.read(), matchers=[match_energies]))

assert result
assert result == {"energies": {"total force_eval": -251.687390311050706}}
assert result == {
"energies": {
"core_hamiltonian": 138.60026809219576,
"electronic_entropic": -1.256404e-08,
"fermi": 0.20382509931778,
"hartree": 343.3040142273272,
"overlap_core": 3.0088e-10,
"self_core": -656.5115154010256,
"total": -251.63989121633358,
"xc": -77.03265812258856,
"total force_eval": -251.687390311050706,
}
}
128 changes: 117 additions & 11 deletions tests/test_optimization.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,35 +13,79 @@
"num_occ_orb": 1,
"num_mol_orb": 1,
"num_orb_func": 10,
"nsteps": 8,
"force_eval_energy": -1.16096052,
"inner_scf": {
"nsteps": 8,
"overlap_core_energy": 6.0512959e-07,
"self_core_energy": -2.82094792,
"core_hamiltonian_energy": 1.07883184,
"hartree_energy": 1.30392538,
"xc_energy": -0.72277042,
"total_energy": -1.16096052,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 2,
"num_occ_orb": 1,
"num_mol_orb": 1,
"num_orb_func": 10,
"nsteps": 4,
"force_eval_energy": -1.16118019,
"inner_scf": {
"nsteps": 4,
"overlap_core_energy": 0.0000012,
"self_core_energy": -2.82094792,
"core_hamiltonian_energy": 1.08370608,
"hartree_energy": 1.30183473,
"xc_energy": -0.72577425,
"total_energy": -1.16118019,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 2,
"num_occ_orb": 1,
"num_mol_orb": 1,
"num_orb_func": 10,
"nsteps": 4,
"force_eval_energy": -1.16118537,
"inner_scf": {
"nsteps": 4,
"overlap_core_energy": 1.07320985e-06,
"self_core_energy": -2.82094792,
"core_hamiltonian_energy": 1.08302262,
"hartree_energy": 1.30210002,
"xc_energy": -0.72536117,
"total_energy": -1.16118537,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 2,
"num_occ_orb": 1,
"num_mol_orb": 1,
"num_orb_func": 10,
"nsteps": 4,
"force_eval_energy": -1.16118537,
"inner_scf": {
"nsteps": 4,
"overlap_core_energy": 1.07320985e-06,
"self_core_energy": -2.82094792,
"core_hamiltonian_energy": 1.08302280,
"hartree_energy": 1.30210002,
"xc_energy": -0.72536117,
"total_energy": -1.16118537,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
],
},
Expand All @@ -53,44 +97,99 @@
"num_occ_orb": 16,
"num_mol_orb": 16,
"num_orb_func": 55,
"nsteps": 10,
"force_eval_energy": -65.98706397,
"inner_scf": {
"nsteps": 10,
"overlap_core_energy": 0.0,
"self_core_energy": -138.89582544,
"core_hamiltonian_energy": 27.22305865,
"hartree_energy": 56.70225641,
"xc_energy": -11.01655328,
"total_energy": -65.98706366,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 32,
"num_occ_orb": 16,
"num_mol_orb": 16,
"num_orb_func": 55,
"nsteps": 9,
"force_eval_energy": -65.98706506,
"inner_scf": {
"nsteps": 9,
"overlap_core_energy": 0.0,
"self_core_energy": -138.89582544,
"core_hamiltonian_energy": 27.22275912,
"hartree_energy": 56.70220023,
"xc_energy": -11.01647284,
"total_energy": -65.98733894,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 32,
"num_occ_orb": 16,
"num_mol_orb": 16,
"num_orb_func": 55,
"nsteps": 5,
"force_eval_energy": -65.98710767,
"inner_scf": {
"nsteps": 5,
"overlap_core_energy": 0.0,
"self_core_energy": -138.89582544,
"core_hamiltonian_energy": 27.20749466,
"hartree_energy": 56.71104696,
"xc_energy": -11.00902476,
"total_energy": -65.98630859,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 32,
"num_occ_orb": 16,
"num_mol_orb": 16,
"num_orb_func": 55,
"nsteps": 4,
"force_eval_energy": -65.98710770,
"inner_scf": {
"nsteps": 4,
"overlap_core_energy": 0.0,
"self_core_energy": -138.89582544,
"core_hamiltonian_energy": 27.20670707,
"hartree_energy": 56.71083178,
"xc_energy": -11.00884866,
"total_energy": -65.98713525,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
{
"nspin": 1,
"nelec": 32,
"num_occ_orb": 16,
"num_mol_orb": 16,
"num_orb_func": 55,
"nsteps": 1,
"force_eval_energy": -65.98710770,
"inner_scf": {
"nsteps": 1,
"overlap_core_energy": 0.0,
"self_core_energy": -138.89582544,
"core_hamiltonian_energy": 27.20676302,
"hartree_energy": 56.71080149,
"xc_energy": -11.00884677,
"total_energy": -65.98710770,
"electronic_entropic_energy": None,
"fermi_energy": None,
"dispersion_energy": None,
},
},
],
"cell_infos": [
Expand Down Expand Up @@ -149,12 +248,19 @@ def test_optimization(opt_type):
for step_scf, opt_ref in zip(opt_steps_scf + [result.levels[0].sublevels[1]], ref["opt"]):
assert step_scf.converged
assert step_scf.sublevels[0].converged
inner_scf_ref = opt_ref.pop("inner_scf")
inner_scf = step_scf.sublevels[0]
for keyw, val in inner_scf_ref.items():
if val is None:
assert getattr(inner_scf, keyw) is None
elif isinstance(val, float):
assert abs(float(getattr(inner_scf, keyw).magnitude) - val) < 1.0e-7
else:
assert getattr(inner_scf, keyw) == val
assert abs(float(step_scf.force_eval_energy.magnitude) - opt_ref.pop("force_eval_energy")) < 1.0e-7
assert step_scf.sublevels[0].nsteps == opt_ref.pop("nsteps")
for keyw, val in opt_ref.items():
assert getattr(step_scf, keyw) == val
if opt_type == "cell":
print(result.levels[0].cell_infos)
dev-zero marked this conversation as resolved.
Show resolved Hide resolved
assert len(result.levels[0].cell_infos) == len(ref["cell_infos"])
for info, ref_info in zip(result.levels[0].cell_infos, ref["cell_infos"]):
assert abs(float(info.volume.magnitude) - ref_info.pop("volume")) < 1.0e-7
Expand Down