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

1""" 

2This module provides the ClusterSpace class. 

3""" 

4 

5import copy 

6import itertools 

7import pickle 

8import tarfile 

9import tempfile 

10from collections import OrderedDict 

11from math import log10, floor 

12from typing import List, Union, Tuple 

13 

14import numpy as np 

15import spglib 

16 

17from _icet import ClusterSpace as _ClusterSpace 

18from ase import Atoms 

19from ase.io import read as ase_read 

20from ase.io import write as ase_write 

21from icet.core.orbit_list import OrbitList 

22from icet.core.structure import Structure 

23from icet.core.sublattices import Sublattices 

24from icet.tools.geometry import (ase_atoms_to_spglib_cell, 

25 get_occupied_primitive_structure, 

26 get_position_from_lattice_site) 

27 

28 

29class ClusterSpace(_ClusterSpace): 

30 """This class provides functionality for generating and maintaining 

31 cluster spaces. 

32 

33 **Note:** In icet all :class:`ase.Atoms` objects must have 

34 periodic boundary conditions. When carrying out cluster expansions 

35 for surfaces and nanoparticles it is therefore recommended to 

36 surround the structure with vacuum and use periodic boundary 

37 conditions. This can be done using e.g., :func:`ase.Atoms.center`. 

38 

39 Parameters 

40 ---------- 

41 structure : ase.Atoms 

42 atomic structure 

43 cutoffs : list(float) 

44 cutoff radii per order that define the cluster space 

45 

46 Cutoffs are specified in units of Angstrom and refer to the 

47 longest distance between two atoms in the cluster. The first 

48 element refers to pairs, the second to triplets, the third 

49 to quadruplets, and so on. ``cutoffs=[7.0, 4.5]`` thus implies 

50 that all pairs distanced 7 A or less will be included, 

51 as well as all triplets among which the longest distance is no 

52 longer than 4.5 A. 

53 chemical_symbols : list(str) or list(list(str)) 

54 list of chemical symbols, each of which must map to an element 

55 of the periodic table 

56 

57 If a list of chemical symbols is provided, all sites on the 

58 lattice will have the same allowed occupations as the input 

59 list. 

60 

61 If a list of list of chemical symbols is provided then the 

62 outer list must be the same length as the `structure` object and 

63 ``chemical_symbols[i]`` will correspond to the allowed species 

64 on lattice site ``i``. 

65 symprec : float 

66 tolerance imposed when analyzing the symmetry using spglib 

67 position_tolerance : float 

68 tolerance applied when comparing positions in Cartesian coordinates 

69 

70 Examples 

71 -------- 

72 The following snippets illustrate several common situations:: 

73 

74 >>> from ase.build import bulk 

75 >>> from ase.io import read 

76 >>> from icet import ClusterSpace 

77 

78 >>> # AgPd alloy with pairs up to 7.0 A and triplets up to 4.5 A 

79 >>> prim = bulk('Ag') 

80 >>> cs = ClusterSpace(structure=prim, cutoffs=[7.0, 4.5], 

81 ... chemical_symbols=[['Ag', 'Pd']]) 

82 >>> print(cs) 

83 

84 >>> # (Mg,Zn)O alloy on rocksalt lattice with pairs up to 8.0 A 

85 >>> prim = bulk('MgO', crystalstructure='rocksalt', a=6.0) 

86 >>> cs = ClusterSpace(structure=prim, cutoffs=[8.0], 

87 ... chemical_symbols=[['Mg', 'Zn'], ['O']]) 

88 >>> print(cs) 

89 

90 >>> # (Ga,Al)(As,Sb) alloy with pairs, triplets, and quadruplets 

91 >>> prim = bulk('GaAs', crystalstructure='zincblende', a=6.5) 

92 >>> cs = ClusterSpace(structure=prim, cutoffs=[7.0, 6.0, 5.0], 

93 ... chemical_symbols=[['Ga', 'Al'], ['As', 'Sb']]) 

94 >>> print(cs) 

95 

96 >>> # PdCuAu alloy with pairs and triplets 

97 >>> prim = bulk('Pd') 

98 >>> cs = ClusterSpace(structure=prim, cutoffs=[7.0, 5.0], 

99 ... chemical_symbols=[['Au', 'Cu', 'Pd']]) 

100 >>> print(cs) 

101 

102 """ 

103 

104 def __init__(self, 

105 structure: Atoms, 

106 cutoffs: List[float], 

107 chemical_symbols: Union[List[str], List[List[str]]], 

108 symprec: float = 1e-5, 

109 position_tolerance: float = None) -> None: 

110 

111 if not isinstance(structure, Atoms): 111 ↛ 112line 111 didn't jump to line 112, because the condition on line 111 was never true

112 raise TypeError('Input configuration must be an ASE Atoms object' 

113 ', not type {}'.format(type(structure))) 

114 if not all(structure.pbc): 

115 raise ValueError('Input structure must have periodic boundary conditions') 

116 if symprec <= 0: 

117 raise ValueError('symprec must be a positive number') 

118 

119 self._config = {'symprec': symprec} 

120 self._cutoffs = cutoffs.copy() 

121 self._input_structure = structure.copy() 

122 self._input_chemical_symbols = copy.deepcopy(chemical_symbols) 

123 chemical_symbols = self._get_chemical_symbols() 

124 

125 self._pruning_history = [] 

126 

127 # set up primitive 

128 occupied_primitive, primitive_chemical_symbols = get_occupied_primitive_structure( 

129 self._input_structure, chemical_symbols, symprec=self.symprec) 

130 self._primitive_chemical_symbols = primitive_chemical_symbols 

131 assert len(occupied_primitive) == len(primitive_chemical_symbols) 

132 

133 # derived tolerances 

134 if position_tolerance is None: 

135 self._config['position_tolerance'] = symprec 

136 else: 

137 if position_tolerance <= 0: 

138 raise ValueError('position_tolerance must be a positive number') 

139 self._config['position_tolerance'] = position_tolerance 

140 effective_box_size = abs(np.linalg.det(occupied_primitive.cell)) ** (1 / 3) 

141 tol = self.position_tolerance / effective_box_size 

142 tol = min(tol, self._config['position_tolerance'] / 5) 

143 self._config['fractional_position_tolerance'] = round(tol, -int(floor(log10(abs(tol))))) 

144 

145 # set up orbit list 

146 self._orbit_list = OrbitList( 

147 structure=occupied_primitive, 

148 cutoffs=self._cutoffs, 

149 symprec=self.symprec, 

150 position_tolerance=self.position_tolerance, 

151 fractional_position_tolerance=self.fractional_position_tolerance) 

152 self._orbit_list.remove_inactive_orbits(primitive_chemical_symbols) 

153 

154 # call (base) C++ constructor 

155 _ClusterSpace.__init__(self, 

156 chemical_symbols=primitive_chemical_symbols, 

157 orbit_list=self._orbit_list, 

158 position_tolerance=self.position_tolerance, 

159 fractional_position_tolerance=self.fractional_position_tolerance) 

160 

161 def _get_chemical_symbols(self): 

162 """ Returns chemical symbols using input structure and 

163 chemical symbols. Carries out multiple sanity checks. """ 

164 

165 # setup chemical symbols as List[List[str]] 

166 if all(isinstance(i, str) for i in self._input_chemical_symbols): 

167 chemical_symbols = [ 

168 self._input_chemical_symbols] * len(self._input_structure) 

169 elif not all(isinstance(i, list) for i in self._input_chemical_symbols): 

170 raise TypeError("chemical_symbols must be List[str] or List[List[str]], not {}".format( 

171 type(self._input_chemical_symbols))) 

172 elif len(self._input_chemical_symbols) != len(self._input_structure): 

173 msg = 'chemical_symbols must have same length as structure. ' 

174 msg += 'len(chemical_symbols) = {}, len(structure)= {}'.format( 

175 len(self._input_chemical_symbols), len(self._input_structure)) 

176 raise ValueError(msg) 

177 else: 

178 chemical_symbols = copy.deepcopy(self._input_chemical_symbols) 

179 

180 for i, symbols in enumerate(chemical_symbols): 

181 if len(symbols) != len(set(symbols)): 

182 raise ValueError( 

183 'Found duplicates of allowed chemical symbols on site {}.' 

184 ' allowed species on site {}= {}'.format(i, i, symbols)) 

185 

186 if len([tuple(sorted(s)) for s in chemical_symbols if len(s) > 1]) == 0: 

187 raise ValueError('No active sites found') 

188 

189 return chemical_symbols 

190 

191 def _get_chemical_symbol_representation(self): 

192 """Returns a str version of the chemical symbols that is 

193 easier on the eyes. 

194 """ 

195 sublattices = self.get_sublattices(self.primitive_structure) 

196 nice_str = [] 

197 for sublattice in sublattices.active_sublattices: 

198 sublattice_symbol = sublattice.symbol 

199 

200 nice_str.append('{} (sublattice {})'.format( 

201 list(sublattice.chemical_symbols), sublattice_symbol)) 

202 return ', '.join(nice_str) 

203 

204 def _get_string_representation(self, 

205 print_threshold: int = None, 

206 print_minimum: int = 10) -> str: 

207 """ 

208 String representation of the cluster space that provides an overview of 

209 the orbits (order, radius, multiplicity etc) that constitute the space. 

210 

211 Parameters 

212 ---------- 

213 print_threshold 

214 if the number of orbits exceeds this number print dots 

215 print_minimum 

216 number of lines printed from the top and the bottom of the orbit 

217 list if `print_threshold` is exceeded 

218 

219 Returns 

220 ------- 

221 multi-line string 

222 string representation of the cluster space. 

223 """ 

224 

225 def repr_orbit(orbit, header=False): 

226 formats = {'order': '{:2}', 

227 'radius': '{:8.4f}', 

228 'multiplicity': '{:4}', 

229 'index': '{:4}', 

230 'orbit_index': '{:4}', 

231 'multi_component_vector': '{:}', 

232 'sublattices': '{:}'} 

233 s = [] 

234 for name, value in orbit.items(): 

235 str_repr = formats[name].format(value) 

236 n = max(len(name), len(str_repr)) 

237 if header: 

238 s += ['{s:^{n}}'.format(s=name, n=n)] 

239 else: 

240 s += ['{s:^{n}}'.format(s=str_repr, n=n)] 

241 return ' | '.join(s) 

242 

243 # basic information 

244 # (use largest orbit to obtain maximum line length) 

245 prototype_orbit = self.orbit_data[-1] 

246 width = len(repr_orbit(prototype_orbit)) 

247 s = [] # type: List 

248 s += ['{s:=^{n}}'.format(s=' Cluster Space ', n=width)] 

249 s += [' {:38} : {}'.format('space group', self.space_group)] 

250 s += [' {:38} : {}' 

251 .format('chemical species', self._get_chemical_symbol_representation())] 

252 s += [' {:38} : {}'.format('cutoffs', ' ' 

253 .join(['{:.4f}'.format(co) for co in self._cutoffs]))] 

254 s += [' {:38} : {}'.format('total number of parameters', len(self))] 

255 t = ['{}= {}'.format(k, c) 

256 for k, c in self.get_number_of_orbits_by_order().items()] 

257 s += [' {:38} : {}'.format('number of parameters by order', ' '.join(t))] 

258 for key, value in sorted(self._config.items()): 

259 s += [' {:38} : {}'.format(key, value)] 

260 

261 # table header 

262 s += [''.center(width, '-')] 

263 s += [repr_orbit(prototype_orbit, header=True)] 

264 s += [''.center(width, '-')] 

265 

266 # table body 

267 index = 0 

268 orbit_list_info = self.orbit_data 

269 while index < len(orbit_list_info): 

270 if (print_threshold is not None and 

271 len(self) > print_threshold and 

272 index >= print_minimum and 

273 index <= len(self) - print_minimum): 

274 index = len(self) - print_minimum 

275 s += [' ...'] 

276 s += [repr_orbit(orbit_list_info[index])] 

277 index += 1 

278 s += [''.center(width, '=')] 

279 

280 return '\n'.join(s) 

281 

282 def __repr__(self) -> str: 

283 """ String representation. """ 

284 return self._get_string_representation(print_threshold=50) 

285 

286 def print_overview(self, 

287 print_threshold: int = None, 

288 print_minimum: int = 10) -> None: 

289 """ 

290 Print an overview of the cluster space in terms of the orbits (order, 

291 radius, multiplicity etc). 

292 

293 Parameters 

294 ---------- 

295 print_threshold 

296 if the number of orbits exceeds this number print dots 

297 print_minimum 

298 number of lines printed from the top and the bottom of the orbit 

299 list if `print_threshold` is exceeded 

300 """ 

301 print(self._get_string_representation(print_threshold=print_threshold, 

302 print_minimum=print_minimum)) 

303 

304 @property 

305 def symprec(self) -> float: 

306 """ tolerance imposed when analyzing the symmetry using spglib """ 

307 return self._config['symprec'] 

308 

309 @property 

310 def position_tolerance(self) -> float: 

311 """ tolerance applied when comparing positions in Cartesian coordinates """ 

312 return self._config['position_tolerance'] 

313 

314 @property 

315 def fractional_position_tolerance(self) -> float: 

316 """ tolerance applied when comparing positions in fractional coordinates """ 

317 return self._config['fractional_position_tolerance'] 

318 

319 @property 

320 def space_group(self) -> str: 

321 """ space group of the primitive structure in international notion (via spglib) """ 

322 structure_as_tuple = ase_atoms_to_spglib_cell(self.primitive_structure) 

323 return spglib.get_spacegroup(structure_as_tuple, symprec=self._config['symprec']) 

324 

325 @property 

326 def orbit_data(self) -> List[OrderedDict]: 

327 """ 

328 list of orbits with information regarding 

329 order, radius, multiplicity etc 

330 """ 

331 data = [] 

332 zerolet = OrderedDict([('index', 0), 

333 ('order', 0), 

334 ('radius', 0), 

335 ('multiplicity', 1), 

336 ('orbit_index', -1), 

337 ('multi_component_vector', '.'), 

338 ('sublattices', '.')]) 

339 sublattices = self.get_sublattices(self.primitive_structure) 

340 data.append(zerolet) 

341 index = 1 

342 while index < len(self): 

343 multi_component_vectors_by_orbit = self.get_multi_component_vectors_by_orbit(index) 

344 orbit_index = multi_component_vectors_by_orbit[0] 

345 mc_vector = multi_component_vectors_by_orbit[1] 

346 orbit = self.get_orbit(orbit_index) 

347 repr_sites = orbit.sites_of_representative_cluster 

348 orbit_sublattices = '-'.join( 

349 [sublattices[sublattices.get_sublattice_index(ls.index)].symbol 

350 for ls in repr_sites]) 

351 local_Mi = self.get_number_of_allowed_species_by_site( 

352 self._get_primitive_structure(), orbit.sites_of_representative_cluster) 

353 mc_vectors = orbit.get_mc_vectors(local_Mi) 

354 mc_permutations = self.get_multi_component_vector_permutations( 

355 mc_vectors, orbit_index) 

356 mc_index = mc_vectors.index(mc_vector) 

357 mc_permutations_multiplicity = len(mc_permutations[mc_index]) 

358 cluster = self.get_orbit(orbit_index).representative_cluster 

359 

360 multiplicity = len(self.get_orbit( 

361 orbit_index).equivalent_clusters) 

362 record = OrderedDict([('index', index), 

363 ('order', cluster.order), 

364 ('radius', cluster.radius), 

365 ('multiplicity', multiplicity * 

366 mc_permutations_multiplicity), 

367 ('orbit_index', orbit_index)]) 

368 record['multi_component_vector'] = mc_vector 

369 record['sublattices'] = orbit_sublattices 

370 data.append(record) 

371 index += 1 

372 return data 

373 

374 def get_number_of_orbits_by_order(self) -> OrderedDict: 

375 """ 

376 Returns the number of orbits by order. 

377 

378 Returns 

379 ------- 

380 an ordered dictionary where keys and values represent order and number 

381 of orbits, respectively 

382 """ 

383 count_orbits = {} # type: dict[int, int] 

384 for orbit in self.orbit_data: 

385 k = orbit['order'] 

386 count_orbits[k] = count_orbits.get(k, 0) + 1 

387 return OrderedDict(sorted(count_orbits.items())) 

388 

389 def get_cluster_vector(self, structure: Atoms) -> np.ndarray: 

390 """ 

391 Returns the cluster vector for a structure. 

392 

393 Parameters 

394 ---------- 

395 structure 

396 atomic configuration 

397 

398 Returns 

399 ------- 

400 the cluster vector 

401 """ 

402 if not isinstance(structure, Atoms): 402 ↛ 403line 402 didn't jump to line 403, because the condition on line 402 was never true

403 raise TypeError('Input structure must be an ASE Atoms object') 

404 

405 try: 

406 cv = _ClusterSpace.get_cluster_vector( 

407 self, 

408 structure=Structure.from_atoms(structure), 

409 fractional_position_tolerance=self.fractional_position_tolerance) 

410 except Exception as e: 

411 self.assert_structure_compatibility(structure) 

412 raise(e) 

413 return cv 

414 

415 def get_coordinates_of_representative_cluster(self, orbit_index: int) -> List[Tuple[float]]: 

416 """ 

417 Returns the positions of atoms in the selected orbit 

418 

419 Parameters 

420 ---------- 

421 orbit_index 

422 index of the orbit from which to calculate the positions of the atoms 

423 

424 Returns 

425 ------- 

426 list of positions of atoms in the selected orbit 

427 

428 """ 

429 # Raise exception if chosen orbit index not in current list of orbit indices 

430 if not (orbit_index in range(len(self._orbit_list))): 

431 raise ValueError('The input orbit index is not in the list of possible values.') 

432 

433 lattice_sites = self._orbit_list.get_orbit(orbit_index).sites_of_representative_cluster 

434 positions = [] 

435 

436 for site in lattice_sites: 

437 pos = get_position_from_lattice_site(structure=self.primitive_structure, 

438 lattice_site=site) 

439 positions.append(pos) 

440 

441 return positions 

442 

443 def _prune_orbit_list(self, indices: List[int]) -> None: 

444 """ 

445 Prunes the internal orbit list 

446 

447 Parameters 

448 ---------- 

449 indices 

450 indices to all orbits to be removed 

451 """ 

452 size_before = len(self._orbit_list) 

453 

454 self._prune_orbit_list_cpp(indices) 

455 for index in sorted(indices, reverse=True): 

456 self._orbit_list.remove_orbit(index) 

457 self._compute_multi_component_vectors() 

458 

459 size_after = len(self._orbit_list) 

460 assert size_before - len(indices) == size_after 

461 self._pruning_history.append(indices) 

462 

463 @property 

464 def primitive_structure(self) -> Atoms: 

465 """ Primitive structure on which cluster space is based """ 

466 structure = self._get_primitive_structure().to_atoms() 

467 # Decorate with the "real" symbols (instead of H, He, Li etc) 

468 for atom, symbols in zip(structure, self._primitive_chemical_symbols): 

469 atom.symbol = min(symbols) 

470 return structure 

471 

472 @property 

473 def chemical_symbols(self) -> List[List[str]]: 

474 """ Species identified by their chemical symbols """ 

475 return self._primitive_chemical_symbols.copy() 

476 

477 @property 

478 def cutoffs(self) -> List[float]: 

479 """ 

480 Cutoffs for different n-body clusters. The cutoff radius (in 

481 Angstroms) defines the largest interatomic distance in a 

482 cluster. 

483 """ 

484 return self._cutoffs 

485 

486 @property 

487 def orbit_list(self): 

488 """Orbit list that defines the cluster in the cluster space""" 

489 return self._orbit_list 

490 

491 def get_possible_orbit_occupations(self, orbit_index: int) \ 

492 -> List[List[str]]: 

493 """Returns possible occupation of the orbit. 

494 

495 Parameters 

496 ---------- 

497 orbit_index 

498 """ 

499 orbit = self.orbit_list.orbits[orbit_index] 

500 

501 indices = [ 

502 lattice_site.index for lattice_site in orbit.sites_of_representative_cluster] 

503 

504 allowed_species = [self.chemical_symbols[index] for index in indices] 

505 

506 return list(itertools.product(*allowed_species)) 

507 

508 def get_sublattices(self, structure: Atoms) -> Sublattices: 

509 """ 

510 Returns the sublattices of the input structure. 

511 

512 Parameters 

513 ---------- 

514 structure 

515 structure the sublattices are based on 

516 """ 

517 sl = Sublattices(self.chemical_symbols, 

518 self.primitive_structure, 

519 structure, 

520 fractional_position_tolerance=self.fractional_position_tolerance) 

521 return sl 

522 

523 def assert_structure_compatibility(self, structure: Atoms, vol_tol: float = 1e-5) -> None: 

524 """ Raises error if structure is not compatible with ClusterSpace. 

525 

526 Todo 

527 ---- 

528 Add check for if structure is relaxed. 

529 

530 Parameters 

531 ---------- 

532 structure 

533 structure to check if compatible with ClusterSpace 

534 """ 

535 # check volume 

536 prim = self.primitive_structure 

537 vol1 = prim.get_volume() / len(prim) 

538 vol2 = structure.get_volume() / len(structure) 

539 if abs(vol1 - vol2) > vol_tol: 

540 raise ValueError('Volume per atom of structure does not match the volume of ' 

541 'ClusterSpace.primitive_structure') 

542 

543 # check occupations 

544 sublattices = self.get_sublattices(structure) 

545 sublattices.assert_occupation_is_allowed(structure.get_chemical_symbols()) 

546 

547 # check pbc 

548 if not all(structure.pbc): 

549 raise ValueError('Input structure must have periodic boundary conditions') 

550 

551 def is_supercell_self_interacting(self, structure: Atoms) -> bool: 

552 """ 

553 Checks whether an structure has self-interactions via periodic 

554 boundary conditions. 

555 

556 Parameters 

557 ---------- 

558 structure 

559 structure to be tested 

560 

561 Returns 

562 ------- 

563 bool 

564 If True, the structure contains self-interactions via periodic 

565 boundary conditions, otherwise False. 

566 """ 

567 ol = self.orbit_list.get_supercell_orbit_list( 

568 structure=structure, 

569 fractional_position_tolerance=self.fractional_position_tolerance) 

570 orbit_indices = set() 

571 for orbit in ol.orbits: 

572 for sites in orbit.equivalent_clusters: 

573 indices = tuple(sorted([site.index for site in sites])) 

574 if indices in orbit_indices: 

575 return True 

576 else: 

577 orbit_indices.add(indices) 

578 return False 

579 

580 def write(self, filename: str) -> None: 

581 """ 

582 Saves cluster space to a file. 

583 

584 Parameters 

585 --------- 

586 filename 

587 name of file to which to write 

588 """ 

589 

590 with tarfile.open(name=filename, mode='w') as tar_file: 

591 

592 # write items 

593 items = dict(cutoffs=self._cutoffs, 

594 chemical_symbols=self._input_chemical_symbols, 

595 pruning_history=self._pruning_history, 

596 symprec=self.symprec, 

597 position_tolerance=self.position_tolerance) 

598 temp_file = tempfile.TemporaryFile() 

599 pickle.dump(items, temp_file) 

600 temp_file.seek(0) 

601 tar_info = tar_file.gettarinfo(arcname='items', fileobj=temp_file) 

602 tar_file.addfile(tar_info, temp_file) 

603 temp_file.close() 

604 

605 # write structure 

606 temp_file = tempfile.NamedTemporaryFile() 

607 ase_write(temp_file.name, self._input_structure, format='json') 

608 temp_file.seek(0) 

609 tar_info = tar_file.gettarinfo(arcname='atoms', fileobj=temp_file) 

610 tar_file.addfile(tar_info, temp_file) 

611 temp_file.close() 

612 

613 @staticmethod 

614 def read(filename: str): 

615 """ 

616 Reads cluster space from filename. 

617 

618 Parameters 

619 --------- 

620 filename 

621 name of file from which to read cluster space 

622 """ 

623 if isinstance(filename, str): 

624 tar_file = tarfile.open(mode='r', name=filename) 

625 else: 

626 tar_file = tarfile.open(mode='r', fileobj=filename) 

627 

628 # read items 

629 items = pickle.load(tar_file.extractfile('items')) 

630 

631 # read structure 

632 temp_file = tempfile.NamedTemporaryFile() 

633 temp_file.write(tar_file.extractfile('atoms').read()) 

634 temp_file.seek(0) 

635 structure = ase_read(temp_file.name, format='json') 

636 

637 tar_file.close() 

638 

639 # ensure backward compatibility 

640 if 'symprec' not in items: # pragma: no cover 

641 items['symprec'] = 1e-5 

642 if 'position_tolerance' not in items: # pragma: no cover 

643 items['position_tolerance'] = items['symprec'] 

644 

645 cs = ClusterSpace(structure=structure, 

646 cutoffs=items['cutoffs'], 

647 chemical_symbols=items['chemical_symbols'], 

648 symprec=items['symprec'], 

649 position_tolerance=items['position_tolerance']) 

650 for indices in items['pruning_history']: 

651 cs._prune_orbit_list(indices) 

652 return cs 

653 

654 def copy(self): 

655 """ Returns copy of ClusterSpace instance. """ 

656 cs_copy = ClusterSpace(structure=self._input_structure, 

657 cutoffs=self.cutoffs, 

658 chemical_symbols=self._input_chemical_symbols, 

659 symprec=self.symprec, 

660 position_tolerance=self.position_tolerance) 

661 for indices in self._pruning_history: 

662 cs_copy._prune_orbit_list(indices) 

663 return cs_copy