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"""Definition of the canonical annealing class.""" 

2 

3from ase import Atoms 

4from ase.data import chemical_symbols 

5from ase.units import kB 

6from typing import Dict, List, Any 

7 

8from .. import DataContainer 

9from ..calculators import ClusterExpansionCalculator 

10from .thermodynamic_base_ensemble import ThermodynamicBaseEnsemble 

11from .semi_grand_canonical_ensemble import get_chemical_potentials 

12from .canonical_annealing import available_cooling_functions 

13 

14 

15class SGCAnnealing(ThermodynamicBaseEnsemble): 

16 """Instances of this class allow one to carry out simulated annealing 

17 in the semi grand canonical ensemble, i.e. the temperature is varied in 

18 pre-defined fashion while the chemical potential is kept fixed. See 

19 :class:`mchammer.ensembles.SemiGrandCanonicalEnsemble` for more information 

20 about the ensemble. 

21 

22 The temperature control scheme is selected via the 

23 ``cooling_function`` keyword argument, while the initial and final 

24 temperature are set via the ``T_start`` and ``T_stop`` arguments. 

25 Several pre-defined temperature control schemes are available 

26 including `linear` and `exponential`. In the latter case the 

27 temperature varies logarithmatically as a function of the MC step, 

28 emulating the exponential temperature dependence of the atomic 

29 exchange rate encountered in many materials. It is also possible 

30 to provide a user defined cooling function via the keyword 

31 argument. This function must comply with the following function 

32 header:: 

33 

34 def cooling_function(step, T_start, T_stop, n_steps): 

35 T = ... # compute temperature 

36 return T 

37 

38 Here ``step`` refers to the current MC trial step. 

39 

40 Parameters 

41 ---------- 

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

43 atomic configuration to be used in the Monte Carlo simulation; 

44 also defines the initial occupation vector 

45 chemical_potentials : Dict[str, float] 

46 chemical potential for each species :math:`\\mu_i`; the key 

47 denotes the species, the value specifies the chemical potential in 

48 units that are consistent with the underlying cluster expansion 

49 calculator : :class:`ClusterExpansionCalculator 

50 <mchammer.calculators.ClusterExpansionCalculator>` 

51 calculator to be used for calculating the potential changes 

52 that enter the evaluation of the Metropolis criterion 

53 T_start : float 

54 temperature from which the annealing is started 

55 T_stop : float 

56 final temperature for annealing 

57 n_steps : int 

58 number of steps to take in the annealing simulation 

59 cooling_function : str/function 

60 to use the predefined cooling functions provide a string 

61 `linear` or `exponential`, otherwise provide a function 

62 boltzmann_constant : float 

63 Boltzmann constant :math:`k_B` in appropriate 

64 units, i.e. units that are consistent 

65 with the underlying cluster expansion 

66 and the temperature units [default: eV/K] 

67 user_tag : str 

68 human-readable tag for ensemble [default: None] 

69 random_seed : int 

70 seed for the random number generator used in the Monte Carlo 

71 simulation 

72 dc_filename : str 

73 name of file the data container associated with the ensemble 

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

75 data container will be appended, and the file will be 

76 updated/overwritten 

77 data_container_write_period : float 

78 period in units of seconds at which the data container is 

79 written to file; writing periodically to file provides both 

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

81 the data [default: 600 s] 

82 ensemble_data_write_interval : int 

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

84 includes for example the current value of the calculator 

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

86 such as temperature or the number of atoms of different species 

87 trajectory_write_interval : int 

88 interval at which the current occupation vector of the atomic 

89 configuration is written to the data container. 

90 sublattice_probabilities : List[float] 

91 probability for picking a sublattice when doing a random swap. 

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

93 sum up to 1. 

94 

95 """ 

96 

97 def __init__(self, 

98 structure: Atoms, 

99 calculator: ClusterExpansionCalculator, 

100 T_start: float, 

101 T_stop: float, 

102 n_steps: int, 

103 chemical_potentials: Dict[str, float], 

104 cooling_function: str = 'exponential', 

105 boltzmann_constant: float = kB, 

106 user_tag: str = None, 

107 random_seed: int = None, 

108 dc_filename: str = None, 

109 data_container: str = None, 

110 data_container_write_period: float = 600, 

111 ensemble_data_write_interval: int = None, 

112 trajectory_write_interval: int = None, 

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

114 

115 self._ensemble_parameters = dict(n_steps=n_steps) # type: Dict[str, Any] 

116 

117 self._chemical_potentials = get_chemical_potentials(chemical_potentials) 

118 

119 # add chemical potentials to ensemble parameters 

120 self._chemical_potentials = get_chemical_potentials(chemical_potentials) 

121 for atnum, chempot in self.chemical_potentials.items(): 

122 mu_sym = 'mu_{}'.format(chemical_symbols[atnum]) 

123 self._ensemble_parameters[mu_sym] = chempot 

124 

125 super().__init__( 

126 structure=structure, 

127 calculator=calculator, 

128 user_tag=user_tag, 

129 random_seed=random_seed, 

130 dc_filename=dc_filename, 

131 data_container=data_container, 

132 data_container_class=DataContainer, 

133 data_container_write_period=data_container_write_period, 

134 ensemble_data_write_interval=ensemble_data_write_interval, 

135 trajectory_write_interval=trajectory_write_interval, 

136 boltzmann_constant=boltzmann_constant) 

137 

138 self._temperature = T_start 

139 self._T_start = T_start 

140 self._T_stop = T_stop 

141 self._n_steps = n_steps 

142 

143 self._ground_state_candidate = self.configuration.structure 

144 self._ground_state_candidate_potential = calculator.calculate_total( 

145 occupations=self.configuration.occupations) 

146 

147 # setup cooling function 

148 if isinstance(cooling_function, str): 

149 available = sorted(available_cooling_functions.keys()) 

150 if cooling_function not in available: 150 ↛ 151line 150 didn't jump to line 151, because the condition on line 150 was never true

151 raise ValueError( 

152 'Select from the available cooling_functions {}'.format(available)) 

153 self._cooling_function = available_cooling_functions[cooling_function] 

154 elif callable(cooling_function): 154 ↛ 157line 154 didn't jump to line 157, because the condition on line 154 was never false

155 self._cooling_function = cooling_function 

156 else: 

157 raise TypeError('cooling_function must be either str or a function') 

158 

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

160 self._flip_sublattice_probabilities = self._get_flip_sublattice_probabilities() 

161 else: 

162 self._flip_sublattice_probabilities = sublattice_probabilities 

163 

164 @property 

165 def chemical_potentials(self) -> Dict[int, float]: 

166 """ 

167 chemical potentials :math:`\\mu_i` (see parameters section above) 

168 """ 

169 return self._chemical_potentials 

170 

171 @property 

172 def temperature(self) -> float: 

173 """ Current temperature """ 

174 return self._temperature 

175 

176 @property 

177 def estimated_ground_state(self): 

178 """ Structure with lowest observed potential during run """ 

179 return self._ground_state_candidate.copy() 

180 

181 @property 

182 def estimated_ground_state_potential(self): 

183 """ Lowest observed potential during run """ 

184 return self._ground_state_candidate_potential 

185 

186 def run(self): 

187 """ Runs the annealing. """ 

188 if self.step >= self._n_steps: 

189 raise Exception('Annealing has already finished') 

190 super().run(self._n_steps - self.step) 

191 

192 def _do_trial_step(self) -> int: 

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

194 self._temperature = self._cooling_function( 

195 self.step, self._T_start, self._T_stop, self._n_steps) 

196 sublattice_index = self.get_random_sublattice_index(self._flip_sublattice_probabilities) 

197 return self.do_sgc_flip( 

198 sublattice_index=sublattice_index, chemical_potentials=self.chemical_potentials) 

199 

200 def _get_ensemble_data(self) -> Dict: 

201 """Returns the data associated with the ensemble. For the 

202 CanonicalAnnealing this specifically includes the temperature. 

203 """ 

204 data = super()._get_ensemble_data() 

205 

206 # species counts 

207 data.update(self._get_species_counts()) 

208 

209 data['temperature'] = self.temperature 

210 if data['potential'] < self._ground_state_candidate_potential: 

211 self._ground_state_candidate_potential = data['potential'] 

212 self._ground_state_candidate = self.configuration.structure 

213 return data