Force constant models

ForceConstantPotential

class hiphive.ForceConstantPotential(cs, parameters)[source]

A finalized force constant model. Can produce force constants for any structure compatible with the structure for which the model was set up.

Parameters
get_force_constants(atoms)[source]

Return the force constants of a compatible structure.

Parameters

atoms (ase.Atoms) – input structure

Returns

force constants

Return type

ForceConstants

orbit_data

list of dictionaries containing detailed information for each orbit, e.g. cluster radius and force constant

Type

list

parameters
primitive_structure

atomic structure

Type

ase.Atoms

static read(f)[source]

Reads a force constant potential from file.

Parameters

f (str or file object) – name of input file (str) or stream to load from (file object)

Returns

the original object as stored in the file

Return type

ForceConstantPotential

symprec
write(f)[source]

Writes a force constant potential to file.

Parameters

f (str or file object) – name of input file (str) or stream to write to (file object)

ForceConstantCalculator

class hiphive.calculators.ForceConstantCalculator(fcs)[source]

This class provides an ASE calculator that can be used in conjunction with integrators and optimizers with the atomic simulation environment (ASE). To initialize an object of this class one must provide the ideal atomic configuration along with a compatible force constant model.

Parameters

fcs (ForceConstants) – the force constants instance must include the atomic configuration

calculate(atoms=None, properties=['energy'], system_changes=['positions', 'numbers', 'cell', 'pbc', 'initial_charges', 'initial_magmoms'])[source]

Do the calculation.

properties: list of str

List of what needs to be calculated. Can be any combination of ‘energy’, ‘forces’, ‘stress’, ‘dipole’, ‘charges’, ‘magmom’ and ‘magmoms’.

system_changes: list of str

List of what has changed since last calculation. Can be any combination of these six: ‘positions’, ‘numbers’, ‘cell’, ‘pbc’, ‘initial_charges’ and ‘initial_magmoms’.

Subclasses need to implement this, but can ignore properties and system_changes if they want. Calculated properties should be inserted into results dictionary like shown in this dummy example:

self.results = {'energy': 0.0,
                'forces': np.zeros((len(atoms), 3)),
                'stress': np.zeros(6),
                'dipole': np.zeros(3),
                'charges': np.zeros(len(atoms)),
                'magmom': 0.0,
                'magmoms': np.zeros(len(atoms))}

The subclass implementation should first call this implementation to set the atoms attribute.

compute_energy_and_forces()[source]

Compute energy and forces.

Returns

energy and forces

Return type

float, list(list(float))

implemented_properties = ['energy', 'forces']

ForceConstants

This module provides functionality for storing and handling of force constants.

class hiphive.force_constants.ForceConstants(supercell)[source]

Base class for force constants

assert_acoustic_sum_rules(order=None, tol=1e-06)[source]

Asserts that acoustic sum rules are enforced for force constants.

Parameters
  • order (int) – specifies which order to check, if None all are checked

  • tol (float) – numeric tolerance for checking sum rules

Raises

AssertionError – if acoustic sum rules are not enforced

clusters

sorted list of clusters

Type

list

compute_gamma_frequencies()[source]

Computes gamma frequencies using the second-order force constants.

Returns

gamma_frequencies – Gamma frequencies in THz

Return type

numpy.ndarray

classmethod from_arrays(supercell, fc2_array=None, fc3_array=None)[source]

Constructs FCs from numpy arrays.

One or both of fc2_array and fc3_array must not be None

Parameters
  • supercell (ase.Atoms) – supercell structure

  • fc2_array (np.ndarray) – second-order force constant in phonopy format, i.e. must have shape (N, N, 3, 3)

  • fc3_array (np.ndarray) – third-order force constant in phonopy format, i.e. must have shape (N, N, N, 3, 3, 3)

classmethod from_dense_dict(fc_dict, supercell)[source]

All permutations of clusters that are not zero must be listed, if label symmetries are fullfilled will return a SortedForceConstants

Parameters
  • fc_dict (dict) – keys corresponding to clusters and values to the force constants

  • supercell (ase.Atoms) –

classmethod from_sparse_dict(fc_dict, supercell)[source]

Assumes label symmetries, meaning only one cluster for each permuation should be included

Parameters
  • fc_dict (dict) – keys corresponding to clusters and values to the force constants

  • supercell (ase.Atoms) –

get_fc_array(order, format='phonopy')[source]

Returns force constants in array format for specified order.

Parameters
  • order (int) – force constants for this order will be returned

  • format (str) – specify which format (shape) the NumPy array should have, possible values are phonopy and ase

Returns

NumPy array with shape (N,)*order + (3,)*order where N is the number of atoms

Return type

numpy.ndarray

get_fc_dict(order=None)[source]

Returns force constant dictionary for one specific order.

The returned dict may be sparse or may be dense depending on the underlying force constants.

Parameters

order (int) – force constants returned for this order

Returns

dictionary with keys corresponding to clusters and values to respective force constant

Return type

dict

n_atoms

number of atoms

Type

int

print_force_constant(cluster)[source]

Prints force constants for a cluster in a nice format.

Parameters

cluster (tuple(int)) – sites belonging to the cluster

classmethod read(fname)[source]

Reads ForceConstants from file.

Parameters

fname (str) – name of file from which to read

classmethod read_phono3py(supercell, fname)[source]

Reads force constants from a phono3py calculation.

Parameters
  • supercell (ase.Atoms) – supercell structure (SPOSCAR)

  • fname (str) – name of third-order force constant file

classmethod read_phonopy(supercell, fname, format=None)[source]

Reads force constants from a phonopy calculation.

Parameters
  • supercell (ase.Atoms) – supercell structure (SPOSCAR)

  • fname (str) – name of second-order force constant file

  • format (str) – format for second-order force constants

classmethod read_shengBTE(supercell, fname, prim)[source]

Reads third order force constants from a shengBTE calculation.

shengBTE force constants will be mapped onto a supercell.

Parameters
  • supercell (str) – supercell structure

  • fname (str) – name of third-order force constant file

  • prim (ase.Atoms) – primitive configuration (must be equivalent to structure used in the shengBTE calculation)

supercell

supercell associated with force constants

Type

ase.Atoms

write(fname)[source]

Writes entire ForceConstants object to file.

Parameters

fname (str) – name of file to which to write

write_to_phono3py(fname)[source]

Writes force constants in phono3py format.

Parameters

fname (str) – name of file to which to write third-order force constant

write_to_phonopy(fname, format=None)[source]

Writes force constants in phonopy format.

Parameters
  • fname (str) – name of file to which to write second-order force constant

  • format (str) – format for second-order force constants

write_to_shengBTE(fname, prim, **kwargs)[source]

Writes third order force constants in shengBTE format.

Parameters
  • fname (str) – name of file to which to write third-order force constant

  • prim (ase.Atoms) – primitive configuration (must be equivalent to structure used in the shengBTE calculation)

class hiphive.force_constants.RawForceConstants(fc_dict, supercell)[source]

Force constants without label symmetries.

Parameters
  • fc_dict (dict) – keys corresponding to clusters and values to the force constants, should contain all clusters with nonzero force constants

  • supercell (ase.Atoms) –

orders

orders for which force constants exist

Type

list

class hiphive.force_constants.SortedForceConstants(fc_dict, supercell)[source]

Force constants with label symmetries.

Parameters
  • fc_dict (dict) – keys corresponding to clusters and values to the force constants, should only contain sorted clusters

  • supercell (ase.Atoms) –

orders

orders for which force constants exist

Type

list

hiphive.force_constants.array_to_dense_dict(fc_array, fc_tol=1e-10)[source]

Constructs a dense dict from an fc array in phonopy format.

Force constants with norm smaller than fc_tol will be considered zero and therefore not included in the fc_dict.

Parameters
  • fc_array (np.ndarray) – force constant array in phonopy format

  • fc_tol (float) – tolerance for considering force constants zero or not

hiphive.force_constants.check_label_symmetries(fc_dict)[source]

Checks label symmetries for dense fc dict.

Todo

tol, which one to use etc

Parameters

fc_dict (dict) – keys corresponding to clusters and values to the force constants

hiphive.force_constants.dense_dict_to_sparse_dict(fc_dict)[source]

Converts dense dict to sparse dict.

This does not check if label symmetry is True, but rather will just keep the sorted clusters and their force constants.

Parameters

fc_dict (dict) – keys corresponding to clusters and values to the force constants

hiphive.force_constants.symbolize_force_constant(fc, tol=1e-10)[source]

Carries out a symbolic symmetrization of a force constant tensor.

Parameters
  • fc (numpy.ndarray) – force constant tensor

  • tol (float) – tolerance used to decide whether two elements are identical

Returns

symbolic representation of force constant matrix

Return type

numpy.ndarray

Constraints

hiphive.enforce_rotational_sum_rules(cs, parameters, sum_rules, **kwargs)[source]

Enforces rotational sum rules by projecting parameters.

Note

The interface to this function might change in future releases.

Parameters
  • cs (ClusterSpace) – the underlying cluster space

  • parameters (numpy.ndarray) – parameters to be constrained

  • sum_rules (list(str)) – type of sum rules to enforce; possible values: ‘Huang’, ‘Born-Huang’

  • ridge_alpha (float) – hyperparameter to the ridge regression algorithm; keyword argument passed to the optimizer; larger values specify stronger regularization, i.e. less correction but higher stability [default: 1e-6]

  • iterations (int) – number of iterations to run the projection since each step projects the solution down to each nullspace in serial; keyword argument passed to the optimizer [default: 10]

Returns

constrained parameters

Return type

numpy.ndarray

Examples

The rotational sum rules can be enforced to the parameters before constructing a force constant potential as illustrated by the following snippet:

cs = ClusterSpace(reference_structure, cutoffs)
sc = StructureContainer(cs)
# add structures to structure container
opt = Optimizer(sc.get_fit_data())
opt.train()
new_params = enforce_rotational_sum_rules(cs, opt.parameters,
    sum_rules=['Huang', 'Born-Huang'])
fcp = ForceConstantPotential(cs, new_params)