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 

3import numpy as np 

4 

5from ase import Atoms 

6from ase.units import kB 

7from typing import Dict, List 

8 

9from .. import DataContainer 

10from ..calculators import ClusterExpansionCalculator 

11from .thermodynamic_base_ensemble import ThermodynamicBaseEnsemble 

12 

13 

14class CanonicalAnnealing(ThermodynamicBaseEnsemble): 

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

16 in the canonical ensemble, i.e. the temperature is varied in 

17 pre-defined fashion while the composition is kept fixed. See 

18 :class:`mchammer.ensembles.CanonicalEnsemble` for more information 

19 about the standard canonical ensemble. 

20 

21 The canonical annealing ensemble can be useful, for example, for 

22 finding ground states or generating low energy configurations. 

23 

24 The temperature control scheme is selected via the 

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

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

27 Several pre-defined temperature control schemes are available 

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

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

30 emulating the exponential temperature dependence of the atomic 

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

32 to provide a user defined cooling function via the keyword 

33 argument. This function must comply with the following function 

34 header:: 

35 

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

37 T = ... # compute temperature 

38 return T 

39 

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

41 

42 Parameters 

43 ---------- 

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

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

46 also defines the initial occupation vector 

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

48 calculator to be used for calculating the potential changes 

49 that enter the evaluation of the Metropolis criterion 

50 T_start : float 

51 temperature from which the annealing is started 

52 T_stop : float 

53 final temperature for annealing 

54 n_steps : int 

55 number of steps to take in the annealing simulation 

56 cooling_function : str/function 

57 to use the predefined cooling functions provide a string 

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

59 boltzmann_constant : float 

60 Boltzmann constant :math:`k_B` in appropriate 

61 units, i.e. units that are consistent 

62 with the underlying cluster expansion 

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

64 user_tag : str 

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

66 random_seed : int 

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

68 simulation 

69 dc_filename : str 

70 name of file the data container associated with the ensemble 

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

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

73 updated/overwritten 

74 data_container_write_period : float 

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

76 written to file; writing periodically to file provides both 

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

78 the data [default: 600 s] 

79 ensemble_data_write_interval : int 

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

81 includes for example the current value of the calculator 

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

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

84 trajectory_write_interval : int 

85 interval at which the current occupation vector of the atomic 

86 configuration is written to the data container. 

87 sublattice_probabilities : List[float] 

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

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

90 sum up to 1. 

91 

92 """ 

93 

94 def __init__(self, 

95 structure: Atoms, 

96 calculator: ClusterExpansionCalculator, 

97 T_start: float, 

98 T_stop: float, 

99 n_steps: int, 

100 cooling_function: str = 'exponential', 

101 user_tag: str = None, 

102 boltzmann_constant: float = kB, 

103 random_seed: int = None, 

104 dc_filename: str = None, 

105 data_container_write_period: float = 600, 

106 ensemble_data_write_interval: int = None, 

107 trajectory_write_interval: int = None, 

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

109 

110 self._ensemble_parameters = dict(n_steps=n_steps) 

111 

112 # add species count to ensemble parameters 

113 for sl in calculator.sublattices: 

114 for symbol in sl.chemical_symbols: 

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

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

117 self._ensemble_parameters[key] = count 

118 

119 super().__init__( 

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

121 random_seed=random_seed, 

122 dc_filename=dc_filename, 

123 data_container_class=DataContainer, 

124 data_container_write_period=data_container_write_period, 

125 ensemble_data_write_interval=ensemble_data_write_interval, 

126 trajectory_write_interval=trajectory_write_interval, 

127 boltzmann_constant=boltzmann_constant) 

128 

129 self._temperature = T_start 

130 self._T_start = T_start 

131 self._T_stop = T_stop 

132 self._n_steps = n_steps 

133 

134 self._ground_state_candidate = self.configuration.structure 

135 self._ground_state_candidate_potential = calculator.calculate_total( 

136 occupations=self.configuration.occupations) 

137 

138 # setup cooling function 

139 if isinstance(cooling_function, str): 

140 available = sorted(available_cooling_functions.keys()) 

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

142 raise ValueError( 

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

144 self._cooling_function = available_cooling_functions[cooling_function] 

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

146 self._cooling_function = cooling_function 

147 else: 

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

149 

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

151 self._swap_sublattice_probabilities = self._get_swap_sublattice_probabilities() 

152 else: 

153 self._swap_sublattice_probabilities = sublattice_probabilities 

154 

155 @property 

156 def temperature(self) -> float: 

157 """ Current temperature """ 

158 return self._temperature 

159 

160 @property 

161 def T_start(self) -> float: 

162 """ Starting temperature """ 

163 return self._T_start 

164 

165 @property 

166 def T_stop(self) -> float: 

167 """ Starting temperature """ 

168 return self._T_stop 

169 

170 @property 

171 def n_steps(self) -> int: 

172 """ Number of steps to carry out """ 

173 return self._n_steps 

174 

175 @property 

176 def estimated_ground_state(self): 

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

178 return self._ground_state_candidate.copy() 

179 

180 @property 

181 def estimated_ground_state_potential(self): 

182 """ Lowest observed potential during run """ 

183 return self._ground_state_candidate_potential 

184 

185 def run(self): 

186 """ Runs the annealing. """ 

187 if self.step >= self.n_steps: 

188 raise Exception('Annealing has already finished') 

189 super().run(self.n_steps - self.step) 

190 

191 def _do_trial_step(self) -> int: 

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

193 self._temperature = self._cooling_function( 

194 self.step, self.T_start, self.T_stop, self.n_steps) 

195 sublattice_index = self.get_random_sublattice_index(self._swap_sublattice_probabilities) 

196 return self.do_canonical_swap(sublattice_index=sublattice_index) 

197 

198 def _get_ensemble_data(self) -> Dict: 

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

200 CanonicalAnnealing this specifically includes the temperature. 

201 """ 

202 data = super()._get_ensemble_data() 

203 data['temperature'] = self.temperature 

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

205 self._ground_state_candidate_potential = data['potential'] 

206 self._ground_state_candidate = self.configuration.structure 

207 return data 

208 

209 

210def _cooling_linear(step, T_start, T_stop, n_steps): 

211 return T_start + (T_stop-T_start) * step / (n_steps - 1) 

212 

213 

214def _cooling_exponential(step, T_start, T_stop, n_steps): 

215 return T_start - (T_start - T_stop) * np.log(step+1) / np.log(n_steps) 

216 

217 

218available_cooling_functions = dict(linear=_cooling_linear, exponential=_cooling_exponential)