Langevin

class Langevin(initial_velocity=None, time_step=None, reservoir_temperature=None, heating_rate=None, friction=None, random_seed=None)

NVT Langevin integrator class which implements the Langevin thermostat.

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.0 * 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
  • 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: 0.0 * Kelvin s**-1
  • thermostat_timescale (PhysicalQuantity of type time) – The time constant for Berendsen temperature coupling.
    Default: 100.0 * fs
  • friction (PhysicalQuantity of type 1/time) – The friction constant in the Langevin equation.
    Default: 0.01 * fs**-1
  • random_seed (int) – The seed for the random generator providing the stochastic forces. Must be between 0 and 2**32.
    Default: The default random seed
kineticEnergy(configuration)
Parameters:configuration (DistributedConfiguration) – The current configuration to calculate the kinetic energy of.
Returns:The kinetic energy of the current configuration.
Return type:PhysicalQuantity of type energy
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 Langevin thermostat:

# 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 MD method
method = Langevin(
    time_step=1*femtoSecond,
    reservoir_temperature=300*Kelvin,
    friction=0.01*femtoSecond**-1,
    heating_rate=0*Kelvin/picoSecond,
    initial_velocity=None
)

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

Notes

  • The Langevin thermostat can be used to simulate a canonical NVT-ensemble, by integrating the Langevin equation

    \[m_i \frac{d^2 \bf{r}_i}{dt^2} = -m_i \gamma \bf{v}_i + \bf{F}_i + \bf{\Gamma}_i \, ,\]

    which adds a friction term with a friction constant \(\gamma\) and a stochastic random force \(\Gamma\) to Newton’s equations of motion. QuantumATK implements the algorithm proposed by Goga et al [GRdV+12].

  • In contrast to other NVT methods, such as NVTNoseHoover or NVTBerendsen, the dynamics of the system can be altered significantly by the Langevin thermostat, and the evolution of the system is no longer deterministic. It should therefore primarily be used to sample configurations from the canonical ensemble, rather than to calculate dynamical properties.

  • You can specify one or more thermostats acting on tagged sub-groups of atoms of the configuration, by giving a list of (tag-name, temperature)-tuples instead of a single, global reservoir_temperature value.

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

  • Setting the reservoir_temperature to 0K, will switch off the stochastic force, resulting in damped dynamics, which can be used to coarsely relax a configuration.

[GRdV+12]N. Goga, A. J. Rzepiela, A. H. de Vries, S. J. Marrink, and H. J. C. Berendsen. Efficient algorithms for langevin and dpd dynamics. Journal of Chemical Theory and Computation, 8(10):3637–3649, 2012. PMID: 26593009. doi:10.1021/ct3000876.