choderalab / openmmtools

A batteries-included toolkit for the GPU-accelerated OpenMM molecular simulation engine.
http://openmmtools.readthedocs.io
MIT License
235 stars 76 forks source link

Is there an "intended" way to rerun a simulation at new alchemical states? #475

Open mrshirts opened 4 years ago

mrshirts commented 4 years ago

I'm interested in, after running a simulation, generating the energies at a number of new energy functions given the configurations I sampled. openmmtools handles running lots of new alchemical states at the same time as the simulation with the multistate functionality, but in my case I don't know which energy functions to test until after.

I was thinking something like DummyIntegrator was designed for this use case, but it doesn't seem to include a way to move to the next configuration on disk.

Would it still make the most sense to doing an explicit loop over loading new coordinates, calling the DummyIntegrator, and then reading off the energies each step?

Any suggestions about best way to do this with OpenMM tools (or an example!) would be much appreciated.

andrrizzi commented 4 years ago

Would it still make the most sense to doing an explicit loop over loading new coordinates, calling the DummyIntegrator, and then reading off the energies each step?

I'm not aware of facilities for performing this task. I would suggest doing the loop, either with OpenMM or OpenMMTools. With the latter, it would probably look something like this (untested code):

sampler_state = states.SamplerState(...)
thermodynamic_state = states.ThermodynamicState(...)
any_integrator = openmm.VerletIntegrator(...)
context = thermodynamic_state.create_context(any_integrator)

for positions, box_vectors in trajectory:
    sampler_state.positions = positions
    sampler_state.box_vectors = box_vectors

    # This sends the coordinates to the GPU. It's an expensive step.
    sampler_state.apply_to_context(context)

    # Note that this also computes the pV work if thermodynamic_state.pressure != None
    reduced_potential = thermodynamic_state.reduced_potential(context)

I don't see particular advantages in using openmmtools over openmm in this case anyway.

mrshirts commented 4 years ago

OK, I will try to do this to see how it works! Seems like a useful thing to have, since one does not always know the unsampled states one desires before running the simulation.

andrrizzi commented 4 years ago

I just realized I could have misunderstood your scenario. If you have to evaluate each sample over many thermodynamic states, OpenMMTools does offer a speedup if you are using the AlchemicalStates as they will recompute for the same sample only the parts of the Hamiltonian that change from one state to the other. If this is the case, I can put together an example.

mrshirts commented 4 years ago

We do want to evaluate each sample over multiple states. It would probably be easiest to treat them as single perturbations of different parameters, (maybe we want to slightly perturb bonds/angles/and torsions, all in the same pass).

Is it possible to lump together a number of alchemical states that are related by different lambdas (i.e. I want to vary the bonds lengths +/- 10%, the bond force constants by +/- 50%, and the torsional potentials by +/- 50%, combinatorially). Or should I just loop through a bunch of thermodynamic states with openmm?

Any example would be much appreciated.

andrrizzi commented 4 years ago

I just remembered I wrote a section of the advanced tutorial covering exactly this case. Let me know if it helps: https://openmmtools.readthedocs.io/en/0.18.1/devtutorial.html#computing-the-reduced-potential-of-one-configuration-at-multiple-thermodynamic-states

mrshirts commented 4 years ago

Great, I will take a look!