-
Notifications
You must be signed in to change notification settings - Fork 20
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
Bug hunt fixes #225
Bug hunt fixes #225
Conversation
Codecov ReportBase: 92.44% // Head: 92.52% // Increases project coverage by
Additional details and impacted files@@ Coverage Diff @@
## main #225 +/- ##
==========================================
+ Coverage 92.44% 92.52% +0.08%
==========================================
Files 57 57
Lines 4170 4270 +100
==========================================
+ Hits 3855 3951 +96
- Misses 315 319 +4
Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here. ☔ View full report at Codecov. |
@@ -115,8 +116,7 @@ class creation of LambdaProtocol. | |||
context, context_integrator = context_cache.get_context( | |||
compound_thermostate_copy) | |||
# TODO: move to our own MD tasks | |||
# feptasks.minimize(compound_thermodynamic_state_copy, | |||
# sampler_state) | |||
#minimize(compound_thermostate_copy, sampler_state) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The original perses code does the minimization here, however in my opinion it's a bit problematic to have it hidden here. It would be better for this to be a) controlled, b) easy to spot in the protocol.
def _minimize_replica(self, replica_id, tolerance, max_iterations): | ||
"""Minimize the specified replica. | ||
""" | ||
|
||
# Retrieve thermodynamic and sampler states. | ||
thermodynamic_state_id = self._replica_thermodynamic_states[replica_id] | ||
thermodynamic_state = self._thermodynamic_states[thermodynamic_state_id] | ||
sampler_state = self._sampler_states[replica_id] | ||
|
||
# Use the FIRE minimizer | ||
integrator = FIREMinimizationIntegrator(tolerance=tolerance) | ||
|
||
# Get context and bound integrator from energy_context_cache | ||
context, integrator = self.energy_context_cache.get_context(thermodynamic_state, integrator) | ||
# inform of platform used in current context | ||
logger.debug(f"{type(integrator).__name__}: Minimize using {context.getPlatform().getName()} platform.") | ||
|
||
# Set initial positions and box vectors. | ||
sampler_state.apply_to_context(context) | ||
|
||
# Compute the initial energy of the system for logging. | ||
initial_energy = thermodynamic_state.reduced_potential(context) | ||
logger.debug('Replica {}/{}: initial energy {:8.3f}kT'.format( | ||
replica_id + 1, self.n_replicas, initial_energy)) | ||
|
||
# Minimize energy. | ||
try: | ||
if max_iterations == 0: | ||
logger.debug('Using FIRE: tolerance {} minimizing to convergence'.format(tolerance)) | ||
while integrator.getGlobalVariableByName('converged') < 1: | ||
integrator.step(50) | ||
else: | ||
logger.debug('Using FIRE: tolerance {} max_iterations {}'.format(tolerance, max_iterations)) | ||
integrator.step(max_iterations) | ||
except Exception as e: | ||
if str(e) == 'Particle coordinate is nan': | ||
logger.debug('NaN encountered in FIRE minimizer; falling back to L-BFGS after resetting positions') | ||
sampler_state.apply_to_context(context) | ||
openmm.LocalEnergyMinimizer.minimize(context, tolerance, max_iterations) | ||
else: | ||
raise e | ||
|
||
# Get the minimized positions. | ||
sampler_state.update_from_context(context) | ||
|
||
# Compute the final energy of the system for logging. | ||
final_energy = thermodynamic_state.reduced_potential(sampler_state) | ||
|
||
# If energy > 0 kT and on a GPU device attempt to use CPU L-BFGS minimizer | ||
if final_energy > 0 and context.getPlatform().getName() in ['CUDA', 'OpenCL']: | ||
logger.debug(f'Positive final FIRE minimizer energy {final_energy}; falling back to CPU L-BFGS') | ||
sampler_state.apply_to_context(context) | ||
integrator = openmm.VerletIntegrator(1.0) | ||
cpu_platform = openmm.Platform.getPlatformByName('CPU') | ||
context = thermodynamic_state.create_context(integrator, cpu_platform) | ||
sampler_state.apply_to_context(context, ignore_velocities=True) | ||
openmm.LocalEnergyMinimizer.minimize(context, tolerance, max_iterations) | ||
|
||
# Get the minimized positions | ||
sampler_state.update_from_context(context) | ||
|
||
# Get the final energy | ||
final_energy = thermodynamic_state.reduced_potential(sampler_state) | ||
|
||
logger.debug('Replica {}/{}: final energy {:8.3f}kT'.format( | ||
replica_id + 1, self.n_replicas, final_energy)) | ||
|
||
# Clean up the integrator | ||
del context | ||
|
||
# Return minimized positions. | ||
return sampler_state.positions |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This probably can go away for now, but we should try to push something like this to openmmtools. Ideally we need a) optional platform choice for minimisation (that isn't in the context cache), and b) fallback to CPU L-BFGS when GPU optimization doesn't get anywhere.
Also - I think here we want the positions to be reset to their pre GPU minimised values before the CPU L-BFGS fallback. I've had some interesting failures yesterday that indicate that the failed GPU minimisation might take things to a bad configurational state.
@@ -2179,6 +2187,86 @@ def _handle_old_new_exceptions(self): | |||
sigma_new, epsilon_new, 0, 1] | |||
) | |||
|
|||
def _impose_virtual_bonds(self): |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This should be reasonably easy to test. One option is to make this a static method and get it to return the bond force?
b629f87
to
6b3092e
Compare
test for constraint detection added |
@RiesBen can this gathering stuff not be in a separate PR? It doesn't seem like part of the original bug issue? |
Virtual bond addition doesn't seem to be working for inter-chain things, for some reason in one of my system I seem to be getting a bond between a monomer protein and water O.o |
@@ -143,6 +148,79 @@ class creation of LambdaProtocol. | |||
self.create(thermodynamic_states=thermodynamic_state_list, | |||
sampler_states=sampler_state_list, storage=reporter) | |||
|
|||
def _minimize_replica(self, replica_id, tolerance, max_iterations): |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
As discussed, let's remove this for now and move it to a separate PR
still only works for H constraints mind
making room for non-H constraint checking
added test for constraint resolution failure
no protocol_unit_success yet
e643a3b
to
1d8ea72
Compare
Work done in this PR:
To-do: