Skip to content

Add PolarMACE electrostatic ML/MM embedding via graph-longrange - #169

Draft
CheukHinHoJerry wants to merge 22 commits into
openmm:mainfrom
CheukHinHoJerry:feat/graph-longrange-polarmace-embedding
Draft

CheukHinHoJerry wants to merge 22 commits into
openmm:mainfrom
CheukHinHoJerry:feat/graph-longrange-polarmace-embedding

Conversation

@CheukHinHoJerry

@CheukHinHoJerry CheukHinHoJerry commented Sep 17, 2026

Copy link
Copy Markdown

Draft — opened for early feedback on the approach after asking claude to refactor the implementation, not as a finished PR. Trying to adapt to some convention from EMLE too. Two questions:

  • Suggestions on tests that I should be adding?
  • There are a huge portion of the code here spent on checking validity of the MM forcefield but I saw EMLE implementation not doing all these yet. How should I go with this?

Summary

Adds PolarMACE electrostatic ML/MM embedding to the MACE potential, using the reusable external-source APIs from graph-longrange. The MM field enters the model's descriptor and polarises the ML density inside a single forward pass, and the MM back-reaction forces are returned and scattered into the full force array.

Dependency

Requires the external_field branch of graph-longrange, which provides the external-source energy and feature blocks:

pip install 'git+https://github.com/WillBaldwin0/graph_electrostatics.git@external_field'

devtools/requirements/mace-aimnet-torchmd.txt is pinned to that branch which supports electrostatic embedding by allowing PolarMACE taking external point charges input.

How it fits the existing structure

MACEPotentialImpl.getSupportedEmbeddings() advertises "electrostatic", so MLPotential.createMixedSystem() routes that name here and leaves "mechanical" to the generic embedding plugin. addForces() still handles both, since the generic plugin calls it.

Minimal example on how to use the code

from openmm import unit
from openmm.app import PDBFile, ForceField, Simulation, PME, HBonds
from openmm import LangevinMiddleIntegrator
from openmmml import MLPotential

# load system
pdb = PDBFile("system.pdb")
# MM forcefield setup
ff  = ForceField("amber14-all.xml", "amber14/tip3pfb.xml")
mm  = ff.createSystem(pdb.topology, nonbondedMethod=PME, constraints=HBonds)

# The ML subset — the QM region.
mlAtoms = list(range(120))

# using mace potential, can also be some trained MACE-POLAR models
potential = MLPotential("mace-polar-1-small")

# Only PolarMACE-class models advertise it; plain MACE returns [].
assert "electrostatic" in potential.getSupportedEmbeddings()

mixed = potential.createMixedSystem(
    pdb.topology, mm, mlAtoms,
    embedding="electrostatic",
    charge=0,            # ML-region total charge
    multiplicity=1,
)

# standard NVT dynamics
sim = Simulation(pdb.topology, mixed,
                 LangevinMiddleIntegrator(300*unit.kelvin, 1/unit.picosecond,
                                          1*unit.femtosecond))
sim.context.setPositions(pdb.positions)
sim.step(1000)

CheukHinHoJerry and others added 22 commits September 10, 2026 12:08
…polarmace-embedding

# Conflicts:
#	doc/userguide.md
#	openmmml/models/macepotential.py
#	test/TestMACEPotential.py
#	test/TestMechanicalEmbedding.py
The MACE backend requirements still pinned graph_electrostatics v0.4.0, which
predates the external-source work: installing from it yields a graph_longrange
without external_source_energy / external_source_features, so PolarMACE
electrostatic embedding fails at import.

Pin the external_field branch instead, and name it (with the install command) in
the ImportError raised when the blocks are missing.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011TXm1UZZ1hYRkFHZgVFwAN
…rostatic path

Upstream openmm#167 restricts PythonForce to the ML atoms via setParticles(), so the
callback state no longer carries MM coordinates. PolarMACE electrostatic
embedding needs them in the same forward pass (the MM field polarises the ML
density) and returns MM back-reaction forces, so that path keeps the explicit
index slice/scatter and does not call setParticles(). The plain ML path uses
upstream's restriction unchanged; with indices=None and mmInfo=None
_computeMACE is behaviourally identical to upstream.

Same pattern as the EMLE embedding PR (openmm#159): a full-system PythonForce that
indexes ML and MM regions itself.

Verified on OpenMM 8.6.1 (isolated env) with graph_longrange@external_field:
TestPolarMACEElectrostaticEmbedding + TestElectrostaticEmbeddingPME +
TestMechanicalEmbedding: 104 passed, 6 skipped.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
No behaviour change. Verified by fingerprinting the mixed System built before and
after (every NonbondedForce parameter and exception, CustomNonbondedForce
exclusions and parameters, plus energy and forces from a Reference Context at
fixed positions, PBC and non-PBC, sharing one saved PolarMACE model): the two
fingerprints are byte-identical, energies -84.29790229105993 and
-56.26526271299811 kJ/mol in both.

- The `embedding != "electrostatic"` guard is unreachable: MLPotential routes an
  embedding to this method only when it appears in getSupportedEmbeddings(), and
  MACE never advertises "mechanical", so that case already falls back to the
  generic embedding plugin. Kept as an explicit internal invariant and commented
  as such rather than left looking like user-facing validation. addForces() still
  handles embedding="mechanical", because the generic plugin calls it.

- Model support was validated twice, in createMixedSystem and again in
  addForces() via _should_use_mm_embedding(). The helper is now the single owner
  of that message and is called early, so an unsupported model is rejected before
  any force-field surgery.

- Dropped dead work in the NonbondedForce block: the loop zeroing chargeProd for
  every exception touching an ML atom was followed by an unconditional
  addException(..., replace=True) over all ML-ML pairs, which overwrote half of
  it. It now handles ML-MM exceptions only.

- New _customNonbondedChargeIndex() resolves the named charge parameter and
  raises if absent, used both by the early validation and by the force loop, so
  the message has one source and no force indices are carried across systems.

- Copied upstream's comment onto setExceptionsUsePeriodicBoundaryConditions() and
  ordered the gates cheap-first: interpolation, model, force-field structure,
  then surgery. Precision parsing moved ahead of model loading so an invalid
  value fails without loading the model.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011TXm1UZZ1hYRkFHZgVFwAN
- Revert the precision-parsing reorder from the previous commit: the block is
  byte-identical to how it was before the cleanup, parsed after the model is
  loaded as upstream does.

- Restore the sentence this branch had dropped, verbatim from upstream:
  "According to the MACE documentation, 'single' precision is recommended for MD
  (faster but less accurate), while 'double' precision is recommended for
  geometry optimization." It was replaced by a shorter paraphrase in a8b89e1;
  the diff should not silently reword upstream documentation.

Fingerprint of the mixed System (all NonbondedForce parameters and exceptions,
CustomNonbondedForce exclusions and parameters, plus energy and forces from a
Reference Context at fixed positions, PBC and non-PBC) is still byte-identical to
the pre-cleanup baseline.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011TXm1UZZ1hYRkFHZgVFwAN
…ents

Style alignment with openmmml/embeddings/emleembedding.py (upstream PR openmm#159), which is the
reference for how an ML/MM embedding is written here. No behaviour change: the mixed System
fingerprint (every NonbondedForce parameter and exception, CustomNonbondedForce exclusions and
parameters, plus energy and forces from a Reference Context at fixed positions, PBC and non-PBC,
against one saved PolarMACE model) is byte-identical to before, energies -84.29790229105993 and
-56.26526271299811 kJ/mol. Tests: 104 passed, 6 skipped.

- camelCase for identifiers we own: modelDevice, positionsFull, includedAtoms, useMMEmbedding,
  mlAtoms/mmAtoms/mmCharges/mmPositions/mmForces, and the private helpers _supportsMMEmbedding,
  _prepareExternalSources, _enablePolarMACEExternalSources. External API names are untouched and
  must stay snake_case: the MACE input/result dict keys ("mm_positions", "mm_charges",
  "mm_source_batch", "external_field", "node_attrs", "fermi_level", "mm_forces"), the compute_*
  and requires_grad kwargs, and the mace_off/mace_mp/mace_polar model-family strings.

- The mmInfo dict no longer crosses into the callback. _prepareMMEmbedding is unpacked at the call
  site and the values are bound by name, as emleembedding.py does with mlIndices/mmIndices/
  mmCharges, so _computeMACE's signature documents its own inputs:

      def _computeMACE(state, model, ptr, nodeAttrs, batch, pbc, returnEnergyType, charge,
                       multiplicity, periodic, mlIndices=None, mmIndices=None, mmCharges=None)

  The electrostatic path is selected by `mmIndices is not None` exactly where it previously used
  `mmInfo is not None`, and still does not call setParticles().

- _shouldUseMMEmbedding both raised and returned a flag, and one call site discarded the flag, so
  it read as a no-op. It is now _validateMMEmbedding, which only raises; the flag is derived where
  it is used (`useMMEmbedding = embedding == "electrostatic"`). Its "unsupported mode" message now
  names what MACE implements and what it delegates.

- Wrapped the four error messages this branch added that ran to 122-187 characters, and split the
  _computeMACE signature over two lines. Upstream's own long lines and the KNOWN_MODELS URL table
  are left alone.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011TXm1UZZ1hYRkFHZgVFwAN
Removes the two test files this branch added (TestPolarMACEElectrostaticEmbedding.py,
TestElectrostaticEmbeddingPME.py) and restores test/TestMACEPotential.py and
test/TestMechanicalEmbedding.py to upstream. The test/ tree is now byte-identical to upstream and
the PR is source-only: macepotential.py plus the requirements pin.

Tests will be added back once their scope is agreed; the previous set was ~1200 lines across four
files, most of it single-assertion checks on the force-field surgery at a granularity the rest of
this repository does not use.

This also fixes a real breakage: TestMACEPotential.py imported _should_use_mm_embedding, which the
camelCase pass renamed, so that file failed at collection. It was not in the set I had been running,
so the "104 passed" figure reported earlier never covered it. Restoring the upstream file removes
the stale import; 129 tests now collect from the two upstream files.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011TXm1UZZ1hYRkFHZgVFwAN

This branch has not been deployed

No deployments
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.

1 participant