-
Notifications
You must be signed in to change notification settings - Fork 12
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
Add support for replica exchange and REST2 #273
Merged
Merged
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
… bugs might still be lurking somewhere
… also renamed the core alchemical angle/dihedral restraint c++ files
Now fixed. |
…_2' into feature_repex [ci skip]
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
This PR adds the low-level changes to
sire
required to support hamiltonian replica exchange and REST2 simulations usingsomd2
. The PR modifies some existingsire
functionality, adds some new functionality, and provides several fixes and enhancements.Modified functionality:
In order to facilitate replica exchange (repex) it was necessary to modify
dynamics.run()
such that samples (trajectory and energy frames) were not always recorded on exit. When performing repex we usedynamics.run()
to perform many short cycles, only saving energies after each. Performance was massively impacted by the existing code where trajectory frames are always saved. The change closes [BUG] Frames and energies are always saved after dynamics.run() #256, with the sampling frequencies for energies and frames now decoupled from the run method, i.e. the are only ever recorded at a multiple of the frequency. This is validated in this test. Note that this doesn't affect the ability to run a short dynamics plot and commit changes to the system, since the coordinates are updated from the current OpenMM state when.commit()
is called, so is unrelated to trajectory sampling, etc.The handling of
NaN
errors has been tweaked so that it's possible to catch and fix via minimisation when an exception is raised on any step within thedynamics.run()
, rather than just on the first step, as long as no sampling has been performed (energy or trajectory frames). This improves the robustness of the dynamics blocks withinsomd2
and allows us to make use ofauto_fix_minimise=True
for all simulation types. Previously we ran things within our owntry/except
block with custom error handling.New functionality:
REST2 is supported by the new
rest2_scale
andrest2_selection
kwargs to the.dynamics
method. The first defines the scale factor to use for REST2 simulations, which specifies the temperature of the REST2 region relative to the rest of the system. I chose temperature, since this is what is commonly reported in the literature, although some codes use beta instead. I'm not sure if it's easier to think about making the REST2 region hotter or the interactions softer. I'm happy to use the inverse if preferred. (Internally the scale factor is inverted before scaling force constants, etc.) Therest2_selection
is used to specify additional atoms to include in the REST2 region via asire
selection string. By default, all atoms in perturbable molecules are included. The selection provides a way of specifying other relevant atoms, e.g. ones in the protein binding site. When perturbable molecules are specified, then only those perturbable atoms are placed in the REST2 region, i.e. atoms that aren't specified are unaffected. This provides a way of selecting a sub-region of a perturbable molecule if desired.Modifications to dihedral and nonbonded terms are handled in two ways. Firstly we deal with terms in the perturbable region via the lambda lever code. To do so, we create a mask for atoms in the perturbable molecule that are affected by the REST2 modifications, which is determined via the selection described above. Secondly, we need to keep track of which torsion terms are proper dihedrals, since impropers are not modified. With this in place, we can then apply the REST2 scaling to the required terms when
LambdaLever::setLambda()
is called. Terms corresponding to atoms that are not within perturbable molecules are handled separately using the OpenMM Python API via the approach used by meld. The terms that are modified are pre-computed when instantiating theSOMMContext
and are modified when callingSOMMContext.set_lambda()
when arest2_selection
containing non-perturbable atoms has been specified. There is a fairly extensive unit test to make sure the appropriate terms are scaled correctly. This covers: 1) A perturbable molecule; 2) A selection of atoms within a perturbable molecule, and 3) A non-perturbable molecule.When calling
dynamics.run()
a user can now pass an optionalrest2_scale_factors
kwarg. This specifies the REST2 scale factor to use at each of thelambda_windows
that are specified for energy evaluation. We only validate that there are the correct number of entries, not that the scale factors are valid, e.g. they should be 1 at the end states, etc. This validation is done at thesomd2
level.Fixes:
I've fixed an issue with the Sire OpenMM minimiser which meant that exceptions where thrown prior to the thread state being restored, which would cause the program to crash. This closes Minimisation timeout exception not always caught somd2#63.
I've added a small update to the
sire.qm
code to search for theQMForce
by name before returning. (It is always namedQMForce
.) Currently it is always the first force that is added to the system so we can search by index. However, I'd rather make this change now so that we don't run into any surprises if the ordering changes in future. (This change is still covered by an existing unit test.)Closes Sire::IO::SDF can't handle a stereo value of 3 #270 by adding handling for double bond stereo info from SDF files.
Closes [BUG] Positive formal_charge isn't converted to SDF charge field correctly #274 by adding handling for positive formal charge on SDF write.
Checklist
devel
into this branch before issuing this pull request (e.g. by runninggit pull origin devel
): [y]@mb2055: Adding you since you are currently running simulations so can comment on any issues you encounter.
@chryswoods: Adding you so that you are aware of changes to the existing functionality. Let me know if any of these changes impact anything in ways that I might not have anticipated.