optimalcontrol

Python/Rust package for NMR spin dynamics implementing analytical ROPE/CROP theories, numerical GRAPE optimisation, and the Seedless band/restraint front-end for isolated spin-1/2 pulse design, with a Spinach-compatible Python API, ensemble support, and paper-reproduction examples. GRAPE propagation, exact coherent gradients, Bloch ensemble profiles, and the analytic spin-1/2 Seedless kernel run in a parallel native Rust extension.

Installation

pip install optimalcontrol-nmr

The import name is still optimalcontrol — only the PyPI distribution is named optimalcontrol-nmr (the bare name is blocked by an existing similar project).

Python 3.10 or newer is required. Environments without a compatible prebuilt wheel also need a stable Rust toolchain so pip can compile the native extension from the sdist.

Getting started

New to the package? Follow the step-by-step beginner manual in docs/user_manual.md: install, design a pulse by two different routes (GRAPE and Seedless), read the result, verify it against an independent Bloch model, and write a shape file. The topic-specific docs/guide_*.md files cover each subsystem in depth.

MCP server (use from Claude and other LLM clients)

The package ships an MCP server exposing the Seedless route, so an LLM client can design and verify pulses conversationally.

Easiest — install as a Claude Code plugin (this repository is its own plugin marketplace; requires uv, no pip install needed — the server is fetched from PyPI automatically):

/plugin marketplace add deepnmr/optimalcontrol
/plugin install optimalcontrol@optimalcontrol

Or register the MCP server directly:

pip install "optimalcontrol-nmr[mcp]"
claude mcp add optimalcontrol -- optimalcontrol-mcp   # Claude Code

No standing server is involved — the client launches optimalcontrol-mcp on demand and talks to it over stdio. Two tools are exposed:

• design_seedless_pulse — declarative band/restraint input; runs a multi-seed phase-only optimisation, verifies the winner against the independent Bloch model, and writes a Bruker .shape file plus a design JSON. Returns file paths and summary numbers only.
• bloch_offset_profile — reloads a design JSON and computes the mx/my/mz offset response with the independent Bloch model.

Any MCP-capable client works; for others, register the same optimalcontrol-mcp command in that client's server configuration.

Source references

• JMR 2003 (ROPE): Unterbeck & Glaser, Journal of Magnetic Resonance 160 (2003) 88–101 — analytical optimal control for heteronuclear transfer under relaxation.
• PNAS 2003 (CROP): Unterbeck & Glaser, Proc. Natl. Acad. Sci. USA 100 (2003) 5172–5177 — cross-correlated relaxation-optimised pulses.
• Nat. Commun. 2025 (Seedless): Buchanan et al., Nature Communications 16, 7276 (2025) — "Seedless: on-the-fly pulse calculation for NMR experiments"; the band/restraint formalism implemented by optimalcontrol.ocseed.
• Spinach: https://spindynamics.org/wiki/index.php?title=Main_Page — MATLAB spin dynamics library whose grape_xy / control struct API this package mirrors.

Local development commands

Install a stable Rust toolchain (rustc and cargo) first. On macOS with Homebrew:

brew install rust

# Install in editable mode with dev and MCP dependencies
python3 -m pip install -e ".[dev,mcp]"

# Lint
ruff check optimalcontrol/ tests/

# Typecheck
mypy optimalcontrol

# Test
python3 -m pytest

Set OPTIMALCONTROL_DISABLE_RUST=1 to run the NumPy/SciPy fallback for numerical comparisons. Normal installations build and use the Rust extension automatically. GitHub Actions runs the full native suite, targeted fallback tests (including MCP), Ruff, mypy, and Clippy on Ubuntu with Python 3.12 and stable Rust for pull requests and pushes to master. See the CI workflow for the exact commands; Rust accelerator tests only run in the native suite.

Performance

On an Apple Silicon development machine, the native path reduced the full example regression runtime from 20.08 s to 2.91 s (6.9x). The focused 72-slice GRAPE benchmark improved single-member gradient evaluation by 2.2x and five-member ensemble gradient evaluation by 3.4x. Run python benchmarks/bench_grape_hotpath.py and pytest tests/test_examples.py to measure the local machine.

Seedless per-step suppression now uses one cumulative adjoint sweep, making its gradient linear in pulse length while preserving the hold objective at every prefix. Run .venv/bin/python benchmarks/bench_seedless.py for 32-, 128-, and 512-step cases; add OPTIMALCONTROL_DISABLE_RUST=1 to measure the NumPy fallback. The GRAPE benchmark also includes density-matrix gradients, which now use batched matrix products to avoid repeated tensor-contraction planning.

Measured before/after on Apple M4 (Python 3.14.6, NumPy 2.5.0, SciPy 1.18.0; median of five timing runs, seconds per objective-plus-gradient call):

Case	Before	After	Speedup
Seedless suppression, Rust, 512 steps, 21 members	0.01562	0.0004694	33.3x
Seedless suppression, NumPy, 512 steps, 21 members	1.333	0.008141	163.8x
Density GRAPE, NumPy, 72 steps, 2x2 state	0.003084	0.0006324	4.9x

These are focused workload measurements, not end-to-end optimizer speedups. In these comparisons, objective values matched exactly and the maximum absolute gradient difference was 1.2e-16 or less. Timings vary with hardware and workload.

Symmetric methyl 180 pulse with water preservation

python -m examples.methyl_water_binary_symmetric_180 writes a Bruker shape and a dense-grid diagnostic plot for a 1.2 GHz proton spectrometer with the carrier at water (4.7 ppm). The 1.740 ms pulse covers methyl protons from -3 to 3 ppm, preserves water Iz, uses only 0/180 degree phase, is exactly time symmetric, and is capped at 10 kHz. It implements the inversion-quality requirement motivated by Kay's methyl-HMQC refocusing-artifact analysis (JBNMR 2019).

The cached candidate was validated at 2401 methyl offsets and 9 water offsets. Worst fidelities are 0.999186 (Ix -> Ix), 0.999128 (-Iy -> Iy), 0.999098 (Iz -> -Iz), and 0.999892 for water Iz -> Iz; the worst predicted inner artifact is 0.09027% of the central line. It is the shortest passing point on the tested local duration grid: 1.740 ms passed while 1.735 ms failed. This is a numerical grid-search result, not a proof of the continuous global minimum. The recorded duration audit is in examples/expected/methyl_water_binary_symmetric_180_duration_search.csv.

Versioning policy

This package follows semantic versioning: MAJOR.MINOR.PATCH.

• PATCH releases contain backwards-compatible fixes only and must not introduce breaking API, file-format, or numerical-contract changes.
• MINOR releases may add features and deprecate existing APIs; deprecations must emit warnings before removal.
• MAJOR releases may remove deprecated APIs or introduce intentional breaking changes, with migration notes recorded in CHANGELOG.md.