Data containers#

DataContainer#

class mchammer.DataContainer(structure, ensemble_parameters, metadata={})[source]#

Data container for storing information concerned with Monte Carlo simulations performed with mchammer.

This is the data container the thermodynamic ensembles write. On top of the storage and access the base class provides, it analyzes scalar observables through get_average() and analyze_data(), which estimates the correlation length and the statistical error. The section Data containers describes the container in more detail.

Parameters:

  • structure (Atoms) – Reference atomic structure associated with the data container.

  • ensemble_parameters (dict) – Parameters associated with the underlying ensemble.

  • metadata (dict) – Metadata associated with the data container.

Example

The following snippet runs a short simulation in the canonical ensemble and analyzes the potential along the trajectory. The parameters of the cluster expansion are made up to keep the example self-contained:

>>> from ase.build import bulk
>>> from icet import ClusterExpansion, ClusterSpace
>>> from mchammer.calculators import ClusterExpansionCalculator
>>> from mchammer.ensembles import CanonicalEnsemble

>>> # prepare cluster expansion
>>> prim = bulk('Au')
>>> cs = ClusterSpace(prim, cutoffs=[4.3], chemical_symbols=['Ag', 'Au'])
>>> ce = ClusterExpansion(cs, [0, 0, 0.1, -0.02])

>>> # prepare initial configuration
>>> structure = prim.repeat(3)
>>> for k in range(5):
...     structure[k].symbol = 'Ag'

>>> # set up and run the simulation
>>> calculator = ClusterExpansionCalculator(structure, ce)
>>> mc = CanonicalEnsemble(structure=structure, calculator=calculator,
...                        temperature=600)
>>> mc.run(len(structure) * 50)

>>> # analyze the potential, skipping the first steps as equilibration
>>> dc = mc.data_container
>>> print(dc.data)
>>> mean = dc.get_average('potential', start=len(structure) * 10)
>>> summary = dc.analyze_data('potential', start=len(structure) * 10)
>>> print(summary['error_estimate'])
analyze_data(tag, start=None, max_lag=None)[source]#

Returns a statistical analysis of a scalar observable.

Parameters:

  • tag (str) – Name of the observable to analyze.

  • start (int | None) – Smallest trial step to include. By default all records are included.

  • max_lag (int | None) – Maximum lag between two points in the data series when computing the autocorrelation function. By default the length of the data series minus one.

Returns:

The mean, the standard deviation, the correlation length in trial steps, and the error estimate of the mean at 95 % confidence, under the keys mean, std, correlation_length, and error_estimate. The correlation length and the error estimate are nan when the autocorrelation function does not decay within the data series.

Return type:

dict

Raises:

  • ValueError – If the observable is not in the data container.

  • ValueError – If the observable is not scalar.

  • ValueError – If the observations are not evenly spaced in trial steps.

append(mctrial, record)#

Appends a record to the data container.

Parameters:

  • mctrial (int) – Monte Carlo trial step at which the observations were made.

  • record (dict[str, int | float | list]) – Observations as tag-value pairs.

Raises:

  • TypeError – If the input parameters have the wrong type.

  • ValueError – If mctrial is smaller than the step of the last record.

Return type:

None

apply_observer(observer)#

Applies an observer to the stored trajectory and adds its output to the records.

The observer is run for every record that holds an occupation vector, so its interval is ignored. This makes it possible to observe a quantity after the simulation that was not observed while it ran.

Parameters:

observer (BaseObserver) – Observer to apply.

Return type:

None

property data: DataFrame#

All records as a pandas.DataFrame with one column per observable and the trial step in the column mctrial. The occupation vectors are left out.

property ensemble_parameters: dict#

Parameters of the ensemble that produced the data (copy).

get(*tags, start=0)#

Returns the accumulated data for the requested observables.

The configurations stored in the data container are returned as a list of Atoms objects when 'trajectory' is one of the tags. Only records that hold every requested observable are returned, so requesting observables with different intervals together returns the steps at which all of them were recorded.

Parameters:

  • tags (str) – Names of the requested observables. The trial step is available as 'mctrial'.

  • start (int) – Smallest trial step to include. By default all records are included.

Returns:

The data of the requested observable, or a tuple with one entry per tag when several tags are given.

Return type:

ndarray | list[Atoms] | tuple[ndarray, list[Atoms]]

Raises:

  • TypeError – If tags is empty.

  • ValueError – If an observable is requested that is not in the data container.

Examples

Below the get() method is illustrated, but first we require a data container.

>>> from ase.build import bulk
>>> from icet import ClusterExpansion, ClusterSpace
>>> from mchammer.calculators import ClusterExpansionCalculator
>>> from mchammer.ensembles import CanonicalEnsemble

>>> # prepare cluster expansion
>>> prim = bulk('Au')
>>> cs = ClusterSpace(prim, cutoffs=[4.3], chemical_symbols=['Ag', 'Au'])
>>> ce = ClusterExpansion(cs, [0, 0, 0.1, -0.02])

>>> # prepare initial configuration
>>> structure = prim.repeat(3)
>>> for k in range(5):
...     structure[k].symbol = 'Ag'

>>> # set up and run MC simulation
>>> calc = ClusterExpansionCalculator(structure, ce)
>>> mc = CanonicalEnsemble(structure=structure, calculator=calc,
...                        temperature=600,
...                        dc_filename='myrun_canonical.dc')
>>> mc.run(100)  # carry out 100 trial swaps

We can now access the data container by reading it from file by using the read() method. For the purpose of this example, however, we access the data container associated with the ensemble directly.

>>> dc = mc.data_container

The following lines illustrate how to use the get() method for extracting data from the data container.

>>> # obtain all values of the potential represented by
>>> # the cluster expansion along the trajectory
>>> p = dc.get('potential')

>>> import matplotlib.pyplot as plt
>>> # as above but this time the MC trial step is included as well
>>> s, p = dc.get('mctrial', 'potential')
>>> _ = plt.plot(s, p)
>>> plt.show(block=False)

>>> # obtain configurations along the trajectory along with
>>> # their potential
>>> p, confs = dc.get('potential', 'trajectory')
get_average(tag, start=None)[source]#

Returns the average of a scalar observable.

Parameters:

  • tag (str) – Name of the observable to average.

  • start (int | None) – Smallest trial step to include. By default all records are included.

Raises:

  • ValueError – If the observable is not in the data container.

  • ValueError – If the observable is not scalar.

Return type:

float

property metadata: dict#

Metadata associated with the data container.

property observables: list[str]#

Names of the observables in the data container.

classmethod read(infile)#

Reads a data container from file.

Parameters:

infile (str | BinaryIO | TextIO) – File from which to read, as a file name or an open file object.

Raises:

  • FileNotFoundError – If infile is not found.

  • KeyError – If infile is a tar archive that lacks one of the members atoms, reference_data and runtime_data.

  • TypeError – If infile is not a tar archive.

Return type:

BaseDataContainer

write(outfile)#

Writes the data container to file.

The file is a tar archive that holds the reference structure, the records, and the parameters, metadata, and last state as JSON.

Parameters:

outfile (bytes | str) – File to which to write.

Return type:

None

WangLandauDataContainer#

class mchammer.WangLandauDataContainer(structure, ensemble_parameters, metadata={})[source]#

Data container for storing information concerned with Wang-Landau simulations performed with mchammer.

On top of the storage and access the base class provides, it records the fill factor, the histogram, and the entropy of the Wang-Landau algorithm, and it gives access to the density of states and to observables averaged over it. The functions listed under Analysis functions combine the containers of several simulations that cover different energy ranges.

Parameters:

  • structure (Atoms) – Reference atomic structure associated with the data container.

  • ensemble_parameters (dict) – Parameters associated with the underlying ensemble.

  • metadata (dict) – Metadata associated with the data container.

append(mctrial, record)#

Appends a record to the data container.

Parameters:

  • mctrial (int) – Monte Carlo trial step at which the observations were made.

  • record (dict[str, int | float | list]) – Observations as tag-value pairs.

Raises:

  • TypeError – If the input parameters have the wrong type.

  • ValueError – If mctrial is smaller than the step of the last record.

Return type:

None

apply_observer(observer)#

Applies an observer to the stored trajectory and adds its output to the records.

The observer is run for every record that holds an occupation vector, so its interval is ignored. This makes it possible to observe a quantity after the simulation that was not observed while it ran.

Parameters:

observer (BaseObserver) – Observer to apply.

Return type:

None

property data: DataFrame#

All records as a pandas.DataFrame with one column per observable and the trial step in the column mctrial. The occupation vectors are left out.

property ensemble_parameters: dict#

Parameters of the ensemble that produced the data (copy).

property fill_factor: float#

Final value of the fill factor in the Wang-Landau algorithm.

property fill_factor_history: DataFrame#

Evolution of the fill factor in the Wang-Landau algorithm.

get(*tags, fill_factor_limit=None)[source]#

Returns the accumulated data for the requested observables.

The configurations stored in the data container are returned as a list of Atoms objects when 'trajectory' is one of the tags.

Parameters:

  • tags (str) – Names of the requested observables. The trial step is available as 'mctrial'.

  • fill_factor_limit (float | None) – Only return data recorded up to the trial step at which the fill factor first reached this value. By default all data is returned.

Returns:

The data of the requested observable, or a tuple with one entry per tag when several tags are given.

Return type:

ndarray | list[Atoms] | tuple[ndarray, list[Atoms]]

Raises:

  • TypeError – If tags is empty.

  • ValueError – If an observable is requested that is not in the data container.

Examples

Below the get() method is illustrated, but first we require a data container.

>>> from ase import Atoms
>>> from icet import ClusterExpansion, ClusterSpace
>>> from mchammer.calculators import ClusterExpansionCalculator
>>> from mchammer.ensembles import WangLandauEnsemble

>>> # prepare cluster expansion
>>> prim = Atoms('Au', positions=[[0, 0, 0]], cell=[1, 1, 10], pbc=True)
>>> cs = ClusterSpace(prim, cutoffs=[1.1], chemical_symbols=['Ag', 'Au'])
>>> ce = ClusterExpansion(cs, [0, 0, 2])

>>> # prepare initial configuration
>>> structure = prim.repeat((4, 4, 1))
>>> for k in range(8):
...     structure[k].symbol = 'Ag'

>>> # set up and run Wang-Landau simulation
>>> calculator = ClusterExpansionCalculator(structure, ce)
>>> mc = WangLandauEnsemble(structure=structure,
...                         calculator=calculator,
...                         energy_spacing=1,
...                         dc_filename='ising_2d_run.dc',
...                         fill_factor_limit=0.3)
>>> # N.B.: in practice one requires more steps
>>> mc.run(number_of_trial_steps=len(structure) * 3000)

We can now access the data container by reading it from file by using the read() method. For the purpose of this example, however, we access the data container associated with the ensemble directly.

>>> dc = mc.data_container

The following lines illustrate how to use the get() method for extracting data from the data container.

>>> # obtain all values of the potential represented by
>>> # the cluster expansion and the MC trial step along the
>>> # trajectory
>>> import matplotlib.pyplot as plt
>>> s, p = dc.get('mctrial', 'potential')
>>> _ = plt.plot(s, p)

>>> # as above but this time only included data recorded up to
>>> # the point when the fill factor reached below 0.6
>>> s, p = dc.get('mctrial', 'potential', fill_factor_limit=0.6)
>>> _ = plt.plot(s, p)
>>> plt.show(block=False)

>>> # obtain configurations along the trajectory along with
>>> # their potential
>>> p, confs = dc.get('potential', 'trajectory')
get_entropy(fill_factor_limit=None)[source]#

Returns the (relative) entropy from this data container accumulated during a Wang-Landau simulation. Returns None if the data container does not contain the required information.

Parameters:

fill_factor_limit (float | None) – Return the entropy recorded up to the point when the specified fill factor limit was reached, or None if the entropy history is empty or the last fill factor is above the limit. Otherwise return the entropy for the last state.

Return type:

DataFrame

get_histogram()[source]#

Returns the histogram from this data container accumulated since the last update of the fill factor. Returns None if the data container does not contain the required information.

Return type:

DataFrame

property metadata: dict#

Metadata associated with the data container.

property observables: list[str]#

Names of the observables in the data container.

classmethod read(infile)[source]#

Reads data container from file.

Parameters:

infile (str | BinaryIO | TextIO) – file from which to read

Raises:

  • FileNotFoundError – if file is not found (str)

  • ValueError – if file is of incorrect type (not a tarball)

Return type:

WangLandauDataContainer

write(outfile)#

Writes the data container to file.

The file is a tar archive that holds the reference structure, the records, and the parameters, metadata, and last state as JSON.

Parameters:

outfile (bytes | str) – File to which to write.

Return type:

None

Analysis functions#

mchammer.data_containers.get_average_observables_wl(dcs, temperatures, observables=None, boltzmann_constant=8.617330337217213e-05, fill_factor_limit=None)[source]#

Returns the average and the standard deviation of the energy from a Wang-Landau simulation for the temperatures specified. If the observables keyword argument is specified the function will also return the mean and standard deviation of the specified observables.

Parameters:

  • dcs (WangLandauDataContainer | dict[Any, WangLandauDataContainer]) – Data container(s) from which to extract density of states as well as observables.

  • temperatures (list[float]) – Temperatures at which to compute the averages.

  • observables (list[str] | None) – Observables for which to compute averages. The observables must refer to fields in the data container.

  • boltzmann_constant (float) – Boltzmann constant \(k_B\) in appropriate units, i.e., units that are consistent with the underlying cluster expansion and the temperature units. By default eV/K.

  • fill_factor_limit (float | None) – Use data recorded up to the point when the specified fill factor limit was reached when computing averages. Otherwise use data for the last state.

Raises:

  • ValueError – If the data container(s) do(es) not contain entropy data from Wang-Landau simulation.

  • ValueError – If data container(s) do(es) not contain requested observable.

Return type:

DataFrame

mchammer.data_containers.get_average_cluster_vectors_wl(dcs, cluster_space, temperatures, boltzmann_constant=8.617330337217213e-05, fill_factor_limit=None)[source]#

Returns the average cluster vectors from a Wang-Landau simulation for the temperatures specified.

Parameters:

  • dcs (WangLandauDataContainer | dict) – Data container(s), from which to extract density of states as well as observables.

  • cluster_space (ClusterSpace) – Cluster space to use for calculation of cluster vectors.

  • temperatures (list[float]) – Temperatures at which to compute the averages.

  • boltzmann_constant (float) – Boltzmann constant \(k_B\) in appropriate units, i.e., units that are consistent with the underlying cluster expansion and the temperature units. By default eV/K.

  • fill_factor_limit (float | None) – Use data recorded up to the point when the specified fill factor limit was reached when computing the average cluster vectors. Otherwise use data for the last state.

Raises:

ValueError – If the data container(s) do(es) not contain entropy data from Wang-Landau simulation.

Return type:

DataFrame

mchammer.data_containers.get_density_of_states_wl(dcs, fill_factor_limit=None)[source]#

Returns a pandas DataFrame with the total density of states from a Wang-Landau simulation. If a dict of data containers is provided the function also returns a dictionary that contains the standard deviation between the entropy of neighboring data containers in the overlap region. These errors should be small compared to the variation of the entropy across each bin.

The function can handle both a single data container and a dict thereof. In the latter case the data containers must cover a contiguous energy range and must at least partially overlap.

Parameters:

  • dcs (WangLandauDataContainer | dict[Any, WangLandauDataContainer]) – Data container(s), from which to extract the density of states.

  • fill_factor_limit (float | None) – Calculate the density of states using the entropy recorded up to the point when the specified fill factor limit was reached. Otherwise return the density of states for the last state.

Raises:

  • TypeError – If dcs does not correspond to not a single (dictionary) of data container(s) from which the entropy can retrieved.

  • ValueError – If the data container does not contain entropy information.

  • ValueError – If a fill factor limit has been provided and the data container either does not contain information about the entropy history or if the last fill factor is higher than the specified limit.

  • ValueError – If multiple data containers are provided and there are inconsistencies with regard to basic simulation parameters such as system size or energy spacing.

  • ValueError – If multiple data containers are provided and there is at least one energy region without overlap.

Return type:

tuple[DataFrame, dict]

mchammer.free_energy_tools.get_free_energy_thermodynamic_integration(dc, cluster_space, forward, max_temperature=inf, sublattice_probabilities=None, boltzmann_constant=8.617330337217213e-05)[source]#

Returns the free energy calculated via thermodynamic integration using the ThermodynamicIntegrationEnsemble.

The temperature dependence of the free energy can be extracted from the thermodynamic integration as

\[F(T) = \frac{F_0(\lambda)}{\lambda} + \frac{T_0}{\lambda} S_\text{B}\]

where \(S_\text{B}\) is the Boltzmann entropy,

\[T = \frac{T_0}{\lambda}\]

and

\[F_0(\lambda) = \int_0^\lambda \left\langle\frac{\mathrm{d}H(\lambda)} {\mathrm{d}\lambda}\right\rangle_{H} \mathrm{d}\lambda\]
Parameters:

  • dc (DataContainer) – Data container from the thermodynamic integration simulation.

  • cluster_space (ClusterSpace) – Cluster space used to construct the cluster expansion used for the simulation.

  • forward (bool) – If True the thermodynamic integration was carried out from the reference system to the cluster expansion, otherwise in the opposite direction.

  • max_temperature (float) – Largest temperature to extract from the thermodynamic integration.

  • sublattice_probababilites – Sublattice probabilties that were provided to the thermodynamic integration simulation.

  • boltzmann_constant (float) – Boltzmann constant in the units used for the thermodynamic integration simulation.

Return type:

(ndarray, ndarray)

mchammer.free_energy_tools.get_free_energy_temperature_integration(dc, cluster_space, forward, temperature_reference, free_energy_reference=None, sublattice_probabilities=None, max_temperature=inf, boltzmann_constant=8.617330337217213e-05)[source]#

Returns the free energy calculated using temperature integration and the corresponding temperature.

\[\frac{A(T_{2})}{T_{2}} = \frac{A(T_{1})}{T_{1}} - \int_{T_{1}}^{T_{2}}\frac{U(T)}{T^2}\mathrm{d}T\]
Parameters:

  • dc (DataContainer) – Data container from a canonical annealing simulation. The first (last for forward=False) temperature in the data container has to be at the same temperature as temperature_reference.

  • cluster_space (ClusterSpace) – Cluster space used to construct the cluster expansion that was used in the simulations.

  • forward (bool) – If True the canonical annealing simulation was carried out from high to low temperature, otherwise the opposite is assumed.

  • temperature_reference (float) – Temperature at which free_energy_reference was calculated.

  • free_energy_reference (float | None) – Reference free energy. By default the free energy at temperature_reference is approximated by \(-T S_B\), where \(S_B\) is the ideal mixing entropy.

  • sublattice_probabilities (list[float] | None) – Sublattice probabilities that were provided to the canonical annealing simulation.

  • max_temperature (float) – Largest temperature to extract from the temperature integration.

  • boltzmann_constant (float) – Boltzmann constant in the units used for the thermodynamic integration simulation.

Return type:

(ndarray, ndarray)

mchammer.data_analysis.analyze_data(data, max_lag=None)[source]#

Returns a statistical analysis of a data series.

Parameters:

  • data (ndarray) – Data series to analyze.

  • max_lag (int | None) – Maximum lag between two data points used for computing the autocorrelation function. By default the length of the data series minus one.

Returns:

The mean, the standard deviation, the correlation length, and the error estimate of the mean at 95 % confidence, under the keys mean, std, correlation_length, and error_estimate. The correlation length and the error estimate are nan when the autocorrelation function does not decay within the data series.

Return type:

dict

mchammer.data_analysis.get_autocorrelation_function(data, max_lag=None)[source]#

Returns the autocorrelation function of a data series.

The autocorrelation function is computed using pandas.Series.autocorr.

Parameters:

  • data (ndarray) – Data series for which to compute the autocorrelation function.

  • max_lag (int | None) – Maximum lag between two data points. By default the length of the data series minus one.

Returns:

The autocorrelation for each lag from zero to max_lag minus one.

Return type:

ndarray

Raises:

ValueError – If max_lag is smaller than one or not smaller than the length of the data series.

mchammer.data_analysis.get_correlation_length(data)[source]#

Returns an estimate of the correlation length of a data series.

The correlation length is the first lag at which the autocorrelation function drops below \(\exp(-2)\). If the autocorrelation function never drops below that value within the data series, the function returns None.

Parameters:

data (ndarray) – Data series for which to estimate the correlation length.

Return type:

int | None

mchammer.data_analysis.get_error_estimate(data, confidence=0.95)[source]#

Returns an estimate of the standard error of the mean of a data series at the given confidence level via

\[\mathrm{error} = t_\mathrm{factor} * \mathrm{std}(\mathrm{data}) / \sqrt{N_s}\]

where \(t_\mathrm{factor}\) is the factor corresponding to the confidence level and \(N_s\) is the number of independent measurements, which is the length of the data series divided by the correlation length.

If the correlation length cannot be estimated because the autocorrelation function does not decay within the data series, the function returns None.

Parameters:

  • data (ndarray) – Data series for which to estimate the error.

  • confidence (float) – Confidence level of the error estimate. By default 0.95.

Return type:

float | None