Skip to content

Generate GROMACS mdp run schedules at build time - #67

Merged
maxscheurer merged 5 commits into
feature/issue-61-per-component-parametrization-specsfrom
feature/issue-62-generate-mdp-run-schedules
Nov 10, 2026
Merged

maxscheurer merged 5 commits into
feature/issue-61-per-component-parametrization-specsfrom
feature/issue-62-generate-mdp-run-schedules

Conversation

@maxscheurer

@maxscheurer maxscheurer commented Nov 6, 2026 •

Copy link
Copy Markdown
Collaborator

Closes #62

Replaces the package-shipped GROMACS .mdp templates and run_schedules.yaml. Runs are now configured through a public SimulationSettings Pydantic contract that complements BuildInput:

  • Settings: temperature, pressure, timestep, EM limits, stage lengths, output intervals, and per-stage GROMACS overrides.
  • Hash: only the state point (temperature, pressure) enters it. A longer production run keeps the same hash.

The mdps are rendered from plain dicts at build time, and existing build directories can be re-rendered without a rebuild. The file names and the keys that orchestration reads stay the same.

The nonbonded block follows the force-field family for every system type (CHARMM force-switch 1.0→1.2 nm, Sage potential-switch 0.8→0.9 nm). lnp builds get the mixedbox schedule. Removing the template directory also removes the force-include that causes #55.

Stacked on #65 (base feature/issue-61-per-component-parametrization-specs), which is itself stacked on #63 and #59.

The implementation plan (revision 2) is posted as a comment below.

@maxscheurer

maxscheurer commented Nov 6, 2026 •

Copy link
Copy Markdown
Collaborator Author

Implementation Plan (revision 2)

Revised 2026-11-06 after discussion with the user. Changes from the first version:

  • The mdp files are rendered by a thin, dict-based layer instead of a typed model of mdp blocks.
  • A public SimulationSettings contract makes runs configurable.
  • Temperature and pressure are part of the run identity. No other setting is.
  • Existing build directories can be re-rendered from new settings without a rebuild.

Stack: develop ← PR 59 ← PR 63 ← PR 65 (feature/issue-61-per-component-parametrization-specs) ← this PR.

  • The issue says "base develop, after PR 59". PR 65, however, renamed the run-schedule selector to nonbonded_scheme and added spec.nonbonded_scheme, and this work builds on both.
  • The branch is therefore cut from PR 65's head (571dedc) and rebased as the lower PRs merge.

Decisions settled in planning (with the user)

  1. The nonbonded block follows the force-field family for every system type. This was decided after a check of primary sources (see "Nonbonded evidence"). It brings in scope the "related finding" that CGenFF and SMIRNOFF mixedbox/bilayer builds currently run with plain 1.0 nm cutoffs.

  2. lnp gets the mixedbox schedule. lnp builds become simulable; today simulate skips them as "build incomplete".

  3. SimulationSettings is a public Pydantic contract that complements BuildInput.

    • BuildInput describes what is simulated. SimulationSettings describes how it is run, in physical units.
    • The v1 settings are temperature, pressure, timestep, EM limits, NVT/NPT/production lengths and output intervals.
    • A per-stage escape hatch accepts raw GROMACS keys.
  4. Identity rule (in code): the hash covers the BuildInput plus the state point of the settings (temperature, pressure) and nothing else.

    • A different temperature is a different data point: a new build directory and a new run-DB row.
    • A longer production run, another timestep, other output intervals, EM limits and overrides never change the hash.
  5. Coupling runs one way. The settings are resolved and validated against the build:

    • defaults depend on the simulation type (bilayer production defaults to 1 µs);
    • some parts are derived from the build: the nonbonded block from the family, semi-isotropic coupling from bilayer, and the Protein Non-Protein groups from the protein types;
    • overrides may not touch keys that the build or a typed field owns.

    The build does not read the settings. The OpenMM preparation stages keep 300 K; see Out of scope.

  6. Supplying the settings.

    • A protocol file or object per build call (--simulation-settings protocol.yaml), applied to every system.
    • An optional inline simulation_settings: section in the input YAML, and dotted CSV columns.
    • Precedence: inline row, then protocol file, then per-type defaults.
    • Inputs with the same hash but different resolved settings are rejected as duplicates.
  7. Re-render without a rebuild.

    • The settings of existing build directories can be changed from a settings file.
    • The mdps are rewritten only for stages without outputs, or after simulate --clean.
    • Temperature/pressure changes are refused, because they would change the identity.
  8. Lean rendering, and no mdp package (see "Package survey").

    • The mdps are rendered from plain dicts: a base dict per stage, overlays per simulation type and nonbonded scheme, the values from the settings, and the overrides.
    • A writer of a few lines turns the dicts into files.
    • There is no mdp object model and there are no block classes. grompp is the validator.
  9. Clean cut (project convention). Deleted with no compatibility layer:

    • RunScheduleManager, with get_settings, list_* and get_run_schedules_path (none of these has a caller);
    • run_schedules.yaml, including the dead settings: block (num_gpus/num_cores/walltime);
    • the template_dir/selector fallback;
    • all 16 template files.
  10. Regression check: a one-off report posted on this PR, plus permanent contract tests. This follows the project's test convention. The only exact values asserted are the decided ones: the nonbonded table, the defaults, and the protein schedules the issue names.

Nonbonded evidence (primary sources, retrieved during planning)

  • CHARMM (C36/C36m, CGenFF, CHARMM TIP3P)
    • Every CHARMM .prm in the SILCSBIO 2026.1 install contains NONBONDED ... vfswitch ... ctofnb 12.0 ctonnb 10.0, i.e. a force switch from 10 to 12 Å.
    • CHARMM-GUI (Lee et al., JCTC 2016) says the cutoff scheme "should be maintained when different simulation programs are used", for the whole C36 force field, small molecules included.
    • The GROMACS manual's CHARMM36 block is Force-switch from 1.0 to 1.2 nm, PME, DispCorr = no. Dispersion correction "should be applied in the case of lipid monolayers, but not bilayers".
    • The C36 lipids were parametrized without a long-range correction (Klauda 2010).
    • CGenFF's pure-liquid validation, however, used force-switching from 10 to 12 Å with a long-range LJ correction (Vanommeslaeghe, JCC 2010).
  • OpenFF / SMIRNOFF
    • The .offxml files encode a 9 Å cutoff and a 1 Å vdW switch width.
    • The spec "assumes that the long-range dispersion handling is isotropic", and OpenFF regards these settings as part of the force field.
    • The Sage 2.3.0 paper used 0.9 nm, a 1 Å LJ switch, a dispersion correction and PME.
    • OpenMM's switch is a potential switch, so the GROMACS equivalent is Potential-switch from 0.8 nm with EnerPres.
    • There is no distinction by system type. Sage is not validated for lipids (OpenFF roadmap).
  • Amber
    • Plain truncation at 8–10 Å, PME, and the default continuum correction (vdwmeth = 1).
    • The GROMACS Amber-port validation (Miletić 2026) used 1.0 nm with potential-shift and EnerPres.
    • The existing sage block for Amber proteins is within Amber practice and matches openmmforcefields' Sage + Amber example.

Package survey (acceptance criterion 4)

Package License What it offers Verdict
GromacsWrapper 0.9.2 GPL-3+ MDP dict class for editing templates GPL; we would still ship templates
mdpeditor 1.2.5 GPL-3+ composition of "blocks", including CHARMM/Amber blocks GPL; last release 2022
asyncmd 0.5.1 GPL-3+ MDP dict for read/modify/write GPL
gromacs-py 2.0.3, BioSimSpace GPL-2 / GPL-3 templates or protocol objects GPL; BioSimSpace is heavy
OpenFF Interchange MDConfig (already a dependency) MIT the nonbonded block of a single-point-energy mdp cannot express CHARMM Force-switch; hard-codes disp-corr = Ener
biobb_gromacs 5.3.2 Apache-2.0 mdp_preset() and write_mdp()/create_mdp() see below

mdfactory is MIT-licensed, so the GPL tools are excluded. biobb_gromacs is the only permissive option:

  • Its whole mdp layer is about 110 lines of dict assignments, with generic presets that do not depend on the force field.
  • On top of that sit a 10-line writer and a merge of about 15 lines.
  • It would pull in biobb_common==5.3.1 (exact pin), biopython, requests, jsonschema and rocrate, to save about 25 lines.

mdfactory therefore keeps the same shape in-house without the dependency. The force-field-aware protocol data has to live in mdfactory in any case.

Problem analysis

  • Static files. RunScheduleManager picks a yaml entry and byte-copies 16 static, near-duplicate mdp files (482 key = value lines over only 45 distinct keys). Users cannot configure runs except by hand-editing the copies, which a rebuild overwrites.
  • lnp has a yaml entry but no directory. The builder copies nothing and logs "Run schedule files copied."
  • Silent grompp skip. _finalize_protein_build runs grompp only if em_mdp.is_file(), so a missing em.mdp skips validation silently.
  • Consumers read the files by name and key. None of them changes in this PR:
    • StageSpec.mdp_file fixes the file names, _missing_build_files requires them, and clean_simulation_outputs removes <stem>_rescue_t*.mdp/mdout.mdp.
    • apply_rescue_tier reads emstep/dt/nsteps through parse_mdp (lowercase, - treated as _, first match wins). A missing key only logs a warning, and whether the output is an integer depends on "." not in current_val.
    • _extract_expected_frames_from_mdp computes int(nsteps) // int(nstxout_compressed). A value written as 5e7 returns None, and frame validation silently degrades.
    • StageSpec.maxwarn (EM/NVT 1, NPT/Production 2) absorbs today's grompp warnings, so a generator that adds a warning breaks simulate.
  • Obsolete keys. grompp reports title and ns_type as obsolete.
  • Packaging. pyproject.toml force-includes mdfactory/run_schedules/gromacs, which is the cause of issue 55.
  • The hash and the build directory. The hash is the SHA-1 of BuildInput.model_dump_json(). It names the build directory and the run-DB row, and the build directory holds <hash>.yaml, which is model_dump(mode="json"). Settings carried inside the BuildInput document therefore travel through Parsl, the build directory and the analysis layer without a new serialization path. Only the hash has to exclude the settings that are not part of the identity.

Target settings (deliverable contract)

SimulationSettings v1. Field names are illustrative; the defaults reproduce today's protocol.

Field Unit Default Identity
temperature K 300 yes
pressure bar 1.0 yes
timestep fs 2 (at most 2: H-bond constraints, no HMR) no
minimization.max_steps / minimization.force_tolerance – / kJ mol⁻¹ nm⁻¹ 50,000 / 1000 no
nvt.length / nvt.output_interval ps 100 / 1 no
npt.length / npt.output_interval ps 100 / 1 no
production.length ns 100 (bilayer: 1000) no
production.output_interval ps 10 (xtc, energies, log) no
gromacs_overrides.{minimization,nvt,npt,production} raw mdp keys none no
  • Validation:

    • values are positive;
    • lengths and intervals are integer multiples of the timestep;
    • the production length is a multiple of its output interval, so frame accounting is exact;
    • overrides use known stage names and no owned keys.
  • Owned keys, rejected in overrides:

    • keys set from typed fields: nsteps, dt, emtol, ref-t, gen-temp, ref-p, nstxout, nstvout, nstfout, nstenergy, nstlog, nstxout-compressed;
    • keys derived from the build: the nonbonded block, tc-grps, pcoupltype, compressibility.

    Everything else can be overridden, for example pcoupl = C-rescale, tau-p, nstcalcenergy or define. The overrides are applied last.

simulation_type: bilayer
parametrization: charmm
system: { ... }
simulation_settings:            # optional; omitted fields take per-type defaults
  temperature: 310.0            # K, part of the hash
  production: {length: 500.0, output_interval: 20.0}   # ns / ps, not part of the hash
  gromacs_overrides:
    npt: {pcoupl: C-rescale}

Nonbonded block, derived from the family and identical in all four stages:

Scheme Families Electrostatics vdW DispCorr
charmm charmm PME, rcoulomb = 1.2 vdwtype = cutoff, rvdw = 1.2, vdw-modifier = Force-switch, rvdw-switch = 1.0 EnerPres for mixedbox and lnp; no for bilayer, proteinbox and protein_mixedbox
sage openff, amber PME, rcoulomb = 0.9 vdwtype = cutoff, rvdw = 0.9, vdw-modifier = Potential-switch, rvdw-switch = 0.8 EnerPres

Intended changes against today's default output:

  • charmm mixedbox/bilayer: from 1.0 nm potential-shift to 1.2 nm Force-switch.
  • openff mixedbox/bilayer: from 1.0 nm to 0.9 nm Potential-switch. openff bilayer DispCorr changes from no to EnerPres.
  • lnp: files for the first time.
  • EM: carries the full block.
  • Protein types: unchanged.

Derived from the simulation type (unchanged; lnp equals mixedbox):

mixedbox, lnp bilayer proteinbox, protein_mixedbox
tc-grps System System Protein Non-Protein
NPT barostat Berendsen, isotropic, tau_p = 2.0 Berendsen, semi-isotropic, tau_p = 5.0 Berendsen, isotropic, tau_p = 2.0
Production barostat Berendsen, isotropic, tau_p = 2.0 Berendsen, semi-isotropic, tau_p = 5.0 Parrinello-Rahman, isotropic, tau_p = 5.0

Constants as today:

  • EM: steepest descent, emstep = 0.01.
  • NVT/NPT: define = -DPOSRES.
  • Thermostat: V-rescale, tau_t = 0.1.
  • Pressure coupling: compressibility 4.5e-5.
  • Constraints: H-bonds with LINCS.
  • Neighbour search: Verlet.
  • PME: order 4, fourierspacing = 0.16.
  • Production: xtc only (compressed-x-grps = System, no trr), refcoord-scaling = com.
  • NVT generates velocities; NPT and production continue.

Deliverables

  1. SimulationSettings model (new: mdfactory/models/simulation_settings.py)
    • Frozen, extra="forbid", with the fields, defaults and validation above.
    • resolve(simulation_type) fills the per-type defaults.
    • An explicit list of identity fields.
    • NumPy docstrings with units.
  2. BuildInput integration (models/input.py)
    • An optional simulation_settings field, resolved in the existing after-validator next to resolve_parametrization.
    • hash covers the BuildInput with the settings reduced to their identity fields.
    • metadata, flatten_system_parameters and the dry-run line show temperature and pressure.
    • <hash>.yaml in the build directory records the full resolved settings.
  3. Rendering (mdfactory/run_schedules/)
    • Base dicts per stage, overlays per simulation type and nonbonded scheme, settings converted to exact step counts, overrides applied last, and a writer of a few lines.
    • A header naming the simulation type, the scheme, the state point and the mdfactory version.
    • Integers rendered as plain integers, dt/emstep as plain decimals, and no obsolete keys.
    • File names equal StageSpec.mdp_file. A test checks this; the build path does not import mdfactory.orchestration, to avoid a Parsl/TUI import cycle.
  4. Builder wiring (mdfactory/build.py)
    • build_mixedbox, build_bilayer, build_lnp and _finalize_protein_build call the writer with the BuildInput; the settings are inside it.
    • The protein grompp check on em.mdp runs unconditionally.
  5. Supplying the settings (cli.py, orchestration/build.py, prepare.py)
    • --simulation-settings on mdfactory build and prepare-build, and build_systems(..., simulation_settings=...) in Python.
    • Deep-merge under the inline values.
    • The duplicate check (same hash, different resolved settings) runs before any build starts.
    • Dotted CSV columns simulation_settings.* go through the existing generic mapping.
  6. Re-render API and CLI
    • A typed API: apply_simulation_settings(build_dirs, settings, dry_run).
    • A thin Cyclopts subcommand (working name mdfactory apply-settings).
    • It merges the new settings over each directory's recorded settings and validates them.
    • It refuses identity changes and changes to stages that already have outputs.
    • It rewrites the mdps of the other stages and the directory's <hash>.yaml, and reports per directory.
  7. Removals
    • RunScheduleManager, run_schedules.yaml and mdfactory/run_schedules/gromacs/ (16 files).
    • The pyproject.toml force-include entry, with its comment. This removes the cause of issue 55.
  8. Regression report (one-off, posted on this PR)
    • A key-by-key comparison of generated files (default settings) against the static files at 571dedc, for mixedbox, bilayer, proteinbox and protein_mixedbox with each scheme.
    • Keys are normalized; enumerations are compared case-insensitively, numbers numerically, and multi-valued fields token-wise.
    • Every difference is classified as one of: an intended nonbonded change, a dropped obsolete key, the EM nonbonded block, or a key written at its GROMACS default. Anything else is a bug.
    • On the cluster e2e builds, a gmx grompp -po effective-parameter diff for every stage.
  9. Docs
    • A new user-guide page, "Simulation settings": the contract and its defaults, the identity rule, the protocol file, the inline/CSV forms, overrides and owned keys, and re-rendering.
    • yaml-format.mdx, csv-format.mdx and cli-reference.mdx.
    • configuration.mdx: the run-schedule section points to the new page, and the reference to config_templates/run_schedules.yaml is removed.
    • parametrization.mdx: the nonbonded scheme, with the "planned" callout replaced.
    • adding-parametrization.mdx: a new scheme means a new overlay.
    • architecture.mdx, troubleshooting.mdx (small boxes with the 1.2 nm cutoff), SKILLS.md.
    • One example YAML with simulation_settings.

Acceptance criteria (mapped to the issue, plus the contract)

  • No .mdp files, run_schedules.yaml or template directories remain under mdfactory/run_schedules/, and RunScheduleManager is gone.

  • Without the force-include, hatch build produces a wheel with no run-schedule .mdp.

  • All five builders write em/nvt/npt/md.mdp from code at build time, and lnp gets the mixedbox schedule.

  • The nonbonded settings follow the table and are identical in all four stages.

  • The CHARMM and the Amber/Sage protein schedules are reproduced with default settings, with no effective difference from today's files beyond listed zero-effect differences.

  • SimulationSettings is a public Pydantic model. The inline YAML section, dotted CSV columns, the protocol file and the Python object all produce the same resolved settings, with precedence row > file > defaults.

  • The hash changes with temperature or pressure and with no other setting. Omitted settings and spelled-out defaults give the same hash.

  • Inputs that differ only in settings outside the identity are rejected as duplicates.

  • The settings reach the mdps:

    • temperature into every ref-t and into gen-temp;
    • pressure into ref-p (both values for semi-isotropic coupling);
    • lengths and intervals into exact step counts;
    • overrides applied last.

    Overrides on owned keys are rejected.

  • Re-rendering rewrites stages without outputs, refuses stages with outputs and temperature/pressure changes, updates <hash>.yaml, and keeps the hash.

  • The regression report lists and classifies every difference, and the grompp -po comparison agrees.

  • The unchanged orchestration code works on the generated files: rescue tiers, expected frames, the required-file check, clean, and build-time grompp. Every stage passes grompp within its StageSpec.maxwarn.

  • The package decision is documented (above), and the docs describe the settings contract and the generated schedules.

Files to create or modify

  • New:
    • mdfactory/models/simulation_settings.py
    • mdfactory/run_schedules/gromacs.py: dicts, overlays, render, writer
    • mdfactory/run_schedules/__init__.py (rewritten): writer and re-render API
  • mdfactory/models/input.py: the field, resolution and identity hash.
  • mdfactory/build.py: writer call sites and unconditional grompp.
  • mdfactory/cli.py: --simulation-settings and apply-settings.
  • mdfactory/orchestration/build.py: merge, duplicate check and dry-run line.
  • mdfactory/prepare.py: only if the dotted mapping needs a tweak.
  • mdfactory/analysis/utils.py: flatten_system_parameters.
  • mdfactory/models/parametrization.py: docstrings of NonbondedScheme and nonbonded_scheme.
  • pyproject.toml: force-include.
  • Deleted: run_schedules.yaml and mdfactory/run_schedules/gromacs/**.
  • Tests:
    • mdfactory/tests/test_simulation_settings.py (new)
    • mdfactory/tests/test_run_schedules.py (new)
    • test_protein_mixedbox.py, test_amber_gmx_build.py, test_build.py
    • the CLI and orchestration build tests
  • Docs:
    • docs/content/docs/user-guide/{simulation-settings (new),yaml-format,csv-format,cli-reference,configuration,parametrization,troubleshooting}.mdx plus meta.json
    • docs/content/docs/developer-guide/{adding-parametrization,architecture}.mdx
    • SKILLS.md and examples/

Testing approach

New tests have no "build" in their names, so CI runs them with pytest -k "not build" mdfactory.

  • test_simulation_settings.py
    • Per-type defaults; for example, bilayer production is 1 µs.
    • The identity rule:
      • production length, timestep, intervals, EM limits and overrides leave the hash unchanged;
      • temperature and pressure change it;
      • omitted settings give the same hash as spelled-out defaults;
      • the list of identity fields is pinned, so a new field forces a decision.
    • Validation errors: non-positive values, values that are not multiples of the timestep, a timestep above 2 fs, an unknown stage, and an owned key in an override.
    • Every supply path (YAML section, dotted CSV columns, protocol file, Python object) resolves to the same result, and the precedence order holds.
    • The duplicate rejection.
  • test_run_schedules.py
    • For every simulation type × scheme combination (10 in total):
      • exactly four files, named as in STAGE_REGISTRY;
      • every file parses with parse_mdp, with no duplicate normalized key;
      • the header is present;
      • the writer overwrites existing files.
    • Orchestration contract, using the real functions:
      • apply_rescue_tier halves emstep/dt and doubles nsteps, with integers kept as integers;
      • _extract_expected_frames_from_mdp returns length / output_interval;
      • _missing_build_files is empty.
    • Settings propagation:
      • temperature reaches every ref-t and gen-temp;
      • pressure reaches ref-p, including both semi-isotropic values;
      • step counts are exact;
      • overrides end up in the files.
    • The nonbonded table, asserted exactly.
    • Protocol invariants per simulation type, asserted as rules:
      • semi-isotropic coupling for bilayer;
      • the tc-grps choice and the matching counts;
      • -DPOSRES only in NVT/NPT;
      • Parrinello-Rahman in protein production;
      • lnp identical to mixedbox apart from the header.
    • Sage cross-check: the sage cutoff and switch distance equal the vdW handler of DEFAULT_OPENFF_FORCEFIELD and of the default water .offxml, loaded with openff-toolkit.
  • Re-render (temporary build directories):
    • stages without outputs are rewritten;
    • a stage with outputs is refused unless it was cleaned;
    • a temperature/pressure change is refused;
    • <hash>.yaml is updated and the hash kept;
    • dry_run writes nothing;
    • the CLI wrapper delegates to the API.
  • CLI and orchestration:
    • --simulation-settings in the dry-run path, for both build and prepare-build;
    • the dry-run line shows temperature and pressure;
    • the duplicate error is raised before any submission.
  • test_protein_mixedbox.py:
    • the template-path tests are removed;
    • test_finalize_installs_* checks the scheme's schedule and the unconditional validate_with_grompp, using a monkeypatched recorder.
  • Opt-in real GROMACS (module names contain "build"; GROMACS ≥ 2026):
    • test_amber_gmx_build.py uses the writer and runs grompp on all four stages within each stage's maxwarn, passing -r for restrained stages. It covers:
      • the protein schedules, on the peptide system;
      • the mixedbox, bilayer and lnp schedules, on a water box of at least 3 nm;
      • one run with non-default settings and an override.
    • test_build.py::test_build_mixedbox/test_build_bilayer also check the four files and the charmm block.
  • Cluster e2e (cpu/gpu partitions):
    • Rebuild the PR 65 matrix with the new code.
    • Check the four files and run grompp on every stage.
    • Run EM → NVT with mdfactory simulate; lnp runs on gpu.
    • Run one temperature series of two rows (two hashes) and one re-render that changes the production length (same hash).
    • Produce the grompp -po comparison for the regression report.
    • Post the results on this PR.
  • Full suite: pytest -k "not build" mdfactory on a compute node, plus ruff 0.9.2 check and format.

Risks and open questions

  • The identity rule lives in code. A field added to SimulationSettings later must be classified explicitly. A test pins the identity fields, and the docstring states the rule.
  • Every hash changes once more, because the resolved temperature and pressure enter the identity. PR 65 already changes every hash, and both ship together.
  • Records after a re-render. The directory's <hash>.yaml and the mdp headers are authoritative. The run-DB input_data keeps the settings from build time, and this is documented. The DB records identity; the directory records the protocol.
  • Re-render is not "extend". Changing a stage that already has outputs (for example, extending a finished production run) is refused. Extension semantics (tpr extension, checkpoints) are a follow-up.
  • Overrides can produce invalid or poor mdps. grompp catches invalid ones at simulate; for protein types the build also checks EM. Owned keys are protected, so overrides cannot contradict the force field, the coupling geometry or the typed settings.
  • The protocol changes on purpose. New mixedbox/bilayer/lnp builds are not comparable with earlier 1.0 nm campaigns.
  • Minimum box size. The 1.2 nm CHARMM cutoff raises the minimum box edge to about 2.6 nm. The repository examples are far larger, and the troubleshooting docs mention it.
  • Silent drift and number formatting are covered by the regression report and the contract tests.
  • grompp warnings versus maxwarn are checked per stage in the opt-in test and the e2e runs.
  • Sage for lipids is not validated by OpenFF; the docs say so.
  • lnp under the mixedbox schedule uses DispCorr = EnerPres despite its lipid/water interfaces. This is the user's decision and is documented.
  • Stack churn. PRs 59, 63 and 65 are open, and this PR is rebased as they merge. nonbonded_scheme is extended, never renamed.
  • PR 56 (the issue 55 fix, base develop) edits the same force-include line and adds a .gitignore exception plus a test that expects 8 templates in the wheel.
    • If it merges first, the rebase removes those parts.
    • If this PR lands first, PR 56 is superseded and issue 55 can be closed.

Out of scope (follow-ups)

  • The OpenMM preparation stages, which keep their hard-coded 1.0 nm cutoff and 300 K. A follow-up could use the settings' temperature, since it is part of the identity.
  • Extending finished runs, and several protocols per build directory.
  • Typed fields for ensemble algorithms (thermostat, barostat, restraints). Meanwhile they are reachable through the overrides.
  • Modernizing the protocol defaults (C-rescale, refcoord-scaling in NPT, LJ-PME) and a dedicated lnp protocol.
  • Run-DB columns for the state point.
  • The protein system-type docs (issue 66).

Plan created by mach6 (revision 2)

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Progress Update

Plan revision 2 is implemented in full. The commit below adds the following:

  • SimulationSettings (mdfactory/models/simulation_settings.py) is a public, frozen Pydantic contract.
    • It covers the fields, defaults and validation of the plan.
    • IDENTITY_FIELDS = {temperature, pressure}. BuildInput.hash excludes every other setting.
    • Overrides may not set keys that mdfactory owns. Key matching ignores case, - and _, as GROMACS does, and also catches the old name nstxtcout.
  • Rendering (mdfactory/run_schedules/gromacs.py) builds each stage from a base dict.
    • Overlays add the coupling of the simulation type (COUPLING) and the nonbonded block of the scheme (NONBONDED, DISPCORR).
    • The settings become exact step counts, and the per-stage overrides are applied last.
    • A small writer adds a header with the type, the scheme, the state point and the mdfactory version. lnp uses the mixedbox schedule.
  • Builders
    • All five builders call write_mdp_files(inp).
    • The protein em.mdp grompp check now always runs.
  • Removed: RunScheduleManager, run_schedules.yaml, the 16 templates and the force-include that causes issue 55.
  • Supplying settings
    • An inline simulation_settings: section.
    • Dotted CSV columns; blank cells fall back to the protocol file or the defaults.
    • --simulation-settings on build and prepare-build, and build_systems(simulation_settings=...).
    • Precedence: row, then file, then defaults. Inputs with the same hash are rejected before anything is built.
  • Re-render: mdfactory.run_schedules.apply_simulation_settings(build_dirs, settings, dry_run) and mdfactory apply-settings.
  • Docs: a new "Simulation Settings" page, the other pages from the plan, SKILLS.md, and the example examples/mixedbox/water_ethanol_smirnoff_310K.yaml.

Decisions made while implementing (beyond the plan text)

  • Build record: run_build_from_dict now writes <hash>.yaml. Single-YAML builds had no record before, so the re-render could not find their settings.
  • BuildInput objects: build_systems(simulation_settings=...) raises when it gets BuildInput objects, whose settings are already resolved. It merges only into dict inputs.
  • Identical inputs: inputs that are fully identical are rejected as duplicates, as are inputs that differ only in settings outside the identity.
  • Summary YAML: --simulation-settings together with a summary YAML is rejected. Those systems were prepared with their settings; apply-settings changes them.
  • Refused re-render: the refusal lists the outputs of the blocking stage. simulate --clean reruns the stage right away, so it cannot clear outputs ahead of a re-render.

Bug found during the e2e runs and fixed: <hash>.yaml was written with yaml.safe_dump(sort_keys=True). Sorting reordered dict fields such as LNP input_fractions, so the reloaded record had a different hash from its directory. This already happened for LNP prepare-build records before this PR. A regression test is added.

Commit: 4392963

E2E results and regression report: generated run schedules (issue 62), 2026-11-06

Uncommitted working tree on top of 4ea82ba. Jobs: mdfactory_e2e/scratch62/jobs/job62_{cpu,gpu,lnp2,lnp_sim2}.sh. Check script: scratch62/check_dirs.py. Static files for comparison: a worktree at 571dedc.

1. Regression report: generated files versus the static files at 571dedc (default settings)

Keys are compared the way GROMACS reads them: case, - and _ are ignored. Enumerations are compared case-insensitively, numbers numerically, and multi-valued fields token by token. The comparison covers mixedbox, bilayer, proteinbox and protein_mixedbox with each nonbonded scheme, all four stages, and 159 differing keys. Each difference falls into one of four classes:

Class Count Where
Intended nonbonded change 71 Only mixedbox and bilayer: charmm 1.0 nm → 1.2 nm with Force-switch from 1.0; openff 1.0 nm → 0.9 nm with Potential-switch from 0.8; openff bilayer DispCorr no → EnerPres
EM nonbonded block 16 mixedbox/bilayer em.mdp now carries vdwtype, vdw-modifier, rvdw-switch and DispCorr
Dropped obsolete key 56 title and ns_type
Key written at its GROMACS default 16 nstfout = 0 in NVT/NPT
Other (would be a bug) 0

The protein schedules (charmm, and amber/sage) differ from the old files only in dropped obsolete keys and nstfout = 0. Full table: scratch62/regression_static.md.

2. grompp -po effective-parameter diff on the built systems (gmx 2026.3)

For each build and stage, grompp ran on the same structure and topology twice, once with the static file and once with the generated file. The effective mdout parameters were then compared. The realized random seeds gen-seed and ld-seed were ignored.

  • proteinbox and protein_mixedbox, charmm and amber: no effective difference in any stage.
  • mixedbox and bilayer, charmm and openff: only the nonbonded keys (rcoulomb, rvdw, rvdw-switch, vdw-modifier, DispCorr), as in the table above.
  • Builds with non-default settings (the 310 K example and the T-series after re-rendering) additionally differ only in the keys those settings own: ref-t, gen-temp, nsteps, the nst* intervals, and the pcoupl = C-rescale override.

Warnings per stage were counted with each stage's input structure as simulate uses it. They stay within StageSpec.maxwarn everywhere and equal the counts for the static files:

Stage Small-molecule types (mixedbox, bilayer, lnp) Protein types
EM 1 / 1: atom names from the CGenFF/OpenFF topology differ from system.pdb (pre-existing) 0 / 1
NVT 1 / 1: -DPOSRES is not used by the topology 0 / 1
NPT 2 / 2: Berendsen, plus -DPOSRES 2 / 2: Berendsen, plus absolute position restraints with pressure coupling
Production 1 / 2: Berendsen 0 / 2

With C-rescale overriding the barostat (the 310 K example), NPT has 1 warning out of 2. Reports: scratch62/report_{cpu,gpu,lnp}.md.

3. Builds and simulations

System Family (scheme) Env Build 4 mdps + record EM → NVT
mixedbox water/ethanol charmm (charmm) dev, cpu OK OK OK (gmx 2026.3, cpu)
mixedbox water/ethanol openff (sage) dev, cpu OK OK OK
mixedbox water/ethanol, water_ethanol_smirnoff_310K.yaml (310 K, 20 ns, C-rescale NPT override) openff (sage) dev, cpu OK OK OK
bilayer POPC charmm (charmm) dev, cpu OK OK OK
bilayer POPC openff (sage) dev, cpu OK OK OK
lnp ALC-0315 openff (sage) dev, gpu OK OK OK (--slurm, gpu, gmx-cuda 2026.3)
lnp ALC-0315 charmm (charmm) dev, gpu OK on the second attempt* OK OK (--slurm, gpu)
proteinbox lysozyme charmm (charmm) dev + gmx_mpi 2024.3, gpu OK (build-time grompp) OK OK (gmx 2026.3, cpu)
protein_mixedbox lysozyme + ethanol charmm (charmm) dev + gmx_mpi 2024.3, gpu OK (build-time grompp) OK OK
proteinbox lysozyme amber (sage) gmx-cuda 2026.3, gpu OK (build-time grompp) OK OK
protein_mixedbox lysozyme + ethanol amber (sage) gmx-cuda 2026.3, gpu OK (build-time grompp) OK OK

* The first charmm LNP attempt failed before any mdp was written: OpenMM reported Energy is NaN in create_sphere during the OpenMM preparation. This is the known random blow-up (issue 58 class). The retry built cleanly.

LNP builds now get the mixedbox schedule and simulate; before this change, simulate skipped them.

4. Temperature series and re-render

  • Temperature series. A CSV with simulation_settings.temperature = 300 and 310 gave two hashes. They equal the hashes of the 300 K example and of the 310 K example, which uses other production settings and an override. Only the state point enters the identity.
  • Production length. After EM and NVT had run, apply-settings longer.yaml (production: {length: 200}) updated both directories. npt.mdp, md.mdp and <hash>.yaml were rewritten: md.mdp now has nsteps = 100000000, and the hashes are unchanged. nvt.mdp, which has outputs, was not touched.
  • Temperature change. apply-settings hotter.yaml (320 K) was refused for both directories with exit status 1 and wrote nothing.
  • NVT length. apply-settings shorter_nvt.yaml (NVT 50 ps) was refused with exit status 1, because NVT already has outputs. The refusal lists nvt.tpr, nvt.cpt, nvt.log, nvt.edr, nvt.gro.

5. Tests

  • Unit tests: pytest -k "not build" mdfactory gives 3406 passed and 1 skipped.
  • Orchestration and CLI modules: -k "not build" deselects test_orchestration_build.py and test_orchestration_cli.py entirely, so they ran separately. 37 passed and 2 failed. The two failures (test_shutdown_parsl_*) also fail at 571dedc when run inside a SLURM job.
  • Opt-in real-GROMACS tests (test_amber_gmx_build.py, gmx 2026.3): 7 passed. They run grompp on all four stages within each stage's maxwarn, passing -r for restrained stages. The cases are the amber protein schedule on the neutralized peptide, the mixedbox, bilayer and lnp schedules on a 3.2 nm water box with both schemes, and one run with non-default settings and overrides.
  • test_build.py::test_build_mixedbox/test_build_bilayer (CGenFF): these tests now also check the four files and the charmm block. The bilayer ion count (n_ions == 6) is flaky on both 571dedc and this branch, because the water count after the OpenMM compression varies close to the salt-count rounding threshold. In three runs on this branch it was 8, 8 and 6, and it passed on 571dedc.
  • Wheel: a hatchling wheel built without the force-include contains no run-schedule .mdp.

Found during e2e and fixed

<hash>.yaml was written with yaml.safe_dump's default sort_keys=True. That reordered dict fields such as LNP input_fractions, so a reloaded record had a different hash from its directory. The bug predates this change: it affected prepare-build records of LNPs. Records are now written with sort_keys=False, and a regression test covers it.


Progress tracked by mach6

@maxscheurer
maxscheurer added this pull request to stack #64 November 6, 2026 16:39
@maxscheurer
maxscheurer marked this pull request as ready for review November 6, 2026 16:39
@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Unverified Review Candidates — Pending Assessment

Review round: 1
Reviewed commit: 4392963

These are unverified candidates. Severity reflects reviewer confidence; do not treat any item as a merge blocker until the assessment comment is posted.

Critical

None.

Important

1. Integer-valued GROMACS overrides are rendered as 100.0 when their CSV column has blank cells (code-reviewer 80, error-auditor 85)
mdfactory/prepare.py (df_to_build_input_models) → mdfactory/run_schedules/gromacs.py (_apply_overrides, format_mdp_value). If any cell in a simulation_settings.gromacs_overrides.<stage>.<key> column is blank, pandas makes the column float64. Both reviewers reproduced nstcalcenergy = 100.0 in the rendered npt.mdp. GROMACS reads integer keys with strtol, so grompp would reject the value. Small-molecule builds have no build-time grompp check, so the failure would only show up in simulate. Suggested fix: convert integral floats in overrides to int, plus a regression test.

2. check-csv and push --csv ignore the protocol file, so they hash to the wrong directories (error-auditor 85)
mdfactory/cli.py (check_csv), mdfactory/utils/push.py (load_models_from_csv). Both call load_build_input_csv(csv) with no settings. The state point is now part of the hash. A CSV built with --simulation-settings protocol.yaml (for example temperature: 310) is therefore re-hashed at 300 K. check-csv misses systems that already exist, and push --csv reports "No folder found for hash". Suggested fix: pass --simulation-settings through these commands, or warn.

3. One bad build record or a write error aborts the whole apply-settings batch, and the writes are not atomic (error-auditor 80)
mdfactory/run_schedules/__init__.py (_load_record, _try_update).

  • _load_record catches only TypeError/ValueError. yaml.YAMLError, OSError and UnicodeDecodeError escape, as do OSErrors from write_mdp_files/write_build_record.
  • _try_update only catches _Refused, so the batch stops with a traceback. Earlier directories are already rewritten, and no per-directory summary is printed.
  • The mdps are written before the record. A failure in between leaves new mdps with the old record.

Suggestions

4. Test gap: the section-level deep merge across layers is never exercised (test-reviewer 85). No precedence or re-render test has two layers that set different fields of the same section. Replacing the deep merge in merge_simulation_settings with a shallow one would pass every test. The user-visible effect would be that a partial production: in a row or in apply-settings resets the length to the type default.

5. Test gap: builder wiring of write_mdp_files is only asserted in CI-deselected test_build.py, and lnp not at all (test-reviewer 85). If the call were dropped from build_mixedbox/build_bilayer/build_lnp, CI would still pass.

6. Test gap: run_build_from_dict writing <hash>.yaml is untested (test-reviewer 85). Every test mocks it. Re-render depends on this record for single-YAML and Parsl builds.

7. Test gap: owned-key rejection is spot-checked (test-reviewer 82). The tests cover 6 of the 13 settings-owned keys and 4 of the 10 build-owned keys. Nothing ties the owned set to the keys the renderer derives. Dropping ref-p or gen-temp from the set would let an override change the state point without changing the hash.

8. Test gap: prepare-build duplicate rejection is untested (test-reviewer 80). Only _build_from_csv and build_systems are covered. Since _prepare_system_directories now always writes the record, a missing check would let one row silently overwrite another row's record.

9. Test gap: re-render refusal paths beyond state point and outputs (test-reviewer 80). The untested paths are invalid merged settings, a hash mismatch or corrupt record, and a mixed batch where one directory is refused.

10. Test gap: grompp acceptance and the non-protein sequential-YAML build run only outside CI (test-reviewer 80). Grompp acceptance is only checked by the opt-in GROMACS-2026 test. The sequential YAML build is only asserted in the deselected test_orchestration_cli.py.

11. Simplify: SimulationSettings.identity() has no caller (simplifier 92). mdfactory/models/simulation_settings.py.

12. Simplify: render_stage repeats the nsteps/dt/interval logic in both branches (simplifier 82). mdfactory/run_schedules/gromacs.py.

13. Simplify: apply_settings counts statuses by hand; collections.Counter would do (simplifier 80). mdfactory/cli.py.

Strengths

  • The identity rule is small and enforced in code: IDENTITY_FIELDS drives BuildInput.hash, and a test pins it.
  • Override key matching follows GROMACS (case-insensitive, - and _ treated alike, the old name nstxtcout caught). Owned keys are protected.
  • The renderer is a lean dict layer: base, then overlays, then settings, then overrides. Step counts are exact, and integers and decimals are written the way the unchanged rescue and frame code needs.
  • Precedence is implemented once and reused by every supply path. Duplicates are rejected before any submission.
  • The clean cut is complete: no RunScheduleManager, templates, yaml or force-include remain, and no compatibility layer was added.
  • The regression report is thorough: 159 differences, all classified, 0 unexplained, and the grompp -po comparison agrees. The sort_keys hash bug found during e2e is fixed and has a regression test.
  • The completeness check found every deliverable and acceptance criterion of plan revision 2 implemented.

Agents run: code-reviewer, error-auditor, test-reviewer, completeness-checker, simplifier


Reviewed by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Review Assessment

Findings comment: #67 (comment)

Round 1 was assessed by the independent assessor. Its results were cross-checked against the code.

Classifications

Finding Classification Reasoning
1. Integer overrides rendered as 100.0 when the CSV column has blanks merge blocker Factual: Reproduced. A CSV column simulation_settings.gromacs_overrides.npt.nstcalcenergy with values 100 and blank is read as float64. row.dropna() (prepare.py:184) and MdpValue (simulation_settings.py:81) keep it as a float. format_mdp_value (gromacs.py:277-278) then writes nstcalcenergy = 100.0. With every cell filled, the value is written as 100. GROMACS 2026.3's libgromacs contains the get_eint error "Right hand side '%s' for parameter '%s' in parameter file is not an integer value", so grompp rejects the line. Scope: Dotted CSV columns are an approved supply path (deliverable 5), and blank-cell fallback was added and documented by this PR. Acceptance criterion 6 requires that all supply paths produce the same resolved settings; here YAML gives 100 and CSV gives 100.0. Practical: A user runs a CSV campaign with an integer override on some rows only. The build succeeds, because small-molecule builds run no build-time grompp. simulate then fails at that stage on every affected row, possibly after EM/NVT compute has been spent.
2. check-csv / push --csv ignore the protocol file useful follow-up Factual: True. cli.py:865 and push.py:125 hash rows without settings. The divergence is new in this PR, but it only appears when a protocol file sets temperature or pressure. Scope: The plan does not cover push or check-csv. The limitation is documented with a workaround (simulation-settings.mdx:102-107). Practical: The failure is loud ("No folder found for hash"), nothing wrong is uploaded, and prepare-build re-checks the DB with the correct settings.
3. apply-settings batch aborts on a corrupt record or OSError; writes not atomic deferred Factual: True (run_schedules/__init__.py:86,89,103,128-129,136). Scope: Refusals are reported per directory as planned. Corrupt-YAML and filesystem errors are not named in the plan. Practical: This needs a hand-corrupted record or a disk/quota failure. The batch stops with a loud traceback, and rerunning recovers.
4. Test gap: cross-layer deep merge within one section deferred Factual: The code is correct (simulation_settings.py:349). Only the test design leaves a shallow-merge regression undetected. Practical: No current defect; this would harden against a future change.
5. Test gap: builder wiring only in deselected test_build.py, lnp none deferred Factual: True. Scope: The plan placed these checks in test_build.py, and lnp was verified in the e2e runs. Practical: All five builders call write_mdp_files.
6. Test gap: run_build_from_dict record write deferred Factual: True, the function is mocked everywhere. Practical: If this regressed, apply-settings would refuse loudly with "no build record".
7. Test gap: owned-key rejection spot-checked deferred Factual: True. The owned sets are complete per the plan today. Practical: This only guards future edits. A test tying the owned set to the rendered keys would be useful hardening.
8. Test gap: prepare-build duplicate rejection deferred Factual: The check is present (cli.py:127), and the shared helper is tested through other callers. Practical: No current defect.
9. Test gap: further re-render refusal paths deferred Factual: Partly true; the missing-record case is tested. Practical: No current defect.
10. grompp acceptance and sequential YAML build outside CI discarded observation Scope: The opt-in GROMACS tests plus e2e runs are what the approved plan specifies, and all of them passed (7 opt-in tests, maxwarn respected, grompp -po diff agrees).
11. Unused SimulationSettings.identity() nitpick Factual: No caller (simulation_settings.py:304). Practical: No effect, but removing it fits the lean-scope convention.
12. Duplicated logic in render_stage nitpick Practical: Style only.
13. Hand-counted statuses vs Counter nitpick Practical: Style only.

Action Plan

  1. Finding 1: write override values that are whole-number floats as int, for example in a GromacsOverrides validator or in the CSV mapping. Real-valued GROMACS keys accept an integer, and dt/emstep are written outside the override path, so the rescue-tier decimal-point contract is unaffected. Add a restrained regression test: a two-row CSV with an integer override in one row and a blank cell in the other must render nstcalcenergy = 100.

Assessment by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Progress Update

This fixes finding 1, the only merge blocker from the round-1 assessment.

The problem. pandas reads a numeric CSV column with any blank cell as float. An integer override such as simulation_settings.gromacs_overrides.npt.nstcalcenergy = 100 was therefore written to the mdp as 100.0. grompp rejects that for integer parameters.

Changes:

  • mdfactory/prepare.py: a new helper, _row_cells, replaces row.dropna() in df_to_build_input_models. For simulation_settings.gromacs_overrides.* columns it converts whole-number floats back to int.
    • The fix sits at the CSV boundary, so YAML and Python inputs are unchanged.
    • A whole-number value written with a decimal point in a CSV (for example tau-p = 2.0) is now written as 2. GROMACS reads 2 for real-valued keys as well.
    • The rescue tiers are unaffected: dt is owned by mdfactory, and an emstep override written as 1 still halves to 0.5.
  • mdfactory/tests/test_simulation_settings.py: new test test_csv_integer_override_with_blank_cells_stays_integer.
    • It uses a two-row CSV: the override is set in one row and blank in the other.
    • The CSV row resolves to the same settings as the YAML form, the blank row has no override, and the mdp gets 100.
    • The test fails without the fix.
  • docs/content/docs/user-guide/csv-format.mdx: one sentence on the integer handling.

Checks (cpu partition):

  • test_simulation_settings.py and test_run_schedules.py: 122 passed.
  • pytest -k "not build" mdfactory: 3407 passed, 1 skipped.
  • ruff check and format: clean.

Commit: 0c6191e


Progress tracked by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Unverified Review Candidates — Pending Assessment

Review round: 2
Reviewed commit: 0c6191e

These are unverified candidates. Severity reflects reviewer confidence; do not treat any item as a merge blocker until the assessment comment is posted.

Round-1 merge blocker (finding 1 of round 1): fixed. _row_cells in mdfactory/prepare.py turns whole-number floats in simulation_settings.gromacs_overrides.* columns back into int. Three reviewers checked it independently:

  • The code-reviewer loaded a two-row CSV with a blank cell and got nstcalcenergy = 100 in the rendered md.mdp.
  • The error-auditor confirmed that 0.5 stays a float, and that strings and bools pass through unchanged.
  • The completeness-checker confirmed the new regression test test_csv_integer_override_with_blank_cells_stays_integer. Its final assertion, the rendered value "100", fails without the fix.

Critical

None.

Important

1. CSV override cells such as None, NA, n/a or null are read as blank and silently dropped (error-auditor 82)
mdfactory/prepare.py:268 calls pd.read_csv with pandas' default NA handling, and _row_cells calls row.dropna().

  • Trigger: a CSV cell simulation_settings.gromacs_overrides.nvt.comm-mode = None. None is a legal GROMACS value for comm-mode.
  • What happens: pandas reads the cell as NaN, and the override disappears with no error or log line. nvt.mdp keeps the default comm-mode = Linear.
  • YAML comparison: the YAML path keeps the string None, so the CSV and YAML paths resolve differently.
  • Suggested fix: read the CSV with keep_default_na=False, na_values=[""], or document lowercase none as the workaround. Add a regression test.

2. A misspelled override key is accepted, and grompp's warning for it is absorbed by the stage's maxwarn (error-auditor 80)
mdfactory/models/simulation_settings.py (GromacsOverrides._check_keys) checks only the syntax of override keys and the owned-key list. An unknown key, for example nstcalcenrgy, is written to the mdp as is.

  • grompp: GROMACS 2026.3 grompp reports Unknown left-hand 'nstcalcenrgy' as a warning. With -maxwarn 1 it continues.
  • Budget: the stage budgets in orchestration/stages.py are EM/NVT 1 and NPT/production 2. A stage that has a spare warning therefore runs, and the override has no effect.
  • Build-time check: the protein build-time check covers only em.mdp, with -maxwarn 0. A typo in an EM override fails the build, while the same typo in NVT passes.
  • Docs: simulation-settings.mdx says "grompp validates the result when the stage runs". That is not true for unknown keys.
  • Suggested fix: treat "Unknown left-hand" as fatal, or check keys against a list of known mdp keys. At minimum, correct the docs sentence.

Suggestions

3. Test gap: the CSV integer fix is only tested on a whole-number value (test-reviewer 80)
No test covers a fractional override, such as tau_p = 4.5, in a column that has blank cells. If the value.is_integer() guard were dropped, 4.5 would be truncated to 4, and the new test would still pass.

4. Simplify: the row-label uniqueness check is repeated (simplifier 80)
mdfactory/cli.py (prepare_build, _build_from_csv) contains the line _check_unique_or_exit(models, [f"row {i}" for i in df.index]) twice. A small helper would keep the labelling in one place.

5. Simplify: Optional[PositiveFloat] in simulation_settings.py (simplifier 80)
The rest of the PR uses X | None, and this is the only use of Optional in the module.

6. Simplify: the "str/Path means load the settings file" branch is duplicated (simplifier 80)
The branch appears in both orchestration/build.py (build_systems) and run_schedules/__init__.py (apply_simulation_settings). It could live in _as_mapping in simulation_settings.py. The explicit rejection of settings together with BuildInput objects in build_systems must stay.

Strengths

  • The round-1 blocker is fixed at the CSV boundary with a targeted regression test, and the YAML and Python paths are unchanged.
  • The code-reviewer exercised apply_simulation_settings end to end on a real mixedbox directory. It checked a production-only re-render, a timestep change, a refused temperature change, an empty update, refusal versus skip for a stage with outputs, and replacement of a differently spelled override key.
  • The orchestration consumers keep working on the generated values:
    • The rescue tiers read emstep, dt and nsteps, and dt stays a decimal after halving.
    • The expected-frame count reads plain integers.
    • The file names match STAGE_REGISTRY.
  • The build record round-trips to the same hash.
  • The completeness check found every deliverable and acceptance criterion of plan revision 2 implemented, with the scope additions disclosed.

Agents run: code-reviewer, error-auditor, test-reviewer, completeness-checker, simplifier


Reviewed by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Review Assessment

Findings comment: #67 (comment)

Round 2 was assessed by the independent assessor. The parent reproduced the pandas behavior behind finding 1 in the dev environment. The round-1 merge blocker is fixed (_row_cells, commit 0c6191e).

Classifications

Finding Classification Reasoning
1. CSV override cells None/NA/n/a/null read as blank and silently dropped merge blocker Factual: Reproduced with pandas 2.3.3. load_build_input_csv calls pd.read_csv(csv_path) (prepare.py:268) with pandas' default missing-value strings, and _row_cells calls row.dropna(). A cell simulation_settings.gromacs_overrides.nvt.comm-mode = None resolves to nvt={}, with no error or log line, so the stage runs with GROMACS' default comm-mode = Linear. YAML comm-mode: None loads as the string 'None'. Scope: Plan revision 2 requires that the YAML section, dotted CSV columns, the protocol file and the Python object "all produce the same resolved settings", and the plan names the overrides as the route to ensemble algorithms. The silent drop is new in this PR. The PR added simulation_settings. to OPTIONAL_FIELD_PREFIXES (prepare.py:27-34); before that, the same cell in a required column raised loudly. Practical: A CSV campaign sets a per-row override whose value is on pandas' missing-value list. The realistic case is comm-mode = None, spelled as in the GROMACS manual. The build succeeds, and the record and mdp files simply lack the override, so the run silently uses a different protocol than requested. The docs promise a fallback only for blank cells, and grompp sees a valid file. The trigger is narrow, but the user's input is lost silently and the parity criterion is broken directly.
2. Misspelled override key accepted; grompp's warning absorbed by maxwarn useful follow-up Factual: Partly true. _check_keys (simulation_settings.py:192-215) checks only key syntax, owned keys and duplicate spellings. Whether grompp's "Unknown left-hand" warning is absorbed depends on the stage's spare budget. With the posted per-stage counts, a typo fails in small-molecule EM, NVT and NPT, and in protein NPT; protein EM is caught at build time. A typo is absorbed in small-molecule production, protein NVT and protein production. Scope: The plan names this risk and makes grompp the validator, with "no mdp object model". A known-key list goes against that design. YAML and CSV behave the same, and the maxwarn budgets are older than this PR. Practical: A user has to mistype a key in a stage that has spare budget, and the warning still appears in the stage log. That is no worse than a typo in a hand-edited mdp before this PR. Worth doing later: soften the sentence at simulation-settings.mdx:124, or treat "Unknown left-hand" as fatal.
3. Test gap: fractional CSV override in a column with blanks deferred Factual: True, the new test uses only 100. The code is correct: the value.is_integer() guard keeps 4.5 a float. Scope: Not required by the plan. Practical: It only guards against a future edit that drops the guard. It could be folded into the finding 1 test as one extra column.
4. Repeated row-label uniqueness call nitpick Factual: Two identical lines (cli.py:127, cli.py:338). Practical: Style only.
5. Optional[...] vs X | None nitpick Factual: The only use in the module (simulation_settings.py:23,164). Practical: No runtime effect.
6. Duplicated "str/Path means load the settings file" branch nitpick Factual: True (orchestration/build.py:80-81, run_schedules/__init__.py:171-172). Practical: Two short call sites that behave the same.

Action Plan

  1. Finding 1: in load_build_input_csv (mdfactory/prepare.py:268), read with pd.read_csv(csv_path, keep_default_na=False, na_values=[""]), so that only truly empty cells count as blank.
    • The parent checked three things: None stays the string 'None', empty cells still become NaN (so the required-column check and the blank-cell fallback still work), and numeric columns stay numeric.
    • This is the single CSV entry point, so build, prepare-build, check-csv and push all get the fix.
    • Add a regression test: a CSV cell ...gromacs_overrides.nvt.comm-mode = None resolves to the same settings as YAML {'nvt': {'comm-mode': 'None'}}, and nvt.mdp contains comm-mode = None.
    • Optionally add a sentence to csv-format.mdx saying that only empty cells count as blank.

Assessment by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Progress Update

This fixes finding 1, the only merge blocker from the round-2 assessment.

The problem. pandas reads the cell text None, NA, null, n/a and similar as missing. A CSV override such as simulation_settings.gromacs_overrides.nvt.comm-mode = None was therefore treated as a blank cell and dropped without any message. The YAML path keeps the string None, so the two supply paths resolved differently.

Changes:

  • mdfactory/prepare.py: load_build_input_csv now reads with pd.read_csv(csv_path, keep_default_na=False, na_values=[""]).
    • Only an empty cell counts as blank, so the required-column check and the blank-cell fallback still work.
    • Numeric columns stay numeric.
    • The setting applies to every column. A required cell containing text such as NA now reaches validation as text, instead of failing early as a NaN.
    • This function is the single CSV entry point, so build, prepare-build, check-csv and push --csv all get the fix.
  • mdfactory/tests/test_simulation_settings.py: new test test_csv_none_override_is_a_value_not_a_blank.
    • The CSV row with None resolves to the same settings as the YAML form {'nvt': {'comm-mode': 'None'}}, and the rendered NVT stage has comm-mode = None.
    • A blank cell in the other row still falls back to the default.
    • The test fails without the fix.
  • docs/content/docs/user-guide/csv-format.mdx: one sentence saying that only an empty cell counts as blank.

Checks (cpu partition):

  • pytest -k "not build" mdfactory: 3408 passed, 1 skipped.
  • test_orchestration_cli.py and test_orchestration_build.py: 37 passed, 2 failed. The two failures are the test_shutdown_parsl_* tests. They fail inside a SLURM job at 571dedc too.
  • ruff check and format (pre-commit): clean.

Finding 2 (useful follow-up), finding 3 (deferred) and the nitpicks in findings 4-6 are not addressed.

Commit: 542972b


Progress tracked by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Unverified Review Candidates — Pending Assessment

Review round: 3
Reviewed commit: 542972b

These are unverified candidates. Severity reflects reviewer confidence; do not treat any item as a merge blocker until the assessment comment is posted.

Round-2 merge blocker (finding 1 of round 2): fixed. load_build_input_csv now reads with keep_default_na=False, na_values=[""], so a CSV cell None stays the string "None" and renders as comm-mode = None, while a truly blank cell still falls back. The code-reviewer and error-auditor ran the CSV loader on probe files and confirmed this. The test-reviewer confirmed that test_csv_none_override_is_a_value_not_a_blank fails without the change.

Round-1 merge blocker (finding 1 of round 1): still fixed. _row_cells restores 100 for whole-number floats in override columns. A fractional value such as tau-p = 2.5 stays a float.

Side effects of the global CSV NA change (checked by three agents): NA, null and None in a required or typed cell used to abort the file with the NaN message. They now fail per row in pydantic. Nothing was found that silently loses a value or produces a wrong result. No repository CSV contains an NA-like token.

Critical

None.

Important

None.

Suggestions

1. Simplify: SimulationSettings.identity() has no caller (simplifier 90). mdfactory/models/simulation_settings.py. BuildInput.hash uses IDENTITY_FIELDS directly. (Seen in round 1.)

2. Simplify: render_stage repeats the nsteps/dt logic in the production and equilibration branches (simplifier 80). mdfactory/run_schedules/gromacs.py. (Seen in round 1.)

3. Simplify: apply_simulation_settings accepts a Path, but only one test uses that branch (simplifier 80). mdfactory/run_schedules/__init__.py. The CLI already loads the file and passes a mapping.

Strengths

  • Both earlier blockers are fixed at the CSV boundary with targeted regression tests that fail without the fix.
  • The completeness check found every deliverable and acceptance criterion of plan revision 2 implemented, with the scope additions disclosed.
  • No repo CSV, protein column or push/check-csv consumer is broken by the stricter NA handling; failures are loud and per row.
  • The test-reviewer found no important missing regression test across identity/hash, supply-path precedence, override validation, duplicate rejection, rendering per type × scheme, the orchestration contract, and the re-render paths.

Agents run: code-reviewer, error-auditor, test-reviewer, completeness-checker, simplifier


Reviewed by mach6

@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Review Assessment

#67 (comment)

Round 3 was assessed by the independent assessor, the developer's advocate and the devil's advocate. The independent assessor also ran test_simulation_settings.py and test_run_schedules.py on the cpu partition: 123 passed. Both earlier merge blockers were verified fixed (_row_cells in 0c6191e; keep_default_na=False, na_values=[""] in 542972b), and each has a regression test that fails without its fix. No agent found an additional concrete merge-blocking problem in the full PR.

Classifications

Finding Classification Reasoning
1. Unused SimulationSettings.identity() deferred Factual: True, no caller (simulation_settings.py:304); BuildInput.hash uses IDENTITY_FIELDS directly. Scope: The plan asks for an explicit list of identity fields, which is IDENTITY_FIELDS and is pinned by a test. Practical: No behavior is affected. Removing the method fits the lean-scope convention and can be done when the file is next touched. The developer's advocate agreed it is not worth its own round.
2. Duplicated nsteps/dt logic in render_stage nitpick Factual: True, two lines repeated. Scope: Style only. Practical: Output is correct. Restructuring the function that writes every mdp would make the regression report and e2e evidence stale for no user benefit.
3. apply_simulation_settings accepts a Path nitpick Factual: True, only a test uses the branch. Scope: The method is the documented public Python API, and build_systems takes the same union type. Practical: No harm. Removing it would make the two entry points inconsistent.

Previously assessed non-blockers are unchanged by the round-3 delta: check-csv/push --csv ignoring the protocol file (useful follow-up), apply-settings batch abort on a corrupt record (deferred), a misspelled override key absorbed by maxwarn together with the overclaiming docs sentence (useful follow-up), and the listed test gaps (deferred).

The devil's advocate found no missing acceptance evidence for the two fixed blockers or for the criterion that all supply paths resolve identically. It needs no new test.

Action Plan

No merge blockers. Nothing to fix before merge.


Assessment by mach6

@maxscheurer
maxscheurer force-pushed the feature/issue-62-generate-mdp-run-schedules branch from 542972b to 6631062 Compare November 10, 2026 09:06
@maxscheurer

Copy link
Copy Markdown
Collaborator Author

Progress Update

Documentation and examples for the whole stack (protein types, Amber route, per-component parametrization, SimulationSettings).

Docs

  • New tutorial user-guide/tutorial-protein-mixedbox.mdx: proteinbox vs protein_mixedbox, Amber and CHARMM routes, simulation_settings, CSV campaign, SLURM executor config, analysis, apply-settings. Linked from the sidebar, index, quick-start, system-types and the mixed-box tutorial.
  • Stale pages fixed: yaml-format (protein tabs, ProteinSpecies table), system-types, csv-format (protein and sizing columns), architecture, workflows, index, quick-start, README, SKILLS, resname limit (3 characters).
  • running-on-hpc: GROMACS 2026 via pixi gmx/gmx-cuda on workers, SimulationSettings with SLURM, max_rescue. cli-reference: config environment. troubleshooting: version gate, pdb2gmx and box-size errors.

Examples

  • examples/pipeline/: CSV, shared protocol, GPU SLURM config (pixi gmx-cuda), sbatch wrapper, single-system YAML with inline settings, README.
  • examples/simulation_settings/: full protocol, short test protocol, apply-settings input, CSV with dotted simulation_settings.* columns.
  • New inputs: multi-chain insulin (proteinbox and merged protein_mixedbox), packed ions, explicit OpenFF family, bilayer at 310 K, fetch_pdbs.sh, examples/README.md.

Verification

  • Example parse tests pass; all changed MDX pages compile.
  • The five new build inputs built on a GPU node.
  • The full pipeline (CSV, protocol, SLURM build, simulate, analysis) ran for all four systems using the short test protocol.

Found while testing, not fixed

  • mdfactory simulate exits 0 when all stages fail.
  • A per-stage gres: null override for EM on the GPU partition fails in GROMACS (GPU visible, no rank count), so the shipped config has no EM override.
  • analysis run <dir> needs prod.xtc; the pipeline uses the summary YAML.

Commit: eed1131


Progress tracked by mach6

@gregorweiss gregorweiss left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

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

please ship

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants