openmmtools.testsystems.DiatomicFluid

class openmmtools.testsystems.DiatomicFluid(nmolecules=250, K=Quantity(value=424.0, unit=kilocalorie/(angstrom**2*mole)), r0=Quantity(value=1.383, unit=angstrom), m1=Quantity(value=14.01, unit=dalton), m2=Quantity(value=14.01, unit=dalton), epsilon=Quantity(value=0.17, unit=kilocalorie/mole), sigma=Quantity(value=1.824, unit=angstrom), charge=Quantity(value=0.0, unit=elementary charge), reduced_density=0.05, switch_width=Quantity(value=0.5, unit=angstrom), cutoff=None, constraint=False, dispersion_correction=True, **kwargs)[source]

Create a diatomic fluid.

Parameters:
nmoleculesint, optional, default=250

Number of molecules.

Kopenmm.unit.Quantity, optional, default=290.1 * unit.kilocalories_per_mole / unit.angstrom**2

harmonic bond potential. default is GAFF c-c bond

r0openmm.unit.Quantity, optional, default=1.550 * unit.amu

bond length. Default is Amber GAFF c-c bond.

constraintbool, default=False

if True, the bond length will be constrained

m1openmm.unit.Quantity, optional, default=12.01 * unit.amu

particle1 mass

m2openmm.unit.Quantity, optional, default=12.01 * unit.amu

particle2 mass

epsilonopenmm.unit.Quantity, optional, default=0.1700 * unit.kilocalories_per_mole

particle Lennard-Jones well depth

sigmaopenmm.unit.Quantity, optional, default=1.8240 * unit.angstroms

particle Lennard-Jones sigma

chargeopenmm.unit.Quantity, optional, default=0.0 * unit.elementary_charge

charge to place on atomic centers to create a dipole

reduced_densityfloat, optional, default=0.05

Reduced density (density * sigma**3); default is appropriate for gas

cutoffopenmm.unit.Quantity, optional, default=None

if specified, the specified cutoff will be used; otherwise, 3.0 * sigma will be used

switch_widthopenmm.unit.Quantity with units compatible with angstroms, optional, default=0.2*unit.angstroms

switching function is turned on at cutoff - switch_width If None, no switch will be applied (e.g. hard cutoff).

dispersion_correctionbool, optional, default=True

if True, will use analytical dispersion correction (if not using switching function)

Notes

The natural period of a harmonic oscillator is T = sqrt(m/K), so you will want to use an integration timestep smaller than ~ T/10.

Examples

Create an uncharged Diatomic fluid.

>>> diatom = DiatomicFluid()
>>> system, positions = diatom.system, diatom.positions

Create a dipolar fluid.

>>> diatom = DiatomicFluid(charge=1.0*unit.elementary_charge)
>>> system, positions = diatom.system, diatom.positions

Create a Diatomic fluid with constraints instead of harmonic bonds

>>> diatom = DiatomicFluid(constraint=True)
>>> system, positions = diatom.system, diatom.positions

Specify a different system size.

>>> diatom = DiatomicFluid(constraint=True, nmolecules=200)
>>> system, positions = diatom.system, diatom.positions
Attributes:
analytical_properties

A list of available analytical properties, accessible via ‘get_propertyname(thermodynamic_state)’ calls.

mdtraj_topology

The mdtraj.Topology object corresponding to the test system (read-only).

name

The name of the test system.

positions

The openmm.unit.Quantity object containing the particle positions, with units compatible with openmm.unit.nanometers.

system

The openmm.System object corresponding to the test system.

topology

The openmm.app.Topology object corresponding to the test system.

Methods

get_potential_expectation(state)

Return the expectation of the potential energy, computed analytically or numerically.

reduced_potential_expectation(...)

Calculate the expected potential energy in state_sampled_from, divided by kB * T in state_evaluated_in.

serialize()

Return the System and positions in serialized XML form.

__init__(nmolecules=250, K=Quantity(value=424.0, unit=kilocalorie/(angstrom**2*mole)), r0=Quantity(value=1.383, unit=angstrom), m1=Quantity(value=14.01, unit=dalton), m2=Quantity(value=14.01, unit=dalton), epsilon=Quantity(value=0.17, unit=kilocalorie/mole), sigma=Quantity(value=1.824, unit=angstrom), charge=Quantity(value=0.0, unit=elementary charge), reduced_density=0.05, switch_width=Quantity(value=0.5, unit=angstrom), cutoff=None, constraint=False, dispersion_correction=True, **kwargs)[source]

Abstract base class for test system.

Parameters:

Methods

__init__([nmolecules, K, r0, m1, m2, ...])

Abstract base class for test system.

get_potential_expectation(state)

Return the expectation of the potential energy, computed analytically or numerically.

reduced_potential_expectation(...)

Calculate the expected potential energy in state_sampled_from, divided by kB * T in state_evaluated_in.

serialize()

Return the System and positions in serialized XML form.

Attributes

analytical_properties

A list of available analytical properties, accessible via 'get_propertyname(thermodynamic_state)' calls.

mdtraj_topology

The mdtraj.Topology object corresponding to the test system (read-only).

name

The name of the test system.

positions

The openmm.unit.Quantity object containing the particle positions, with units compatible with openmm.unit.nanometers.

system

The openmm.System object corresponding to the test system.

topology

The openmm.app.Topology object corresponding to the test system.