Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

"""Definition of the canonical ensemble class.""" 

 

import numpy as np 

 

from ase import Atoms 

from ase.units import kB 

from typing import List 

 

from .. import DataContainer 

from ..calculators.base_calculator import BaseCalculator 

from .thermodynamic_base_ensemble import ThermodynamicBaseEnsemble 

 

 

class CanonicalEnsemble(ThermodynamicBaseEnsemble): 

"""Instances of this class allow one to simulate systems in the 

canonical ensemble (:math:`N_iVT`), i.e. at constant temperature 

(:math:`T`), number of atoms of each species (:math:`N_i`), and 

volume (:math:`V`). 

 

The probability for a particular state in the canonical ensemble is 

proportional to the well-known Boltzmann factor, 

 

.. math:: 

 

\\rho_{\\text{C}} \\propto \\exp [ - E / k_B T ]. 

 

Since the concentrations or equivalently the number of atoms of each 

species is held fixed in the canonical ensemble, a trial step must 

conserve the concentrations. This is accomplished by randomly picking two 

unlike atoms and swapping their identities. The swap is accepted with 

probability 

 

.. math:: 

 

P = \\min \\{ 1, \\, \\exp [ - \\Delta E / k_B T ] \\}, 

 

where :math:`\\Delta E` is the change in potential energy caused by the 

swap. 

 

The canonical ensemble provides an ideal framework for studying the 

properties of a system at a specific concentration. Properties such as 

potential energy or phenomena such as chemical ordering at a specific 

temperature can conveniently be studied by simulating at that temperature. 

The canonical ensemble is also a convenient tool for "optimizing" a 

system, i.e., finding its lowest energy chemical ordering. In practice, 

this is usually achieved by simulated annealing, i.e. the system is 

equilibrated at a high temperature, after which the temperature is 

continuously lowered until the acceptance probability is almost zero. In a 

well-behaved system, the chemical ordering at that point corresponds to a 

low-energy structure, possibly the global minimum at that particular 

concentration. 

 

Parameters 

---------- 

structure : :class:`Atoms <ase.Atoms>` 

atomic configuration to be used in the Monte Carlo simulation; 

also defines the initial occupation vector 

calculator : :class:`BaseCalculator <mchammer.calculators.ClusterExpansionCalculator>` 

calculator to be used for calculating the potential changes 

that enter the evaluation of the Metropolis criterion 

temperature : float 

temperature :math:`T` in appropriate units [commonly Kelvin] 

boltzmann_constant : float 

Boltzmann constant :math:`k_B` in appropriate 

units, i.e. units that are consistent 

with the underlying cluster expansion 

and the temperature units [default: eV/K] 

user_tag : str 

human-readable tag for ensemble [default: None] 

data_container : str 

name of file the data container associated with the ensemble 

will be written to; if the file exists it will be read, the 

data container will be appended, and the file will be 

updated/overwritten 

random_seed : int 

seed for the random number generator used in the Monte Carlo 

simulation 

ensemble_data_write_interval : int 

interval at which data is written to the data container; this 

includes for example the current value of the calculator 

(i.e. usually the energy) as well as ensembles specific fields 

such as temperature or the number of atoms of different species 

data_container_write_period : float 

period in units of seconds at which the data container is 

written to file; writing periodically to file provides both 

a way to examine the progress of the simulation and to back up 

the data [default: np.inf] 

trajectory_write_interval : int 

interval at which the current occupation vector of the atomic 

configuration is written to the data container. 

sublattice_probabilities : List[float] 

probability for picking a sublattice when doing a random swap. 

This should be as long as the number of sublattices and should 

sum up to 1. 

 

 

Example 

------- 

The following snippet illustrate how to carry out a simple Monte Carlo 

simulation in the canonical ensemble. Here, the parameters of the cluster 

expansion are set to emulate a simple Ising model in order to obtain an 

example that can be run without modification. In practice, one should of 

course use a proper cluster expansion:: 

 

from ase.build import bulk 

from icet import ClusterExpansion, ClusterSpace 

from mchammer.calculators import ClusterExpansionCalculator 

from mchammer.ensembles import CanonicalEnsemble 

 

# prepare cluster expansion 

# the setup emulates a second nearest-neighbor (NN) Ising model 

# (zerolet and singlet ECIs are zero; only first and second neighbor 

# pairs are included) 

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, 

data_container='myrun_canonical.dc') 

mc.run(100) # carry out 100 trial swaps 

""" 

 

def __init__(self, structure: Atoms, calculator: BaseCalculator, 

temperature: float, user_tag: str = None, 

boltzmann_constant: float = kB, 

data_container: DataContainer = None, random_seed: int = None, 

data_container_write_period: float = np.inf, 

ensemble_data_write_interval: int = None, 

trajectory_write_interval: int = None, 

sublattice_probabilities: List[float] = None) -> None: 

 

self._ensemble_parameters = dict(temperature=temperature) 

 

# add species count to ensemble parameters 

symbols = set([symbol for sub in calculator.sublattices 

for symbol in sub.chemical_symbols]) 

for symbol in symbols: 

key = 'n_atoms_{}'.format(symbol) 

count = structure.get_chemical_symbols().count(symbol) 

self._ensemble_parameters[key] = count 

 

super().__init__( 

structure=structure, calculator=calculator, user_tag=user_tag, 

data_container=data_container, 

random_seed=random_seed, 

data_container_write_period=data_container_write_period, 

ensemble_data_write_interval=ensemble_data_write_interval, 

trajectory_write_interval=trajectory_write_interval, 

boltzmann_constant=boltzmann_constant) 

 

159 ↛ 162line 159 didn't jump to line 162, because the condition on line 159 was never false if sublattice_probabilities is None: 

self._swap_sublattice_probabilities = self._get_swap_sublattice_probabilities() 

else: 

self._swap_sublattice_probabilities = sublattice_probabilities 

 

@property 

def temperature(self) -> float: 

""" Current temperature """ 

return self._ensemble_parameters['temperature'] 

 

def _do_trial_step(self): 

""" Carries out one Monte Carlo trial step. """ 

sublattice_index = self.get_random_sublattice_index(self._swap_sublattice_probabilities) 

self.do_canonical_swap(sublattice_index=sublattice_index)