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