Skip to content

Address confusion due to similarity of Simulation and Simulator names #65

Closed

Description

@dwhswenson

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?

Activity

  1. jhprinz commented on Nov 11, 2014

    @jhprinz
    Contributor

    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:

    1. Start with a topology that described our molecule in a phycal way
    2. 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.
    3. Add a configuration and a momentum which must be compatible. with the topology/model
    4. 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...
    5. 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:

    1. 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)
    2. 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.
    3. 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.
    4. Simulator() or Simulation() the overall object that knows all of the above and can generate new Snapshots() by using the integrator

    Questions:

    1. 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.
    2. 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
    3. Each Snapshot() should contain a reference to the Model used!

    Is this overkill? I imagine then in the final code something like

    simulation = OpenMMSimulation()
    simulation.model = OpenMMSimulation.ForceFieldModel(
        topology=...pdb, 
        nonbondedMethod=PME, 
        nonbondedCutoff=1*nanometer, 
        constraints=HBonds, 
        forcefield_solute="amber96.xml", 
        forcefield_solvent="tip3p.xml"
    )
    simulation.dynamics = OpenMMSimulation.VVVRIntegrator(
        temperature = 300.0 * kelvin, 
        collision_rate=1.0 / picoseconds, 
        timestep=2.0 * femtoseconds
    )
    simulation.snapshot = Snapshot(
        configuration=Configuration(0),
        velocities=model.randomvelocities(300.0 * kelvin)
    )
    simulation.run(1000)
    snapshot_new = simulation.snapshot
    momentum_new = simulation.momentum
    

    I guess I stop here to keep the discussion running...

  2. jhprinz commented on Nov 11, 2014

    @jhprinz
    Contributor

    Maybe an example for the toy_simulation

    simulation = ToySimulation()
    simulation.model = ToySimulation.1DPotentialModel(
        potential = TwoWellGaussPotential(parameters)
    )
    simulation.dynamics = ToySimulation.LangevinBAOABIntegrator(
        dt = 2.0 * femtoseconds,
        temperature = 300.0 * kelvin,
        gamma = 10.0 / picoseconds
    )
    simulation.snapshot = Snapshot(
        coordinates=np.array([[0,0,0]]),   # if we still assume 3 spatial coordinates although only x is used
        momentum=simulation.model.randomvelocities(300.0 * kelvin)
    )
    simulation.run(1000)
    
  3. dwhswenson commented on Nov 11, 2014

    @dwhswenson
    MemberAuthor

    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 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.

  4. jhprinz commented on Nov 12, 2014

    @jhprinz
    Contributor

    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.

  5. jhprinz commented on Nov 12, 2014

    @jhprinz
    Contributor

    I think I move the discussion about Calculation to a new issue...

  6. dwhswenson commented on Nov 13, 2014

    @dwhswenson
    MemberAuthor

    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.

  7. jhprinz commented on Nov 14, 2014

    @jhprinz
    Contributor

    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

    sf.model = OpenMM.FFModel(...)
    sf.dynamics = OpenMM.VerletIntergrator(...)
    

    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.

  8. kyleabeauchamp commented on Nov 14, 2014

    @kyleabeauchamp

    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:

    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:

    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

    sf.model = OpenMM.FFModel(...)
    sf.dynamics = OpenMM.VerletIntergrator(...)

    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).

  9. jhprinz commented on Nov 14, 2014

    @jhprinz
    Contributor

    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!

  10. dwhswenson commented on Nov 14, 2014

    @dwhswenson
    MemberAuthor

    Simulation/Simulator

    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.)

  11. kyleabeauchamp commented on Nov 14, 2014

    @kyleabeauchamp

    Sound reasonable, I'm not sugggesting that you actually reuse the stuff, just wanted to point it out for reference.

  12. jhprinz commented on Nov 14, 2014

    @jhprinz
    Contributor

    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 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!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions