diff --git a/README.md b/README.md index 4160515..5fd541f 100644 --- a/README.md +++ b/README.md @@ -70,6 +70,10 @@ every detected occurrence, every refinement with its band and saddle check, and `reactions.sqlite3` index under `output_root//`. Nothing is overwritten, and a resubmitted job continues each unfinished trajectory from its last checkpoint. +Refining one event stops after `max_hours_per_barrier` hours, 4 by default, and records +`time_limit` with no barrier, so that one slow band cannot hold up a trajectory for days +([time limit per barrier](docs/campaigns.md#time-limit-per-barrier)). MD has no time limit. + ## Models MACE, UMA, AIMNet2, and ANI-1xnr are built in. `reactionflow models` lists their models, heads, and diff --git a/docs/campaigns.md b/docs/campaigns.md index f4feb78..cd29dee 100644 --- a/docs/campaigns.md +++ b/docs/campaigns.md @@ -247,6 +247,31 @@ NVT barriers are potential energies; NPT barriers are enthalpies. `--json` print including every barrier value, the counts of each pathway, frequency and connectivity status, trajectory ID, and one example occurrence per class, for scripted analysis. +## Time limit per barrier + +Refining one event stops after `max_hours_per_barrier` hours, 4 by default. Set it at the top of the +campaign to a positive number, or to `null` for no limit. The event then ends as `time_limit`, with +no barrier and its band as far as it got, and the trajectory moves on. MD has no time limit. The +steps of an event refined in steps share one limit: steps already refined keep their results, and +those not started end as `time_limit` with no images. + +The clock counts the time a refinement has run, not the time a job waits in the queue. It carries +across resubmissions, because an interrupted refinement starts over when it resumes, and a job +shorter than the refinement would otherwise repeat it forever. So work lost to an interruption counts +against the limit. The count is saved once a minute, so each resubmission can lose up to a minute of +it, and a job that refines for under a minute adds nothing. It is checked before each model +evaluation, so a single evaluation that never returns is not stopped. + +The setting is not part of a trajectory's contract, so you can change it when you resubmit, and the +new limit applies to the event being refined and to later ones. An event that already ended as +`time_limit` stays that way, because results are never overwritten. + +Four hours comes from the file times and job logs of 288-atom acetonitrile runs at high pressure +with UMA on one GPU, on older commits. The two refinements that reached a connected saddle, both at +700 K, took 1.1 and 1.8 hours. None of the nine longest, 8.5 to 26.7 hours at 2000 K or hotter, +reached a saddle: four did not converge and five had no imaginary mode. More atoms, more images, a +slower model, or a CPU can need more, and `status` counts the refinements that hit the limit. + ## Perlmutter mapping The bundled ANI environment is prepared by diff --git a/docs/how-it-works.md b/docs/how-it-works.md index 14f5ebd..6779fdd 100644 --- a/docs/how-it-works.md +++ b/docs/how-it-works.md @@ -61,11 +61,14 @@ first stage that fails, with a status that says why. | NEB or SSNEB | `neb_failed` | | Climbing image | `ci_neb_failed` | | Saddle check | `saddle_validated`, `no_saddle`, `saddle_rejected`, or `saddle_unchecked` | +| Any stage, past the time limit | `time_limit` | `collapsed` means both endpoints relaxed to the same bond topology. `unresolved` means the change could not be mapped onto the relaxed endpoints, because a bond stayed between the two thresholds, the change disappeared, or relaxation changed bonds elsewhere in the cell. `stepwise` means the -event was refined in steps, each with its own outcome. `failed` means an unexpected error. The +event was refined in steps, each with its own outcome. `time_limit` means the refinement ran past +`max_hours_per_barrier`, 4 hours by default ([time limit per +barrier](campaigns.md#time-limit-per-barrier)). `failed` means an unexpected error. The saddle check's four verdicts are described in [pathway refinement](pathway-refinement.md#saddle-check). In every case the result is saved, the checkpoint is restored, and MD continues. @@ -231,7 +234,8 @@ directory. To continue after a job ends, submit the same command again. Unfinished trajectories resume from their last checkpoint and finished ones are left as they are. Each output directory is bound to -its structure, settings, and model, and a changed campaign is rejected instead of resumed. +its structure, trajectory conditions, reaction settings, and model, and a change to those is +rejected instead of resumed. [`examples/perlmutter/run-campaign.sbatch`](../examples/perlmutter/run-campaign.sbatch) is a job script for Perlmutter's four-GPU nodes. Thirty-two trajectories on eight nodes are submitted like diff --git a/docs/pathway-refinement.md b/docs/pathway-refinement.md index 36fce0f..aad5c7d 100644 --- a/docs/pathway-refinement.md +++ b/docs/pathway-refinement.md @@ -70,8 +70,10 @@ once with that lease, so each NEB, CI-NEB, or SSNEB step evaluates only the movi `PathwayOutcome.status` is `unresolved`, `collapsed`, `relaxation_failed`, `neb_failed`, `ci_neb_failed`, or `failed` when refinement stops before the band converges. A converged band ends -with the saddle check's verdict below. Outcomes retain calculator-free endpoint or NEB snapshots as -far as the attempt progressed. Every converged band keeps its image energies, and only +with the saddle check's verdict below. `time_limit` can come at any stage, the saddle check +included, and keeps the images so far but no energies or verdict; a later step that it closes before +the step starts has no images. Outcomes retain calculator-free endpoint or NEB snapshots as far as +the attempt progressed. Every band with a verdict keeps its image energies, and only `saddle_validated` adds a forward barrier in eV: potential energy for NVT, enthalpy for NPT. No status stops MD from resuming. diff --git a/src/reactionflow/campaign.py b/src/reactionflow/campaign.py index e467ecc..44b5ca5 100644 --- a/src/reactionflow/campaign.py +++ b/src/reactionflow/campaign.py @@ -12,7 +12,7 @@ from typing import Any from .pathway import _RETIRED_KEYS -from .run import ReactionRunConfig +from .run import DEFAULT_MAX_HOURS_PER_BARRIER, ReactionRunConfig _TRAJECTORY_ID = re.compile(r"[A-Za-z0-9][A-Za-z0-9._-]*\Z") @@ -167,6 +167,7 @@ class CampaignConfig: adapter_profiles: Mapping[str, AdapterSpec] trajectory_profiles: tuple[str, ...] trajectory_structures: tuple[Path, ...] + max_hours_per_barrier: float | None = DEFAULT_MAX_HOURS_PER_BARRIER require_gpu: bool = True @classmethod @@ -186,6 +187,7 @@ def load(cls, path: str | Path) -> CampaignConfig: "trajectories", "require_gpu", "adapter_profiles", + "max_hours_per_barrier", } if unknown: raise ValueError(f"unknown campaign keys: {sorted(unknown)}") @@ -250,6 +252,11 @@ def structure_path(raw: object, name: str) -> Path: identifiers = [trajectory.id for trajectory in trajectories] if len(set(identifiers)) != len(identifiers): raise ValueError("trajectory IDs must be unique") + max_hours = value.get("max_hours_per_barrier", DEFAULT_MAX_HOURS_PER_BARRIER) + if max_hours is not None: + max_hours = _finite_float(max_hours, "campaign.max_hours_per_barrier", minimum=0.0) + if max_hours == 0: + raise ValueError("campaign.max_hours_per_barrier must be positive, or null") require_gpu = value.get("require_gpu", True) if type(require_gpu) is not bool: raise TypeError("campaign.require_gpu must be a boolean") @@ -261,6 +268,7 @@ def structure_path(raw: object, name: str) -> Path: adapter_profiles=adapter_profiles, trajectory_profiles=tuple(trajectory_profiles), trajectory_structures=tuple(trajectory_structures), + max_hours_per_barrier=max_hours, require_gpu=require_gpu, ) diff --git a/src/reactionflow/cli.py b/src/reactionflow/cli.py index 324a292..44e2722 100644 --- a/src/reactionflow/cli.py +++ b/src/reactionflow/cli.py @@ -177,6 +177,7 @@ def run_selected_trajectory( atoms, runtime_provider=adapter, pathway_calculator_provider=adapter.calculator, + max_hours_per_barrier=campaign.max_hours_per_barrier, total_steps=trajectory.total_steps, pressure_GPa=trajectory.pressure_GPa, ) diff --git a/src/reactionflow/pathway.py b/src/reactionflow/pathway.py index 5fe06ba..35029ee 100644 --- a/src/reactionflow/pathway.py +++ b/src/reactionflow/pathway.py @@ -45,6 +45,14 @@ class _EndpointUnresolved(ValueError): pass +class BarrierTimeLimit(BaseException): + """Ends a refinement as ``time_limit``. + + Not an Exception: the Hessian and the connectivity check record any Exception as their own + failure and carry on, which would report the wrong status and keep running the model. + """ + + @dataclass(frozen=True, slots=True) class PathwayConfig: """Small set of controls for endpoint relaxation, CI-NEB, and the saddle check. @@ -774,6 +782,8 @@ def relaxed(frame: int, stage: str) -> tuple[Atoms, bool]: if made is not None and made > (ends[-1][0] if ends else state): ends.append((made, frame)) ends.append((changed, last)) + except BarrierTimeLimit as exc: + return PathwayOutcome("time_limit", images=_snapshots(images), message=str(exc)) except _EndpointUnresolved as exc: return PathwayOutcome("unresolved", images=_snapshots(images), message=str(exc)) except Exception as exc: @@ -803,7 +813,9 @@ def relaxed(frame: int, stage: str) -> tuple[Atoms, bool]: reactant_bonds=candidate.reactant_bonds ^ before, product_bonds=candidate.reactant_bonds ^ made, ) - if end is None: + if steps and steps[-1].outcome.status == "time_limit": + outcome = replace(steps[-1].outcome, images=()) + elif end is None: outcome = PathwayOutcome("unresolved", message="the step before ended elsewhere") else: outcome, end = _refine_step( @@ -959,6 +971,8 @@ def _refine_step( finally: for image in images: image.calc = None + except BarrierTimeLimit as exc: + return PathwayOutcome("time_limit", images=_snapshots(images), message=str(exc)), None except _EndpointUnresolved as exc: return PathwayOutcome("unresolved", images=_snapshots(images), message=str(exc)), end except Exception as exc: diff --git a/src/reactionflow/run.py b/src/reactionflow/run.py index f4f7c9a..640aee1 100644 --- a/src/reactionflow/run.py +++ b/src/reactionflow/run.py @@ -4,9 +4,10 @@ import json import shutil +import time import warnings from collections.abc import Iterator, Mapping -from contextlib import contextmanager +from contextlib import contextmanager, suppress from dataclasses import asdict, dataclass, field from numbers import Integral from pathlib import Path @@ -21,13 +22,23 @@ from ._version import __version__ from .candidates import ReactionCandidate, ReactionTracker from .detection import BondChangeDetector, BondDetectorConfig, canonical_copy, transport_copy -from .pathway import CalculatorProvider, PathwayConfig, PathwayOutcome, refine_pathway +from .pathway import ( + BarrierTimeLimit, + CalculatorProvider, + PathwayConfig, + PathwayOutcome, + refine_pathway, +) from .restart import ExactRestartSnapshot from .runtime import ExactDynamicsRuntime, ExactRuntimeProvider from .store import OccurrenceRecord, OccurrenceStore _PHASES = {"new", "running", "refining", "completed", "failed"} +# By file times, refinements that reached a connected saddle took 1.1 and 1.8 hours (288-atom +# acetonitrile, SSNEB, UMA on one GPU), and none of the nine longest, from 8.5 h, reached a saddle. +DEFAULT_MAX_HOURS_PER_BARRIER = 4.0 + @dataclass(frozen=True, slots=True) class ReactionRunConfig: @@ -142,6 +153,8 @@ def __init__(self, root: Path, config: ReactionRunConfig) -> None: self._global_frame = 0 self._pending: list[str] = [] self._failure: dict[str, str] | None = None + self._barrier_seconds = 0.0 + self._barrier_clock = 0.0 self._initial: Atoms | None = None self._detector: BondChangeDetector | None = None self._tracker: ReactionTracker | None = None @@ -182,6 +195,7 @@ def open(cls, root: str | Path) -> ReactionRun: run._pending = list(map(str, value["pending_pathway_ids"])) run._created_with = value.get("created_with_reactionflow") run._failure = value.get("failure") + run._barrier_seconds = float(value.get("barrier_seconds", 0.0)) if value.get("detector_state") is not None: run._detector = BondChangeDetector.from_state( value["detector_state"], @@ -226,6 +240,7 @@ def _write_state(self) -> None: "global_step": self._global_step, "global_frame": self._global_frame, "pending_pathway_ids": list(self._pending), + "barrier_seconds": self._barrier_seconds, "failure": self._failure, "config": self.config.to_dict(), "detector_state": None if self._detector is None else self._detector.export_state(), @@ -485,6 +500,45 @@ def _publish_outcome(self, occurrence_id: str, outcome: PathwayOutcome) -> None: shutil.rmtree(temporary, ignore_errors=True) raise + def _time_limited( + self, provider: CalculatorProvider, max_hours: float | None + ) -> CalculatorProvider: + """A provider that ends a refinement as ``time_limit`` once it has run for ``max_hours``. + + The time is checked before each model evaluation and added up in ``state.json``, because an + interrupted refinement starts over when it resumes: a job shorter than the refinement would + otherwise repeat it forever. + """ + + def check() -> None: + before = self._barrier_seconds + now = time.monotonic() + self._barrier_seconds += now - self._barrier_clock + self._barrier_clock = now + if max_hours is not None and self._barrier_seconds >= max_hours * 3600: + raise BarrierTimeLimit(f"stopped at the limit of {max_hours:g} h per barrier") + if self._barrier_seconds // 60 > before // 60: # save once a minute + with suppress(OSError): # this runs inside the refinement, whose result it would be + self._write_state() + + @contextmanager + def lease(stage: str) -> Iterator[Any]: + check() + with provider(stage) as calculator: + calculate = calculator.calculate + + def checked(*args: Any, **kwargs: Any) -> Any: + check() + return calculate(*args, **kwargs) + + calculator.calculate = checked + try: + yield calculator + finally: + calculator.calculate = calculate + + return lease + def _refine_pending( self, calculator_provider: CalculatorProvider, @@ -499,6 +553,7 @@ def _refine_pending( # Published before an interruption; make its name durable and move on. sync_directory(self.pathways) else: + self._barrier_clock = time.monotonic() candidate = self.occurrences.load(occurrence_id) outcome = refine_pathway( candidate, @@ -513,6 +568,7 @@ def _refine_pending( ) self._publish_outcome(occurrence_id, outcome) self._pending.pop(0) + self._barrier_seconds = 0.0 self._write_state() def _complete(self) -> None: @@ -544,15 +600,22 @@ def run( *, runtime_provider: ExactRuntimeProvider, pathway_calculator_provider: CalculatorProvider, + max_hours_per_barrier: float | None = DEFAULT_MAX_HOURS_PER_BARRIER, total_steps: int, pressure_GPa: float | None = None, ) -> RunSummary: - """Run or resume live detection, serial NEB/CI-NEB, and exact MD to a step target.""" + """Run or resume live detection, serial NEB/CI-NEB, and exact MD to a step target. + + Refining one event stops after ``max_hours_per_barrier`` hours, or never with ``None``. + """ if isinstance(total_steps, bool) or not isinstance(total_steps, Integral): raise ValueError("total_steps must be a non-negative integer") if total_steps < self._global_step: raise ValueError("total_steps cannot precede the durable run counter") + pathway_calculator_provider = self._time_limited( + pathway_calculator_provider, max_hours_per_barrier + ) if self._phase == "new": if atoms is None: raise ValueError("initial atoms are required for a new run") diff --git a/src/reactionflow/status.py b/src/reactionflow/status.py index 89a4480..d5092fd 100644 --- a/src/reactionflow/status.py +++ b/src/reactionflow/status.py @@ -313,6 +313,11 @@ def format_status(report: dict[str, Any]) -> str: f"{verdicts['stepwise']} event{'' if verdicts['stepwise'] == 1 else 's'} refined in" " steps: each step counts as an event of its own class." ) + if verdicts["time_limit"]: + lines.append( + f"{verdicts['time_limit']} refinement{'' if verdicts['time_limit'] == 1 else 's'}" + " stopped at max_hours_per_barrier, with no barrier." + ) radius = report["radius"] lines.append( f"A class is a bond change with the atoms within {radius} bond{'' if radius == 1 else 's'}" diff --git a/tests/test_campaign.py b/tests/test_campaign.py index f4aa041..8a78344 100644 --- a/tests/test_campaign.py +++ b/tests/test_campaign.py @@ -348,3 +348,14 @@ def test_perlmutter_template_maps_four_local_cpu_monitors_to_four_gpus() -> None assert "#SBATCH --cpus-per-task=32" in script assert 'srun --cpu-bind=cores reactionflow run "$campaign"' in script assert "--ntasks=32" not in script + + +def test_barrier_time_limit_defaults_to_four_hours_and_null_turns_it_off(tmp_path) -> None: + path = _campaign(tmp_path) + assert CampaignConfig.load(path).max_hours_per_barrier == 4.0 + value = json.loads(path.read_text()) + path.write_text(json.dumps({**value, "max_hours_per_barrier": None})) + assert CampaignConfig.load(path).max_hours_per_barrier is None + path.write_text(json.dumps({**value, "max_hours_per_barrier": 0})) + with pytest.raises(ValueError, match="max_hours_per_barrier"): + CampaignConfig.load(path) diff --git a/tests/test_run.py b/tests/test_run.py index 081fa8f..6fba32e 100644 --- a/tests/test_run.py +++ b/tests/test_run.py @@ -5,6 +5,7 @@ import shutil from contextlib import contextmanager from pathlib import Path +from types import SimpleNamespace from typing import ClassVar import numpy as np @@ -22,6 +23,7 @@ ReactionRun, ReactionRunConfig, assign_atom_ids, + pathway, ) from reactionflow.campaign import CampaignConfig from reactionflow.cli import main @@ -429,6 +431,83 @@ def fail_frequencies(*_args, **_kwargs): assert leases.live == 0 and leases.max_live == 1 +@pytest.mark.parametrize("stage", ["_minimize", "_hessian", "_check_connectivity"]) +def test_time_limit_ends_a_refinement_at_any_stage_and_md_continues( + tmp_path, monkeypatch, stage +) -> None: + # The Hessian and connectivity check record any Exception as their own failure and carry on, so + # a limit they caught would end as saddle_unchecked. The clock jumps when the stage starts. + clock = {"seconds": 0.0} + original = getattr(pathway, stage) + + def late(*args, **kwargs): + clock["seconds"] = 5 * 3600.0 + return original(*args, **kwargs) + + monkeypatch.setattr(pathway, stage, late) + monkeypatch.setattr( + "reactionflow.run.time", SimpleNamespace(monotonic=lambda: clock["seconds"]) + ) + leases = LeaseCounter() + run = ReactionRun.create(tmp_path, config=run_config()) + + summary = run.run( + pair(0.6), + runtime_provider=ScriptedRuntimeProvider(leases), + pathway_calculator_provider=leases, + total_steps=18, + max_hours_per_barrier=4.0, + ) + + assert (summary.phase, summary.global_step, run.failure) == ("completed", 18, None) + assert leases.live == 0 and leases.max_live == 1 + result = json.loads(next((tmp_path / "pathways").glob("*/result.json")).read_text()) + assert (result["status"], result["barrier_eV"]) == ("time_limit", None) + assert "4 h" in result["message"] + assert json.loads((tmp_path / "state.json").read_text())["barrier_seconds"] == 0.0 + + +def test_time_limit_counts_across_jobs_too_short_to_finish_a_barrier(tmp_path, monkeypatch) -> None: + # An interrupted refinement starts over when it resumes, so without a count that outlives each + # 70 s job, this refinement would repeat forever under a 100 s limit. + clock = {"seconds": 0.0, "job_ends": 70.0} + + class EndsWithTheJob(PairDoubleWell): + def calculate(self, *args, **kwargs): + clock["seconds"] += 10.0 + if clock["seconds"] >= clock["job_ends"]: + raise KeyboardInterrupt("the job ended") + super().calculate(*args, **kwargs) + + @contextmanager + def calculators(_stage): + yield EndsWithTheJob() + + def advance(run, atoms=None): + clock["job_ends"] = clock["seconds"] + 70.0 + return run.run( + atoms, + runtime_provider=ScriptedRuntimeProvider(LeaseCounter()), + pathway_calculator_provider=calculators, + total_steps=18, + max_hours_per_barrier=100 / 3600, + ) + + monkeypatch.setattr( + "reactionflow.run.time", SimpleNamespace(monotonic=lambda: clock["seconds"]) + ) + with pytest.raises(KeyboardInterrupt, match="the job ended"): + advance(ReactionRun.create(tmp_path, config=run_config()), pair(0.6)) + state = json.loads((tmp_path / "state.json").read_text()) + assert state["phase"] == "refining" and 60.0 <= state["barrier_seconds"] <= 70.0 + + summary = advance(ReactionRun.open(tmp_path)) + + assert (summary.phase, summary.pathways) == ("completed", 1) + result = json.loads(next((tmp_path / "pathways").glob("*/result.json")).read_text()) + assert result["status"] == "time_limit" + + def test_reopen_refines_the_reverse_reaction_too_and_serializes_leases( tmp_path, monkeypatch ) -> None: