mdtraj.reporters.HDF5Reporter

class mdtraj.reporters.HDF5Reporter(file, reportInterval, coordinates=True, time=True, cell=True, potentialEnergy=True, kineticEnergy=True, temperature=True, velocities=False, atomSubset=None)

HDF5Reporter stores a molecular dynamics trajectory in the HDF5 format.

This object supports saving all kinds of information from the simulation – more than any other trajectory format. In addition to all of the options, the topology of the system will also (of course) be stored in the file. All of the information is compressed, so the size of the file is not much different than DCD, despite the added flexibility.

Parameters:

file : str, or HDF5TrajectoryFile

Either an open HDF5TrajecoryFile object to write to, or a string specifying the filename of a new HDF5 file to save the trajectory to.

reportInterval : int

The interval (in time steps) at which to write frames.

coordinates : bool

Whether to write the coordinates to the file.

time : bool

Whether to write the current time to the file.

cell : bool

Whether to write the current unit cell dimensions to the file.

potentialEnergy : bool

Whether to write the potential energy to the file.

kineticEnergy : bool

Whether to write the kinetic energy to the file.

temperature : bool

Whether to write the instantaneous temperature to the file.

velocities : bool

Whether to write the velocities to the file.

atomSubset : array_like, default=None

Only write a subset of the atoms, with these (zero based) indices to the file. If None, all of the atoms will be written to disk.

Notes

If you use the atomSubset option to write only a subset of the atoms to disk, the kineticEnergy, potentialEnergy, and temperature fields will not change. They will still refer to the energy and temperature of the whole system, and are not “subsetted” to only include the energy of your subsystem.

Examples

>>> simulation = Simulation(topology, system, integrator)
>>> h5_reporter = HDF5Reporter('traj.h5', 100)
>>> simulation.reporters.append(h5_reporter)
>>> simulation.step(10000)
>>> traj = mdtraj.trajectory.load('traj.lh5')

Attributes

backend

Methods

close() Close the underlying trajectory file
describeNextReport(simulation) Get information about the next report this object will generate.
report(simulation, state) Generate a report.
__init__(file, reportInterval, coordinates=True, time=True, cell=True, potentialEnergy=True, kineticEnergy=True, temperature=True, velocities=False, atomSubset=None)

Create a HDF5Reporter.

Methods

__init__(file, reportInterval[, …]) Create a HDF5Reporter.
close() Close the underlying trajectory file
describeNextReport(simulation) Get information about the next report this object will generate.
report(simulation, state) Generate a report.
close()

Close the underlying trajectory file

describeNextReport(simulation)

Get information about the next report this object will generate.

Parameters:

simulation : simtk.openmm.app.Simulation

The Simulation to generate a report for

Returns:

report_description : tuple

A five element tuple. The first element is the number of steps until the next report. The remaining elements specify whether that report will require positions, velocities, forces, and energies respectively.

report(simulation, state)

Generate a report.

Parameters:

simulation : simtk.openmm.app.Simulation

The Simulation to generate a report for

state : simtk.openmm.State

The current state of the simulation