Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

vdW section should support different methods for periodic and non-periodic systems #51

Closed
mattwthompson opened this issue Mar 17, 2023 · 9 comments · Fixed by #53
Closed
Labels

Comments

@mattwthompson
Copy link
Member

When creating an openmm.System using a gas-phase system,

  • the prescribed vdW "method" method is "cutoff"
  • the Topology has no box vectors
  • the resulting openmm.NonbondedForce uses NoCutoff
>>> import openmm
>>> from openff.toolkit import ForceField, Molecule
>>> from openff.interchange import Interchange
>>> from openff.units import unit
>>>
>>> force_field = ForceField('openff-2.0.0.offxml')
>>> topology = Molecule.from_smiles("CCO").to_topology()
>>> assert topology.box_vectors is None
>>> interchange = Interchange.from_smirnoff(force_field, topology)
>>> assert interchange.box is None
>>> system = interchange.to_openmm(combine_nonbonded_forces=True)
>>>
>>> for force in system.getForces():
...     if isinstance(force, openmm.NonbondedForce):
...         assert force.getNonbondedMethod() == openmm.NonbondedForce.NoCutoff
...
>>> # (runs with no errors)
>>> print(force_field['vdW'].method)
cutoff

This seems obviously incorrect. Using a cutoff is not the same as not using a cutoff. (SMIRNOFF does not fully describe what "cutoff" means, but whatever it means, it probably does not mean "use NoCutoff.) I've brought this up in the past to mixed reception.

This behavior has been in the toolkit since the beginning and I believe that, as a result, all force fields have been fit with this wart in place. What should the behavior here be @lilyminium? Follow the specification or track the existing reference implementation? Below the divider, I will ramble on with more context, but it can be ignored if there's a clear answer to this question.


Splitting out non-bonded forces is not a solution to the problem of mixing PME electrostatics and cutoff vdW interactions. OpenMM will throw an exception, but Interchange already pre-empts this and forbids this combination from existing:

>>> system = interchange.to_openmm(combine_nonbonded_forces=False)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/Users/mattthompson/mambaforge/envs/smirnoff-plugins-test/lib/python3.9/site-packages/openff/interchange/components/interchange.py", line 349, in to_openmm
    return to_openmm_(self, combine_nonbonded_forces=combine_nonbonded_forces)
  File "/Users/mattthompson/mambaforge/envs/smirnoff-plugins-test/lib/python3.9/site-packages/openff/interchange/interop/openmm/__init__.py", line 89, in to_openmm
    particle_map = _process_nonbonded_forces(
  File "/Users/mattthompson/mambaforge/envs/smirnoff-plugins-test/lib/python3.9/site-packages/openff/interchange/interop/openmm/_nonbonded.py", line 108, in _process_nonbonded_forces
    _func(
  File "/Users/mattthompson/mambaforge/envs/smirnoff-plugins-test/lib/python3.9/site-packages/openff/interchange/interop/openmm/_nonbonded.py", line 743, in _create_multiple_nonbonded_forces
    raise UnsupportedCutoffMethodError(
openff.interchange.exceptions.UnsupportedCutoffMethodError: When using `openmm.CustomNonbondedForce`, vdW and electrostatics cutoff methods must agree on whether or not periodic boundary conditions should be used.OpenMM will throw an error. Found vdw method 1, and electrostatics method 0,

Arguably this is a quirk or OpenMM and maybe there's a workaround I'm unaware of. The fact that vanilla NonbondedForce usage combines PME electrostatics and cutoff vdW interactions is awfully convenient for that case, but doesn't work for other cases, such as when CustomNonbondedForce must be used.

Modifying the non-bonded methods to differ from what is specified in the force field is, in my opinion, a non-starter.

It's possible that using no-cutoff interactions in the gas phase is functionally equivalent to using cutoff interactions in the condensed phase provided an isotropic, assuming that the relevant assumptions are valid and the error introduced by the approximation is sufficiently small. It certainly looks good in papers, but I have no clue if that matches what's observed in the wild. This would mean there's functionally little or no difference between gas-phase no-cutoff calculations used in (valence) fitting and how OpenFF force fields are implemented in the condensed phase. Still, having these quirks cancel each other out is separate from how confusing this is to implement.

OpenMM's NoCutoff method might not have clear analogs in other packages. I think Amber supports it in some sense but the documentation is nearly impossible to search. GROMACS does not properly support gas-phase simulations in versions after 2019 and nobody is actively working on it. I'm not exactly sure about LAMMPS, but there are enough pair styles implemented for me to guess that it can be done. I have no idea about CHARMM, Desmond, etc. nor do I track how it's implemented in the novel ML-based engines.

The number of combinations to consider grows pretty rapidly the more you look. Periodicity, combinations of electrostatics and vdW methods, their implementations in different engines, the use of non-LJ functional forms, maybe the use of implicit solvent, ... the list goes on and is unlikely to be trimmed down without removing sections of SMIRNOFF. The two use cases I think we should be concerned about are using a main-line OpenFF force field (PME + cutoff vdW) with a topology that is either periodic or not, but I want to highlight the background of complexity that is can be introduced when making even a single change.

@jthorton
Copy link

jthorton commented Jun 1, 2023

Hi @mattwthompson I have been looking into this as I am now trying to fit double exponential force fields to gas phase QM dimer interaction energies. I believe that in order to match the QM reference I would need my MM system to be fully interacting with NoCutoff electrostatic and vdW terms similar to the workaround we had to do in smirnoff-plugins. I have considered changing the vdW method to NoCutoff but this becomes an issue when doing joint training to physical properties and dimer energies as I want to use a cutoff in the periodic simulations.

I also believe that it is a similar situation with the valence fitting as we are trying to reproduce the gas phase QM energies of the molecules and so the MM model should be fully interacting with no cutoffs. I think then it's best to have interchange match the old behaviour.

@mattwthompson
Copy link
Member Author

mattwthompson commented Jun 6, 2023

I understand the value in matching the previous behavior, but I also am not sure it was correct. This is a difficult pair of options to choose between, I think.

What do you think @lilyminium?

@lilyminium
Copy link
Contributor

Ugh, this is tricky. So let me get this straight:

  1. without a box, combine_nonbonded_forces=True: vdw = "cutoff", electrostatics = "pme" results in vdw="NoCutoff" and electrostatics="NoCutoff"
  2. without a box, combine_nonbonded_forces=False: error, because OpenMM requires cutoff agreement
  3. with a box, combine_nonbonded_forces=True: vdw = "cutoff", electrostatics = "pme" results in vdw="pme" and electrostatics="pme"
  4. with a box, combine_nonbonded_forces=False: vdw = "cutoff", electrostatics = "pme" results in vdw="CutoffPeriodic" and electrostatics="pme"

From a technical/user-friendly standpoint, this is pretty suboptimal -- it really seems that literally following the spec should be done with CutoffNonPeriodic. However, IMO the points raised in #29 are still valid. I'm very keen to avoid suddenly giving different results for the same code, but we should be clearer to users what's going on. Could the same solution there apply to the vdW terms, i.e. specifying the periodic/non-periodic methods?

As to whether the non-periodic method should be NoCutoff at all -- that's a bit harder. It is odd that OpenMM can't support CutoffNonPeriodic and NoCutoff in the same system. But, I think, in general you would want your vacuum simulations to mimic that behaviour -- especially if it's trained that way -- and to the best of my knowledge you can emulate NoCutoff in both AMBER and GROMACS with big boxes and cutoffs.

@jthorton
Copy link

jthorton commented Jun 6, 2023

+1 for being able to specify the periodic/non-periodic methods this should be explicit and flexible.

@mattwthompson
Copy link
Member Author

Here's a more concise summary of what happens than I wrote earlier:

Periodic: True, Combine: True
	[('NonbondedForce', 'PME')]
Periodic: True, Combine: False
	[('CustomNonbondedForce', 'CutoffPeriodic'), ('NonbondedForce', 'PME')]
Periodic: False, Combine: True
	[('NonbondedForce', 'NoCutoff')]
Periodic: False, Combine: False
	When using `openmm.CustomNonbondedForce`, vdW and electrostatics cutoff methods must agree on whether or not periodic boundary conditions should be used.OpenMM will throw an error. Found vdw method 1, and electrostatics method 0, 

(1 is CutoffNonperiodic, 0 is NoCutoff)

I think - please correct me if I'm wrong! - Josh is asking for the last case to be NoCutoff and NoCutoff, mimicking what happens in the combine_nonbonded_forces=True case. It's really unlucky that this can't be used with these SMIRNOFF plugins. My key objection here is that SMIRNOFF doesn't allow for no-cutoff vdW interactions; I'm really hesitant to have a <vdW ... method="cutoff"> force field produce vdW interactions that aren't cut off. I'm really not a fan of the Periodic: False, Combine: True producing this; I think @j-wags talked me out of throwing an exception there. I see now the contradiction of how this is handled differently with combine_nonbonded_forces=False, which I don't think is totally intentional. In any case, if SMIRNOFF were to split out the periodic and non-periodic cases analogously to what resulted from OFF-EP 0005(b), this would no longer be true.

It appears we're in agreement (Lily suggested, Josh offered support, I also support) an OFF-EP that splits out vdW method into two different options, one for periodic and another for non-periodic systems. This would make the force fields a tiny bit more verbose but help me out a ton, make the OpenMM workflows a bit cleaner, and hopefully either improve or not impact the GROMACS/Amber workflows. If I'm understanding everything accurately, I think this is a good proposal - should we move this to the standards tracker?

@lilyminium lilyminium transferred this issue from openforcefield/openff-interchange Jun 6, 2023
@lilyminium
Copy link
Contributor

Thanks @mattwthompson for summarising and clarifying this. I'm in agreement and have moved the issue now :)

@mattwthompson mattwthompson added SMIRNOFF and removed help wanted Extra attention is needed labels Jun 13, 2023
@mattwthompson mattwthompson changed the title openmm.NonbondedForce does not respect vdWHandler.cutoff in the gas phase vdW section should support different methods for periodic and non-periodic systems Jun 13, 2023
@mattwthompson
Copy link
Member Author

This is largely a duplicate of #7 which is celebrating its second birthday today

@mattwthompson
Copy link
Member Author

The relevant changes I'll be proposing (keep an eye out for the EP):

  • potential stays the same (continuing to encode 12-6 LJ or an alternative)
    • (This is a somewhat unfortunate language difference between the Electrostatics section, since "potential" has a different meaning here.)
  • method="cutoff" becomes split into periodic_method and nonperiodic_method.
    • Supported values would be periodic_method="cutoff" and nonperiodic_method="no-cutoff" (perhaps spelled differently)
    • In the future, but not in this EP, there might be something like periodic_method="PME" (again, perhaps spelled differently)
  • Implied default values, and the implied up-conversion with existing SMIRNOFF force fields, would map from method="cutoff" to periodic_method="cutoff" and nonperiodic_method="no-cutoff"

@jchodera
Copy link
Member

jchodera commented Aug 1, 2023

@mattwthompson : This sounds like a good solution to the immediate problem.

By analogy to the Electrostatics tag, we could later consider adding support for something like implicit solvent (e.g., GBSA) or other coase-grained representation.

Note that, if OpenMM was ever creating a NonbondedForce with CutoffPeriodic for non-periodic systems without setting the solvent dielectric to 1, this would have resulted in using a functional form for the electrostatics that dampens the Coulomb force.

We should be mindful that we still have to address the ambiguity in the long-range vdW handling at a future date. I have a proposal we can review.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging a pull request may close this issue.

4 participants