Source code for hiphive.force_constant_potential

This module introduces the ForceConstantPotential object which acts as the
finalized force constant model.

import pickle
import tarfile
import numpy as np

from collections import Counter
from typing import Any, Dict, List
from .core.atoms import Atoms
from .force_constant_model import ForceConstantModel
from .core.orbits import Orbit
from .core.orbits import OrientationFamily
from .core.tensors import rotation_to_cart_coord, rotate_tensor
from .input_output.read_write_files import (add_items_to_tarfile_pickle,
from .input_output.logging_tools import logger

logger = logger.getChild('ForceConstantPotential')

[docs]class ForceConstantPotential: """ A finalized force constant model. Can produce force constants for any structure compatible with the structure for which the model was set up. Parameters ---------- cs : ClusterSpace The cluster space the model is based upon parameters : numpy.ndarray The fitted paramteres metadata : dict metadata dictionary, will be pickled when object is written to file """ def __init__(self, cs, parameters, metadata=None): self._prim = cs.primitive_structure.copy() self.cluster_list = cs.cluster_list.copy() self.atom_list = cs.atom_list.copy() self.orbits = [] self.spacegroup = cs.spacegroup self._config = cs._config self.cs_summary = cs.summary # add metadata if metadata is None: metadata = dict() self._metadata = metadata self._add_default_metadata() # Extract the eigentensors from the cluster space and use the paramters # to construct the finalized force constants parameters = cs._map_parameters(parameters) p = 0 for orb in cs.orbits: new_orbit = Orbit() fc = np.zeros(orb.eigentensors[0].shape) for et, a in zip(orb.eigentensors, parameters[p:]): fc += et * a new_orbit.force_constant = fc new_orbit.order = orb.order new_orbit.radius = orb.radius new_orbit.maximum_distance = orb.maximum_distance for of in orb.orientation_families: new_of = OrientationFamily() new_of.cluster_indices = of.cluster_indices.copy() sym_ind = of.symmetry_index R = rotation_to_cart_coord(cs.rotation_matrices[sym_ind], self.primitive_structure.cell) fc = rotate_tensor(new_orbit.force_constant, R.T) perm = cs.permutations[of.permutation_indices[0]] new_of.force_constant = fc.transpose(perm) new_orbit.orientation_families.append(new_of) self.orbits.append(new_orbit) p += len(orb.eigentensors) @property def symprec(self): return self._config['symprec']
[docs] def write(self, filename): """ Writes a ForceConstantPotential to file. Parameters ---------- filename : str name of file to write ForceConstantPotential to """ # Create a tar archive if isinstance(filename, str): tar_file =, mode='w') else: raise ValueError('filename must be str') # objects with custom write add_list_to_tarfile_custom(tar_file, self.orbits, 'orbits') # prim with its builtin write/read functions. Note prim is a hiphive.core.Atoms object items_custom = {'_prim': self._prim} add_items_to_tarfile_custom(tar_file, items_custom) # Attributes in pickle format pickle_attributes = ['_config', '_metadata', 'spacegroup', 'atom_list', 'cluster_list', 'cs_summary'] items_pickle = dict() for attribute in pickle_attributes: items_pickle[attribute] = self.__getattribute__(attribute) add_items_to_tarfile_pickle(tar_file, items_pickle, 'attributes') # Done! tar_file.close()
[docs] @staticmethod def read(filename): """ Reads a ForceConstantPotentialfrom file. Parameters ---------- filename : str name of input file to load ForceConstantPotential from Returns ------- ForceConstantPotential the original object as stored in the file """ # Try usage of old read format, Remove in hiphive 1.0 try: old_fcp = ForceConstantPotential._read_old(filename) logger.warning('This fcp was written with a version <1.0. Please rewrite it.') return old_fcp except Exception: pass # Instantiate empty cs obj. fcp = ForceConstantPotential.__new__(ForceConstantPotential) # Load from file on disk if type(filename) is str: tar_file ='r', name=filename) else: raise ValueError('filename must be str') # Attributes with custom read fileobj = tar_file.extractfile('_prim') fcp._prim = fcp.orbits = read_list_custom(tar_file, 'orbits', # Attributes attributes = read_items_pickle(tar_file, 'attributes') for name, value in attributes.items(): fcp.__setattr__(name, value) # Done! tar_file.close() return fcp
# TODO: Remove? def _write_old(self, f): """Writes a force constant potential to file using old format. Parameters ---------- f : str or file object name of input file (str) or stream to write to (file object) """ if isinstance(f, str): with open(f, 'wb') as fobj: pickle.dump(self, fobj) else: try: pickle.dump(self, f) except Exception: raise Exception('Failed writing to file.') @staticmethod def _read_old(f): """Reads a force constant potential from file in old format. Parameters ---------- f : str or file object name of input file (str) or stream to load from (file object) Returns ------- ForceConstantPotential the original object as stored in the file """ if isinstance(f, str): with open(f, 'rb') as fobj: # This allows for reading FCPs with ASE-3.17 and 3.18 fcp = pickle.load(fobj) _prim = fcp._prim if hasattr(_prim, '_cell'): # 3.17 cell = _prim._cell else: # 3.18 cell = _prim.cell[:] # assume PBC True (as it has to be True in hiphive) pbc = [True, True, True] new_prim = Atoms( symbols=_prim.symbols, positions=_prim.positions, cell=cell, pbc=pbc) fcp._prim = new_prim return fcp else: return pickle.load(f) @property def metadata(self): """ dict : metadata associated with force constant potential """ return self._metadata @property def primitive_structure(self): """ ase.Atoms : atomic structure """ return self._prim.copy() @property def orbit_data(self) -> List[Dict[str, Any]]: """list of dictionaries containing detailed information for each orbit, e.g. cluster radius and force constant """ data = [] for orbit_index, orbit in enumerate(self.orbits): d = {} d['index'] = orbit_index d['order'] = orbit.order d['radius'] = orbit.radius d['maximum_distance'] = orbit.maximum_distance d['n_clusters'] = len(orbit.orientation_families) types = [] for atom_ind in self.cluster_list[orbit.prototype_index]: types.append(self.primitive_structure.numbers[ self.atom_list[atom_ind].site]) d['prototype_cluster'] = self.cluster_list[orbit.prototype_index] d['prototype_atom_types'] = types d['geometrical_order'] = len(set(d['prototype_cluster'])) d['force_constant'] = orbit.force_constant d['force_constant_norm'] = np.linalg.norm(orbit.force_constant) data.append(d) return data
[docs] def print_tables(self): """ Prints information concerning the underlying cluster space to stdout,including, e.g., the number of cluster, orbits, and parameters by order and number of bodies. """ self.cs_summary.print_tables()
[docs] def get_force_constants(self, atoms): """ Return the force constants of a compatible structure. Parameters ---------- atoms : ase.Atoms input structure Returns ------- ForceConstants force constants """ return ForceConstantModel(atoms, self).get_force_constants()
def __str__(self): orbits = self.orbit_data orbit_counts = Counter([orbit['order'] for orbit in orbits]) cluster_counts = Counter() for orbit in orbits: cluster_counts[orbit['order']] += orbit['n_clusters'] n = 54 s = [] s.append(' ForceConstantPotential '.center(n, '=')) s.append(f'Spacegroup {self.spacegroup}') s.append(f'Cell:\n{self.primitive_structure.cell}') s.append(f'Basis:\n{self.primitive_structure.basis}') s.append(f'Numbers: {self.primitive_structure.numbers}') s.append(f'Cutoff matrix:\n{self.cs_summary._cutoff_matrix}') for order in sorted(orbit_counts.keys()): n_dofs = self.cs_summary.ndofs_by_order[order] s.append(f'Order {order}, #orbits {orbit_counts[order]}, #cluster ' f'{cluster_counts[order]}, #parameters {n_dofs}') s.append(f'Total number of orbits: {len(orbits)}') s.append(f'total number of clusters: {sum(cluster_counts.values())}') s.append(f'total number of parameters: {sum(self.cs_summary.ndofs_by_order.values())}') s.append(''.center(n, '=')) return '\n'.join(s) def __repr__(self): return 'ForceConstantPotential(ClusterSpace({!r}, ...), [...])'.format( self.primitive_structure) def _add_default_metadata(self): """Adds default metadata to metadata dict.""" import getpass import socket from datetime import datetime from . import __version__ as hiphive_version self._metadata['date_created'] ='%Y-%m-%dT%H:%M:%S') self._metadata['username'] = getpass.getuser() self._metadata['hostname'] = socket.gethostname() self._metadata['hiphive_version'] = hiphive_version