diff --git a/t3/pdep/parser.py b/t3/pdep/parser.py index abc4260b..127bb81b 100644 --- a/t3/pdep/parser.py +++ b/t3/pdep/parser.py @@ -304,6 +304,17 @@ class PDepNetwork: the PES diagram by where its energy came from. A species with no thermo comment is simply absent from this dict. Default ``{}`` so every existing construction site keeps working unchanged. + ts_labels_with_statmech (frozenset): The transition state labels whose + ``transitionState(...)`` call declares statmech data + (a non-empty ``modes`` list, or the by-reference + ``transitionState(label, path)`` form) -- see + ``_transition_state_declares_statmech``. A label + absent from this set names a TS whose conformer has + no modes, for which ``Reaction.can_tst()`` is False + and every k(E) comes from the inverse Laplace + transform. Default empty for the same + construction-site-compatibility reason as + ``species_structures``. """ network_id: str path: str @@ -318,6 +329,7 @@ class PDepNetwork: source_hash: str | None = None species_structures: dict = field(default_factory=dict) species_thermo_comments: dict = field(default_factory=dict) + ts_labels_with_statmech: frozenset = frozenset() def expected_net_reaction_count(self) -> int: """ @@ -801,6 +813,39 @@ def parse_pdep_network_file(path: str, require_reactions: bool = True) -> PDepNe source_hash=hash_bytes(data), require_reactions=require_reactions) +def _transition_state_declares_statmech(kwargs: dict) -> bool: + """ + Whether one ``transitionState(...)`` call's keywords declare statmech data for the TS. + + Statmech reaches a transition state through exactly two Arkane DSL spellings (see the + ``_call_keywords`` docstring): the inline form's ``modes = [HarmonicOscillator(...), ...]`` + list, and the by-reference form ``transitionState('TS2', 'qm/TS2.py')`` (which + ``_call_keywords`` maps to a ``path`` keyword) that loads a statmech file. An E0-only block -- + ``label``/``E0``/``spinMultiplicity``/``opticalIsomers``, the shape every un-QM'd transition + state takes in an explored reduced network -- has neither, leaves the TS's + ``conformer.modes`` empty, and therefore makes ``Reaction.can_tst()`` False + (``rmgpy/reaction.py``), forcing the master equation onto the inverse Laplace transform for + every k(E) of its channel. That distinction is what + ``t3.pdep.pes_rounds.e0_sensitivity_is_measurable`` classifies on, so it is detected here, + at the parse, from the file's own AST -- never by executing the file. + + The ``modes`` list is inspected structurally (a non-empty ``ast.List``/``ast.Tuple``), not + literal-evaluated: its elements are calls (``HarmonicOscillator(...)``), which + ``ast.literal_eval`` cannot resolve. An empty ``modes = []`` counts as NO statmech -- it + leaves ``conformer.modes`` empty exactly as omitting the keyword does. + + Args: + kwargs (dict): The call's keyword name -> AST value node mapping, from ``_call_keywords``. + + Returns: + bool: Whether the call declares statmech data. + """ + if 'path' in kwargs: + return True + modes_node = kwargs.get('modes') + return isinstance(modes_node, (ast.List, ast.Tuple)) and len(modes_node.elts) > 0 + + def parse_pdep_network_text(text: str, network_id: str, path: str = '', source_hash: str | None = None, require_reactions: bool = True) -> PDepNetwork: @@ -846,6 +891,7 @@ def parse_pdep_network_text(text: str, network_id: str, path: str = '', species_structures = dict() species_thermo_comments = dict() transition_state_labels = list() + ts_labels_with_statmech = list() path_reactions = list() network_label = None isomers = tuple() @@ -879,6 +925,8 @@ def parse_pdep_network_text(text: str, network_id: str, path: str = '', action="omit it from the transition state labels") if label is not None: transition_state_labels.append(label) + if _transition_state_declares_statmech(kwargs): + ts_labels_with_statmech.append(label) elif call_name == 'reaction': path_reactions.append(_parse_reaction(kwargs, path=path)) @@ -937,6 +985,7 @@ def parse_pdep_network_text(text: str, network_id: str, path: str = '', source_hash=source_hash, species_structures=species_structures, species_thermo_comments=species_thermo_comments, + ts_labels_with_statmech=frozenset(ts_labels_with_statmech), ) diff --git a/t3/pdep/pes_loop.py b/t3/pdep/pes_loop.py index f7a5ca54..0d41ed18 100644 --- a/t3/pdep/pes_loop.py +++ b/t3/pdep/pes_loop.py @@ -51,18 +51,21 @@ import logging import os -from dataclasses import dataclass +from dataclasses import asdict, dataclass from pathlib import Path +import yaml + from t3.pdep.api import explore_pdep_network from t3.pdep.capture import captured_qm_artifact_path from t3.pdep.diagram import draw_pes_diagram from t3.pdep.join import arc_ts_label from t3.pdep.explorer.config import PDepExplorerConfig from t3.pdep.explorer.result import EXPLORATION_STATUS_SUCCEEDED -from t3.pdep.parser import parse_pdep_network_file +from t3.pdep.parser import parse_pdep_network_file, to_json_safe from t3.pdep.pes_qm import adopt_prior_qm -from t3.pdep.pes_rounds import (RoundPaths, adoption_channel_keys_by_ts_label, +from t3.pdep.pes_rounds import (SKIP_BELOW_FLOOR, SKIP_NO_EVIDENCE, SKIP_UNMEASURABLE, + RoundPaths, adoption_channel_keys_by_ts_label, attach_sensitivity_evidence, channel_keys_by_ts_label, hybrid_network_path, round_paths, split_qm_candidates) @@ -82,6 +85,17 @@ # that is barrierless end to end, or whose every barriered channel measured a # response below the floor, is not a loop that "converged" after doing # work, it never had any to do. +# - 'unmeasurable' nothing could be queued, but at least one remaining candidate's E0 +# sensitivity was structurally UNCOMPUTABLE rather than measured +# (t3.pdep.pes_rounds.e0_sensitivity_is_measurable: no statmech modes on a +# bimolecular channel, so under ILT the master equation cannot respond to +# that E0 at all -- its reported coefficient is a structural zero, not a +# measurement). This stop point and the queueing decisions are IDENTICAL to +# the 'converged'/'no_candidates' floor stop; only the claim differs, because +# "converged" was not demonstrated -- the screen was blind to those channels. +# What a run SHOULD do about them (rank on barrier height, seed provisional +# statmech, ...) is an open decision (I-030/I-031) that this status +# deliberately reports rather than takes. # - 'stalled' a round ran the QM runner and it returned no newly-converged TS labels, and # ``config.termination.stop_when_no_new_ts`` says that is a reason to stop now # rather than spend the rest of the round budget on a runner that is not making @@ -100,16 +114,31 @@ PES_LOOP_CONVERGED = 'converged' PES_LOOP_MAX_ROUNDS = 'max_rounds' PES_LOOP_NO_CANDIDATES = 'no_candidates' +PES_LOOP_UNMEASURABLE = 'unmeasurable' PES_LOOP_STALLED = 'stalled' PES_LOOP_DIAGRAM_ONLY = 'diagram_only' PES_LOOP_FAILED = 'failed' +# The durable per-round record file, written into each round's own directory +# (``RoundPaths.root``) by ``_record_round``: the ``RoundRecord`` -- including every skipped +# entry's machine-readable classification and NUMERIC coefficient/delta_ln_k fields -- as plain +# YAML. Before this file existed the record was assembled, returned inside ``PESLoopResult``, and +# dropped, while termination messages directed the reader to "the round record" -- a pointer to an +# artifact that was never written (I-031). +ROUND_RECORD_FILENAME = 'round_record.yml' + @dataclass(frozen=True) class RoundRecord: """ What one round of the PES exploration loop did. + Every record is also PERSISTED, at the moment it is made, to its round's own directory as + ``round_record.yml`` (see ``_record_round``): each skipped entry carries its machine-readable + ``classification`` and its numeric ``coefficient``/``delta_ln_k`` as YAML fields, so both a + human and a script can read back what the round measured -- and what it structurally could + not measure -- without regexing prose. + Attributes: index (int): The zero-based round number. network_path (str): The network file this round explored from (before this round's own @@ -248,6 +277,53 @@ def _draw_round_diagram(explored_network_path: str, diagram_path: str, logger) - return diagram_path +def round_record_path(paths: RoundPaths) -> str: + """ + Where one round's durable ``RoundRecord`` YAML lives. + + Args: + paths (RoundPaths): The round's directory layout. + + Returns: + str: The absolute path of the round's ``round_record.yml``. + """ + return os.path.join(paths.root, ROUND_RECORD_FILENAME) + + +def _record_round(rounds: list, record: RoundRecord, paths: RoundPaths, logger) -> None: + """ + Append one round's record to the in-memory list AND persist it to the round's directory. + + The write is best-effort in exactly the sense ``_draw_round_diagram`` is: the record describes + the round's result, it is not part of it, so a full disk or a read-only directory logs a loud + warning rather than flipping a round that genuinely ran into a failed loop. The in-memory + append is unconditional -- ``PESLoopResult.rounds`` never silently loses a round. + + The YAML is written with ``yaml.safe_dump`` over ``to_json_safe(asdict(record))``, so every + skipped entry's ``classification`` and NUMERIC ``coefficient``/``delta_ln_k`` fields land as + real YAML scalars a script can read directly -- a consumer never has to regex the prose + ``reason`` to recover a number (I-031). ``safe_dump`` also guarantees no ``!!python/...`` tag + can ever land in the file, so it reloads with a plain ``yaml.safe_load``. + + Args: + rounds (list): The loop's accumulating ``RoundRecord`` list, appended in place. + record (RoundRecord): The record to append and persist. + paths (RoundPaths): The round's directory layout (the YAML goes to ``paths.root``). + logger: A T3 ``Logger``, or ``None``. + """ + rounds.append(record) + record_path = round_record_path(paths) + try: + with open(record_path, 'w') as f: + yaml.safe_dump(to_json_safe(asdict(record)), f, sort_keys=False) + except Exception as e: + message = (f'PES loop: could not write the round {record.index} record to ' + f'{record_path!r} ({type(e).__name__}): {e}') + _logger.warning(message) + if logger is not None: + logger.warning(message) + + def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, adopted_ts_labels: frozenset | None = None, logger=None) -> PESLoopResult: """ @@ -432,9 +508,10 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, f"{current_network_path!r} does not exist. The previous round's qm_runner " "must write it (see t3.pdep.pes_rounds.hybrid_network_path) before returning.") prior = rounds[-1] if rounds else None - rounds.append(RoundRecord(index=round_index, network_path=current_network_path, - diagram_path=None, queued_ts_labels=(), skipped=(), - status=PES_LOOP_FAILED, reason=reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=current_network_path, + diagram_path=None, queued_ts_labels=(), skipped=(), + status=PES_LOOP_FAILED, reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_FAILED, reason=reason, final_network_path=prior.network_path if prior else None, final_diagram_path=prior.diagram_path if prior else None) @@ -454,9 +531,10 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, reason = '; '.join(result.reasons) if result.reasons else \ f"exploration ended with status {result.status!r} and no stated reason." prior = rounds[-1] if rounds else None - rounds.append(RoundRecord(index=round_index, network_path=current_network_path, - diagram_path=None, queued_ts_labels=(), skipped=(), - status=PES_LOOP_FAILED, reason=reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=current_network_path, + diagram_path=None, queued_ts_labels=(), skipped=(), + status=PES_LOOP_FAILED, reason=reason), + paths, logger) # The same rule as the missing-hybrid branch above: this round failed, but rounds # 0..N-1 explored real networks and drew real diagrams, and those remain the best # result this run has. Reporting None here throws them away and makes the caller (e.g. @@ -488,10 +566,11 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, f'structures, or structurally duplicated channels). QM computed for them ' f'could not be carried across rounds, so every round would re-submit the ' f'same transition states to ARC; refusing to spend the budget on that.') - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=(), - skipped=split_qm_candidates(network, frozenset()).skipped, - status=PES_LOOP_FAILED, reason=reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=(), + skipped=split_qm_candidates(network, frozenset()).skipped, + status=PES_LOOP_FAILED, reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_FAILED, reason=reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -508,9 +587,10 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, round0_has_adopted = round_index == 0 and prior_adopted_artifacts if not split.candidates and not round0_has_adopted: status = PES_LOOP_NO_CANDIDATES if round_index == 0 else PES_LOOP_CONVERGED - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=(), - skipped=split.skipped, status=status)) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=(), + skipped=split.skipped, status=status), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=status, reason='', final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -521,13 +601,14 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, # justify and no capture to evidence), so scope='sensitive' cannot rank here: the # record simply reports the first `limit` candidates in network-file order. offered = split.candidates[:config.qm.max_transition_states_per_round] - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, - queued_ts_labels=tuple(candidate.ts_label - for candidate in offered), - skipped=split.skipped, status=PES_LOOP_DIAGRAM_ONLY, - reason='no qm_runner configured: explored and drew the ' - 'diagram only, nothing was computed.')) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, + queued_ts_labels=tuple(candidate.ts_label + for candidate in offered), + skipped=split.skipped, status=PES_LOOP_DIAGRAM_ONLY, + reason='no qm_runner configured: explored and drew the ' + 'diagram only, nothing was computed.'), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_DIAGRAM_ONLY, reason=rounds[-1].reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -553,10 +634,11 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, f'transition state must carry real sensitivity evidence ' f'(t3.pdep.capture refuses a captured artifact without it), so this ' f'round cannot queue QM.') - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=(), - skipped=split.skipped, status=PES_LOOP_FAILED, - reason=reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=(), + skipped=split.skipped, status=PES_LOOP_FAILED, + reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_FAILED, reason=reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -570,31 +652,81 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, # fault in the evidence pipeline, not a legitimate "nothing worth computing". reason = (f'round {round_index}: the master-equation sensitivity analysis for ' f'network {network.network_id!r} measured no finite sensitivity ' - f'evidence for any QM candidate, so nothing can be queued (see the ' - f"round record's skipped entries for the per-candidate reasons).") - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=(), - skipped=split.skipped, status=PES_LOOP_FAILED, - reason=reason)) + f'evidence for any QM candidate, so nothing can be queued ' + f'(per-candidate reasons: {round_record_path(paths)}).') + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=(), + skipped=split.skipped, status=PES_LOOP_FAILED, + reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_FAILED, reason=reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) - # Every remaining candidate measured a real response BELOW the floor: a measured - # "nothing here is worth the QM spend", which is the same honest terminal outcome - # as having no candidates at all -- each skip carries its measured value in the - # round record. Never 'failed' (nothing malfunctioned) and never a queue of - # below-floor transition states (their manifests would record a coefficient that - # justified nothing). The reason clause is what separates this, at the summary - # level, from "everything was already computed". - status = PES_LOOP_NO_CANDIDATES if round_index == 0 else PES_LOOP_CONVERGED - reason = (f'round {round_index}: every remaining QM candidate of network ' - f'{network.network_id!r} measured a ln(k) response below the ' - f'min_delta_ln_k floor ({config.qm.min_delta_ln_k:.3e}); see the round ' - f"record's skipped entries for the measured values.") - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=(), - skipped=split.skipped, status=status, reason=reason)) + # Nothing survived the evidence screen, so the loop stops HERE either way (the + # same stop point, and the same queueing decisions, as before the unmeasurable + # classification existed -- I-031 changes what is claimed, never what is queued). + # What is claimed depends on WHAT the screen actually established, per the + # classifications attach_sensitivity_evidence just stamped: + # + # * every remaining candidate MEASURED a real response below the floor -- a + # measured "nothing here is worth the QM spend": 'converged' (or + # 'no_candidates' at round 0), exactly as before. Never 'failed' (nothing + # malfunctioned) and never a queue of below-floor transition states (their + # manifests would record a coefficient that justified nothing). + # * at least one remaining candidate was structurally UNMEASURABLE + # (SKIP_UNMEASURABLE): its reported value is a structural zero, not a + # measurement, so "converged" was not demonstrated -- the screen was blind to + # that channel. Reported as 'unmeasurable', loudly, with the count; what to DO + # about such channels is an open decision this loop refuses to take silently. + # + # The reason also names WHICH criterion fired -- indistinguishable stops (floor + # vs. exhausted budget in the same round) are exactly what made the r002 run's + # "converged" unfalsifiable from the outside. + evidence_skips = { + classification: tuple(s for s in split.skipped + if s.classification == classification) + for classification in (SKIP_BELOW_FLOOR, SKIP_UNMEASURABLE, SKIP_NO_EVIDENCE)} + clauses = [] + if evidence_skips[SKIP_BELOW_FLOOR]: + labels = ', '.join(s.ts_label for s in evidence_skips[SKIP_BELOW_FLOOR]) + clauses.append(f'{len(evidence_skips[SKIP_BELOW_FLOOR])} ({labels}) measured ' + f'a ln(k) response below the min_delta_ln_k floor ' + f'({config.qm.min_delta_ln_k:.3e})') + if evidence_skips[SKIP_UNMEASURABLE]: + labels = ', '.join(s.ts_label for s in evidence_skips[SKIP_UNMEASURABLE]) + clauses.append(f'{len(evidence_skips[SKIP_UNMEASURABLE])} ({labels}) are ' + f'structurally UNMEASURABLE -- no statmech modes on a ' + f'bimolecular channel, so the master equation cannot respond ' + f'to their E0 at all and their reported values are structural ' + f'zeros, not measurements') + if evidence_skips[SKIP_NO_EVIDENCE]: + labels = ', '.join(s.ts_label for s in evidence_skips[SKIP_NO_EVIDENCE]) + clauses.append(f'{len(evidence_skips[SKIP_NO_EVIDENCE])} ({labels}) had no ' + f'finite sensitivity row') + remaining = sum(len(entries) for entries in evidence_skips.values()) + budget_note = (f' This was also the last budgeted round (max_rounds: ' + f'{max_rounds}); the criterion above, not the exhausted budget, ' + f'is what stopped the loop.') if round_index == max_rounds - 1 else '' + if evidence_skips[SKIP_UNMEASURABLE]: + status = PES_LOOP_UNMEASURABLE + reason = (f'round {round_index}: none of the {remaining} remaining QM ' + f'candidate(s) of network {network.network_id!r} can be queued: ' + f'{"; ".join(clauses)}. Criterion: unmeasurable candidates remain ' + f'({len(evidence_skips[SKIP_UNMEASURABLE])} of {remaining}), so ' + f'this stop is NOT a demonstrated convergence. Per-candidate ' + f'values: {round_record_path(paths)}.{budget_note}') + else: + status = PES_LOOP_NO_CANDIDATES if round_index == 0 else PES_LOOP_CONVERGED + reason = (f'round {round_index}: every remaining QM candidate of network ' + f'{network.network_id!r} fell to the evidence screen: ' + f'{"; ".join(clauses)}. Criterion: the min_delta_ln_k floor. ' + f'Per-candidate measured values: ' + f'{round_record_path(paths)}.{budget_note}') + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=(), + skipped=split.skipped, status=status, reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=status, reason=reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -651,9 +783,11 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, if config.termination.stop_when_no_new_ts and not made_progress: reason = ('qm_runner returned no newly-converged transition states and ' 'termination.stop_when_no_new_ts is set: continuing would not make progress.') - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=queued_ts_labels, - skipped=split.skipped, status=PES_LOOP_STALLED, reason=reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=queued_ts_labels, + skipped=split.skipped, status=PES_LOOP_STALLED, + reason=reason), + paths, logger) return PESLoopResult(rounds=tuple(rounds), status=PES_LOOP_STALLED, reason=reason, final_network_path=explored_network_path, final_diagram_path=diagram_path) @@ -669,10 +803,11 @@ def run_pes_loop(config: PESLoopConfig, project_directory: str, qm_runner=None, 'qm_runner converged no new transition states this round; ' 're-exploring the same network next round rather than advancing to ' 'a hybrid file that was never written.') - rounds.append(RoundRecord(index=round_index, network_path=explored_network_path, - diagram_path=diagram_path, queued_ts_labels=queued_ts_labels, - skipped=split.skipped, status='continuing', - reason=no_progress_reason)) + _record_round(rounds, RoundRecord(index=round_index, network_path=explored_network_path, + diagram_path=diagram_path, queued_ts_labels=queued_ts_labels, + skipped=split.skipped, status='continuing', + reason=no_progress_reason), + paths, logger) # The next round explores the QM-informed surface: a real qm_runner (Task 6) writes this # round's hybrid network input file (t3.pdep.hybrid) to hybrid_network_path(paths, diff --git a/t3/pdep/pes_rounds.py b/t3/pdep/pes_rounds.py index 7d6c07b4..5148660b 100644 --- a/t3/pdep/pes_rounds.py +++ b/t3/pdep/pes_rounds.py @@ -37,6 +37,63 @@ # loop converts the same handful of structures once per round per consumer. _CANONICAL_STRUCTURE_CACHE: dict = dict() +# The closed vocabulary for a SkippedChannel's ``classification`` -- WHY a channel was not sent to +# ARC, as a machine-readable code alongside the prose ``reason`` (a consumer must never have to +# regex a sentence to recover a fact). The first three come from ``split_qm_candidates``; the last +# three from ``attach_sensitivity_evidence``. ``SKIP_UNMEASURABLE`` is deliberately distinct from +# ``SKIP_BELOW_FLOOR``: a below-floor skip is a MEASURED "not worth the spend", while an +# unmeasurable skip is a channel whose E0 sensitivity the master-equation model structurally +# cannot compute at all (see ``e0_sensitivity_is_measurable``) -- rendering the latter as a +# measurement is a false statement about the physics (I-030/I-031). +SKIP_NO_TRANSITION_STATE = 'no_transition_state' +SKIP_ALREADY_COMPUTED = 'already_computed' +SKIP_BARRIERLESS = 'barrierless' +SKIP_NO_EVIDENCE = 'no_evidence' +SKIP_BELOW_FLOOR = 'below_floor' +SKIP_UNMEASURABLE = 'unmeasurable' + + +def e0_sensitivity_is_measurable(path_reaction: PDepPathReaction, + ts_declares_statmech: bool) -> bool: + """ + Whether a master-equation dln(k)/dE0 sensitivity is structurally computable for this path + reaction's transition state. + + The screen this loop runs perturbs a transition state's ``E0`` and measures the ln(k) + response. That response is structurally zero -- not small, uncomputable -- for a specific + class of candidates (diagnosed and control-verified on the r002 CHO2 run, I-030): + + * A TS with no statmech modes makes ``Reaction.can_tst()`` False (it is literally + ``len(self.transition_state.conformer.modes) > 0``, ``rmgpy/reaction.py``), so RRKM is + never used and every k(E) comes from the inverse Laplace transform + (``rmgpy/pdep/reaction.pyx``). Under ILT the TS ``E0`` enters ONLY as a discrete threshold + gate (``if e_list[r] > E0``), never as a continuous partition-function input. + * For a **bimolecular** channel (an association/dissociation) that gate never binds: the rate + is fixed by the bimolecular asymptote, the detailed-balance clamp, and HPL renormalization, + so perturbing the saddle's ``E0`` is an exact no-op -- the control experiment left every + k(T,P) bit-identical, and the SA reports 0.0 (or ~1e-18 solver roundoff). + * For a **unimolecular** channel the gate CAN bind (a high saddle truncates the k(E) + integrand), so the sensitivity is a real measurement even without modes -- r002 round 0's + isomerization saddle measured a −12% ln(k) response this way and was rightly queued. + + Hence the predicate: measurable iff the TS declares statmech (RRKM engages and ``E0`` is a + continuous input) OR the channel is unimolecular on both sides. This is the STRUCTURAL test, + chosen over inspecting the SA result, because it does not depend on a structural zero + surviving floating-point arithmetic (the r002 entrance channels measured ~1e-18 in one round + and exactly 0.0 in the next -- no epsilon separates those from a genuinely small response). + + Args: + path_reaction (PDepPathReaction): The parsed path reaction. + ts_declares_statmech (bool): Whether the network file declares statmech data for this + reaction's TS (``PDepNetwork.ts_labels_with_statmech``). + + Returns: + bool: Whether the ME E0 sensitivity is structurally computable for this candidate. + """ + if ts_declares_statmech: + return True + return len(path_reaction.reactants) == 1 and len(path_reaction.products) == 1 + @dataclass(frozen=True) class QMCandidate: @@ -56,12 +113,19 @@ class QMCandidate: delta_ln_k (float | None): The corresponding dimensionless rate response, ``abs(coefficient) * perturbation`` -- same convention as ``t3.pdep.selector.SensitiveTransitionState``. + e0_sensitivity_measurable (bool): Whether the ME E0 sensitivity is structurally + computable for this candidate (``e0_sensitivity_is_measurable``), stamped by + ``split_qm_candidates`` from the network file's own declarations. Classification + only: it never gates queueing (an unmeasurable candidate is skipped by exactly the + same evidence/floor tests as before -- it is just no longer recorded as having + "measured" its structural zero). """ path_reaction: PDepPathReaction ts_label: str family: str | None coefficient: float | None = None delta_ln_k: float | None = None + e0_sensitivity_measurable: bool = True @dataclass(frozen=True) @@ -71,10 +135,25 @@ class SkippedChannel: Attributes: label (str): The path reaction's label. - reason (str): Why it was skipped, for the log and the round record. + reason (str): Why it was skipped, for the log and the round record -- prose, for a human. + ts_label (str | None): The network-local transition state label, when the channel declares + one (``None`` for a ``SKIP_NO_TRANSITION_STATE`` skip). + classification (str): The machine-readable skip class, one of the ``SKIP_*`` constants. + ``''`` only for records predating this field. + coefficient (float | None): The signed dln(k)/dE0 value (mol/J) the round's ME SA + reported for this TS, when the skip happened at the evidence stage and the SA + reported one. For a ``SKIP_UNMEASURABLE`` skip this is the reported STRUCTURAL ZERO, + carried verbatim (never clamped or substituted) -- the ``classification`` is what + says it is not a measurement. + delta_ln_k (float | None): The corresponding dimensionless ``abs(coefficient) * + perturbation`` response, same provenance and caveat as ``coefficient``. """ label: str reason: str + ts_label: str | None = None + classification: str = '' + coefficient: float | None = None + delta_ln_k: float | None = None @dataclass(frozen=True) @@ -109,20 +188,25 @@ def split_qm_candidates(network: PDepNetwork, computed_ts_labels: frozenset) -> skipped.append(SkippedChannel( label=path_reaction.label, reason=f"'{path_reaction.label}': declares no transition state in the network " - f'file, so there is nothing to compute.')) + f'file, so there is nothing to compute.', + classification=SKIP_NO_TRANSITION_STATE)) continue if ts_label in computed_ts_labels: skipped.append(SkippedChannel( label=path_reaction.label, reason=f"'{path_reaction.label}': transition state {ts_label} already has QM; not " - f'queueing it again.')) + f'queueing it again.', + ts_label=ts_label, classification=SKIP_ALREADY_COMPUTED)) continue verdict = classify_barrierless(path_reaction) if verdict.is_barrierless: - skipped.append(SkippedChannel(label=path_reaction.label, reason=verdict.reason)) + skipped.append(SkippedChannel(label=path_reaction.label, reason=verdict.reason, + ts_label=ts_label, classification=SKIP_BARRIERLESS)) continue - candidates.append(QMCandidate(path_reaction=path_reaction, ts_label=ts_label, - family=verdict.family)) + candidates.append(QMCandidate( + path_reaction=path_reaction, ts_label=ts_label, family=verdict.family, + e0_sensitivity_measurable=e0_sensitivity_is_measurable( + path_reaction, ts_label in network.ts_labels_with_statmech))) return CandidateSplit(candidates=tuple(candidates), skipped=tuple(skipped)) @@ -148,6 +232,14 @@ def attach_sensitivity_evidence(split: CandidateSplit, definitionally unreachable. The floor applies under BOTH ``qm.scope`` values -- 'sensitive' ranks and 'all' does not, but neither may queue below it. + A skipped candidate that is structurally UNMEASURABLE (``QMCandidate.e0_sensitivity_measurable`` + False -- see ``e0_sensitivity_is_measurable``) takes the SAME two skip branches, so queueing is + unchanged, but is classified ``SKIP_UNMEASURABLE`` rather than ``SKIP_NO_EVIDENCE``/ + ``SKIP_BELOW_FLOOR``, with a reason that says the model could not respond to its E0 at all: + its reported value is a structural zero, and recording it as a measurement would be a false + statement about the physics. The reported value itself is still carried verbatim on the + skip's numeric fields, never clamped or substituted. + Args: split (CandidateSplit): The split to stamp, from ``split_qm_candidates``. evidence_by_ts_label (dict): Network-local TS label -> ``(coefficient, delta_ln_k)``, @@ -163,23 +255,60 @@ def attach_sensitivity_evidence(split: CandidateSplit, candidates, skipped = [], list(split.skipped) for candidate in split.candidates: pair = evidence_by_ts_label.get(candidate.ts_label) + # The unmeasurable classification (e0_sensitivity_is_measurable, stamped by + # split_qm_candidates) NEVER moves a candidate between the queue and the skips: the + # evidence/floor tests below are exactly the ones that decided before it existed, so a + # run makes the same QM decisions -- the classification only changes what the record SAYS + # about a skip that was happening anyway (a structural zero is not a measurement). + unmeasurable_clause = ( + 'its E0 sensitivity is structurally uncomputable, not small: its ' + 'transitionState(...) block declares no statmech modes, so Reaction.can_tst() is ' + 'False and every k(E) of its bimolecular channel comes from the inverse Laplace ' + 'transform, where the TS E0 enters only as a threshold gate that never binds ' + '(rmgpy/pdep/reaction.pyx) -- perturbing this E0 is an exact no-op') if pair is None: + if not candidate.e0_sensitivity_measurable: + reason = (f"'{candidate.path_reaction.label}': transition state " + f'{candidate.ts_label} could not measure an E0 sensitivity: ' + f'{unmeasurable_clause}, and this round\'s master-equation sensitivity ' + f'analysis reported no finite row for it; not queueing it rather than ' + f'inventing a number.') + skipped.append(SkippedChannel( + label=candidate.path_reaction.label, reason=reason, + ts_label=candidate.ts_label, classification=SKIP_UNMEASURABLE)) + continue skipped.append(SkippedChannel( label=candidate.path_reaction.label, reason=f"'{candidate.path_reaction.label}': transition state " f'{candidate.ts_label} has no finite sensitivity evidence in this ' f"round's master-equation sensitivity analysis, and a captured artifact " f'must carry the evidence that justified selecting it ' - f'(t3.pdep.capture); not queueing it rather than inventing a number.')) + f'(t3.pdep.capture); not queueing it rather than inventing a number.', + ts_label=candidate.ts_label, classification=SKIP_NO_EVIDENCE)) continue coefficient, delta_ln_k = pair if delta_ln_k < min_delta_ln_k: + if not candidate.e0_sensitivity_measurable: + reason = (f"'{candidate.path_reaction.label}': transition state " + f'{candidate.ts_label} could not measure an E0 sensitivity: ' + f'{unmeasurable_clause}. The SA reported {delta_ln_k:.3e}, which is a ' + f'structural zero, not a measurement; it falls below the ' + f'min_delta_ln_k floor ({min_delta_ln_k:.3e}) and the channel is ' + f'skipped exactly as before, but "below the floor" is not a fact this ' + f'screen established about it.') + skipped.append(SkippedChannel( + label=candidate.path_reaction.label, reason=reason, + ts_label=candidate.ts_label, classification=SKIP_UNMEASURABLE, + coefficient=coefficient, delta_ln_k=delta_ln_k)) + continue skipped.append(SkippedChannel( label=candidate.path_reaction.label, reason=f"'{candidate.path_reaction.label}': transition state " f'{candidate.ts_label} measured a ln(k) response of {delta_ln_k:.3e}, ' f'below the min_delta_ln_k floor ({min_delta_ln_k:.3e}); its leverage on ' - f'the network is too small to justify the QM spend.')) + f'the network is too small to justify the QM spend.', + ts_label=candidate.ts_label, classification=SKIP_BELOW_FLOOR, + coefficient=coefficient, delta_ln_k=delta_ln_k)) continue candidates.append(replace(candidate, coefficient=coefficient, delta_ln_k=delta_ln_k)) return CandidateSplit(candidates=tuple(candidates), skipped=tuple(skipped)) diff --git a/tests/test_pdep/test_parser.py b/tests/test_pdep/test_parser.py index 444c88ac..d8688117 100644 --- a/tests/test_pdep/test_parser.py +++ b/tests/test_pdep/test_parser.py @@ -1431,3 +1431,65 @@ def test_positional_exception_is_narrow_mixed_positional_keyword_refused(self): (call,) = [n.value for n in ast.parse("transitionState('TS2', path='qm/TS2.py')").body] with pytest.raises(ValueError, match='POSITIONAL'): _call_keywords(call, call_name='transitionState') + + +class TestTransitionStateStatmechDetection(object): + """``PDepNetwork.ts_labels_with_statmech``: which ``transitionState(...)`` blocks declare + statmech data (I-031). A TS absent from this set has an empty ``conformer.modes``, so + ``Reaction.can_tst()`` is False and every k(E) of its channel comes from the inverse Laplace + transform -- the fact ``t3.pdep.pes_rounds.e0_sensitivity_is_measurable`` classifies on.""" + + E0_ONLY = ("transitionState(\n" + " label = 'TS1',\n" + " E0 = (21.347, 'kJ/mol'),\n" + " spinMultiplicity = 1,\n" + " opticalIsomers = 1,\n" + ")\n") + WITH_MODES = ("transitionState(\n" + " label = 'TS2',\n" + " E0 = (135.3, 'kJ/mol'),\n" + " modes = [HarmonicOscillator(frequencies=([1000.0, 2000.0], 'cm^-1'))],\n" + " spinMultiplicity = 1,\n" + " opticalIsomers = 1,\n" + ")\n") + BY_REFERENCE = "transitionState('TS3', 'qm/TS3.py')\n" + EMPTY_MODES = ("transitionState(\n" + " label = 'TS4',\n" + " E0 = (10.0, 'kJ/mol'),\n" + " modes = [],\n" + ")\n") + + def _parse(self, ts_blocks: str): + return parse_pdep_network_text(text=ts_blocks, network_id='synthetic_statmech', + require_reactions=False) + + def test_e0_only_block_declares_no_statmech(self): + """The exact shape every un-QM'd TS takes in an explored reduced network (verified on + r002 round 1's network0_reduced.py: E0/spinMultiplicity/opticalIsomers only).""" + network = self._parse(self.E0_ONLY) + assert network.transition_state_labels == ('TS1',) + assert network.ts_labels_with_statmech == frozenset() + + def test_non_empty_modes_list_declares_statmech(self): + """The inline form ``t3.pdep.hybrid`` writes for a QM'd TS. The list elements are calls + (``HarmonicOscillator(...)``), so the detection must be structural, never literal_eval.""" + network = self._parse(self.WITH_MODES) + assert network.ts_labels_with_statmech == frozenset({'TS2'}) + + def test_by_reference_path_form_declares_statmech(self): + """``transitionState(label, path)`` loads a statmech file (arkane/input.py); the TS it + names will have modes once Arkane runs, so it is statmech-bearing.""" + network = self._parse(self.BY_REFERENCE) + assert network.ts_labels_with_statmech == frozenset({'TS3'}) + + def test_empty_modes_list_declares_no_statmech(self): + """``modes = []`` leaves ``conformer.modes`` empty exactly as omitting the keyword does, + and an empty-modes TS still fails ``can_tst()``; counting it would relabel a structurally + unmeasurable channel as measured.""" + network = self._parse(self.EMPTY_MODES) + assert network.ts_labels_with_statmech == frozenset() + + def test_mixed_file_keeps_the_sets_apart(self): + network = self._parse(self.E0_ONLY + self.WITH_MODES + self.BY_REFERENCE + self.EMPTY_MODES) + assert network.transition_state_labels == ('TS1', 'TS2', 'TS3', 'TS4') + assert network.ts_labels_with_statmech == frozenset({'TS2', 'TS3'}) diff --git a/tests/test_pdep/test_pes_loop.py b/tests/test_pdep/test_pes_loop.py index 3b9c90a6..1cfb4ee3 100644 --- a/tests/test_pdep/test_pes_loop.py +++ b/tests/test_pdep/test_pes_loop.py @@ -5,6 +5,7 @@ import os import pytest +import yaml from arc.molecule.molecule import Molecule @@ -15,6 +16,7 @@ from t3.pdep.parser import PDepNetwork, PDepPathReaction from t3.pdep.pes_loop import (PES_LOOP_CONVERGED, PES_LOOP_DIAGRAM_ONLY, PES_LOOP_FAILED, PES_LOOP_MAX_ROUNDS, PES_LOOP_NO_CANDIDATES, PES_LOOP_STALLED, + PES_LOOP_UNMEASURABLE, _build_explorer_config, hybrid_network_path, run_pes_loop) from t3.pdep.pes_rounds import round_paths @@ -1222,3 +1224,192 @@ def _explode(network_path, output_path): assert result.final_diagram_path is None assert result.final_network_path is not None, \ 'the exploration succeeded; only the drawing failed' + + +def _stub_network(monkeypatch, tmp_path, path_reactions, species_structures, + ts_labels_with_statmech=frozenset()): + """Patch the explorer seams with one fixed, caller-built network -- the bimolecular / + statmech-bearing shapes ``_stub_explorer``'s fabricated unimolecular channels cannot + express (I-031).""" + network = PDepNetwork(network_id='network1_1', path='/abs/network1_1.py', label=None, + species_labels=tuple(species_structures), + transition_state_labels=tuple(r.transition_state + for r in path_reactions + if r.transition_state), + path_reactions=tuple(path_reactions), isomers=(), reactant_channels=(), + product_channels=(), species_structures=species_structures, + ts_labels_with_statmech=ts_labels_with_statmech) + calls = [] + + def _fake_explore(*, network_path, config, logger=None): + calls.append(network_path) + output_path = os.path.join(str(tmp_path), f'explored_round{len(calls) - 1}.py') + open(output_path, 'w').close() + return PDepExplorationResult(network_id='network1_1', + status=EXPLORATION_STATUS_SUCCEEDED, + network_paths=(output_path,)) + + def _fake_parse(path, require_reactions=True): + return dataclasses.replace(network, path=path) + + def _fake_draw(network_path, output_path): + open(output_path, 'w').close() + + monkeypatch.setattr('t3.pdep.pes_loop.explore_pdep_network', _fake_explore) + monkeypatch.setattr('t3.pdep.pes_loop.parse_pdep_network_file', _fake_parse) + monkeypatch.setattr('t3.pdep.pes_loop.draw_pes_diagram', _fake_draw) + + +def _entrance_channel_fixture(): + """Two bimolecular association channels plus one unimolecular isomerization, no TS declaring + statmech -- the r002 round-1 shape.""" + adj = {'H': _alkane_adjlist(1), 'CO2': _alkane_adjlist(2), + 'HOCO': _alkane_adjlist(3), 'HCO2': _alkane_adjlist(4)} + + def _reaction(label, ts, reactants, products): + return PDepPathReaction(label=label, reactants=reactants, products=products, + transition_state=ts, kinetics_type='Arrhenius', + kinetics_comment='family: R_Addition_MultipleBond') + + return (_reaction('reaction1', 'TS1', ('HOCO',), ('HCO2',)), + _reaction('reaction2', 'TS2', ('H', 'CO2'), ('HOCO',)), + _reaction('reaction3', 'TS3', ('H', 'CO2'), ('HCO2',))), adj + + +class TestUnmeasurableTermination(object): + """I-031: a stop where the remaining candidates could not MEASURE their E0 sensitivity is + reported as 'unmeasurable', never as convergence; the message names which criterion fired + and how many candidates were blind; queueing is untouched.""" + + @staticmethod + def _never_called_runner(candidates, paths, cfg, network_id, adopted=None): + raise AssertionError('qm_runner must not be called when nothing survives the screen.') + + def test_all_remaining_unmeasurable_is_not_convergence(self, tmp_path, monkeypatch, config): + reactions, adj = _entrance_channel_fixture() + _stub_network(monkeypatch, tmp_path, reactions, adj) + + def _fake_sa(*, network_path, sa_dir, method, timeout=None, logger=None): + return {'TS1': (0.0, 0.0), 'TS2': (1.66e-18, 1.39e-14), 'TS3': (0.0, 0.0)} + + monkeypatch.setattr('t3.pdep.pes_loop.run_round_me_sensitivity', _fake_sa) + result = run_pes_loop(config, project_directory=str(tmp_path), + qm_runner=self._never_called_runner) + assert result.status == PES_LOOP_UNMEASURABLE + assert result.status != PES_LOOP_CONVERGED and result.status != PES_LOOP_NO_CANDIDATES + assert '2 (TS2, TS3) are structurally UNMEASURABLE' in result.reason + assert 'NOT a demonstrated convergence' in result.reason + assert '1 (TS1) measured a ln(k) response below the min_delta_ln_k floor' in result.reason + assert 'round_record.yml' in result.reason + classes = {s.ts_label: s.classification for s in result.rounds[0].skipped} + assert classes == {'TS1': 'below_floor', 'TS2': 'unmeasurable', 'TS3': 'unmeasurable'} + + def test_statmech_bearing_channels_keep_the_converged_claim(self, tmp_path, monkeypatch, + config): + """Mutation control at the loop level: the SAME network with statmech declared for the + bimolecular TSs is measurable, so the same all-below-floor stop IS a measured 'done'.""" + reactions, adj = _entrance_channel_fixture() + _stub_network(monkeypatch, tmp_path, reactions, adj, + ts_labels_with_statmech=frozenset({'TS2', 'TS3'})) + + def _fake_sa(*, network_path, sa_dir, method, timeout=None, logger=None): + return {'TS1': (0.0, 0.0), 'TS2': (1.66e-18, 1.39e-14), 'TS3': (0.0, 0.0)} + + monkeypatch.setattr('t3.pdep.pes_loop.run_round_me_sensitivity', _fake_sa) + result = run_pes_loop(config, project_directory=str(tmp_path), + qm_runner=self._never_called_runner) + assert result.status == PES_LOOP_NO_CANDIDATES + assert 'Criterion: the min_delta_ln_k floor' in result.reason + + def test_the_floor_stop_names_its_criterion_and_the_exhausted_budget(self, tmp_path, + monkeypatch): + """When the floor fires on the LAST budgeted round, 'converged' and 'ran out of rounds' + were previously indistinguishable from outside; the message now says which fired.""" + config = PESLoopConfig(pes={'network': '/abs/network1_1.py', 'source': ['HOCHO'], + 'bath_gas': {'He': 1.0}}, + termination={'max_rounds': 1}) + _stub_explorer(monkeypatch, tmp_path, families=['1,2_Insertion_CO']) + monkeypatch.setattr('t3.pdep.pes_loop.run_round_me_sensitivity', + lambda **kwargs: {'TS0': (0.0, 0.0)}) + result = run_pes_loop(config, project_directory=str(tmp_path), + qm_runner=self._never_called_runner) + assert result.status == PES_LOOP_NO_CANDIDATES + assert 'Criterion: the min_delta_ln_k floor' in result.reason + assert 'last budgeted round (max_rounds: 1)' in result.reason + assert 'not the exhausted budget' in result.reason + + +class TestRoundRecordPersistence(object): + """I-031: the round record is WRITTEN, not just returned -- one YAML per round, in the + round's own directory, with numeric fields a script reads without regexing prose.""" + + def test_terminal_round_record_is_persisted_with_numeric_fields(self, tmp_path, monkeypatch, + config): + reactions, adj = _entrance_channel_fixture() + _stub_network(monkeypatch, tmp_path, reactions, adj) + + def _fake_sa(*, network_path, sa_dir, method, timeout=None, logger=None): + return {'TS1': (0.0, 0.0), 'TS2': (1.66e-18, 1.39e-14), 'TS3': (0.0, 0.0)} + + monkeypatch.setattr('t3.pdep.pes_loop.run_round_me_sensitivity', _fake_sa) + result = run_pes_loop(config, project_directory=str(tmp_path), + qm_runner=TestUnmeasurableTermination._never_called_runner) + record_path = os.path.join(round_paths(str(tmp_path), 0).root, 'round_record.yml') + assert os.path.isfile(record_path) + with open(record_path) as f: + payload = yaml.safe_load(f) + assert payload['index'] == 0 + assert payload['status'] == PES_LOOP_UNMEASURABLE + assert payload['reason'] == result.reason + assert payload['queued_ts_labels'] == [] + by_ts = {entry['ts_label']: entry for entry in payload['skipped']} + assert by_ts['TS2']['classification'] == 'unmeasurable' + assert by_ts['TS2']['coefficient'] == 1.66e-18 + assert by_ts['TS2']['delta_ln_k'] == 1.39e-14 + assert by_ts['TS1']['classification'] == 'below_floor' + assert by_ts['TS1']['coefficient'] == 0.0 + + def test_every_round_gets_its_own_record_file(self, tmp_path, monkeypatch, config): + """A converging two-round loop writes round_0 and round_1 records matching the returned + RoundRecords -- the 'continuing' round is persisted too, not only the terminal one.""" + _stub_explorer(monkeypatch, tmp_path, families=['1,2_Insertion_CO', '1,2_Insertion_CO']) + + def _runner(candidates, paths, cfg, network_id, adopted=None): + _touch_hybrid_file(paths, network_id) + return frozenset(c.ts_label for c in candidates), frozenset(c.ts_label + for c in candidates) + + result = run_pes_loop(config, project_directory=str(tmp_path), qm_runner=_runner) + assert result.status == PES_LOOP_CONVERGED + for record in result.rounds: + record_path = os.path.join(round_paths(str(tmp_path), record.index).root, + 'round_record.yml') + assert os.path.isfile(record_path) + with open(record_path) as f: + payload = yaml.safe_load(f) + assert payload['index'] == record.index + assert payload['status'] == record.status + assert payload['queued_ts_labels'] == list(record.queued_ts_labels) + assert [entry['label'] for entry in payload['skipped']] == \ + [skip.label for skip in record.skipped] + + def test_a_failed_record_write_never_fails_the_round(self, tmp_path, monkeypatch, config, + caplog): + """Best-effort exactly as the diagram: the record describes the result, it is not part + of it.""" + _stub_explorer(monkeypatch, tmp_path, families=['1,2_Insertion_CO']) + + def _runner(candidates, paths, cfg, network_id, adopted=None): + _touch_hybrid_file(paths, network_id) + return frozenset(c.ts_label for c in candidates), frozenset(c.ts_label + for c in candidates) + + def _refuse_dump(*args, **kwargs): + raise OSError('disk full') + + monkeypatch.setattr('t3.pdep.pes_loop.yaml.safe_dump', _refuse_dump) + with caplog.at_level(logging.WARNING, logger='t3.pdep.pes_loop'): + result = run_pes_loop(config, project_directory=str(tmp_path), qm_runner=_runner) + assert result.status == PES_LOOP_CONVERGED + assert len(result.rounds) == 2 + assert 'could not write the round 0 record' in caplog.text diff --git a/tests/test_pdep/test_pes_rounds.py b/tests/test_pdep/test_pes_rounds.py index cf6a96e9..1fa265a6 100644 --- a/tests/test_pdep/test_pes_rounds.py +++ b/tests/test_pdep/test_pes_rounds.py @@ -12,8 +12,11 @@ from t3.pdep import barrierless from t3.pdep.parser import PDepNetwork, PDepPathReaction, parse_pdep_network_file from t3.pdep.pes_rounds import (CandidateSplit, PES_LOOP_DIAGRAM_FILENAME, + SKIP_ALREADY_COMPUTED, SKIP_BARRIERLESS, SKIP_BELOW_FLOOR, + SKIP_NO_EVIDENCE, SKIP_NO_TRANSITION_STATE, SKIP_UNMEASURABLE, adoption_channel_keys_by_ts_label, attach_sensitivity_evidence, - channel_keys_by_ts_label, round_paths, split_qm_candidates, + channel_keys_by_ts_label, e0_sensitivity_is_measurable, + round_paths, split_qm_candidates, structural_channel_key) @@ -435,3 +438,109 @@ def test_a_structurally_unkeyable_channel_is_still_refused(self): network = self._network(['family: H_Abstraction']) network = dataclasses.replace(network, species_structures={'A': 'not an adjacency list'}) assert adoption_channel_keys_by_ts_label(network) == {} + + +def _bimolecular_rxn(label: str, ts: str) -> PDepPathReaction: + return PDepPathReaction(label=label, reactants=('H', 'CO2'), products=('HOCO',), + transition_state=ts, kinetics_type='Arrhenius', + kinetics_comment='family: R_Addition_MultipleBond') + + +class TestE0SensitivityMeasurability(object): + """``e0_sensitivity_is_measurable`` (I-031): the ME E0 sensitivity is structurally + uncomputable -- not small -- exactly for a no-statmech TS on a bimolecular channel (ILT's + threshold gate never binds for an association; verified by control experiment on r002).""" + + def test_no_modes_bimolecular_is_unmeasurable(self): + assert e0_sensitivity_is_measurable(_bimolecular_rxn('r1', 'TS1'), + ts_declares_statmech=False) is False + + def test_no_modes_unimolecular_is_measurable(self): + """The ILT threshold gate CAN bind for a high unimolecular saddle -- r002 round 0's + isomerization saddle measured a real -12% ln(k) response without modes and was queued, + so 'no modes' alone must never classify a channel unmeasurable.""" + assert e0_sensitivity_is_measurable(_rxn('r1', 'family: 1,2_Insertion_CO'), + ts_declares_statmech=False) is True + + def test_statmech_makes_any_channel_measurable(self): + """With modes, can_tst() is True, RRKM engages, and E0 is a continuous input.""" + assert e0_sensitivity_is_measurable(_bimolecular_rxn('r1', 'TS1'), + ts_declares_statmech=True) is True + + def test_split_stamps_the_classification_from_the_network_declarations(self): + network = dataclasses.replace( + _network([_bimolecular_rxn('r1', 'TS1'), _bimolecular_rxn('r2', 'TS2'), + _rxn('r3', 'family: 1,2_Insertion_CO', ts='TS3')]), + ts_labels_with_statmech=frozenset({'TS2'})) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + measurable = {c.ts_label: c.e0_sensitivity_measurable for c in split.candidates} + assert measurable == {'TS1': False, 'TS2': True, 'TS3': True} + + +class TestSkippedChannelClassification(object): + """Every skip carries a machine-readable ``classification`` and, at the evidence stage, the + NUMERIC coefficient/delta_ln_k -- a consumer must never regex the prose ``reason`` (I-031).""" + + def test_split_skips_carry_their_classes(self): + network = _network([_rxn('r1', 'family: 1,2_Insertion_CO', ts='TS_done'), + _rxn('r2', 'in family R_Recombination.'), + PDepPathReaction(label='r3', reactants=('A',), products=('B',), + transition_state=None, kinetics_type='Arrhenius', + kinetics_comment='')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset({'TS_done'})) + classes = {s.label: s.classification for s in split.skipped} + assert classes == {'r1': SKIP_ALREADY_COMPUTED, 'r2': SKIP_BARRIERLESS, + 'r3': SKIP_NO_TRANSITION_STATE} + ts_labels = {s.label: s.ts_label for s in split.skipped} + assert ts_labels == {'r1': 'TS_done', 'r2': 'TS_r2', 'r3': None} + + def test_unmeasurable_below_floor_is_classified_and_carries_the_structural_zero(self): + """The reported value is carried VERBATIM (never clamped or hidden); the classification + -- not a substituted number -- is what says it is not a measurement.""" + network = _network([_bimolecular_rxn('r1', 'TS1')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + split = attach_sensitivity_evidence(split, {'TS1': (1.66e-18, 1.39e-14)}, + min_delta_ln_k=1e-3) + assert split.candidates == () + (skip,) = split.skipped + assert skip.classification == SKIP_UNMEASURABLE + assert skip.coefficient == 1.66e-18 + assert skip.delta_ln_k == 1.39e-14 + assert 'structurally uncomputable' in skip.reason + assert 'not a measurement' in skip.reason + assert 'measured a ln(k) response' not in skip.reason + + def test_measurable_below_floor_keeps_the_measured_classification(self): + network = _network([_rxn('r1', 'family: 1,2_Insertion_CO')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + split = attach_sensitivity_evidence(split, {'TS_r1': (0.0, 0.0)}, min_delta_ln_k=1e-3) + (skip,) = split.skipped + assert skip.classification == SKIP_BELOW_FLOOR + assert skip.coefficient == 0.0 + assert skip.delta_ln_k == 0.0 + assert 'below the min_delta_ln_k floor' in skip.reason + + def test_unmeasurable_with_no_evidence_row_is_still_unmeasurable(self): + network = _network([_bimolecular_rxn('r1', 'TS1')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + split = attach_sensitivity_evidence(split, {}, min_delta_ln_k=1e-3) + (skip,) = split.skipped + assert skip.classification == SKIP_UNMEASURABLE + assert skip.coefficient is None and skip.delta_ln_k is None + + def test_measurable_with_no_evidence_row_is_no_evidence(self): + network = _network([_rxn('r1', 'family: 1,2_Insertion_CO')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + split = attach_sensitivity_evidence(split, {}, min_delta_ln_k=1e-3) + (skip,) = split.skipped + assert skip.classification == SKIP_NO_EVIDENCE + + def test_classification_never_gates_queueing(self): + """THE non-goal pin: an unmeasurable candidate above the floor is queued exactly as + before -- the classification changes what the record says, never what gets queued.""" + network = _network([_bimolecular_rxn('r1', 'TS1')]) + split = split_qm_candidates(network, computed_ts_labels=frozenset()) + split = attach_sensitivity_evidence(split, {'TS1': (-3.0e-5, 0.25)}, + min_delta_ln_k=1e-3) + assert [c.ts_label for c in split.candidates] == ['TS1'] + assert split.skipped == ()