NPTMartynaTobiasKlein

class NPTMartynaTobiasKlein(initial_velocity=None, time_step=None, reservoir_temperature=None, thermostat_timescale=None, reservoir_pressure=None, barostat_timescale=None, coupling_mask=None, heating_rate=None, compression_rate=None, chain_length=None)

The NPTMartynaTobiasKlein integrator class which implements the barostat proposed by Martyna, Tobias, and Klein.

Parameters:
  • initial_velocity (ConfigurationVelocities | ZeroVelocities | MaxwellBoltzmannDistribution) – A class that implements a distribution of initial velocities for the particles in the MD simulation.
    Default: ZeroVelocities
  • time_step (PhysicalQuantity of type time) – The time-step interval used in the MD simulation.
    Default: 1 * fs
  • reservoir_temperature (PhysicalQuantity of type temperature | list) – The reservoir temperature in the simulation. The temperature can be given as a single temperature value for the entire system, or as a list of 2-tuples of str and PhysicalQuantity of type temperature, applying local thermostats to the tagged groups of atoms. E.g. [('group1', 280.0 * Kelvin), ('group2', 320.0 * Kelvin)].
    Default: 300.0 * Kelvin
  • thermostat_timescale (PhysicalQuantity of type time) – The time constant for Berendsen temperature coupling.
    Default: 100.0 * fs
  • reservoir_pressure (PhysicalQuantity of type pressure) – The reservoir pressure in the simulation. The pressure can be given either as a scalar value, as a vector of length 3 denoting the hydrostatic components, as a vector of length 6 in Voigt notation, or as a 3x3 tensor. A scalar value will result in isotropic pressure coupling, whereas for all other representations, each pressure component will be coupled to the barostat independently.
    Default: 1.0 * bar
  • barostat_timescale (PhysicalQuantity of type time) – The time constant for Berendsen pressure coupling.
    Default: 1000.0 * fs
  • coupling_mask (array of bools) – The mask determining which elements of the stress tensor, the barostat should couple to. Can be given as a vector of length 3, denoting the hydrostatic components, as a vector of length 6 in Voigt notation, or as a 3x3 tensor. By default the barostat couples only to the diagonal elements.
    Default: True for diagonal elements, False for shear components.
  • heating_rate (PhysicalQuantity of type temperature/time | None) – The heating rate of the target temperature. A value of None disables the heating of the system.
    Default: None
  • compression_rate (PhysicalQuantity of type pressure/time | None) – The compression rate of the target pressure.
    Default: 0.0 * bar / fs
  • chain_length (int) – The number of subsequent thermostats.
    Default: 3
isotropicCoupling()
Returns:Whether the barostat applies isotropic or anisotropic pressure coupling.
Return type:bool
kineticEnergy(configuration)
Returns:The kinetic energy of the current configuration.
Return type:PhysicalQuantity of type energy
thermostats()
Returns:The list of thermostats.
Return type:list
timeStep()
Returns:The time step.
Return type:PhysicalQuantity of type time

Usage Example

Perform a molecular dynamics run of 50 steps on FCC Si, using the MartynaTobiasKlein barostat with isotropic pressure coupling:

# Set up configuration
bulk_configuration = BulkConfiguration(
    bravais_lattice=FaceCenteredCubic(5.4306*Angstrom),
    elements=[Silicon, Silicon],
    fractional_coordinates=[[0.0, 0.0, 0.0], [0.25, 0.25, 0.25]]
    )

# Set calculator
calculator = TremoloXCalculator(parameters=Tersoff_Si_1988b())
bulk_configuration.setCalculator(calculator)

# Set up a MTK-barostat with isotropic pressure coupling.
method = NPTMartynaTobiasKlein(
    time_step=1*femtoSecond,
    reservoir_temperature=300*Kelvin,
    thermostat_timescale=100*femtoSecond,
    reservoir_pressure=1.0*bar,
    barostat_timescale=500.0*femtoSecond,
    heating_rate=0*Kelvin/picoSecond,
    compression_rate=0.0*bar/picoSecond,
    chain_length=3,
    initial_velocity=None
)

# Run MD simulation
md_trajectory = MolecularDynamics(
    bulk_configuration,
    constraints=[],
    trajectory_filename='trajectory.nc',
    steps=50,
    log_interval=10,
    method=method
)

Apply anisotropic pressure coupling only to the zz-direction of the pressure tensor.

method = NPTMartynaTobiasKlein(
    time_step=1*femtoSecond,
    reservoir_temperature=300*Kelvin,
    thermostat_timescale=100*femtoSecond,
    reservoir_pressure=[0.0, 0.0, 10.0]*GPa,
    barostat_timescale=500.0*femtoSecond,
    heating_rate=0*Kelvin/picoSecond,
    coupling_mask=[False, False, True, False, False, False],
    compression_rate=0.0*bar/picoSecond,
    chain_length=3,
    initial_velocity=None
)

Notes

  • The MartynaTobiasKlein barostat can be used to simulate a canonical NPT-ensemble. Similar to the NVTNoseHoover algorithm, chains of subsequent thermostats and barostats are used to suppress large temperature and pressure oscillations. The equations of motion are taken from Ref. [MTK94], whereas the integration algorithm is implemented as described in Ref [TALopezRendon+06].

  • If a scalar reservoir_pressure value is given, isotropic pressure coupling is applied, i.e. all cell vectors are rescaled by the same factor, and the isotropic pressure \(P=1/3 (P_{xx} + P_{yy} + P_{zz})\) is used.

    If a pressure vector (of length 3 or 6, in Voigt notation), or a 3x3-tensor is given as reservoir_pressure, anisotropic pressure coupling is applied, meaning that each degree of freedom of the cell vectors will be rescaled independently.

You can enable/disable the coupling to selected components of the pressure tensor via the parameter coupling_matrix. This parameter only has an effect with anisotropic pressure coupling.

By default, the barostat couples only to the hydrostatic pressure, i.e. the diagonal elements of the pressure tensor.

  • If a non-zero heating_rate is specified, the reservoir temperature will be changed linearly during the simulation, according to the specified heating rate.

    If a non-zero compression_rate is specified, the reservoir pressure will be changed linearly during the simulation, according to the specified compression rate.

[MTK94]Glenn J. Martyna, Douglas J. Tobias, and Michael L. Klein. Constant pressure molecular dynamics algorithms. J. Chem. Phys., 101(5):4177–4189, 1994. doi:http://dx.doi.org/10.1063/1.467468.
[TALopezRendon+06]Mark E Tuckerman, José Alejandre, Roberto López-Rendón, Andrea L Jochim, and Glenn J Martyna. A liouville-operator derived measure-preserving integrator for molecular dynamics simulations in the isothermal–isobaric ensemble. J. Phys. A: Mathematical and General, 39(19):5629, 2006. URL: http://stacks.iop.org/0305-4470/39/i=19/a=S18.