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

reportIntervalint

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

coordinatesbool

Whether to write the coordinates to the file.

timebool

Whether to write the current time to the file.

cellbool

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

potentialEnergybool

Whether to write the potential energy to the file.

kineticEnergybool

Whether to write the kinetic energy to the file.

temperaturebool

Whether to write the instantaneous temperature to the file.

velocitiesbool

Whether to write the velocities to the file.

atomSubsetarray_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(self)

Close the underlying trajectory file

describeNextReport(self, simulation)

Get information about the next report this object will generate.

report(self, simulation, state)

Generate a report.

__init__(self, file, reportInterval, coordinates=True, time=True, cell=True, potentialEnergy=True, kineticEnergy=True, temperature=True, velocities=False, atomSubset=None)

Create a HDF5Reporter.

Methods

__init__(self, file, reportInterval[, …])

Create a HDF5Reporter.

close(self)

Close the underlying trajectory file

describeNextReport(self, simulation)

Get information about the next report this object will generate.

report(self, simulation, state)

Generate a report.

Attributes

backend