SemiGrandCanonicalTransitionMatrixHook¶
- class SemiGrandCanonicalTransitionMatrixHook(elements, chemical_potential=None, temperature=None, swap_ratio=None, random_seed=None, call_interval=None, pnpd_update_interval=None, measurement_interval=None)¶
A hook function that implements semi-grand canonical Monte Carlo simulation
- Parameters:
elements (Sequence of
PeriodicTableElement
) – The elements being swapped between.chemical_potential (PhysicalQuantity of type energy) – The chemical potential of swapping from the first to the second atom. Default: 0 eV
temperature (PhysicalQuantity of type temperature) – The Monte Carlo selection temperature Default: 300 Kelvin
swap_ratio (float) – The fraction of available atoms that are attempted to be swapped at each call. Default: 0.01
random_seed (int) – The seed for the random number generator used to select the Monte Carlo moves.
call_interval (int) – The number of dynamics steps between each Monte Carlo step. Default: 1
pnpd_update_interval (int) – The number of Monte Carlo steps between each update of the particle number probability distribution. Note that this is based on Monte Carlo steps, not dynamics steps. Default: 1
measurement_interval (int) – The number of dynamics steps between recording the particle number probability distribution. Default: 1
- callInterval()¶
- Returns:
The call interval of this hook function.
- Return type:
int
- chemicalPotential()¶
- Returns:
The chemical potential difference used in calculating swap probabilities.
- Return type:
PhysicalQuantity of type energy
- collectionMatrix()¶
- Returns:
The tri-banded collection in which sampled transition probabilities are collected.
- Return type:
ndarray
- elements()¶
- Returns:
List of the elements being swapped in the Monte Carlo step.
- Return type:
list of
PeriodicTableElement
- measurementInterval()¶
- Returns:
The number of dynamics steps between writing data to the trajectory.
- Return type:
int
- monteCarloQuantity()¶
The value being averaged in the Monte Carlo simulation. In this simulation this is the fraction of the second atom in the configuration.
- Returns:
The ratio of second atoms to the first in the configuration.
- Return type:
float
- particleNumberProbabilityDistribution()¶
- Returns:
The relative probability of being in each composition state.
- Return type:
ndarray
- pnpdUpdateInterval()¶
- Returns:
The number of dynamics steps between measurements of the PNPD.
- Return type:
int
- randomSeed()¶
- Returns:
The random seed used to create the Monte Carlo random number generator.
- Return type:
int
- sampleCount()¶
- Returns:
The number of transition samples from each state.
- Return type:
ndarray
- swapElementCount(configuration)¶
Count the number of atoms that have been swapped to the second element.
- Parameters:
configuration (
BulkConfiguration
) – The current configuration.- Returns:
The number of elements
- Return type:
int
- swapRatio()¶
- Returns:
The fraction of atoms to attempt to swap in each Monte Carlo cycle.
- Return type:
float
- temperature()¶
- Returns:
The temperature of the simulation.
- Return type:
PhysicalQuantity of type temperature
- uniqueString()¶
Return a unique string representing the state of the object.
In this example the semi-grand canonical transition matrix Monte Carlo method is used to study the
composition of a gold-silver alloy. The result of the calculation is a
ConfigurationDataContainer
that contains the particle number probability distribution
(PNPD) for the system. Since these metals are miscible, the PNPD should have a single maximum at
the average composition value of around 25% gold and 75% silver.
A workflow for the calculation can be downloaded here
with a corresponding python script Au_Ag_Transition_Matrix.py
. This simulation may take around several hours to run, given the
number of samples required.
Notes¶
The SemiGrandCanonicalTransitionMatrixHook
defines a transition matrix Monte Carlo
simulation in the semi-grand canonical ensemble. This is a complimentary technique to the
SemiGrandCanonicalMonteCarloHook
method. Both simulations are carried out in the
semi-grand canonical ensemble and can be used to determine the composition of materials at
different chemical potentials. While in normal Monte Carlo methods the property of interest is
directly sampled, in transition matrix Monte Carlo the probability of transitioning from one state
to another is sampled[1]. Here this is the number of each element in the
material. As only one atom is swapped at a time, the probability of being in any state can be
calculated. This is known as the particle number probability density (PNPD). In order to
effectively calculate the PNPD it is necessary to sample all states. This is done by biasing
the selection by the negative change in the PNPD. Since the PNPD is needed to effectively sample
the states, this is calculated and updated during the simulation. Once the converges the simulation
samples evenly across all states.
The PNPD relates straightforwardly to the relative free energy of each state such that:
Here \(\Pi(N; \mu, P, T)\) is the PNPD at given particle number \(N\) and chemical potential, temperature and pressure and \(W\) is similarly the free energy. To get the average system property from the PNPD, that property is averaged over the available states. For instance, the average number of elements in the material \(\left<N\right>\) can be given as:
The PNPD can also be re-weighted to different chemical potential values. This allows exploring the whole range of chemical potentials with one simulation. The PNPD at a different chemical potential can be given as:
Here \(\mu_0\) is the chemical potential that the simulations was originally carried out at, \(\mu\) is the new chemical potential and \(C\) is a normalization constant.
One of the main advantages of transition matrix simulations is that it samples the entire space of compositions, which works well with systems that have multiple stable states. In the case of immiscible metals, there can be a minimum in the free energy close to each pure state, and then a large free energy barrier as the elements mix. As the chemical potential changes the relative energies of the minima change, leading to one being more populated than the other. This situation is very difficult to model with traditional Monte Carlo techniques, but can be explored with transition Matrix techniques. The trade-off here is that transition Matrix Monte Carlo methods require at least an order of magnitude more steps than a traditional Monte Carlo simulation to reach convergence. This is due to the transition Matrix simulation scanning the entire space of configurations rather than just a small amount around a minimum in the free energy.
To practically carry out a transition Matrix Monte Carlo simulation, it is recommended that the
chemical potential should be near the difference in free energy between the two materials. The PNPD
stored in the trajectory is the natural log of the PNPD, as this is numerically much easier to
handle. In addition to the arguments to the arguments from the
SemiGrandCanonicalMonteCarloHook
, the SemiGrandCanonicalTransitionMatrixHook
has two additional arguments that control the frequency of updating and storing the PNPD. Frequent
updating of the PNPD can slow the simulation in some cases. This is controlled by the argument
pnpd_update_interval
, which gives the number of Monte Carlo cycles between PNPD updates.
Similarly measurement_interval
gives the frequency with which the PNPD is written to the
trajectory. This should be done often enough so that convergence in the PNPD can be observed, but
not so much as to take large amounts of space in the trajectory. This frequency is also based on
the number of dynamics steps between recording each PNPD. As well as the PNPD the trajectory also
records the number of samples from each state. This is stored as the measurement sample_count
.
For a converged PNPD each state should be sampled enough times to have an accurate probability of
transitioning from that state.