def initialize(self): """ Initializes the associated storage to index momentums in it """ error_if_no_simtk_unit("KineticContainerStore") super(KineticContainerStore, self).initialize() self.create_variable( 'velocities', 'numpy.float32', dimensions=('n_atoms', 'n_spatial'), description="the velocity of atom 'atom' in dimension " + "'coordinate' of momentum 'momentum'.", chunksizes=('n_atoms', 'n_spatial'), simtk_unit=unit.nanometers / unit.picoseconds)
def initialize(self): super(StaticContainerStore, self).initialize() error_if_no_simtk_unit("StaticContainerStore") self.create_variable( 'coordinates', 'numpy.float32', dimensions=('n_atoms', 'n_spatial'), description="coordinate of atom '{ix[1]}' in dimension " + "'{ix[2]}' of configuration '{ix[0]}'.", chunksizes=('n_atoms', 'n_spatial'), simtk_unit=unit.nanometers) self.create_variable( 'box_vectors', 'numpy.float32', dimensions=('n_spatial', 'n_spatial'), chunksizes=('n_spatial', 'n_spatial'), simtk_unit=unit.nanometers)
def snapshot_from_testsystem(testsystem, simple_topology=False, periodic=True): """ Construct a Snapshot from openmm topology and state objects Parameters ---------- testsystem : openmmtools.Topology The filename of the .pdb file to be used simple_topology : bool if `True` only a simple topology with n_atoms will be created. This cannot be used with complex CVs but loads and stores very fast periodic : bool True (default) if system is periodic; if False, box vectors are None Returns ------- :class:`openpathsampling.engines.Snapshot` the constructed Snapshot """ error_if_no_simtk_unit("snapshot_from_testsystem") u_nm = unit.nanometers u_ps = unit.picoseconds velocities = unit.Quantity(np.zeros(testsystem.positions.shape), u_nm / u_ps) if simple_topology: topology = Topology(*testsystem.positions.shape) else: topology = MDTrajTopology(md.Topology.from_openmm(testsystem.topology)) if periodic: sys_box_vectors = testsystem.system.getDefaultPeriodicBoxVectors() box_vectors = np.array([v / u_nm for v in sys_box_vectors]) * u_nm else: box_vectors = None snapshot = Snapshot.construct( coordinates=testsystem.positions, box_vectors=box_vectors, velocities=velocities, engine=OpenMMToolsTestsystemEngine(topology, testsystem.name) ) return snapshot
def empty_snapshot_from_openmm_topology(topology, simple_topology=False): """ Return an empty snapshot from an openmm.Topology object Velocities will be set to zero. Parameters ---------- topology : openmm.Topology the topology representing the structure and number of atoms simple_topology : bool if `True` only a simple topology with n_atoms will be created. This cannot be used with complex CVs but loads and stores very fast Returns ------- openpathsampling.engines.Snapshot the complete snapshot with zero coordinates and velocities """ error_if_no_simtk_unit("empty_snapshot_from_openmm_topology") u_nm = unit.nanometers u_ps = unit.picoseconds n_atoms = topology.n_atoms if simple_topology: topology = Topology(n_atoms, 3) else: error_if_no_mdtraj("empty_snaphsot_from_openmm_topology") topology = MDTrajTopology(md.Topology.from_openmm(topology)) snapshot = Snapshot.construct( coordinates=unit.Quantity(np.zeros((n_atoms, 3)), u_nm), box_vectors=unit.Quantity(topology.setUnitCellDimensions(), u_nm), velocities=unit.Quantity(np.zeros((n_atoms, 3)), u_nm / u_ps), engine=TopologyEngine(topology) ) return snapshot
def trajectory_from_mdtraj(mdtrajectory, simple_topology=False, velocities=None): """ Construct a Trajectory object from an mdtraj.Trajectory object Parameters ---------- mdtrajectory : mdtraj.Trajectory Input mdtraj.Trajectory simple_topology : bool if `True` only a simple topology with n_atoms will be created. This cannot be used with complex CVs but loads and stores very fast velocities : np.array velocities in units of nm/ps Returns ------- openpathsampling.engines.Trajectory the constructed Trajectory instance """ error_if_no_simtk_unit("trajectory_from_mdtraj") trajectory = Trajectory() u_nm = unit.nanometer u_ps = unit.picosecond vel_unit = u_nm / u_ps if simple_topology: topology = Topology(*mdtrajectory.xyz[0].shape) else: topology = MDTrajTopology(mdtrajectory.topology) if velocities is None: empty_vel = unit.Quantity(np.zeros(mdtrajectory.xyz[0].shape), vel_unit) if mdtrajectory.unitcell_vectors is not None: box_vects = unit.Quantity(mdtrajectory.unitcell_vectors, unit.nanometers) else: box_vects = [None] * len(mdtrajectory) engine = TopologyEngine(topology) for frame_num in range(len(mdtrajectory)): # mdtraj trajectories only have coordinates and box_vectors coord = unit.Quantity(mdtrajectory.xyz[frame_num], u_nm) if velocities is not None: vel = unit.Quantity(velocities[frame_num], vel_unit) else: vel = empty_vel box_v = box_vects[frame_num] statics = Snapshot.StaticContainer(coordinates=coord, box_vectors=box_v, engine=engine) kinetics = Snapshot.KineticContainer(velocities=vel, engine=engine) snap = Snapshot(statics=statics, kinetics=kinetics, engine=engine) trajectory.append(snap) return trajectory