Representing FCs from external sources as a FCP

There exist several codes that allow one to compute force constants and hiphive provides functionality to import some of these formats, including phonopy, phono3py, and thirdorder.py (ShengBTE). It can be convenient to represent these force constants as ForceConstantPotential objects, for example in order to enforce rotational sum rules.

Firstly these force constants have to be imported as a ForceConstants object. Secondly, one has to construct the ClusterSpace on which the force constants will be projected. Here, we recommend using the longest possible cutoffs for the supercell when constructing the ClusterSpace. Then one can obtain the parameters as follows:

from hiphive.utilities import extract_parameters
parameters = extract_parameters(fcs, cs)

Now a ForceConstantPotential can be created by combining the ClusterSpace and the parameters in the usual manner:

fcp = ForceConstantPotential(cs, parameters)

Note that the ForceConstantPotential might not correspond perfectly to underlying (external) force constants. Small errors can arise during the projection due to e.g.

1. If the external force constants do not obey the same crystal symmetries or acoustic sum rules as the ClusterSpace, then there will be some errors in the projection.

2. The force constants outside the cutoff for the ClusterSpace will not be captured in the FCP.

3. Potentially some issues with supercell/force constant folding.

Here, we note that 1. is often not a problem as most programs calculating force constants enforce symmetries and acoustic sum rules. Problems regarding 2. and 3. might be an indication that the supercell used for the force constant calculation is too small, and using a larger supercell will likely resolves these problems.

A minimal example can be found in
tests/integration/fcs_sensing.py

import numpy as np

from ase.build import bulk
from hiphive import ClusterSpace, ForceConstantPotential
from hiphive.utilities import extract_parameters


def test_fcs_sensing():
    tols = [1e-12, 1e-10]
    methods = ['numpy', 'scipy']

    cutoffs = [5.0, 4.0]
    prim = bulk('Al', a=4.05)
    dim = 4

    ideal = prim.repeat(dim)

    cs = ClusterSpace(prim, cutoffs)
    parameters = np.random.random(cs.n_dofs)
    fcp = ForceConstantPotential(cs, parameters)
    fcs = fcp.get_force_constants(ideal)
    for tol, method in zip(tols, methods):
        fitted_parameters = extract_parameters(fcs, cs, lstsq_method=method)
        assert np.linalg.norm(fitted_parameters - parameters) < tol