Merge pull request #629 from choderalab/serialize-on-nan

Serialize System, State, and Integrator on NaN fail during propagation, 0.15.0

Yank/repex.py

@@ -48,7 +48,7 @@
 
 LICENSE
 
-This code is licensed under the latest available version of the GNU General Public License.
+This code is licensed under the latest available version of the MIT License.
 
 """
 

Yank/sampling.py

@@ -718,7 +718,7 @@ def _propagate_replica(self, replica_index):
             if (not np.isnan(u_new)) and ((du <= 0.0) or (np.random.rand() < np.exp(-du))):
                 self.rotation_trials_accepted += 1
                 self.replica_positions[replica_index] = perturbed_positions
-            #print("rotation du = %f (%d)" % (du, self.rotation_trials_accepted))
+            # print("rotation du = %f (%d)" % (du, self.rotation_trials_accepted))
             # Accumulate timing information.
             final_time = time.time()
             elapsed_time = final_time - initial_time
@@ -734,7 +734,7 @@ def _propagate_replica(self, replica_index):
         MAX_NAN_RETRIES = 6
         nan_counter = 0
         completed = False
-        while (not completed):
+        while not completed:
             try:
                 # Set box vectors.
                 box_vectors = self.replica_box_vectors[replica_index]
@@ -744,11 +744,12 @@ def _propagate_replica(self, replica_index):
                 if np.any(np.isnan(positions / unit.angstroms)):
                     raise Exception('Initial particle positions for replica %d before propagation are NaN' % replica_index)
                 # Set positions.
-                positions = self.replica_positions[replica_index]
                 context.setPositions(positions)
                 setpositions_end_time = time.time()
                 # Assign Maxwell-Boltzmann velocities.
                 context.setVelocitiesToTemperature(state.temperature, int(np.random.randint(0, MAX_SEED)))
+                # Get velocities in case we have to record them on NaN
+                velocities = context.getState(getVelocities=True).getVelocities(asNumpy=True)
                 setvelocities_end_time = time.time()
                 # Check if initial potential energy is NaN.
                 if np.isnan(context.getState(getEnergy=True).getPotentialEnergy() / state.kT):
@@ -777,7 +778,25 @@ def _propagate_replica(self, replica_index):
                     nan_counter += 1
                     if nan_counter >= MAX_NAN_RETRIES:
                         raise Exception('Maximum number of NAN retries (%d) exceeded.' % MAX_NAN_RETRIES)
-                    logger.info('NaN detected in replica %d. Retrying (%d / %d).' % (replica_index, nan_counter, MAX_NAN_RETRIES))
+                    logger.info('NaN detected in replica %d. Retrying (%d / %d).' % (replica_index, nan_counter,
+                                                                                     MAX_NAN_RETRIES))
+                    # Serialize to XML
+                    from yank.utils import serialize_openmm_object_to_file
+                    positions = self.replica_positions[replica_index]
+                    box_vectors = self.replica_box_vectors[replica_index]
+                    context.setPeriodicBoxVectors(box_vectors[0, :], box_vectors[1, :], box_vectors[2, :])
+                    context.setPositions(positions)
+                    context.setVelocities(velocities)
+                    openmm_state = context.getState(getPositions=True, getVelocities=True, getEnergy=True,
+                                                    getForces=True, getParameters=True)
+                    prefix = 'iteration%d-replica%d-state%d-nan%d' % (self.iteration, replica_index, state_index,
+                                                                      nan_counter)
+                    serialize_openmm_object_to_file(prefix + '.system.xml', context.getSystem())
+                    serialize_openmm_object_to_file(prefix + '.integrator.xml', context.getIntegrator())
+                    serialize_openmm_object_to_file(prefix + '.state.xml', openmm_state)
+                    logger.info('Serialized initial System, State, and Integrator to %s.*' % prefix)
+                    del openmm_state
+
                 else:
                     # It's not an exception we recognize, so re-raise it
                     raise e

Yank/tests/test_yank.py

@@ -13,7 +13,7 @@
 
 LICENSE
 
-This code is licensed under the latest available version of the GNU General Public License.
+This code is licensed under the latest available version of the MIT License.
 
 """
 

Yank/utils.py

@@ -21,6 +21,7 @@
 import numpy as np
 from simtk import unit
 from schema import Optional, Use
+from simtk.openmm import XmlSerializer
 
 from openmoltools.utils import wraps_py2, unwrap_py2  # Shortcuts for other modules
 
@@ -628,9 +629,25 @@ def _combinations_generator(self, leaf_paths, leaf_vals):
             yield copy.deepcopy(template_tree._d)
 
 
-#========================================================================================
+# ========================================================================================
 # Miscellaneous functions
-#========================================================================================
+# ========================================================================================
+
+def serialize_openmm_object_to_file(filename, openmm_object):
+    """Serialize an OpenMM System, State, or Integrator to specified file.
+    
+    Parameters
+    ----------
+    filename : str
+        The file to write to
+    openmm_object : System, State, or Integrator
+        The OpenMM object to serialize
+    """
+    serialized_object = XmlSerializer.serialize(openmm_object)
+    outfile = open(filename, 'w')
+    outfile.write(serialized_object)
+    outfile.close()
+
 
 def get_data_filename(relative_path):
     """Get the full path to one of the reference files shipped for testing
@@ -641,7 +658,7 @@ def get_data_filename(relative_path):
 
     Parameters
     ----------
-    name : str
+    relative_path : str
         Name of the file to load, with respect to the yank egg folder which
         is typically located at something like
         ~/anaconda/lib/python2.7/site-packages/yank-*.egg/examples/

devtools/README.md

@@ -28,7 +28,7 @@ Checklist for all Updates
 Checklist for Major Revisions (YANK X.Y+1.0)
 ------------------------------------------
 - [ ] Make sure all issues related to the milestone will be closed by this commit or moved to future releases
-- [ ] Update `docs/changelog.rst`
+- [ ] Update `docs/whatsnew.rst`
 - [ ] Update `setup.py` with version number and `ISRELEASED` to `True`
 - [ ] Do the steps for All Upates
 - [ ] Create a new release on GitHub, reference the tag and copy the changes in `docs/whatsnew.rst`

docs/whatsnew.rst

@@ -6,13 +6,23 @@ This section features and improvements of note in each release.
 
 The full release history can be viewed `at the github yank releases page <https://github.com/choderalab/yank/releases>`_.
 
+0.15.0 Backend and Helpful Debugging Build
+------------------------------------------
+- Added support for ``solvent_dsl`` in user defined systems of YAML pages
+- Removed Command Line Interface ability to do ``yank prepare`` and ``yank run``
+- Added ability to overwrite individual YAML commands from command line
+- Added YAML feature to ``extend_simulation`` without modifying YAML files or command line every iteration
+- NaN's generated during simulations serialize system, state, and integrator which can be passed off for debugging to others
+- Backend website updating and pushes improved
+- Improved GROMACS extension file handling
+
 0.14.1 Early Access of 1.0 Release
 ----------------------------------
 - YAML Syntax Structure Frozen. YANK YAML Version 1.0. All YAML scripts from this version will be compatible with future versions until YAML 2.0
   New features may appear in the time meantime, but scripts will be forwards compatible.
 - Initial support for OpenMM XML systems and PDB files
 - Support for separate solvent configurations for the two phases when defined from amber/gromacs/openmm files
-- `clearance` in YAML no mandatory parameter of explicit solvent, but only when molecule setup goes through pipeline
+- ``clearance`` in YAML now mandatory parameter of explicit solvent, but only when molecule setup goes through pipeline
 - Boresch Orientational Restraints fully implemented and documented.
 - Long range anisotropic dispersion correction improved to work on both ends of thermodynamic cycle leg
 - Documentation updated with better algorithms and theory sections.
@@ -47,19 +57,19 @@ The full release history can be viewed `at the github yank releases page <https:
 0.9.0 (development)
 -------------------
 - Changed YAML Syntax
-- New Command `yank analyze extrat-trajectory` to extract data from NetCDF4 file in a common trajectory format.
+- New Command ``yank analyze extrat-trajectory`` to extract data from NetCDF4 file in a common trajectory format.
 - Support for solvation free energy calculations.
 - Automatic detection of MPI.
 - Various bug fixes.
 
 0.8.0 (development)
 -------------------
-- `alchemy` split to a standalone repository
+- ``alchemy`` split to a standalone repository
 - YAML based input files for setting up and running simulations. Uses an AmberTools-based pipeline
 
 0.7.0 (development)
 -------------------
-- Convert to single `Context` Hamiltonian Replica Exchange
+- Convert to single ``Context`` Hamiltonian Replica Exchange
 
 v0.6.1 (development)
 --------------------

docs/yamlpages/yaml.rst

@@ -147,6 +147,7 @@ Detailed Options List
     * :ref:`phase1_path <yaml_systems_user_defined>`
     * :ref:`phase2_path <yaml_systems_user_defined>`
     * :ref:`ligand_dsl <yaml_systems_user_defined>`
+    * :ref:`solvent_dsl` <yaml_systems_user_defined>`
     * :ref:`solvent <yaml_systems_user_defined>`
     * :ref:`gromacs_include_dir <yaml_systems_user_defined>`
 

setup.py

@@ -16,8 +16,8 @@
 DOCLINES = __doc__.split("\n")
 
 ########################
-VERSION = "0.15.0dev0"
-ISRELEASED = False
+VERSION = "0.15.0"
+ISRELEASED = True
 __version__ = VERSION
 ########################
 CLASSIFIERS = """\
@@ -39,6 +39,7 @@
 # Writing version control information to the module
 ################################################################################
 
+
 def git_version():
     # Return the git revision as a string
     # copied from numpy setup.py