Skip to content

0.6.1 — selective REST2: choosing the hot region

Validated on CUDA; the energy-check tolerance was rewritten and half of its replacement is refuted

This page describes v0.6.1, tagged on the commit that carries it. The CUDA evidence exists: tests/test_selective_rest2_cuda.py passed 12 tests with 0 skips under --error-on-skip, and two 8-rung TYK2 ladders ran to completion on four cards (ligand only, and ligand plus 19 pocket sidechains; 5 ns per state, 40 ns aggregate each).

What is NOT clean: gate G2-6, which recomputes the recorded exchange energies against the saved states, is a recorded FAIL on both TYK2 ladders — 0.077/0.093 and 0.106/0.121 kT against a 0.05 kT tolerance that had been calibrated on a system 30x smaller. The failure is the tolerance, not the Hamiltonian: the discrepancy is 3.1e-7 relative where the fixture's was 7.3e-7, and a deliberately wrong Hamiltonian misses by 8,000–15,500 kT.

Its replacement, G2-6b, was registered with a falsifiable prediction and has now been measured on a third system size (10,972 particles, |u| = 64,994 kT, used to calibrate neither of the other two). Its upper bound — the half that gates anything — holds, with 8.6x margin on the absolute discrepancy and 5.0x on the cross-Hamiltonian one, and a deliberately wrong Hamiltonian is caught by 3,478x the bound against a 100x requirement. That is a real validation of the size-aware criterion at a size neither ladder nor fixture defined.

Its lower bound is refuted as written, and that matters more for the next version than for this one: the measured discrepancy came in below the floor the row predicted, and the floor is not being lowered to fit. The finding underneath it is that the statistic itself is wrong — a max over a finite sample is biased low, and the same system gave 0.00563 kT at N=9 and 0.00766 kT at N=64. A future G2-6c must fix the estimator (a stated N, coordinates from a stated stage, a quantile or RMS rather than a max) before any level is set again. Against the validated upper bound, ladder A passes and ladder B's cross-Hamiltonian discrepancy is outside it by 1.16x; both ladders' G2-6 status stays FAIL and their figures stand as measured.

The full gate lane is green on c19988c: fast 2385 passed / 3 skipped, slow 438 / 7, gpu 123 passed / 1 skipped — the single gpu skip is the --cuda-evidence one, so nothing skipped silently, and the MPS daemon was torn down with its directory verified gone.

What the gate does not cover: the TYK2 ladder rows. No test in any lane recomputes a production ladder's exchange energies, so a green gate says nothing about them. Ladder B sits 1.16x outside the G2-6b bound and both ladders stay recorded FAIL.

Nothing is removed or renamed, and no existing configuration changes meaning. A REST2 ladder with no new key resolves exactly as it did in 0.6.0, its scaled states are byte-identical, and a 0.6.0 ladder still resumes.

The hot region can be chosen

Until now a REST2 ladder heated the whole solute: every atom that is not solvent or a counter-ion. For anything larger than a peptide that is wasteful, and often counter-productive — the degrees of freedom you care about are a loop, a few sidechains, or one ligand, and heating the rest costs sampling and can unfold what you wanted kept folded.

Three optional keys in scaler.config choose a region instead:

backbone_scaling_list: ":45,46,59"
sidechain_scaling_list: ":45-50"
ligand_scaling_dict:
  L01:
    mask: ":201"
    torsion_exclusions: auto        # or an md-tools-torsion-exclusions/1 file

With no selector the selection is the legacy full solute. With any selector it is explicit: only the named categories and instances are hot, an omitted category means none, and a region that heats nothing is refused. Selective scaling is explicit-solvent only; a selector on an implicit System is refused, because a generalised-Born energy is not separable per atom.

Masks are a documented subset of AMBER's grammar (md-tools-residue-mask/1): ":45", ":45,46,59", ":45-50", ":45,48-52", numbered by one-based topology residue index across every chain without a reset. Everything else — atom and chain selectors, operators, wildcards, residue names — is refused by name with the supported subset printed.

Two questions, answered separately

Which atoms carry the nonbonded factors and which torsions are scaled are different questions, and a selective region answers them separately (md-tools-selective-rest2/1):

  • Atoms. Backbone atoms of backbone-selected residues, sidechain atoms of sidechain-selected residues, every atom of a selected ligand instance. S–S pairs and exceptions carry (1−τ)², S–E (1−τ), E–E 1.
  • Torsions. Every eligible proper term across a selected central bond, where each bond has a documented owner: φ/ψ to the residue's backbone, any bond touching a sidechain atom to the sidechain (so χ1 is the sidechain's although its quartet reaches backbone N and CA), the peptide bond to the residue that owns ω, a disulfide jointly to both sidechains, and anything else between residues refused rather than guessed.
  • CMAP is scaled only when both of a term's torsions are selected. A map shared with an unselected residue is duplicated, never scaled in place.

The existing exclusions are untouched: ordinary amide ω, aromatic rings, other double bonds and every improper stay unscaled, and the proline-like exception is preserved.

Ligands are selected per copy

mask names one residue, so selecting one copy of a compound never selects its twin. A per-instance torsion_exclusions file (md-tools-torsion-exclusions/1) adds protected bonds and can never un-protect one; it is bound to a parameter package, written in package-local atom names, and its resolved contents are copied into the record, because a path is not provenance.

The compact form L01: <path> is refused for now, with the explicit form printed: the ligand mapping record has no instance alias, and a residue name cannot say which copy is meant.

Two helpers for choosing a region

  • Which residues line a ligand: python -m md_tools.rest2.pocket build/built.pdb --ligand ":201" --cutoff-nm 0.5 prints the criterion it used, a table (index, chain, author residue id, name, minimum distance) and a mask to paste. Heavy atoms only, minimum distance over all pairs, strictly inside the cutoff, minimum image when there is a box, solvent excluded. It resolves nothing at run time.
  • Which sidechain bonds are scalable: md_tools.rest2.sidechains (md-tools-sidechain-rotatability/1) is a per-residue table of every sidechain central bond — each χ, each terminal rotation, and each never-scaled bond with its class. It is checked against the torsion classifier residue by residue, and a disagreement refuses a build.

What a run records

scaler.yaml gains, beside the unchanged solute block:

  • selection — md-tools-solute-selection/2.0: the mode, the masks as written, the residue map, the hot atoms, the owned torsion bonds, the CMAP decisions and the ligand instances;
  • selection_sha256 (the identity) and selection_provenance_sha256 (the whole record);
  • scaler_arguments — exactly what each state was built with, so every state rebuilds from the record alone, including in an exported reference bundle.

Identity is the Hamiltonian, not its provenance. md-tools-hamiltonian-identity/v3 hashes what determines the energy — hot atoms, torsion and CMAP decisions, protected bonds, each ligand instance's residue, package and resolved exclusions — and never mask spelling, labels or file paths. A respelled mask or a relabelled instance resumes; one more hot atom, or one flipped torsion or CMAP decision, is refused. A 0.6.0 (v2) identity is still accepted, but only for a legacy selection and only when every v2 field matches.

In a build-md REST2 configuration the three keys are claims about the saved states, checked against the record and refused on disagreement, never applied — the region is chosen once, by build-top --rest2-scaler.

Two fixes that change behaviour

  • The topology digest is canonical. It hashed bonds in the order the file listed them, so reading and rewriting a structure with a cross-residue CONECT bond flipped it between two values forever. Bonds are now hashed as a sorted set (atoms-in-index-order+sorted-bond-set/1). Selection records written before 0.6.1 still validate through a named compatibility branch.
  • A ladder's rungs are never re-derived at run time. md-run and the generated scripts always read them from saved states; the direct Python API did not. A plan with no rungs is now refused, and a caller that supplies its own must declare rung_source="caller-supplied" with a reason. Every rung's origin — saved-state <i> (verified) or caller-modified — is recorded in restart.json by canonical digest, and joins the ladder's identity. Ladders started before 0.6.1 still resume.

Not in this release

  • Selective implicit-solvent (GB) scaling: refused.
  • Chain qualifiers and the rest of the AMBER mask grammar.
  • Nucleic acids: refused by name, because no backbone/sidechain definition exists for them.
  • Any change to cMD or AIS behaviour.

Known defects in this release

export-reference cannot read a ladder that build-md's own run.sh produced. It looks for <run>/<stage>.log; run.sh writes remd_records/<stage>_prod1.log, so the tutorial's section 6 command stops with REST2.log: cannot be read. Workaround, documented on that page: cp remd_records/REST2_prod1.log REST2.log before exporting. This is not a regression — 0.6.0 shipped the same behaviour — and the honest fix is to resolve the log the way the RUN recorded it rather than guess a filename, which is a change with its own failure modes. It is the first item of 0.6.2, with a test that fails before it and that exercises a ladder produced by build-md rather than one assembled by a fixture.

export-reference also verifies the source structure and the ligand package by digest, beside the build. A run tree that has been moved, or registered as a dataset without them, cannot be exported. The two registered tutorial datasets are in that position.

Evidence

Deterministic, on CPU: analytic S–S/S–E/E–E factors per pair and per exception; the PME nonbonded energy as a²E_SS + aE_SE + E_EE against an independent decomposition; every torsion term against the ownership rule recomputed from the topology, under ff19SB and ff14SB; CMAP shared and mixed maps; τ = 0 identity; byte-identical legacy states against the frozen 0.6.0 scaler; boundary cases; reordered atoms; a mutated exclusion file; end-to-end refusals and resumes; and a selective ladder's recorded exchange energies against its saved states at stored coordinates.

On CUDA: four selective ladders (backbone, sidechain, ligand, combined) with their exchange energies recomputed at stored coordinates, an interrupted and resumed selective ladder, and four protein–ligand ladders on the TYK2 complex using a registered ligand package. Three of those are registered datasets, and the tutorial cites them.

The full test lane on the released tree: 2946 passed, 11 skipped, 0 failed across four lanes against a 2957-test collection, every lane verified against the same commit before and after.