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

1import itertools 

2import random 

3from typing import List, Dict, Tuple, Union 

4import numpy as np 

5from mchammer.ensembles import TargetClusterVectorAnnealing 

6from mchammer.calculators import (TargetVectorCalculator, 

7 compare_cluster_vectors) 

8from icet.tools import (enumerate_structures, 

9 enumerate_supercells) 

10from icet import ClusterSpace 

11from icet.input_output.logging_tools import logger 

12from ase import Atoms 

13from ase.data import chemical_symbols as periodic_table 

14 

15 

16def generate_target_structure_from_supercells(cluster_space: ClusterSpace, 

17 supercells: List[Atoms], 

18 target_concentrations: dict, 

19 target_cluster_vector: List[float], 

20 T_start: float = 5.0, 

21 T_stop: float = 0.001, 

22 n_steps: int = None, 

23 optimality_weight: float = 1.0, 

24 random_seed: int = None, 

25 tol: float = 1e-5) -> Atoms: 

26 """ 

27 Given a ``cluster_space`` and a ``target_cluster_vector`` and one 

28 or more ``supercells``, generate a structure that as closely as 

29 possible matches that cluster vector. 

30 

31 Internally the function uses a simulated annealing algorithm and the 

32 difference between two cluster vectors is calculated with the 

33 measure suggested by A. van de Walle et al. in Calphad **42**, 13-18 

34 (2013) [WalTiwJon13]_ (for more information, see 

35 :class:`mchammer.calculators.TargetVectorCalculator`). 

36 

37 Parameters 

38 ---------- 

39 cluster_space 

40 a cluster space defining the lattice to be occupied 

41 supercells 

42 list of one or more supercells among which an optimal 

43 structure will be searched for 

44 target_concentrations 

45 concentration of each species in the target structure, per 

46 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

47 single sublattice Au-Pd structure, or 

48 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

49 for a system with two sublattices. 

50 The symbols defining sublattices ('A', 'B' etc) can be 

51 found by printing the `cluster_space` 

52 target_cluster_vector 

53 cluster vector that the generated structure should match as 

54 closely as possible 

55 T_start 

56 artificial temperature at which the simulated annealing starts 

57 T_stop 

58 artifical temperature at which the simulated annealing stops 

59 n_steps 

60 total number of Monte Carlo steps in the simulation 

61 optimality_weight 

62 controls weighting :math:`L` of perfect correlations, see 

63 :class:`mchammer.calculators.TargetVectorCalculator` 

64 random_seed 

65 seed for the random number generator used in the 

66 Monte Carlo simulation 

67 tol 

68 Numerical tolerance 

69 """ 

70 target_concentrations = _validate_concentrations(target_concentrations, cluster_space) 

71 

72 calculators = [] 

73 

74 # Loop over all supercells and intialize 

75 # them with a "random" occupation of symbols that 

76 # fulfill the target concentrations 

77 valid_supercells = [] 

78 warning_issued = False 

79 for supercell in supercells: 

80 supercell_copy = supercell.copy() 

81 try: 

82 occupy_structure_randomly(supercell_copy, cluster_space, 

83 target_concentrations) 

84 except ValueError: 

85 if not warning_issued: 

86 logger.warning('At least one supercell was not commensurate with the specified ' 

87 'target concentrations.') 

88 warning_issued = True 

89 continue 

90 valid_supercells.append(supercell_copy) 

91 calculators.append(TargetVectorCalculator(supercell_copy, cluster_space, 

92 target_cluster_vector, 

93 optimality_weight=optimality_weight, 

94 optimality_tol=tol)) 

95 

96 if len(valid_supercells) == 0: 

97 raise ValueError('No supercells that may host the specified ' 

98 'target_concentrations were supplied.') 

99 

100 ens = TargetClusterVectorAnnealing(structure=valid_supercells, calculators=calculators, 

101 T_start=T_start, T_stop=T_stop, 

102 random_seed=random_seed) 

103 return ens.generate_structure(number_of_trial_steps=n_steps) 

104 

105 

106def generate_target_structure(cluster_space: ClusterSpace, 

107 max_size: int, 

108 target_concentrations: dict, 

109 target_cluster_vector: List[float], 

110 include_smaller_cells: bool = True, 

111 pbc: Union[Tuple[bool, bool, bool], Tuple[int, int, int]] = None, 

112 T_start: float = 5.0, 

113 T_stop: float = 0.001, 

114 n_steps: int = None, 

115 optimality_weight: float = 1.0, 

116 random_seed: int = None, 

117 tol: float = 1e-5) -> Atoms: 

118 """ 

119 Given a ``cluster_space`` and a ``target_cluster_vector``, generate 

120 a structure that as closely as possible matches that cluster vector. 

121 The search is performed among all inequivalent supercells shapes up 

122 to a certain size. 

123 

124 Internally the function uses a simulated annealing algorithm and the 

125 difference between two cluster vectors is calculated with the 

126 measure suggested by A. van de Walle et al. in Calphad **42**, 13-18 

127 (2013) [WalTiwJon13]_ (for more information, see 

128 :class:`mchammer.calculators.TargetVectorCalculator`). 

129 

130 Parameters 

131 ---------- 

132 cluster_space 

133 a cluster space defining the lattice to be occupied 

134 max_size 

135 maximum supercell size 

136 target_concentrations 

137 concentration of each species in the target structure, per 

138 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

139 single sublattice Au-Pd structure, or 

140 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

141 for a system with two sublattices. 

142 The symbols defining sublattices ('A', 'B' etc) can be 

143 found by printing the `cluster_space` 

144 target_cluster_vector 

145 cluster vector that the generated structure should match as 

146 closely as possible 

147 include_smaller_cells 

148 if True, search among all supercell sizes including 

149 ``max_size``, else search only among those exactly matching 

150 ``max_size`` 

151 pbc 

152 Periodic boundary conditions for each direction, e.g., 

153 ``(True, True, False)``. The axes are defined by 

154 the cell of ``cluster_space.primitive_structure``. 

155 Default is periodic boundary in all directions. 

156 T_start 

157 artificial temperature at which the simulated annealing starts 

158 T_stop 

159 artifical temperature at which the simulated annealing stops 

160 n_steps 

161 total number of Monte Carlo steps in the simulation 

162 optimality_weight 

163 controls weighting :math:`L` of perfect correlations, see 

164 :class:`mchammer.calculators.TargetVectorCalculator` 

165 random_seed 

166 seed for the random number generator used in the 

167 Monte Carlo simulation 

168 tol 

169 Numerical tolerance 

170 """ 

171 target_concentrations = _validate_concentrations(target_concentrations, cluster_space) 

172 

173 if pbc is None: 

174 pbc = (True, True, True) 

175 

176 prim = cluster_space.primitive_structure 

177 prim.set_pbc(pbc) 

178 

179 supercells = [] 

180 if include_smaller_cells: 

181 sizes = list(range(1, max_size + 1)) 

182 else: 

183 sizes = [max_size] 

184 

185 for size in sizes: 

186 # For efficiency, make a first check that the current size is 

187 # commensurate with all concentrations (example: {'Au': 0.5, 

188 # 'Pd': 0.5} would not be commensurate with a supercell with 3 

189 # atoms). 

190 supercell = cluster_space.primitive_structure.repeat((size, 1, 1)) 

191 if not _concentrations_fit_structure(structure=supercell, 

192 cluster_space=cluster_space, 

193 concentrations=target_concentrations): 

194 continue 

195 

196 # Loop over all inequivalent supercells and intialize 

197 # them with a "random" occupation of symbols that 

198 # fulfill the target concentrations 

199 for supercell in enumerate_supercells(prim, [size]): 

200 supercell.set_pbc(True) 

201 supercells.append(supercell) 

202 

203 return generate_target_structure_from_supercells(cluster_space=cluster_space, 

204 supercells=supercells, 

205 target_concentrations=target_concentrations, 

206 target_cluster_vector=target_cluster_vector, 

207 T_start=T_start, 

208 T_stop=T_stop, 

209 n_steps=n_steps, 

210 optimality_weight=optimality_weight, 

211 random_seed=random_seed, 

212 tol=tol) 

213 

214 

215def generate_sqs_from_supercells(cluster_space: ClusterSpace, 

216 supercells: List[Atoms], 

217 target_concentrations: dict, 

218 T_start: float = 5.0, 

219 T_stop: float = 0.001, 

220 n_steps: int = None, 

221 optimality_weight: float = 1.0, 

222 random_seed: int = None, 

223 tol: float = 1e-5) -> Atoms: 

224 """ 

225 Given a ``cluster_space`` and one or more ``supercells``, generate 

226 a special quasirandom structure (SQS), i.e., a structure that for 

227 the provided supercells size provides the best possible 

228 approximation to a random alloy [ZunWeiFer90]_. 

229 

230 In the present case, this means that the generated structure will 

231 have a cluster vector that as closely as possible matches the 

232 cluster vector of an infintely large randomly occupated supercell. 

233 Internally the function uses a simulated annealing algorithm and the 

234 difference between two cluster vectors is calculated with the 

235 measure suggested by A. van de Walle et al. in Calphad **42**, 13-18 

236 (2013) [WalTiwJon13]_ (for more information, see 

237 :class:`mchammer.calculators.TargetVectorCalculator`). 

238 

239 Parameters 

240 ---------- 

241 cluster_space 

242 a cluster space defining the lattice to be occupated 

243 supercells 

244 list of one or more supercells among which an optimal 

245 structure will be searched for 

246 target_concentrations 

247 concentration of each species in the target structure, per 

248 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

249 single sublattice Au-Pd structure, or 

250 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

251 for a system with two sublattices. 

252 The symbols defining sublattices ('A', 'B' etc) can be 

253 found by printing the `cluster_space` 

254 T_start 

255 artificial temperature at which the simulated annealing starts 

256 T_stop 

257 artifical temperature at which the simulated annealing stops 

258 n_steps 

259 total number of Monte Carlo steps in the simulation 

260 optimality_weight 

261 controls weighting :math:`L` of perfect correlations, see 

262 :class:`mchammer.calculators.TargetVectorCalculator` 

263 random_seed 

264 seed for the random number generator used in the 

265 Monte Carlo simulation 

266 tol 

267 Numerical tolerance 

268 """ 

269 

270 sqs_vector = _get_sqs_cluster_vector(cluster_space=cluster_space, 

271 target_concentrations=target_concentrations) 

272 return generate_target_structure_from_supercells(cluster_space=cluster_space, 

273 supercells=supercells, 

274 target_concentrations=target_concentrations, 

275 target_cluster_vector=sqs_vector, 

276 T_start=T_start, T_stop=T_stop, 

277 n_steps=n_steps, 

278 optimality_weight=optimality_weight, 

279 random_seed=random_seed, 

280 tol=tol) 

281 

282 

283def generate_sqs(cluster_space: ClusterSpace, 

284 max_size: int, 

285 target_concentrations: dict, 

286 include_smaller_cells: bool = True, 

287 pbc: Union[Tuple[bool, bool, bool], Tuple[int, int, int]] = None, 

288 T_start: float = 5.0, 

289 T_stop: float = 0.001, 

290 n_steps: int = None, 

291 optimality_weight: float = 1.0, 

292 random_seed: int = None, 

293 tol: float = 1e-5) -> Atoms: 

294 """ 

295 Given a ``cluster_space``, generate a special quasirandom structure 

296 (SQS), i.e., a structure that for a given supercell size provides 

297 the best possible approximation to a random alloy [ZunWeiFer90]_. 

298 

299 In the present case, this means that the generated structure will 

300 have a cluster vector that as closely as possible matches the 

301 cluster vector of an infintely large randomly occupated supercell. 

302 Internally the function uses a simulated annealing algorithm and the 

303 difference between two cluster vectors is calculated with the 

304 measure suggested by A. van de Walle et al. in Calphad **42**, 13-18 

305 (2013) [WalTiwJon13]_ (for more information, see 

306 :class:`mchammer.calculators.TargetVectorCalculator`). 

307 

308 Parameters 

309 ---------- 

310 cluster_space 

311 a cluster space defining the lattice to be occupated 

312 max_size 

313 maximum supercell size 

314 target_concentrations 

315 concentration of each species in the target structure, per 

316 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

317 single sublattice Au-Pd structure, or 

318 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

319 for a system with two sublattices. 

320 The symbols defining sublattices ('A', 'B' etc) can be 

321 found by printing the `cluster_space` 

322 include_smaller_cells 

323 if True, search among all supercell sizes including 

324 ``max_size``, else search only among those exactly matching 

325 ``max_size`` 

326 pbc 

327 Periodic boundary conditions for each direction, e.g., 

328 ``(True, True, False)``. The axes are defined by 

329 the cell of ``cluster_space.primitive_structure``. 

330 Default is periodic boundary in all directions. 

331 T_start 

332 artificial temperature at which the simulated annealing starts 

333 T_stop 

334 artifical temperature at which the simulated annealing stops 

335 n_steps 

336 total number of Monte Carlo steps in the simulation 

337 optimality_weight 

338 controls weighting :math:`L` of perfect correlations, see 

339 :class:`mchammer.calculators.TargetVectorCalculator` 

340 random_seed 

341 seed for the random number generator used in the 

342 Monte Carlo simulation 

343 tol 

344 Numerical tolerance 

345 """ 

346 

347 sqs_vector = _get_sqs_cluster_vector(cluster_space=cluster_space, 

348 target_concentrations=target_concentrations) 

349 return generate_target_structure(cluster_space=cluster_space, 

350 max_size=max_size, 

351 target_concentrations=target_concentrations, 

352 target_cluster_vector=sqs_vector, 

353 include_smaller_cells=include_smaller_cells, 

354 pbc=pbc, 

355 T_start=T_start, T_stop=T_stop, 

356 n_steps=n_steps, 

357 optimality_weight=optimality_weight, 

358 random_seed=random_seed, 

359 tol=tol) 

360 

361 

362def generate_sqs_by_enumeration(cluster_space: ClusterSpace, 

363 max_size: int, 

364 target_concentrations: dict, 

365 include_smaller_cells: bool = True, 

366 pbc: Union[Tuple[bool, bool, bool], Tuple[int, int, int]] = None, 

367 optimality_weight: float = 1.0, 

368 tol: float = 1e-5) -> Atoms: 

369 """ 

370 Given a ``cluster_space``, generate a special quasirandom structure 

371 (SQS), i.e., a structure that for a given supercell size provides 

372 the best possible approximation to a random alloy [ZunWeiFer90]_. 

373 

374 In the present case, this means that the generated structure will 

375 have a cluster vector that as closely as possible matches the 

376 cluster vector of an infintely large randomly occupied supercell. 

377 Internally the function uses a simulated annealing algorithm and the 

378 difference between two cluster vectors is calculated with the 

379 measure suggested by A. van de Walle et al. in Calphad **42**, 13-18 

380 (2013) [WalTiwJon13]_ (for more information, see 

381 :class:`mchammer.calculators.TargetVectorCalculator`). 

382 

383 This functions generates SQS cells by exhaustive enumeration, which 

384 means that the generated SQS cell is guaranteed to be optimal with 

385 regard to the specified measure and cell size. 

386 

387 Parameters 

388 ---------- 

389 cluster_space 

390 a cluster space defining the lattice to be occupied 

391 max_size 

392 maximum supercell size 

393 target_concentrations 

394 concentration of each species in the target structure, per 

395 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

396 single sublattice Au-Pd structure, or 

397 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

398 for a system with two sublattices. 

399 The symbols defining sublattices ('A', 'B' etc) can be 

400 found by printing the `cluster_space` 

401 include_smaller_cells 

402 if True, search among all supercell sizes including 

403 ``max_size``, else search only among those exactly matching 

404 ``max_size`` 

405 pbc 

406 Periodic boundary conditions for each direction, e.g., 

407 ``(True, True, False)``. The axes are defined by 

408 the cell of ``cluster_space.primitive_structure``. 

409 Default is periodic boundary in all directions. 

410 optimality_weight 

411 controls weighting :math:`L` of perfect correlations, see 

412 :class:`mchammer.calculators.TargetVectorCalculator` 

413 tol 

414 Numerical tolerance 

415 """ 

416 target_concentrations = _validate_concentrations(target_concentrations, cluster_space) 

417 sqs_vector = _get_sqs_cluster_vector(cluster_space=cluster_space, 

418 target_concentrations=target_concentrations) 

419 # Translate concentrations to the format required for concentration 

420 # restricted enumeration 

421 cr = {} # type: Dict[str, tuple] 

422 sublattices = cluster_space.get_sublattices(cluster_space.primitive_structure) 

423 for sl in sublattices: 

424 mult_factor = len(sl.indices) / len(cluster_space.primitive_structure) 

425 if sl.symbol in target_concentrations: 

426 sl_conc = target_concentrations[sl.symbol] 

427 else: 

428 sl_conc = {sl.chemical_symbols[0]: 1.0} 

429 for species, value in sl_conc.items(): 

430 c = value * mult_factor 

431 if species in cr: 

432 cr[species] = (cr[species][0] + c, 

433 cr[species][1] + c) 

434 else: 

435 cr[species] = (c, c) 

436 

437 # Check to be sure... 

438 c_sum = sum(c[0] for c in cr.values()) 

439 assert abs(c_sum - 1) < tol # Should never happen, but... 

440 

441 orbit_data = cluster_space.orbit_data 

442 best_score = 1e9 

443 

444 if include_smaller_cells: 

445 sizes = list(range(1, max_size + 1)) 

446 else: 

447 sizes = [max_size] 

448 

449 # Prepare primitive structure with the right boundary conditions 

450 prim = cluster_space.primitive_structure 

451 if pbc is None: 

452 pbc = (True, True, True) 

453 prim.set_pbc(pbc) 

454 

455 # Enumerate and calculate score for each structuer 

456 for structure in enumerate_structures(prim, 

457 sizes, 

458 cluster_space.chemical_symbols, 

459 concentration_restrictions=cr): 

460 cv = cluster_space.get_cluster_vector(structure) 

461 score = compare_cluster_vectors(cv_1=cv, cv_2=sqs_vector, 

462 orbit_data=orbit_data, 

463 optimality_weight=optimality_weight, 

464 tol=tol) 

465 

466 if score < best_score: 

467 best_score = score 

468 best_structure = structure 

469 return best_structure 

470 

471 

472def occupy_structure_randomly(structure: Atoms, cluster_space: ClusterSpace, 

473 target_concentrations: dict) -> None: 

474 """ 

475 Occupy a structure with quasirandom order but fulfilling 

476 ``target_concentrations``. 

477 

478 Parameters 

479 ---------- 

480 structure 

481 ASE Atoms object that will be occupied randomly 

482 cluster_space 

483 cluster space (needed as it carries information about sublattices) 

484 target_concentrations 

485 concentration of each species in the target structure, per 

486 sublattice (for example ``{'Au': 0.5, 'Pd': 0.5}`` for a 

487 single sublattice Au-Pd structure, or 

488 ``{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}`` 

489 for a system with two sublattices. 

490 The symbols defining sublattices ('A', 'B' etc) can be 

491 found by printing the `cluster_space` 

492 """ 

493 target_concentrations = _validate_concentrations(cluster_space=cluster_space, 

494 concentrations=target_concentrations) 

495 

496 if not _concentrations_fit_structure(structure, cluster_space, target_concentrations): 

497 raise ValueError('Structure with {} atoms cannot accomodate ' 

498 'target concentrations {}'.format(len(structure), 

499 target_concentrations)) 

500 

501 symbols_all = [''] * len(structure) 

502 for sl in cluster_space.get_sublattices(structure): 

503 symbols = [] # type: List[str] # chemical_symbols in one sublattice 

504 chemical_symbols = sl.chemical_symbols 

505 if len(chemical_symbols) == 1: 

506 symbols += [chemical_symbols[0]] * len(sl.indices) 

507 else: 

508 sl_conc = target_concentrations[sl.symbol] 

509 for chemical_symbol in sl.chemical_symbols: 

510 n_symbol = int(round(len(sl.indices) * sl_conc[chemical_symbol])) 

511 symbols += [chemical_symbol] * n_symbol 

512 

513 # Should not happen but you never know 

514 assert len(symbols) == len(sl.indices) 

515 

516 # Shuffle to introduce randomness 

517 random.shuffle(symbols) 

518 

519 # Assign symbols to the right indices 

520 for symbol, lattice_site in zip(symbols, sl.indices): 

521 symbols_all[lattice_site] = symbol 

522 

523 assert symbols_all.count('') == 0 

524 structure.set_chemical_symbols(symbols_all) 

525 

526 

527def _validate_concentrations(concentrations: dict, 

528 cluster_space: ClusterSpace, 

529 tol: float = 1e-5) -> dict: 

530 """ 

531 Validates concentration specification against a cluster space 

532 (raises `ValueError` if they do not match). 

533 

534 Parameters 

535 ---------- 

536 concentrations 

537 concentration specification 

538 cluster_space 

539 cluster space to check against 

540 tol 

541 Numerical tolerance 

542 

543 Returns 

544 ------- 

545 target concentrations 

546 An adapted version of concentrations, which is always a dictionary 

547 even if there is only one sublattice 

548 """ 

549 sls = cluster_space.get_sublattices(cluster_space.primitive_structure) 

550 

551 if not isinstance(list(concentrations.values())[0], dict): 

552 concentrations = {'A': concentrations} 

553 

554 # Ensure concentrations sum to 1 at each sublattice 

555 for sl_conc in concentrations.values(): 

556 conc_sum = sum(list(sl_conc.values())) 

557 if abs(conc_sum - 1.0) > tol: 

558 raise ValueError('Concentrations must sum up ' 

559 'to 1 for each sublattice (not {})'.format(conc_sum)) 

560 

561 # Symbols need to match on each sublattice 

562 for sl in sls: 

563 if sl.symbol not in concentrations: 

564 if len(sl.chemical_symbols) > 1: 

565 raise ValueError('A sublattice ({}: {}) is missing in ' 

566 'target_concentrations'.format(sl.symbol, 

567 list(sl.chemical_symbols))) 

568 else: 

569 sl_conc = concentrations[sl.symbol] 

570 if tuple(sorted(sl.chemical_symbols)) != tuple(sorted(list(sl_conc.keys()))): 

571 raise ValueError('Chemical symbols on a sublattice ({}: {}) are ' 

572 'not the same as those in the specified ' 

573 'concentrations {}'.format(sl.symbol, list(sl.chemical_symbols), 

574 list(sl_conc.keys()))) 

575 

576 return concentrations 

577 

578 

579def _concentrations_fit_structure(structure: Atoms, 

580 cluster_space: ClusterSpace, 

581 concentrations: Dict[str, Dict[str, float]], 

582 tol: float = 1e-5) -> bool: 

583 """ 

584 Check if specified concentrations are commensurate with a 

585 certain supercell (including sublattices) 

586 

587 Parameters 

588 ---------- 

589 structure 

590 atomic configuration to be checked 

591 cluster_space 

592 corresponding cluster space 

593 concentrations 

594 which concentrations, per sublattice, e.g., ``{'A': {'Ag': 0.3, 'Au': 0.7}}`` 

595 tol 

596 numerical tolerance 

597 """ 

598 # Check that concentrations are OK in each sublattice 

599 for sublattice in cluster_space.get_sublattices(structure): 

600 if sublattice.symbol in concentrations: 

601 sl_conc = concentrations[sublattice.symbol] 

602 for conc in sl_conc.values(): 

603 n_symbol = conc * len(sublattice.indices) 

604 if abs(int(round(n_symbol)) - n_symbol) > tol: 

605 return False 

606 return True 

607 

608 

609def _get_sqs_cluster_vector(cluster_space: ClusterSpace, 

610 target_concentrations: Dict[str, Dict[str, float]]) -> np.ndarray: 

611 """ 

612 Get the SQS vector for a certain cluster space and certain 

613 concentration. Here SQS vector refers to the cluster vector of an 

614 infintely large supercell with random occupation. 

615 

616 Parameters 

617 ---------- 

618 cluster_space 

619 the kind of lattice to be occupied 

620 target_concentrations 

621 concentration of each species in the target structure, 

622 per sublattice (for example `{'A': {'Ag': 0.5, 'Pd': 0.5}}`) 

623 """ 

624 target_concentrations = _validate_concentrations(concentrations=target_concentrations, 

625 cluster_space=cluster_space) 

626 

627 sublattice_to_index = {letter: index for index, 

628 letter in enumerate('ABCDEFGHIJKLMNOPQRSTUVWXYZ')} 

629 all_sublattices = cluster_space.get_sublattices( 

630 cluster_space.primitive_structure) 

631 

632 # Make a map from chemical symbol to integer, later on used 

633 # for evaluating cluster functions. 

634 # Internally, icet sorts species according to atomic numbers. 

635 # Also check that each symbol only occurs in one sublattice. 

636 symbol_to_integer_map = {} 

637 found_species = [] # type: List[str] 

638 for sublattice in all_sublattices: 

639 if len(sublattice.chemical_symbols) < 2: 

640 continue 

641 atomic_numbers = [periodic_table.index(sym) for sym in sublattice.chemical_symbols] 

642 for i, species in enumerate(sorted(atomic_numbers)): 

643 found_species.append(species) 

644 symbol_to_integer_map[periodic_table[species]] = i 

645 

646 # Target concentrations refer to all atoms, but probabilities only 

647 # to the sublattice. 

648 probabilities = {} 

649 for sl_conc in target_concentrations.values(): 

650 if len(sl_conc) == 1: 

651 continue 

652 for symbol in sl_conc.keys(): 

653 probabilities[symbol] = sl_conc[symbol] 

654 

655 # For every orbit, calculate average cluster function 

656 cv = [1.0] 

657 for orbit in cluster_space.orbit_data: 

658 if orbit['order'] < 1: 

659 continue 

660 

661 # What sublattices are there in this orbit? 

662 sublattices = [all_sublattices[sublattice_to_index[letter]] 

663 for letter in orbit['sublattices'].split('-')] 

664 

665 # What chemical symbols do these sublattices refer to? 

666 symbol_groups = [sublattice.chemical_symbols for sublattice in sublattices] 

667 

668 # How many allowed species in each of those sublattices? 

669 nbr_of_allowed_species = [len(symbol_group) 

670 for symbol_group in symbol_groups] 

671 

672 # Calculate contribution from every possible combination of 

673 # symbols weighted with their probability 

674 cluster_product_average = 0 

675 for symbols in itertools.product(*symbol_groups): 

676 cluster_product = 1 

677 for i, symbol in enumerate(symbols): 

678 mc_vector_component = orbit['multi_component_vector'][i] 

679 species_i = symbol_to_integer_map[symbol] 

680 prod = cluster_space.evaluate_cluster_function(nbr_of_allowed_species[i], 

681 mc_vector_component, 

682 species_i) 

683 cluster_product *= probabilities[symbol] * prod 

684 cluster_product_average += cluster_product 

685 cv.append(cluster_product_average) 

686 return np.array(cv)