ForceBiasMonteCarlo¶
- class ForceBiasMonteCarlo(temperature=None, max_atom_displacement=None, heating_rate=None, random_seed=None, max_random_attempts=None)¶
Set up the ForceBiasMonteCarlo object.
- Parameters:
temperature (PhysicalQuantity of type temperature | None) – The temperature at which the Monte Carlo simulation should be run. Default: 300*Kelvin
max_atom_displacement (PhysicalQuantity of type length | None) – The maximum distance an atom can move in each Cartesian direction during a single step. Default: 0.1*Angstrom
heating_rate (PhysicalQuantity of type temperature | None) – The change in temperature per step. Default: 0*Kelvin
random_seed (int | None) – The seed for the random generator. Must be between 0 and 2**32. Default: The default random seed
max_random_attempts (int) – The maximum attempts used in the rejection sampling of the probability distribution. Default: 500
- monteCarloStep(configuration, forces, stress, constraints=None)¶
Perform a Monte Carlo step and apply the new positions to the configuration.
- Parameters:
configuration (
MoleculeConfiguration
|BulkConfiguration
|DeviceConfiguration
|SurfaceConfiguration
) – The configuration on which the Monte Carlo step should be performed.forces (PhysicalQuantity of type energy/length) – The atomic forces for the given configuration.
stress (PhysicalQuantity of type pressure) – The current stress in the unit cell.
constraints (list of type BaseConstraint | None) – The list of constraints to apply.
- monteCarloTimeStep()¶
- Returns:
The average time that elapses between each step according to the time-stamped force bias MC (tfMC) formalism or None if not initialized.
- Return type:
PhysicalQuantity of type time | None
- randomSeed()¶
- Returns:
The random seed used to initialize the random number generator.
- Return type:
int
- reservoirTemperature()¶
- Returns:
The current reservoir temperature.
- Return type:
PhysicalQuantity of type temperature
- uniqueString()¶
Return a unique string representing the state of the object.
Usage Example¶
Use the ForceBiasMonteCarlo object to run a TimeStampedForceBiasMonteCarlo simulation:
# -------------------------------------------------------------
# Time-Stamped Force-Bias Monte Carlo
# -------------------------------------------------------------
method = ForceBiasMonteCarlo(
reservoir_temperature=500.0*Kelvin,
max_atom_displacement=0.3*Ang,
)
mc_trajectory = TimeStampedForceBiasMonteCarlo(
bulk_configuration,
constraints=[],
trajectory_filename='tfmc_trajectory.hdf5',
steps=500,
log_interval=50,
method=method
)
Notes¶
The ForceBiasMonteCarlo class is used in the TimeStampedForceBiasMonteCarlo function to run a time-stamped force-bias Monte Carlo (TFMC) simulation [1].
A Monte Carlo step is taken by displacing each Cartesian component
\(\alpha\) of the position of atom i by a distance
\(\zeta_{i,\alpha}\Delta_i\). The maximum displacement magnitude for each
atom is given by
\(\Delta_i = \Delta \left(m_{\textrm{min}}/m_i\right)^{1/4}\)
[2], where \(\Delta\) is specified by the argument
max_atom_displacement
, and \(m_i\) are the atomic masses.
\(\zeta_{i, \alpha}\) is sampled from the interval \([-1, 1]\) with a probability:
where \(\gamma_{i,\alpha} = 0.5 F_{i, \alpha} \Delta/k_B T\),
\(F_{i,\alpha}\) is the force component on atom i, and T is the
temperature at which the configurations are sampled, specified by the argument
reservoir_temperature
.
By adjusting the max_atom_displacement
parameter one can tune accuracy
vs. efficiency of the simulation. A small value results in a more accurate
sampling of the canonical ensemble, whereas a large value increases the efficiency
at which the phase space is sampled. Typically, values in the range
\(0.1 R_{eq} \, - \, 0.3 R_{eq}\) are a reasonable choice, where \(R_{eq}\)
represents an equilibrium bond length in the system.
If a non-zero heating_rate
is specified, the reservoir temperature will be
changed by the given value after each Monte Carlo step, resulting in an
increase or decrease of the temperature during sampling.