You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
{{ message }}
Repository navigation
Address confusion due to similarity of Simulation and Simulator names #65
We have two objects, Simulator and Simulation, which have very similar names and very different purposes. I think that will be terribly confusing for users. In addition, we have a Calculation class whose purpose might also be confused with the other two.
Here is my description of what each object does (@jhprinz, let me know if I'm missing any aspects):
Simulation: interacts directly with the underlying MD engine. It just transfers state information (positions, velocities, etc.) back and forth between MD engine and our code.
Simulator: this one is a bit harder to describe. It mainly involves two functions: (1) setting up the system (including setting up our storage and setting up topologies, force fields, etc. for the MD engine); (2) the generate function, which runs the MD and stops when we reach a stopping condition.
Calculation: think of this as the main() for whatever kind of sampling we're doing: it's the part that distinguishes the overall sampling algorithm, e.g., making TIS different from SATIS or AMS.
Originally the Simulator/Simulation divide was because we needed our own class (Simulator) to interact with OpenMM's Simulation. I think we can clean that up now; whether we also need to change the Calculation name is up to you.
Here's my suggestion for refactoring:
Let specific Simulations (e.g., OpenMMSimulation) inherit from Simulator (which might be renamed to Simulation). For cases like OpenMMSimulation, which currently inherits from openmm.app.Simulation, we could use multiple inheritance, or we could make the openmm.app.Simulation into an instance variable for the OpenMMSimulation. So there would be a superclass with the code to set up the storage and for the generic generate function, and the subclasses would support system set-up and dynamics state information for OpenMM, or for my ToyDynamics, or (perhaps one day) for Gromacs.
From a user perspective, this would mean that instead of having to build both a Simulator and a Simulation, you just need a single object which combines information from both. In fact, the final object is kind of an OpenTIS wrapper around the MD engine's objects, providing a few extra functions and the ability to translate to and from our data structures.
If this seems like a good idea, I can implement it. Any thoughts?
This is a very good important point to address. I also think this is right now very confusing. I guess even to us. Let me first lay out the differentiation between Simulation and Simulator, which is what I originally intended and not what it necessarily has to be:
Simulation() : This is the OpenMM simulation object
Simulator() : This was supposed to be the wrapper for the Simulation that contains information about the system to be simulated and would correspond more to OpenMM Simulation class. This is the part that generates snapshots, momenta, configurations etc, knows the topology (and thus the number of atoms which is necessary for the storage. So far it is also attached to the storage, but that should not be necessary.
So, in theory all we would need is a Simulator() object that knows all about the simulation in an abstract way. It knows how to run a simulation (including setting stopping conditions, etc.), it can give Snapshots/Trajectories.
Maybe it just comes down to almost copying the OpenMM scheme.I think there are in general several steps involved to specify how a simulation works. Usually we want:
Start with a topology that described our molecule in a phycal way
Generate a model that can accurately describe the molecule using the topology (by specifying used force fields, PME, constraints, cut-offs, etc). This effectively creates a computable potential energy function.
Add a configuration and a momentum which must be compatible. with the topology/model
Select an integrator which specifies the differential equation to be solved AND the numeric way on how to do it. This contains thermodynamic variables, etc...
Run a Simulator that uses all these ingredients to produce new configurations/momenta in the specified physical ensemble.
Does this make sense?
A problem that I encountered with OpenMM is that their System object does contain the energy function but no reference to the original topology. So this is kind of missing especially for analytical purposes although not required to do a simulation. Also, a topology is something very special used in Force-Field based energy functions, but we might not need this for simple toy models.
So I think we should have ideally 4 objects:
Model() or System() that contains a describtion of the energy function and stuff that was used to generate it. For modelcules the topology as well as used forcefields, etc... For 2D toy models may only the functional form of V(x). In a mathematical sense this know the dimensionality of x(t) in our differential equation and a way to compute/approximate V(x)
Dynamics() or Integrator() that (again) contains the actual way of doing the integration (e.g. VVVRIntegrator) as well as the information used to generate (Langevin/Brownian, temperature, physical ensemble, thermostat, barostat). In mathematical language this know the form of the differential equation as well as a numerical way to compute steps.
Snapshot() or State() that contains the coordinates and velocities of system. In mathematical terms this is a particular instance of x or x(t) if part of a trajectory. This object consists about a Configuration and a Momentum object.
Simulator() or Simulation() the overall object that knows all of the above and can generate new Snapshots() by using the integrator
Questions:
Where to place the tilmestep? Part of the integrator or part of the simulator? Both make sense although from my description it would be more reasonable to have this as part of the integrator.
Is this generic enough for a topology? If we e.g. want to use mdtraj to output to a known file format we could use the topology in the model or a generic one (without any bonds) depending on the model type. or a Model knows how to generate a mdtraj.Topology object
Each Snapshot() should contain a reference to the Model used!
Is this overkill? I imagine then in the final code something like
I think you're suggesting something very similar to what I mean, although I think you're focusing on a few aspects that I don't think we need to worry about. Different simulation objects don't need to have the same internal structure, as long as they implement the methods we need to interact with them.
Let's talk about two different layers to our code: the path sampling layer, which just functions on trajectories, and the dynamics engine layer, which generates new frames of a trajectory.
I agree that all engines can boil down to your model/dynamics/snapshot model (that's how I implemented both my grad-school C and post-doc C++ codes, although I used the names potential/integrator/system, resp. -- even a file reader can be an "integrator"/"dynamics").
However, I think it might be cleaner and more straightforward for anyone who wants to add a new engine to ask them to implement a few simple methods than to ask them to structure their objects in a certain way, especially since we _need_ the methods to have the right behavior, and _all_ we need is for the methods to have the right behavior.
Now the long version:
Path sampling-like approaches function on trajectories (or, more correctly, ordered [usually time-ordered] series of vectors; because of our MD bias, we call those ordered series "trajectories"). As long as the path sampling layer sees something it can identify as a trajectory, it doesn't care how that trajectory was generated.
Then there's the dynamics engine. This does all the dirty work of generating the trajectories, frame by frame. The internal structure of the dynamics engine can vary widely: as you mention, the topology/force-field split makes sense in molecular mechanics, where you decompose interactions into N-body (for fixed maximum N) interactions. However, for ab initio MD, the "topology" is just the isotope identification, and for simple models where you explicitly specify the entire PES, it's just the masses. For a generic system of differential equations, there's no distinction between topology and force field (since the purpose of separating them is to make it easy to define a generic system of differential equations.)
My suggestion is that Model/Integrator/System all get wrapped up into the Simulation object. How that Simulation object implements those internally can be different for each engine. However, each such object needs to implement an interface consisting of the following functions (we can adjust these, but at least this much behavior needs to available):
init_simulation_with_snapshot(snapshot): takes snapshot, a state of the simulation, and initializes the simulation object in that state
load_snapshot(snaphost): loads the current simulation state into snapshot
generate_next_frame(): generates the next frame according to the simulation, and returns it as a Snapshot object
stop(): tells the simulation to stop running
I also implemented load_configuration and load_velocity to make things easier with the other code we already had, but those four are the core. I'm also fully in favor of switching the first two to a @property called current_snapshot; that requires more refactoring than I wanted to do in #48 without clearing it with you, but I think it is a better solution. (Above, you call it snapshot; clarifying "current" might be worth the extra characters.)
So why those four core functions? Well, you know how they work for the OpenMMSimulation or ToySimulation examples. So let's talk about a GromacsSimulation object. It would be pointless to launch an entire new Gromacs run every saved frame, so the best approach would be to have Gromacs running as a sort of coprocess. The first two functions are, of course, pretty trivial (although initializing the simulation involves writing a file and loading some additional file-related info into the GromacsSimulation). But the last two functions are interesting.
For the GromacsSimulation, the function generate_next_frame() doesn't actually generate the frame. Instead, it waits until a new frame is written to the trr file, and then it reads that frame from the file. And then there's the stop() function; when we generated one frame at a time, stop() didn't do anything. However, the GromacsSimulation doesn't know to stop until we tell it. The stop() function here will stop the Gromacs process, and do whatever cleanup is needed.
(For what it's worth, I think Tsjerk Wassenaar has already implemented most of the stuff we need for the Gromacs get_next_frame() and stop() functions as part of the tipsi.py script he's writing to implement a very user-friendly TIS/TPS interface for Gromacs.)
Note that I've replaced your run(nsteps) function with a get_next_frame() function. That's all we really need. We leave it to the implementer to define a number of steps between frames, because the path sampling layer doesn't need to know about that.
All of this was my thinking in #48 and #52, although I guess I didn't discuss it in detail anywhere.
Another reason not to define the model and dynamics as separate objects is that they are not necessarily independent of each other. In particular, they're not likely to be compatible between different simulation engines. For example, it is meaningless to use an OpenMM-defined model on my LangevinBAOABIntegrator dynamics object: the interfaces for OpenMM's force fields are unknown to that integrator. Enforcing the extra structure is pretty on paper, but can't be used in practice.
In general, I think that we should use Python's duck-typing to enhance our flexibility instead of trying to create super-general structures that work with everything. I say this because I spent a lot of time trying to do the latter: that approach led to messes in my C code that moved me to C++, and then to messes in my C++ code that made me eager to move to Python. (I'm sure I'll eventually screw up in Python, but at least I learned what to avoid from C and C++!)
I think the answer to your questions 1 and 2 are "we don't need to dictate the internal structure of a Simulation object as long as it correctly implements the interface", and for question 3, I think the trick is to make it hard (impossible?) to keep more than one Simulation in a single storage. If the (unique) simulation parameters are stored in the storage, then no problem. And as for the user mixing snapshots from two Simulation objects by incorrect instantiation and loading... well, I'm willing with label that a 10T error/nofix.
Finally, a summary of what initial intent of this issue was:
The above info is pretty much how I implemented things in #48 and #52. I intended this issue to propose a change to combine the Simulator and Simulation objects. As we both mentioned above, the separation between the two was originally because the Simulation was the openmm.app.Simulation object. The Simulator was our own. However, as of #48, the Simulation is our own, but a class that inherits from openmm.app.Simulation. What I'm suggesting here is essentially to have all Simulation classes inherit from something that is currently our Simulator class. Since the Simulator is independent of the dynamics engine, this should be possible. Then we only have one class instead of two.
Great. I like these detailed discussions a lot.
Okay, so we agree that my attempt seems nice and ordered, but might be too overkill and unnecessarily restrictive. I just thought that since OpenMM is using this and we were mainly aiming at OpenMM it seemed natural to make this into a clear seperation of the necessary parts.
I would disagree that the seperation into dynamics and system etc need to be interchangeable between different simulations. Maybe I wasn't clear on that. Right now we are effectively using some seperation like this, at least ToyDynamics has a seperate Integrator and so has OpenMM. Of course each integrator does only make sense in the simulation context. So you are right that we do not need to seperate these.
One point I thought about, when doing this was the way to store this and keep simulations reusable with different parameters. In principle we are saying the same thing, just different parts. You described what simulation should do and I was trying to explain how to set the parameters of simulation in a structured way. But back to the original problem: Restructuring of Simulation
I would do it as you said and have only a Simulation object that allows to run a simulation (whatever that is) that generates snapshots in some way. I really like the idea of an asynchronous way of handling, with just starting and then getting elements as they become available. If we already have this in place why not keep this really asynchronous and add callbacks once a new frame is available. This way we could even run stuff on external clusters and continue once a new frame is available, it would even allow running several processes in parallel. But this might again be a little to much for now.
So, I will sum up Davids idea of functions for Simulation()
init_simulation_with_snapshot(snapshot) : takes snapshot, a state of the simulation, and initializes the simulation object in that state
load_snapshot(snaphost) : loads the current simulation state into snapshot
generate_next_frame() : generates the next frame according to the simulation, and returns it as a Snapshot object
stop() : tells the simulation to stop running
I think I would simplify the names of the functions and treat this like a generator so. Also, if we think of this as a worker then we would just start (which does all the initialization and (in case of Grimaces) starts the simulation. Have a function to generate the next snapshot
start(snapshot) : takes snapshot, a state of the simulation, and initializes the simulation object in that state and start the simulation for an independent runner or waits (like OpenMM)
next() : generates the next frame according to the simulation, and returns it as a Snapshot object
stop() : tells the simulation to stop running. For gromacs this will stop the background process and remove the left-over simulation objects.
for convenience I would add properties, which will change after each call of next()
snapshot : return the current Snapshot() object
configuration : return the current Configuration() object
momentum : return the current Momentum() object
In the future we could also introduce so called futures and make next() non-blocking to make next() return a snapshot object that is empty, but the rest of the code will continue until the real content of the Snapshot is needed. This is similar to the way I implemented the loading of Snapshots. This does not work with OpenMM, but could work with Gromacs and means, we can run several simulation.next() commands that will all run a Gromacs simulation and once we do something with these the python program waits. To allow for multiple simulators we just need to be able to create several Simulation() objects. We could even have remote Simulations using Siegetank. How cool would that be.
Again, David, this idea with implementing a next function is fantastic. So I think this is settled.
Quickly back to my idea with a Dynamics/System separation. Since it is really independent of this PR we might not do this at this point. Although I think we could do this for our OpenMM and ToyModel Approach so encourage a more structured code. Kind of as a good practice rather then a must. The main object is just the simulator, however this is instantiated.
Leaves a last thing to think about. How to serialize the Simulation() object? We can use the jsonify algorithm which will turn a simulation into a json string as long as there are no complex objects in it. Right now we have only strings and some number with units. I will add the capability to the json-netcdf serializer to also accept things with a unit. This means in theory we can also just store simulation objects.
So, enough :) David if you want to start implementing a new Simulation() go ahead.
A few related threads going here. Let me try to first describe how I see the difference/connection between the two major ones:
On one hand, we're talking about the possibility of combining what was previously two objects, a Simulator and a Simulation, into one object by having the generic behavior of a Simulator subclassed by individual Simulations for each dynamics engine. It sounds like we're agreed that this is a good idea. On this part, there's still some discussion to have about the specifics of the API we establish (which has to implemented separately for each dynamics engine).
On the other hand, we're talking about the specifics of the implementation of the engine-specific Simulation subclasses. On that, I think we're also in agreement that we shouldn't enforce a specific internal approach for these (since they can't trade information between them); however, we should demonstrate best practices in our own implementations. Here, I still have some questions about what qualifies as best practices in some specific situations.
So a few things about each of these, and then a few other little threads around the edges of this....
Our API for dynamics engines
When you say, "treat this like a generator", do you mean that we might actually make this into a proper generator? I haven't written anything using generators before. This seems like a perfect example for the use of yield, but I'm not sure whether we actually get much benefit from that, since we can't actually do anything else with our trajectory (identify ensemble) until the whole trajectory has been returned.
If we don't make it a proper generator, I think we might be better using a more explicit name than just next(). I think simulation.next() isn't quite as clear as simulation.next_frame() (just next could as easily refer to getting the entire next trajectory, instead of just the next frame).
The idea having a non-blocking next() is really interesting and definitely appealing when thinking about things like waiting for a coprocess to write a new frame. However, I'm not sure how we can use that idea to our advantage. Isn't the algorithm effectively blocked until we get that frame and can analyze it to determine whether we continue the trajectory?
On the details of the API calls. I was thinking of making the snapshot a true @property, with getter and setter. The setter would do the whole initialization. I don't know whether we need a separate start() function, or whether next_frame() should start any sort of coprocess the first time it runs. Of course, with our current engines, if initialization happens in the snapshot setter, then start(), like stop(), doesn't do anything. So we can wait on that decision until we're actually adding support for Gromacs.
As for the configuration and momentum properties: In general, is there a reason for a user to directly access the Configuration and Momentum classes, or does it make more sense to have the user interact with everything via Snapshots?
Best practices for Simulation subclasses
I agree that, at its heart, all the trajectories we might look at are made by combining 3 elements: a model, a current snapshot, and a method for changing the current snapshot to the next one according to the model. That is, a system of (stochastic) different equations, a state vector, and a method for integrating the diffEqs. And I think I should change the Toy Dynamics to reflect this better (I absorbed the system into the simulation object out of laziness: it works, but it isn't structured in the most pedagogically informative way, and the primary purpose of the Toy Dynamics is to be pedagogical.) This should be done, but not in this PR. I'll add a new issue for it, but don't know when I'll get around to fixing it.
It sounds like you're suggesting that we should make a similar pedagogical point in the OpenMMSimulation object. In principle, we combine the force field and the topology to generate a model; do we want a function OpenMMSimulation.create_model(force_field, topology)? Even though that's how I see the world, I'm not sure whether that's worth implementing. I think it might confuse users, and I think users might be better served by just learning how to set up OpenMM simulations from the OpenMM tutorial. If we change how the user interacts with setting up the OpenMM simulation, then a user who learns to use OpenMM in OpenTIS and then wants to use OpenMM in YANK will find different behavior.
Serializing Simulations
In general, the implementation of storage for each Simulation subclass will have to be up to the implementer. In the same way we use a PDB for the topology, the best thing for Gromacs might just be to read in the appropriate mdp file, which contains all the simulation parameters. Would jsonify handle the current Toy Dynamics easily (it's a pretty simple setup)? And did we already have a way to support openmm.app.Simulation objects? Since our OpenMMSimulation is just a subclass of that with new functions wrapped around it for our API/storage, anything that works on openmm.app.Simulation should work on OpenMMSimulation.
It seems this issue contains lots of major issues, but i will still keep this here since we have almost found a solution it seems. So here my answers:
Simulation and/or Simulator
Yes simulation should be just one class that handles the simulation hence the name
Our API for dynamics engines
I also thought about an iterator/generator (and learned that a generator is just some nice way to construct an iterator), but I haven't really made up my mind, if this can be useful.
To make this clearer let me extend this a little. This idea to able to start a separate job/thread somewhere that generates new samples is really great and so I think we should adapt this as the general way to treat a simulation. This is actually the same in siege tank:
Start : This will setup the simulation and start creating new samples (or not!)
Next/Step : Grab the next produced sample, if they are any present or will be (might involve downloading samples that are already ready like in Siegetank, reading a Gromacs output file and grab the next frame, or running OpenMM for some steps and returning the next snapshot). I think also it should always return an object unless the simluation is or has been stopped or has not been started yet. This new snapshot might (optionally) be lazy if the simulator allows for this
Stop : Tell the simulation to stop: No more samples will be generated, although the simulation could have decided to stop by itself, due to a crash, error or some stop criterion in the simulator.
This scheme is effectively an iterator and so we could implement it as such. The problem is that the simulation should not be the one to determine, when to stop the iteration and so in theory all the nice
iterator stuff would produce infinite trajectories or just make no sense. In a PathMover this would actually make sense. So, a mover could return a RunningTrajectory() iterator that is an iterator instead of a list and we could use this trajectory in other functions, although it is still being simulated. if we request e.g. frame 5 then it would wait until frame 5 is ready or return a snapshot object which is the one at that position, but has no coordinates. This gets really crazy. All this futures stuff is really cool and allows to run and write code independent of when the results are being generated. And so it can be easily run in parallel. One of the tricky parts is do deal with cases where we pronised a result, but the simulation dies and so there will never be a result. So, I think we don't really need this, although it is very cool.
I would keep
start
next_frame or just step
stop
this way we have the most general case and it is still simple. And if we add Gromacs and/or siege tank to it than we do not have to change anything.
Of course you are right, for toy_dynamics and openmm this scheme is not so useful , since the stop is unnecessary, I would not make it a property, although that would work with openmm and your integrators, but for the gromacs idea we should explicitly have a function to start the simulation process (which is potentially more costly than running. We should have snapshot as a property only with an getter which prevents errors. I would add configuration/momentum only for convenience and have a more direct access to the simulation. Might not be necessary and could be added easily later.
The last thing was the question in #73 if we sould
Best practices
You make an important point that users should get to know openmm. That is definitely a good point. My idea of this deperation comes from siegetank where you can have several streams per project. In our terms this would be several simulations per project. Each simulation could be an openmm instance on a cluster or a stream in siegetank or a separate thread running one of your integrators. In Siegetank all of these should have the same model and dynamics. But I agree that these do not have to be wrapped as I suggested previously. I made a simple example, how I thought about this:
# A factory for simulations in siegetank
sf = SiegetankSimulationFramework()
storage.topology = sf.topology
simulator0 = sf.create_simulator('sim00')
simulator1 = sf.create_simulator('sim01')
simulator0.start(initial_snapshot)
simulator1.start(initial_snapshot)
pm0 = ForwardShootMover(selector=UniformSelector, ensemble=Length10, simulator=simulator0)
pm1 = ForwardShootMover(selector=UniformSelector, ensemble=Length20, simulator=simulator1)
new_sample_trajectory1 = pm0.move(Trajectory([snapshot0])))
new_sample_trajectory2 = pm1.move(Trajectory([snapshot0])))
# Save trajectories from different simulations in the same storage!
storage.sample.save(new_sample_trajectory1)
storage.sample.save(new_sample_trajectory2)
If we do want to use my idea then we would use the Model/Dynamics wrapper to create the factory like
So a simulation only needs to provide a topology which is the part needed to set up the storage and everything that is necessary for later structural analysis. (Of course the integrator, temperature, etc. are also important). I agree that we should leave this idea of Model/Dynamics only to ToyModels. We can always add these as factories to make things look more structured later.
Serializing Simulations
I think most of the stuff in a simulator will automatically be stored. It is a little tricky though since it just takes all stuff in __dict__ and turns it into a JSON string. In your toy system case it would depend on when the object is stored, esp for variables like _mass etc. But we can write a special safer that takes care of this for simulations.
So I want to throw in my two cents: John has some pretty general MCMC code
in our repex and openmmtools repositories. Might be worth looking at for
inspiration, so that we don't independently solve the same problems twice.
It seems this issue contains lots of major issues, but i will still keep
this here since we have almost found a solution it seems. So here my
answers:
Simulation and/or Simulator
Yes simulation should be just one class that handles the simulation hence
the name
Our API for dynamics engines
I also thought about an iterator/generator (and learned that a generator
is just some nice way to construct an iterator), but I haven't really made
up my mind, if this can be useful.
To make this clearer let me extend this a little. This idea to able to
start a separate job/thread somewhere that generates new samples is really
great and so I think we should adapt this as the general way to treat a
simulation. This is actually the same in siege tank:
Start : This will setup the simulation and start creating new
samples (or not!)
Next/Step : Grab the next produced sample, if they are any present
or will be (might involve downloading samples that are already ready like
in Siegetank, reading a Gromacs output file and grab the next frame, or
running OpenMM for some steps and returning the next snapshot). I think
also it should always return an object unless the simluation is or has been
stopped or has not been started yet. This new snapshot might (optionally)
be lazy if the simulator allows for this
Stop : Tell the simulation to stop: No more samples will be
generated, although the simulation could have decided to stop by itself,
due to a crash, error or some stop criterion in the simulator.
This scheme is effectively an iterator and so we could implement it as
such. The problem is that the simulation should not be the one to
determine, when to stop the iteration and so in theory all the nice
iterator stuff would produce infinite trajectories or just make no sense.
In a PathMover this would actually make sense. So, a mover could return a
RunningTrajectory() iterator that is an iterator instead of a list and we
could use this trajectory in other functions, although it is still being
simulated. if we request e.g. frame 5 then it would wait until frame 5 is
ready or return a snapshot object which is the one at that position, but
has no coordinates. This gets really crazy. All this futures stuff is
really cool and allows to run and write code independent of when the
results are being generated. And so it can be easily run in parallel. One
of the tricky parts is do deal with cases where we pronised a result, but
the simulation dies and so there will never be a result. So, I think we
don't really need this, although it is very cool.
I would keep
start
next_frame or just step
stop
this way we have the most general case and it is still simple. And if we
add Gromacs and/or siege tank to it than we do not have to change anything.
Of course you are right, for toy_dynamics and openmm this scheme is not so
useful , since the stop is unnecessary, I would not make it a property,
although that would work with openmm and your integrators, but for the
gromacs idea we should explicitly have a function to start the simulation
process (which is potentially more costly than running. We should have
snapshot as a property only with an getter which prevents errors. I would
add configuration/momentum only for convenience and have a more direct
access to the simulation. Might not be necessary and could be added easily
later.
The last thing was the question in #73 #73 if we sould
Best practices
You make an important point that users should get to know openmm. That is
definitely a good point. My idea of this deperation comes from siegetank
where you can have several streams per project. In our terms this would be
several simulations per project. Each simulation could be an openmm
instance on a cluster or a stream in siegetank or a separate thread running
one of your integrators. In Siegetank all of these should have the same
model and dynamics. But I agree that these do not have to be wrapped as I
suggested previously. I made a simple example, how I thought about this:
So a simulation only needs to provide a topology which is the part needed
to set up the storage and everything that is necessary for later structural
analysis. (Of course the integrator, temperature, etc. are also important).
I agree that we should leave this idea of Model/Dynamics only to ToyModels.
We can always add these as factories to make things look more structured
later.
Serializing Simulations
I think most of the stuff in a simulator will automatically be stored. It
is a little tricky though since it just takes all stuff in dict and
turns it into a JSON string. In your toy system case it would depend on
when the object is stored, esp for variables like _mass etc. But we can
write a special safer that takes care of this for simulations.
—
Reply to this email directly or view it on GitHub #65 (comment).
You are totally right, this is why we kept the discussion so far to mainly conceptual questions, which often help me to really understand the parts involved. Right now we only use a very simple Alanine simulation that came from an openmm tutorial and David implemented some integrators with potentials that he has used before in his C-codes to test things and where he knows, what should happen when using TIS. To get more models we should definitely look at openmmtools. I admit that it becomes difficult to keep an overview so thanks for reminding us!
One last question on this subtopic: we've got one object; what's the best name for it?
I feel like maybe DynamicsEngine is a better description for what this has become. The naming history is that OpenMM has a Simulation object, and we needed an object to control that, which got named Simulator. However, this has evolved into an object which acts as the interface between our code and whatever underlying dynamics code the user wants to use to generate trajectories.
This object generates one trajectory at a time. For OpenMM, a single trajectory might be your "MD simulation." For TIS, you'd never say that a single trajectory is your "TIS simulation." So I'm concerned that the name Simulator might be confused with the object we're currently calling a Calculaton.
Our API for dynamics engines
The idea with combining start and next_frame would be that the Simulator object knows whether it is managing a process yet, and can start it if it is not managing one. However, I agree that code is clearer with an explicit start. It sounds like we're pretty close to agreed on an API similar to this:
(@property) current_snapshot : getter returns the current state snapshot; setter sets the snapshot and makes sure everything is prepared to start a trajectory (e.g., for Gromacs this might involve writing the Gromacs input files).
start(snapshot=None) : Starts the trajectory. If snapshot is given, the first step (called in super) is current_snapshot = snapshot. If not, it begins the run or error if no snapshot is loaded. Mainly overridden if the engine is external.
next_frame() : generates the next saved frame of the trajectory, and returns that frame as a snapshot
stop(trajectory) : stops the simulation (mainly for external engines) and cleans up. The argument trajectory is the final trajectory seen by our code (this might be useful during the clean-up stages for the external program).
For now I'm also implementing @propertys called momentum and configuration, although we may consider removing those later since I think they serve no purpose beyond what we already get with snapshot, and having two ways of accessing the same information will probably just confuse people.
Best practices
For now, I think I'll just implement the API without adding the 4-object model on OpenMM stuff, but we can revisit the 4-object model for OpenMM later. If we want to add it, it should be easy to add.
You're also addressing the question of how to run several trajectories (namely, shooting moves) in parallel, and that's worth some discussion. Let's move that to another issue, shall we? This one is crowded.
Serializing
Do you think we should add to the API the rule that all information that the user wants stored from a Simulator needs to be in the Simulator's __dict__? Basically, the user has to flatten any information that needs to be stored into the Simulator object. That way we don't have to try to predict what the user might want. This is done in the OpenMMSimulation, but not yet in the ToySimulation. That would be another thing to do when I implement #73.
Existing stuff in repex/openmmtools
Just glancing at it, the MCMC moves in repex are quite similar in layout to what we've developed as PathMoves and GlobalMoves. Wish I'd seen that earlier! I think this point is about our Calculation object discussion (some of which is also in #24), right @kyleabeauchamp? In terms of importing repex directly, I think there's a problem of domain: John's methods appear to me to be moves in configuration space. Our moves are in path space. The other possible issue is whether the approach in repex is actually restricted to MCMC sampling: several of our algorithm are not Markov chains. We're really doing general moves in path space; sometimes under ensembles which are evolving as a part of the move (e.g., AMS).
Thanks for pointing me to openmmtools: I was thinking of eventually implementing OpenMM version of some of models I actually play with in my main research; openmmtools gives examples that should make it very easy for me to learn all I need in order to do that at some actually do that!
However, I should clarify the differences with ToyDynamics. First and foremost, the point of ToyDynamics was to give a proof-of-concept that we had completely disentangled the path sampling layer from the MD engine layer. So it had to be completely free of OpenMM code, since that was the previously-entangled MD engine! Also, the ToyDynamics is a landscape-based approach to dynamics, whereas OpenMM is interaction-based. I don't think that even the simple 2D models I've shown are separable into pairwise interactions. (I suppose you could describe them as a single 3-body interaction, although it's not clear to me whether OpenMM supports generic 3-body interactions.)
I like the term, this is general enough and far enough from simulation. Actually, if, then I would have preferred Simulator, because the Simulation is more like Trajectory/Sample. So let's use
´DynamicsEngine`
API for DynamicsEngine()
Also agreed. The property setter might be a good compromise, although the setting without calling start might cause problems for some engines (esp external), so there is room for error. If we want the setting and the starting to be separate with might want to use engine.init(snapshot) instead of engine.current_snapshot be settable. I think ir should be something that can be read/requested but it's generation is up to the engine
Best practices
Agreed, I will move the multi-simulation discussion to a new issue
Serializing
No, the other way round. I wanted to say that the fall-back serializer looks at the object.dict and serializes it (while also taking care of objects that are itself storable and simtk.units). Exluded is only the id, which makes no sense, since this is the table row in the netcdf file. This is used for some objects like Ensembles, etc. This allows for certain flexibility in the objects internal structure. For other classes this is not necessary and overkill like configurations, momenta, etc. For these classes the JSON is not used and they have their own storage functions. Also, if we want the object to be immutable then all the "temporary" member variables should be separate and not stored, but so far there is no distinction, like we are using in storing the openmm simulator before.
I think we should use the JSON, but add a class property that can specify which elements in __dict__ should be stored and which are to be excluded.
The second question is about loading these objects. Right now the idea was only a way to store these object in a simple way to provide enough information to reconstruct them or make sense out of them. Not to pickle the exact object. More like a ´repr´ output. Importing necessary modules, etc is not the responsibility of the storage. In my implementation the JSON load constructs a LoadedObject() instance which is filled with the JSON content or the JSON load can be applied to an existing object and will overwrite the existing member with the loaded ones, like
# recreates an object with member functions
my_engine = OpenMMEngine()
storage.engine.load(1, engine_my) # restores inside
# recreates an object without member functions
my_engine = storage.engine.load(1) # creates an LoadedObject()
If it is not using JSON then it always creates the full object. Since a simulation needs to do a lot of initialization we could also do something like
# Create a new engine
my_engine = OpenMMEngine(parameters)
my_engine.initialize()
# Load previously saved engine
my_engine = store.engine.load(1, OpenMMEngine())
my_engine.initialize()
which I think would be the best solution, load/__init__ only creates the object and initialize does all the stuff the make it ready to be used
Existing stuff in repex/openmmtools
My apologies. I understood the comment wrong. I wasn't aware that you were also referring to the ReplicaExchange functionality. Sure, that would have been useful. I guess I am still stuck at single replica schemes :) This is true. We should have a look at that and if not using the code it still might make sense to try to keep the way it is used the same, so that someone familiar with repex can directly use the same for pathsampling. Thanks Kyle!
Activity
jhprinz commentedon Nov 11, 2014
This is a very good important point to address. I also think this is right now very confusing. I guess even to us. Let me first lay out the differentiation between Simulation and Simulator, which is what I originally intended and not what it necessarily has to be:
Simulation(): This is the OpenMM simulation objectSimulator(): This was supposed to be the wrapper for theSimulationthat contains information about the system to be simulated and would correspond more to OpenMM Simulation class. This is the part that generates snapshots, momenta, configurations etc, knows the topology (and thus the number of atoms which is necessary for the storage. So far it is also attached to the storage, but that should not be necessary.So, in theory all we would need is a
Simulator()object that knows all about the simulation in an abstract way. It knows how to run a simulation (including setting stopping conditions, etc.), it can give Snapshots/Trajectories.Maybe it just comes down to almost copying the OpenMM scheme.I think there are in general several steps involved to specify how a simulation works. Usually we want:
Does this make sense?
A problem that I encountered with OpenMM is that their
Systemobject does contain the energy function but no reference to the original topology. So this is kind of missing especially for analytical purposes although not required to do a simulation. Also, a topology is something very special used in Force-Field based energy functions, but we might not need this for simple toy models.So I think we should have ideally 4 objects:
Model()orSystem()that contains a describtion of the energy function and stuff that was used to generate it. For modelcules the topology as well as used forcefields, etc... For 2D toy models may only the functional form ofV(x). In a mathematical sense this know the dimensionality of x(t) in our differential equation and a way to compute/approximateV(x)Dynamics()orIntegrator()that (again) contains the actual way of doing the integration (e.g. VVVRIntegrator) as well as the information used to generate (Langevin/Brownian, temperature, physical ensemble, thermostat, barostat). In mathematical language this know the form of the differential equation as well as a numerical way to compute steps.Snapshot()orState()that contains the coordinates and velocities of system. In mathematical terms this is a particular instance ofxorx(t)if part of a trajectory. This object consists about a Configuration and a Momentum object.Simulator()orSimulation()the overall object that knows all of the above and can generate newSnapshots()by using the integratorQuestions:
Snapshot()should contain a reference to the Model used!Is this overkill? I imagine then in the final code something like
I guess I stop here to keep the discussion running...
jhprinz commentedon Nov 11, 2014
Maybe an example for the toy_simulation
dwhswenson commentedon Nov 11, 2014
First the tl;dr summary:
I think you're suggesting something very similar to what I mean, although I think you're focusing on a few aspects that I don't think we need to worry about. Different simulation objects don't need to have the same internal structure, as long as they implement the methods we need to interact with them.
Let's talk about two different layers to our code: the path sampling layer, which just functions on trajectories, and the dynamics engine layer, which generates new frames of a trajectory.
I agree that all engines can boil down to your model/dynamics/snapshot model (that's how I implemented both my grad-school C and post-doc C++ codes, although I used the names potential/integrator/system, resp. -- even a file reader can be an "integrator"/"dynamics").
However, I think it might be cleaner and more straightforward for anyone who wants to add a new engine to ask them to implement a few simple methods than to ask them to structure their objects in a certain way, especially since we _need_ the methods to have the right behavior, and _all_ we need is for the methods to have the right behavior.
Now the long version:
Path sampling-like approaches function on trajectories (or, more correctly, ordered [usually time-ordered] series of vectors; because of our MD bias, we call those ordered series "trajectories"). As long as the path sampling layer sees something it can identify as a trajectory, it doesn't care how that trajectory was generated.
Then there's the dynamics engine. This does all the dirty work of generating the trajectories, frame by frame. The internal structure of the dynamics engine can vary widely: as you mention, the topology/force-field split makes sense in molecular mechanics, where you decompose interactions into N-body (for fixed maximum N) interactions. However, for ab initio MD, the "topology" is just the isotope identification, and for simple models where you explicitly specify the entire PES, it's just the masses. For a generic system of differential equations, there's no distinction between topology and force field (since the purpose of separating them is to make it easy to define a generic system of differential equations.)
My suggestion is that Model/Integrator/System all get wrapped up into the
Simulationobject. How thatSimulationobject implements those internally can be different for each engine. However, each such object needs to implement an interface consisting of the following functions (we can adjust these, but at least this much behavior needs to available):init_simulation_with_snapshot(snapshot): takessnapshot, a state of the simulation, and initializes the simulation object in that stateload_snapshot(snaphost): loads the current simulation state intosnapshotgenerate_next_frame(): generates the next frame according to the simulation, and returns it as aSnapshotobjectstop(): tells the simulation to stop runningI also implemented
load_configurationandload_velocityto make things easier with the other code we already had, but those four are the core. I'm also fully in favor of switching the first two to a@propertycalledcurrent_snapshot; that requires more refactoring than I wanted to do in #48 without clearing it with you, but I think it is a better solution. (Above, you call itsnapshot; clarifying "current" might be worth the extra characters.)So why those four core functions? Well, you know how they work for the
OpenMMSimulationorToySimulationexamples. So let's talk about aGromacsSimulationobject. It would be pointless to launch an entire new Gromacs run every saved frame, so the best approach would be to have Gromacs running as a sort of coprocess. The first two functions are, of course, pretty trivial (although initializing the simulation involves writing a file and loading some additional file-related info into theGromacsSimulation). But the last two functions are interesting.For the
GromacsSimulation, the functiongenerate_next_frame()doesn't actually generate the frame. Instead, it waits until a new frame is written to the trr file, and then it reads that frame from the file. And then there's thestop()function; when we generated one frame at a time,stop()didn't do anything. However, theGromacsSimulationdoesn't know to stop until we tell it. Thestop()function here will stop the Gromacs process, and do whatever cleanup is needed.(For what it's worth, I think Tsjerk Wassenaar has already implemented most of the stuff we need for the Gromacs
get_next_frame()andstop()functions as part of thetipsi.pyscript he's writing to implement a very user-friendly TIS/TPS interface for Gromacs.)Note that I've replaced your
run(nsteps)function with aget_next_frame()function. That's all we really need. We leave it to the implementer to define a number of steps between frames, because the path sampling layer doesn't need to know about that.All of this was my thinking in #48 and #52, although I guess I didn't discuss it in detail anywhere.
Another reason not to define the model and dynamics as separate objects is that they are not necessarily independent of each other. In particular, they're not likely to be compatible between different simulation engines. For example, it is meaningless to use an OpenMM-defined model on my LangevinBAOABIntegrator dynamics object: the interfaces for OpenMM's force fields are unknown to that integrator. Enforcing the extra structure is pretty on paper, but can't be used in practice.
In general, I think that we should use Python's duck-typing to enhance our flexibility instead of trying to create super-general structures that work with everything. I say this because I spent a lot of time trying to do the latter: that approach led to messes in my C code that moved me to C++, and then to messes in my C++ code that made me eager to move to Python. (I'm sure I'll eventually screw up in Python, but at least I learned what to avoid from C and C++!)
I think the answer to your questions 1 and 2 are "we don't need to dictate the internal structure of a
Simulationobject as long as it correctly implements the interface", and for question 3, I think the trick is to make it hard (impossible?) to keep more than oneSimulationin a single storage. If the (unique) simulation parameters are stored in the storage, then no problem. And as for the user mixing snapshots from twoSimulationobjects by incorrect instantiation and loading... well, I'm willing with label that a 10T error/nofix.Finally, a summary of what initial intent of this issue was:
The above info is pretty much how I implemented things in #48 and #52. I intended this issue to propose a change to combine the
SimulatorandSimulationobjects. As we both mentioned above, the separation between the two was originally because theSimulationwas theopenmm.app.Simulationobject. TheSimulatorwas our own. However, as of #48, theSimulationis our own, but a class that inherits fromopenmm.app.Simulation. What I'm suggesting here is essentially to have allSimulationclasses inherit from something that is currently ourSimulatorclass. Since theSimulatoris independent of the dynamics engine, this should be possible. Then we only have one class instead of two.jhprinz commentedon Nov 12, 2014
Great. I like these detailed discussions a lot.
Okay, so we agree that my attempt seems nice and ordered, but might be too overkill and unnecessarily restrictive. I just thought that since OpenMM is using this and we were mainly aiming at OpenMM it seemed natural to make this into a clear seperation of the necessary parts.
I would disagree that the seperation into dynamics and system etc need to be interchangeable between different simulations. Maybe I wasn't clear on that. Right now we are effectively using some seperation like this, at least ToyDynamics has a seperate Integrator and so has OpenMM. Of course each integrator does only make sense in the simulation context. So you are right that we do not need to seperate these.
One point I thought about, when doing this was the way to store this and keep simulations reusable with different parameters. In principle we are saying the same thing, just different parts. You described what simulation should do and I was trying to explain how to set the parameters of simulation in a structured way. But back to the original problem: Restructuring of
SimulationI would do it as you said and have only a Simulation object that allows to run a simulation (whatever that is) that generates snapshots in some way. I really like the idea of an asynchronous way of handling, with just starting and then getting elements as they become available. If we already have this in place why not keep this really asynchronous and add callbacks once a new frame is available. This way we could even run stuff on external clusters and continue once a new frame is available, it would even allow running several processes in parallel. But this might again be a little to much for now.
So, I will sum up Davids idea of functions for Simulation()
init_simulation_with_snapshot(snapshot): takes snapshot, a state of the simulation, and initializes the simulation object in that stateload_snapshot(snaphost): loads the current simulation state into snapshotgenerate_next_frame(): generates the next frame according to the simulation, and returns it as a Snapshot objectstop(): tells the simulation to stop runningI think I would simplify the names of the functions and treat this like a generator so. Also, if we think of this as a worker then we would just start (which does all the initialization and (in case of Grimaces) starts the simulation. Have a function to generate the next snapshot
start(snapshot): takes snapshot, a state of the simulation, and initializes the simulation object in that state and start the simulation for an independent runner or waits (like OpenMM)next(): generates the next frame according to the simulation, and returns it as a Snapshot objectstop(): tells the simulation to stop running. For gromacs this will stop the background process and remove the left-over simulation objects.for convenience I would add properties, which will change after each call of
next()snapshot: return the currentSnapshot()objectconfiguration: return the currentConfiguration()objectmomentum: return the currentMomentum()objectIn the future we could also introduce so called
futuresand makenext()non-blocking to make next() return a snapshot object that is empty, but the rest of the code will continue until the real content of the Snapshot is needed. This is similar to the way I implemented the loading of Snapshots. This does not work with OpenMM, but could work with Gromacs and means, we can run severalsimulation.next()commands that will all run a Gromacs simulation and once we do something with these the python program waits. To allow for multiple simulators we just need to be able to create severalSimulation()objects. We could even have remote Simulations using Siegetank. How cool would that be.Again, David, this idea with implementing a next function is fantastic. So I think this is settled.
Quickly back to my idea with a Dynamics/System separation. Since it is really independent of this PR we might not do this at this point. Although I think we could do this for our OpenMM and ToyModel Approach so encourage a more structured code. Kind of as a good practice rather then a must. The main object is just the simulator, however this is instantiated.
Leaves a last thing to think about. How to serialize the
Simulation()object? We can use the jsonify algorithm which will turn a simulation into a json string as long as there are no complex objects in it. Right now we have only strings and some number with units. I will add the capability to the json-netcdf serializer to also accept things with a unit. This means in theory we can also just store simulation objects.So, enough :) David if you want to start implementing a new Simulation() go ahead.
jhprinz commentedon Nov 12, 2014
I think I move the discussion about Calculation to a new issue...
dwhswenson commentedon Nov 13, 2014
A few related threads going here. Let me try to first describe how I see the difference/connection between the two major ones:
On one hand, we're talking about the possibility of combining what was previously two objects, a
Simulatorand aSimulation, into one object by having the generic behavior of aSimulatorsubclassed by individualSimulations for each dynamics engine. It sounds like we're agreed that this is a good idea. On this part, there's still some discussion to have about the specifics of the API we establish (which has to implemented separately for each dynamics engine).On the other hand, we're talking about the specifics of the implementation of the engine-specific
Simulationsubclasses. On that, I think we're also in agreement that we shouldn't enforce a specific internal approach for these (since they can't trade information between them); however, we should demonstrate best practices in our own implementations. Here, I still have some questions about what qualifies as best practices in some specific situations.So a few things about each of these, and then a few other little threads around the edges of this....
Our API for dynamics engines
When you say, "treat this like a generator", do you mean that we might actually make this into a proper generator? I haven't written anything using generators before. This seems like a perfect example for the use of
yield, but I'm not sure whether we actually get much benefit from that, since we can't actually do anything else with our trajectory (identify ensemble) until the whole trajectory has been returned.If we don't make it a proper generator, I think we might be better using a more explicit name than just
next(). I thinksimulation.next()isn't quite as clear assimulation.next_frame()(justnextcould as easily refer to getting the entire next trajectory, instead of just the next frame).The idea having a non-blocking
next()is really interesting and definitely appealing when thinking about things like waiting for a coprocess to write a new frame. However, I'm not sure how we can use that idea to our advantage. Isn't the algorithm effectively blocked until we get that frame and can analyze it to determine whether we continue the trajectory?On the details of the API calls. I was thinking of making the
snapshota true@property, with getter and setter. The setter would do the whole initialization. I don't know whether we need a separatestart()function, or whethernext_frame()should start any sort of coprocess the first time it runs. Of course, with our current engines, if initialization happens in the snapshot setter, thenstart(), likestop(), doesn't do anything. So we can wait on that decision until we're actually adding support for Gromacs.As for the
configurationandmomentumproperties: In general, is there a reason for a user to directly access theConfigurationandMomentumclasses, or does it make more sense to have the user interact with everything viaSnapshots?Best practices for
SimulationsubclassesI agree that, at its heart, all the trajectories we might look at are made by combining 3 elements: a model, a current snapshot, and a method for changing the current snapshot to the next one according to the model. That is, a system of (stochastic) different equations, a state vector, and a method for integrating the diffEqs. And I think I should change the Toy Dynamics to reflect this better (I absorbed the system into the simulation object out of laziness: it works, but it isn't structured in the most pedagogically informative way, and the primary purpose of the Toy Dynamics is to be pedagogical.) This should be done, but not in this PR. I'll add a new issue for it, but don't know when I'll get around to fixing it.
It sounds like you're suggesting that we should make a similar pedagogical point in the
OpenMMSimulationobject. In principle, we combine the force field and the topology to generate a model; do we want a functionOpenMMSimulation.create_model(force_field, topology)? Even though that's how I see the world, I'm not sure whether that's worth implementing. I think it might confuse users, and I think users might be better served by just learning how to set up OpenMM simulations from the OpenMM tutorial. If we change how the user interacts with setting up the OpenMM simulation, then a user who learns to use OpenMM in OpenTIS and then wants to use OpenMM in YANK will find different behavior.Serializing
SimulationsIn general, the implementation of storage for each
Simulationsubclass will have to be up to the implementer. In the same way we use a PDB for the topology, the best thing for Gromacs might just be to read in the appropriate mdp file, which contains all the simulation parameters. Would jsonify handle the current Toy Dynamics easily (it's a pretty simple setup)? And did we already have a way to supportopenmm.app.Simulationobjects? Since ourOpenMMSimulationis just a subclass of that with new functions wrapped around it for our API/storage, anything that works onopenmm.app.Simulationshould work onOpenMMSimulation.jhprinz commentedon Nov 14, 2014
It seems this issue contains lots of major issues, but i will still keep this here since we have almost found a solution it seems. So here my answers:
Simulation and/or Simulator
Yes simulation should be just one class that handles the simulation hence the name
Our API for dynamics engines
I also thought about an iterator/generator (and learned that a generator is just some nice way to construct an iterator), but I haven't really made up my mind, if this can be useful.
To make this clearer let me extend this a little. This idea to able to start a separate job/thread somewhere that generates new samples is really great and so I think we should adapt this as the general way to treat a simulation. This is actually the same in siege tank:
Start: This will setup the simulation and start creating new samples (or not!)Next/Step: Grab the next produced sample, if they are any present or will be (might involve downloading samples that are already ready like in Siegetank, reading a Gromacs output file and grab the next frame, or running OpenMM for some steps and returning the next snapshot). I think also it should always return an object unless the simluation is or has been stopped or has not been started yet. This new snapshot might (optionally) be lazy if the simulator allows for thisStop: Tell the simulation to stop: No more samples will be generated, although the simulation could have decided to stop by itself, due to a crash, error or some stop criterion in the simulator.This scheme is effectively an iterator and so we could implement it as such. The problem is that the simulation should not be the one to determine, when to stop the iteration and so in theory all the nice
iterator stuff would produce infinite trajectories or just make no sense. In a PathMover this would actually make sense. So, a mover could return a RunningTrajectory() iterator that is an iterator instead of a list and we could use this trajectory in other functions, although it is still being simulated. if we request e.g. frame 5 then it would wait until frame 5 is ready or return a snapshot object which is the one at that position, but has no coordinates. This gets really crazy. All this futures stuff is really cool and allows to run and write code independent of when the results are being generated. And so it can be easily run in parallel. One of the tricky parts is do deal with cases where we pronised a result, but the simulation dies and so there will never be a result. So, I think we don't really need this, although it is very cool.
I would keep
startnext_frameor juststepstopthis way we have the most general case and it is still simple. And if we add Gromacs and/or siege tank to it than we do not have to change anything.
Of course you are right, for toy_dynamics and openmm this scheme is not so useful , since the
stopis unnecessary, I would not make it a property, although that would work with openmm and your integrators, but for the gromacs idea we should explicitly have a function to start the simulation process (which is potentially more costly than running. We should have snapshot as a property only with an getter which prevents errors. I would add configuration/momentum only for convenience and have a more direct access to the simulation. Might not be necessary and could be added easily later.The last thing was the question in #73 if we sould
Best practices
You make an important point that users should get to know openmm. That is definitely a good point. My idea of this deperation comes from siegetank where you can have several streams per project. In our terms this would be several simulations per project. Each simulation could be an openmm instance on a cluster or a stream in siegetank or a separate thread running one of your integrators. In Siegetank all of these should have the same model and dynamics. But I agree that these do not have to be wrapped as I suggested previously. I made a simple example, how I thought about this:
If we do want to use my idea then we would use the Model/Dynamics wrapper to create the factory like
So a simulation only needs to provide a topology which is the part needed to set up the storage and everything that is necessary for later structural analysis. (Of course the integrator, temperature, etc. are also important). I agree that we should leave this idea of Model/Dynamics only to ToyModels. We can always add these as factories to make things look more structured later.
Serializing Simulations
I think most of the stuff in a simulator will automatically be stored. It is a little tricky though since it just takes all stuff in
__dict__and turns it into a JSON string. In your toy system case it would depend on when the object is stored, esp for variables like_massetc. But we can write a special safer that takes care of this for simulations.kyleabeauchamp commentedon Nov 14, 2014
So I want to throw in my two cents: John has some pretty general MCMC code
in our repex and openmmtools repositories. Might be worth looking at for
inspiration, so that we don't independently solve the same problems twice.
On Thu, Nov 13, 2014 at 10:19 PM, Jan-Hendrik Prinz <
notifications@github.com> wrote:
jhprinz commentedon Nov 14, 2014
You are totally right, this is why we kept the discussion so far to mainly conceptual questions, which often help me to really understand the parts involved. Right now we only use a very simple Alanine simulation that came from an openmm tutorial and David implemented some integrators with potentials that he has used before in his C-codes to test things and where he knows, what should happen when using TIS. To get more models we should definitely look at openmmtools. I admit that it becomes difficult to keep an overview so thanks for reminding us!
dwhswenson commentedon Nov 14, 2014
Simulation/Simulator
One last question on this subtopic: we've got one object; what's the best name for it?
I feel like maybe
DynamicsEngineis a better description for what this has become. The naming history is that OpenMM has aSimulationobject, and we needed an object to control that, which got namedSimulator. However, this has evolved into an object which acts as the interface between our code and whatever underlying dynamics code the user wants to use to generate trajectories.This object generates one trajectory at a time. For OpenMM, a single trajectory might be your "MD simulation." For TIS, you'd never say that a single trajectory is your "TIS simulation." So I'm concerned that the name
Simulatormight be confused with the object we're currently calling aCalculaton.Our API for dynamics engines
The idea with combining
startandnext_framewould be that theSimulatorobject knows whether it is managing a process yet, and can start it if it is not managing one. However, I agree that code is clearer with an explicitstart. It sounds like we're pretty close to agreed on an API similar to this:@property)current_snapshot: getter returns the current state snapshot; setter sets the snapshot and makes sure everything is prepared to start a trajectory (e.g., for Gromacs this might involve writing the Gromacs input files).start(snapshot=None): Starts the trajectory. Ifsnapshotis given, the first step (called in super) iscurrent_snapshot = snapshot. If not, it begins the run or error if no snapshot is loaded. Mainly overridden if the engine is external.next_frame(): generates the next saved frame of the trajectory, and returns that frame as a snapshotstop(trajectory): stops the simulation (mainly for external engines) and cleans up. The argumenttrajectoryis the final trajectory seen by our code (this might be useful during the clean-up stages for the external program).For now I'm also implementing
@propertys calledmomentumandconfiguration, although we may consider removing those later since I think they serve no purpose beyond what we already get withsnapshot, and having two ways of accessing the same information will probably just confuse people.Best practices
For now, I think I'll just implement the API without adding the 4-object model on OpenMM stuff, but we can revisit the 4-object model for OpenMM later. If we want to add it, it should be easy to add.
You're also addressing the question of how to run several trajectories (namely, shooting moves) in parallel, and that's worth some discussion. Let's move that to another issue, shall we? This one is crowded.
Serializing
Do you think we should add to the API the rule that all information that the user wants stored from a
Simulatorneeds to be in theSimulator's__dict__? Basically, the user has to flatten any information that needs to be stored into theSimulatorobject. That way we don't have to try to predict what the user might want. This is done in theOpenMMSimulation, but not yet in theToySimulation. That would be another thing to do when I implement #73.Existing stuff in repex/openmmtools
Just glancing at it, the MCMC moves in repex are quite similar in layout to what we've developed as PathMoves and GlobalMoves. Wish I'd seen that earlier! I think this point is about our
Calculationobject discussion (some of which is also in #24), right @kyleabeauchamp? In terms of importing repex directly, I think there's a problem of domain: John's methods appear to me to be moves in configuration space. Our moves are in path space. The other possible issue is whether the approach in repex is actually restricted to MCMC sampling: several of our algorithm are not Markov chains. We're really doing general moves in path space; sometimes under ensembles which are evolving as a part of the move (e.g., AMS).Thanks for pointing me to openmmtools: I was thinking of eventually implementing OpenMM version of some of models I actually play with in my main research; openmmtools gives examples that should make it very easy for me to learn all I need in order to do that at some actually do that!
However, I should clarify the differences with ToyDynamics. First and foremost, the point of ToyDynamics was to give a proof-of-concept that we had completely disentangled the path sampling layer from the MD engine layer. So it had to be completely free of OpenMM code, since that was the previously-entangled MD engine! Also, the ToyDynamics is a landscape-based approach to dynamics, whereas OpenMM is interaction-based. I don't think that even the simple 2D models I've shown are separable into pairwise interactions. (I suppose you could describe them as a single 3-body interaction, although it's not clear to me whether OpenMM supports generic 3-body interactions.)
kyleabeauchamp commentedon Nov 14, 2014
Sound reasonable, I'm not sugggesting that you actually reuse the stuff, just wanted to point it out for reference.
jhprinz commentedon Nov 14, 2014
Agreed, we are almost ready...
DynamicsEngine()I like the term, this is general enough and far enough from simulation. Actually, if, then I would have preferred Simulator, because the Simulation is more like Trajectory/Sample. So let's use
´DynamicsEngine`
API for
DynamicsEngine()Also agreed. The property setter might be a good compromise, although the setting without calling start might cause problems for some engines (esp external), so there is room for error. If we want the setting and the starting to be separate with might want to use
engine.init(snapshot)instead ofengine.current_snapshotbe settable. I think ir should be something that can be read/requested but it's generation is up to the engineBest practices
Agreed, I will move the multi-simulation discussion to a new issue
Serializing
No, the other way round. I wanted to say that the fall-back serializer looks at the object.dict and serializes it (while also taking care of objects that are itself storable and simtk.units). Exluded is only the id, which makes no sense, since this is the table row in the netcdf file. This is used for some objects like Ensembles, etc. This allows for certain flexibility in the objects internal structure. For other classes this is not necessary and overkill like configurations, momenta, etc. For these classes the JSON is not used and they have their own storage functions. Also, if we want the object to be immutable then all the "temporary" member variables should be separate and not stored, but so far there is no distinction, like we are using in storing the openmm simulator before.
I think we should use the JSON, but add a class property that can specify which elements in
__dict__should be stored and which are to be excluded.The second question is about loading these objects. Right now the idea was only a way to store these object in a simple way to provide enough information to reconstruct them or make sense out of them. Not to pickle the exact object. More like a ´repr´ output. Importing necessary modules, etc is not the responsibility of the storage. In my implementation the JSON load constructs a LoadedObject() instance which is filled with the JSON content or the JSON load can be applied to an existing object and will overwrite the existing member with the loaded ones, like
If it is not using JSON then it always creates the full object. Since a simulation needs to do a lot of initialization we could also do something like
which I think would be the best solution,
load/__init__only creates the object andinitializedoes all the stuff the make it ready to be usedExisting stuff in repex/openmmtools
My apologies. I understood the comment wrong. I wasn't aware that you were also referring to the ReplicaExchange functionality. Sure, that would have been useful. I guess I am still stuck at single replica schemes :) This is true. We should have a look at that and if not using the code it still might make sense to try to keep the way it is used the same, so that someone familiar with repex can directly use the same for pathsampling. Thanks Kyle!