Coverage for mchammer/observers/binary_short_range_order_observer.py: 94%
70 statements
« prev ^ index » next coverage.py v7.5.0, created at 2024-05-06 04:14 +0000
1from collections import namedtuple
2from typing import Dict
4import numpy as np
6from ase import Atoms
7from icet import ClusterSpace
8from mchammer.observers import ClusterCountObserver
9from mchammer.observers.base_observer import BaseObserver
11ClusterCountInfo = namedtuple('ClusterCountInfo', ['counts', 'dc_tags'])
14class BinaryShortRangeOrderObserver(BaseObserver):
15 """
16 This class represents a short range order (SRO) observer for a
17 binary system.
19 Parameters
20 ----------
21 cluster_space
22 Cluster space used for initialization.
23 structure
24 Defines the lattice which the observer will work on.
25 radius
26 The maximum radius for the neigbhor shells considered.
27 interval
28 Observation interval. Defaults to ``None`` meaning that if the
29 observer is used in a Monte Carlo simulations, then the :class:`Ensemble` object
30 will determine the interval.
32 Example
33 -------
34 The following snippet illustrates how to use the short-range order (SRO)
35 observer in a Monte Carlo simulation of a bulk supercell. Here, the
36 parameters of the cluster expansion are set to emulate a simple Ising model
37 in order to obtain an example that can be run without modification. In
38 practice, one should of course use a proper cluster expansion::
40 >>> from ase.build import bulk
41 >>> from icet import ClusterExpansion, ClusterSpace
42 >>> from mchammer.calculators import ClusterExpansionCalculator
43 >>> from mchammer.ensembles import CanonicalEnsemble
44 >>> from mchammer.observers import BinaryShortRangeOrderObserver
46 >>> # prepare cluster expansion
47 >>> # the setup emulates a second nearest-neighbor (NN) Ising model
48 >>> # (zerolet and singlet ECIs are zero; only first and second neighbor
49 >>> # pairs are included)
50 >>> prim = bulk('Au')
51 >>> cs = ClusterSpace(prim, cutoffs=[4.3], chemical_symbols=['Ag', 'Au'])
52 >>> ce = ClusterExpansion(cs, [0, 0, 0.1, -0.02])
54 >>> # prepare initial configuration
55 >>> nAg = 10
56 >>> structure = prim.repeat(3)
57 >>> structure.set_chemical_symbols(nAg * ['Ag'] + (len(structure) - nAg) * ['Au'])
59 >>> # set up MC simulation
60 >>> calc = ClusterExpansionCalculator(structure, ce)
61 >>> mc = CanonicalEnsemble(structure=structure, calculator=calc, temperature=600,
62 ... dc_filename='myrun_sro.dc')
64 >>> # set up observer and attach it to the MC simulation
65 >>> sro = BinaryShortRangeOrderObserver(cs, structure, interval=len(structure),
66 radius=4.3)
67 >>> mc.attach_observer(sro)
69 >>> # run 1000 trial steps
70 >>> mc.run(1000)
72 After having run this snippet one can access the SRO parameters via the
73 data container::
75 >>> print(mc.data_container.data)
76 """
78 def __init__(self, cluster_space, structure: Atoms,
79 radius: float, interval: int = None) -> None:
80 super().__init__(interval=interval, return_type=dict,
81 tag='BinaryShortRangeOrderObserver')
83 self._structure = structure
85 self._cluster_space = ClusterSpace(
86 structure=cluster_space.primitive_structure,
87 cutoffs=[radius],
88 chemical_symbols=cluster_space.chemical_symbols)
89 self._cluster_count_observer = ClusterCountObserver(
90 cluster_space=self._cluster_space, structure=structure,
91 interval=interval)
93 self._sublattices = self._cluster_space.get_sublattices(structure)
94 binary_sublattice_counts = 0
95 for symbols in self._sublattices.allowed_species:
96 if len(symbols) == 2:
97 binary_sublattice_counts += 1
98 self._symbols = sorted(symbols)
99 elif len(symbols) > 2:
100 raise ValueError('Cluster space has more than two allowed'
101 ' species on a sublattice. '
102 'Allowed species: {}'.format(symbols))
103 if binary_sublattice_counts != 1:
104 raise ValueError('Number of binary sublattices must equal one,'
105 ' not {}'.format(binary_sublattice_counts))
107 def get_observable(self, structure: Atoms) -> Dict[str, float]:
108 """Returns the value of the property from a cluster expansion
109 model for a given atomic configurations.
111 Parameters
112 ----------
113 structure
114 Input atomic structure.
115 """
117 df = self._cluster_count_observer.get_cluster_counts(structure)
119 symbol_counts = self._get_atom_count(structure)
120 conc_B = self._get_concentrations(structure)[self._symbols[0]]
122 pair_orbit_indices = set(
123 df.loc[df['order'] == 2]['orbit_index'].tolist())
124 N = symbol_counts[self._symbols[0]] + symbol_counts[self._symbols[1]]
125 sro_parameters = {}
126 for k, orbit_index in enumerate(sorted(pair_orbit_indices)):
127 orbit_df = df.loc[df['orbit_index'] == orbit_index]
128 A_B_pair_count = 0
129 total_count = 0
130 total_A_count = 0
131 for i, row in orbit_df.iterrows():
132 total_count += row.cluster_count
133 if self._symbols[0] in row.occupation:
134 total_A_count += row.cluster_count
135 if self._symbols[0] in row.occupation and \
136 self._symbols[1] in row.occupation:
137 A_B_pair_count += row.cluster_count
139 key = 'sro_{}_{}'.format(self._symbols[0], k+1)
140 Z_tot = symbol_counts[self._symbols[0]] * 2 * total_count / N
141 if conc_B == 1 or Z_tot == 0: 141 ↛ 142line 141 didn't jump to line 142, because the condition on line 141 was never true
142 value = 0
143 else:
144 value = 1 - A_B_pair_count/(Z_tot * (1-conc_B))
145 sro_parameters[key] = value
147 return sro_parameters
149 def _get_concentrations(self, structure: Atoms) -> Dict[str, float]:
150 """Returns concentrations for each species relative its sublattice.
152 Parameters
153 ----------
154 structure
155 The configuration to be analyzed.
156 """
157 occupation = np.array(structure.get_chemical_symbols())
158 concentrations = {}
159 for sublattice in self._sublattices:
160 if len(sublattice.chemical_symbols) == 1: 160 ↛ 161line 160 didn't jump to line 161, because the condition on line 160 was never true
161 continue
162 for symbol in sublattice.chemical_symbols:
163 symbol_count = occupation[sublattice.indices].tolist().count(
164 symbol)
165 concentration = symbol_count / len(sublattice.indices)
166 concentrations[symbol] = concentration
167 return concentrations
169 def _get_atom_count(self, structure: Atoms) -> Dict[str, float]:
170 """Returns atom counts for each species relative its sublattice.
172 Parameters
173 ----------
174 structure
175 The configuration to be analyzed.
176 """
177 occupation = np.array(structure.get_chemical_symbols())
178 counts = {}
179 for sublattice in self._sublattices:
180 if len(sublattice.chemical_symbols) == 1: 180 ↛ 181line 180 didn't jump to line 181, because the condition on line 180 was never true
181 continue
182 for symbol in sublattice.chemical_symbols:
183 symbol_count = occupation[sublattice.indices].tolist().count(
184 symbol)
185 counts[symbol] = symbol_count
186 return counts