Coverage for mchammer/observers/site_occupancy_observer.py: 100%
38 statements
« prev ^ index » next coverage.py v7.5.0, created at 2024-12-26 04:12 +0000
1from ase import Atoms
2from icet import ClusterSpace
3from icet.core.structure import Structure
4from mchammer.observers.base_observer import BaseObserver
5from typing import List, Dict
6import numpy as np
9class SiteOccupancyObserver(BaseObserver):
10 """
11 This class represents a site occupation factor (SOF) observer.
13 A SOF observer allows to compute the site occupation factors along the
14 trajectory sampled by a Monte Carlo (MC) simulation.
16 Parameters
17 ----------
18 cluster_space
19 Cluster space from which the allowed species are extracted.
20 structure
21 Supercell consistent with primitive structure in :attr:`cluster_space`.
22 Used to determine which species are allowed on each site.
23 sites
24 Dictionary containing lists of sites that are to be considered.
25 The keys will be taken as the names of the sites; the indices refer
26 to the primitive structure associated with the cluster space.
27 interval
28 Observation interval. Defaults to ``None`` meaning that if the
29 observer is used in a Monte Carlo simulation, then the
30 :class:`Ensemble` object will determine the interval.
32 Example
33 -------
34 The following snippet illustrates how to use the site occupancy factor (SOF)
35 observer in a Monte Carlo simulation of a surface slab. Here, the SOF
36 observer is used to monitor the concentrations of different species at the
37 surface, the first subsurface layer, and the remaining 'bulk'. A minimal
38 cluster expansion is used with slightly modified surface interactions in
39 order to obtain an example that can be run without much ado. In practice,
40 one should of course use a proper cluster expansion::
42 >>> from ase.build import fcc111
43 >>> from icet import ClusterExpansion, ClusterSpace
44 >>> from mchammer.calculators import ClusterExpansionCalculator
45 >>> from mchammer.ensembles import CanonicalEnsemble
46 >>> from mchammer.observers import SiteOccupancyObserver
48 >>> # prepare reference structure
49 >>> prim = fcc111('Au', size=(1, 1, 10), vacuum=10.0)
50 >>> prim.translate((0.1, 0.1, 0.0))
51 >>> prim.wrap()
52 >>> prim.pbc = True # icet requires pbc in all directions
54 >>> # prepare cluster expansion
55 >>> cs = ClusterSpace(prim, cutoffs=[3.7], chemical_symbols=['Ag', 'Au'])
56 >>> params = [0] + 5 * [0] + 10 * [0.1]
57 >>> params[1] = 0.01
58 >>> params[6] = 0.12
59 >>> ce = ClusterExpansion(cs, params)
60 >>> print(ce)
62 >>> # prepare initial configuration based on a 2x2 supercell
63 >>> structure = prim.repeat((2, 2, 1))
64 >>> for k in range(20):
65 >>> structure[k].symbol = 'Ag'
67 >>> # set up MC simulation
68 >>> calc = ClusterExpansionCalculator(structure, ce)
69 >>> mc = CanonicalEnsemble(structure=structure, calculator=calc, temperature=600,
70 ... dc_filename='myrun_sof.dc')
72 >>> # set up observer and attach it to the MC simulation
73 >>> sites = {'surface': [0, 9], 'subsurface': [1, 8],
74 ... 'bulk': list(range(2, 8))}
75 >>> sof = SiteOccupancyObserver(cs, structure, sites, interval=len(structure))
76 >>> mc.attach_observer(sof)
78 >>> # run 1000 trial steps
79 >>> mc.run(1000)
81 After having run this snippet one can access the SOFs via the data
82 container::
84 >>> print(mc.data_container.data)
85 """
87 def __init__(self, cluster_space: ClusterSpace,
88 structure: Atoms,
89 sites: Dict[str, List[int]],
90 interval: int = None) -> None:
91 super().__init__(interval=interval, return_type=dict, tag='SiteOccupancyObserver')
93 self._sites = {site: sorted(indices)
94 for site, indices in sites.items()}
96 self._set_allowed_species(cluster_space, structure)
98 def _set_allowed_species(self,
99 cluster_space: ClusterSpace,
100 structure: Atoms):
101 """
102 Set the allowed species for the selected sites in the Atoms object
104 Parameters
105 ----------
106 cluster_space
107 Cluster space implicitly defining allowed species.
108 structure
109 Specific supercell (consistent with cluster_space) whose
110 allowed species are to be determined.
111 """
113 primitive_structure = Structure.from_atoms(cluster_space.primitive_structure)
114 chemical_symbols = cluster_space.chemical_symbols
116 if len(chemical_symbols) == 1:
117 # If the allowed species are the same for all sites no loop is
118 # required
119 allowed_species = {site: chemical_symbols[0] for
120 site in self._sites.keys()}
121 else:
122 # Loop over the lattice sites to find the allowed species
123 allowed_species = {}
124 for site, indices in self._sites.items():
125 allowed_species[site] = None
126 positions = structure.positions[np.array(indices)]
127 lattice_sites = primitive_structure.find_lattice_sites_by_positions(
128 positions=positions,
129 fractional_position_tolerance=cluster_space.fractional_position_tolerance)
130 for k, lattice_site in enumerate(lattice_sites):
131 species = chemical_symbols[lattice_site.index]
132 # check that the allowed species are equal for all sites
133 if allowed_species[site] is not None and \
134 species != allowed_species[site]:
135 raise Exception('The allowed species {} for the site'
136 ' with index {} differs from the'
137 ' result {} for the previous index'
138 ' ({})!'.format(species, indices[k],
139 allowed_species[site],
140 indices[k-1]))
141 allowed_species[site] = species
143 self._allowed_species = allowed_species
145 def get_observable(self, structure: Atoms) -> Dict[str, List[float]]:
146 """
147 Returns the site occupation factors for a given atomic configuration.
149 Parameters
150 ----------
151 structure
152 Input atomic structure.
153 """
155 chemical_symbols = np.array(structure.get_chemical_symbols())
156 sofs = {}
157 for site, indices in self._sites.items():
158 counts = {species: 0 for species in self._allowed_species[site]}
159 symbols, sym_counts = np.unique(chemical_symbols[indices],
160 return_counts=True)
161 for sym, count in zip(symbols, sym_counts):
162 counts[sym] += count
164 for species in counts.keys():
165 key = 'sof_{}_{}'.format(site, species)
166 sofs[key] = float(counts[species]) / len(indices)
168 return sofs