Coverage for mchammer/ensembles/sgc_annealing.py: 92%

62 statements  

« prev     ^ index     » next       coverage.py v7.5.0, created at 2024-12-26 04:12 +0000

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 r"""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 :attr:`cooling_function` keyword argument, while the initial and final 

24 temperature are set via the :attr:`T_start` and :attr:`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 :attr:`step` refers to the current MC trial step. 

39 

40 Parameters 

41 ---------- 

42 structure 

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

44 also defines the initial occupation vector. 

45 chemical_potentials 

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 

50 Calculator to be used for calculating the potential changes 

51 that enter the evaluation of the Metropolis criterion. 

52 T_start 

53 Temperature from which the annealing is started. 

54 T_stop 

55 Final temperature for annealing. 

56 n_steps 

57 Number of steps to take in the annealing simulation. 

58 cooling_function 

59 to use the predefined cooling functions provide a string 

60 ``'linear'`` or ``'exponential'``, otherwise provide a function. 

61 boltzmann_constant 

62 Boltzmann constant :math:`k_B` in appropriate 

63 units, i.e. units that are consistent 

64 with the underlying cluster expansion 

65 and the temperature units. Default: eV/K. 

66 user_tag 

67 Human-readable tag for ensemble. Default: ``None``. 

68 random_seed 

69 Seed for the random number generator used in the Monte Carlo simulation. 

70 dc_filename 

71 Name of file the data container associated with the ensemble 

72 will be written to. If the file exists it will be read, the 

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

74 updated/overwritten. 

75 data_container_write_period 

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

77 written to file. Writing periodically to file provides both 

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

79 the data. Default: 600 s. 

80 ensemble_data_write_interval 

81 Interval at which data is written to the data container. This 

82 includes for example the current value of the calculator 

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

84 such as temperature or the number of atoms of different species. 

85 Default: Number of sites in the :attr:`structure`. 

86 trajectory_write_interval 

87 Interval at which the current occupation vector of the atomic 

88 configuration is written to the data container. 

89 Default: Number of sites in the :attr:`structure`. 

90 sublattice_probabilities 

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[str, Any] = dict(n_steps=n_steps) 

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 r""" 

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 

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