Other functions¶
Self-consistent phonons¶
-
hiphive.self_consistent_phonons.
self_consistent_harmonic_model
(atoms_ideal, calc, cs, T, alpha, n_iterations, n_structures, parameters_start=None, fit_kwargs={})[source]¶ Constructs a set of self-consistent second-order force constants that provides the closest match to the potential energy surface at a the specified temperature.
- Parameters
atoms_ideal (ase.Atoms) – ideal structure
calc (ASE calculator object) – calculator to be used as reference potential
cs (ClusterSpace) – clusterspace onto which to project the reference potential
T (float) – temperature in K
alpha (float) – stepsize in optimization algorithm
n_iterations (int) – number of iterations in poor mans
n_structures (int) – number of structures to use when fitting
parameters_start (numpy.ndarray) – parameters from which to start the optimization
fit_kwargs (dict) – kwargs to be used in the fitting process (via Optimizer)
- Returns
sequence of parameter vectors generated while iterating to self-consistency
- Return type
list(numpy.ndarray)
Utilities¶
This module contains various support/utility functions.
-
class
hiphive.utilities.
Shell
(types, distance, count=0)[source]¶ Neighbor Shell class
- Parameters
types (list or tuple) – atomic types for neighbor shell
distance (float) – interatomic distance for neighbor shell
count (int) – number of pairs in the neighbor shell
-
hiphive.utilities.
extract_parameters
(fcs, cs)[source]¶ Extracts parameters from force constants.
TODO: Rename this function with more explanatory name?
This function can be used to extract parameters to create a ForceConstantPotential from a known set of force constants. The return values come from NumPy’s lstsq function.
- Parameters
fcs (
ForceConstants
) – force constantscs (
ClusterSpace
) – cluster space
- Return type
- Returns
x ({(N,), (N, K)} ndarray) – Least-squares solution. If b is two-dimensional, the solutions are in the K columns of x.
residuals ({(1,), (K,), (0,)} ndarray) – Sums of residuals; squared Euclidean 2-norm for each column in
b - a*x
. If the rank of a is < N or M <= N, this is an empty array. If b is 1-dimensional, this is a (1,) shape array. Otherwise the shape is (K,).rank (int) – Rank of matrix a.
s ((min(M, N),) ndarray) – Singular values of a.
-
hiphive.utilities.
find_permutation
(atoms, atoms_ref)[source]¶ Returns the best permutation of atoms for mapping one configuration onto another.
- Parameters
atoms (
Atoms
) – configuration to be permutedatoms_ref (
Atoms
) – configuration onto which to map
Examples
After obtaining the permutation via
p = find_permutation(atoms1, atoms2)
the reordered structureatoms1[p]
will give the closest match toatoms2
.- Return type
List
[int
]
-
hiphive.utilities.
get_displacements
(atoms, atoms_ideal, cell_tol=0.0001)[source]¶ Returns the the smallest possible displacements between a displaced configuration relative to an ideal (reference) configuration.
Notes
assumes periodic boundary conditions in all directions
- Parameters
atoms (
Atoms
) – configuration with displaced atomsatoms_ideal (
Atoms
) – ideal configuration relative to which displacements are computedcell_tol (
float
) – cell tolerance; if cell missmatch more than tol value error is raised
- Return type
-
hiphive.utilities.
get_neighbor_shells
(atoms, cutoff, dist_tol=1e-05)[source]¶ Returns a list of neighbor shells.
Distances are grouped into shells via the following algorithm:
Find smallest atomic distance d_min
Find all pair distances in the range d_min + 1 * dist_tol
Construct a shell from these and pop them from distance list
Go to 1.
- Parameters
atoms (
Atoms
) – configuration used for finding shellscutoff (
float
) – exclude neighbor shells which have a distance larger than this valuedist_tol (
float
) – distance tolerance
- Return type
List
[Shell
]
-
hiphive.utilities.
prepare_structure
(atoms, atoms_ideal, calc=None)[source]¶ Prepare a structure in the format suitable for a
StructureContainer
.- Parameters
atoms (
Atoms
) – input structureatoms_ideal (
Atoms
) – reference structure relative to which displacements are computedcalc (
Optional
[SinglePointCalculator
]) – ASE calculator used for computing forces
- Returns
prepared ASE atoms object with forces and displacements as arrays
- Return type
ASE atoms object
-
hiphive.utilities.
prepare_structures
(structures, atoms_ideal, calc=None)[source]¶ Prepares a set of structures in the format suitable for adding them to a
StructureContainer
.structures should represent a list of supercells with displacements while atoms_ideal should provide the ideal reference structure (without displacements) for the given structures.
The structures that are returned will have their positions reset to the ideal structures. Displacements and forces will be added as arrays to the atoms objects.
If no calculator is provided, then there must be an ASE SinglePointCalculator <ase.calculators.singlepoint> object attached to the structures or the forces should already be attached as arrays to the structures.
If a calculator is provided then it will be used to compute the forces for all structures.
Example
The following example illustrates the use of this function:
db = connect('dft_training_structures.db') training_structures = [row.toatoms() for row in db.select()] training_structures = prepare_structures(training_structures, atoms_ideal) for s in training_structures: sc.add_structure(s)
- Parameters
structures (
List
[Atoms
]) – list of input displaced structuresatoms_ideal (
Atoms
) – reference structure relative to which displacements are computedcalc (
Optional
[SinglePointCalculator
]) – ASE calculator used for computing forces
- Returns
- Return type
list of prepared structures with forces and displacements as arrays
Enforcing rotational sum rules¶
-
hiphive.core.rotational_constraints.
enforce_rotational_sum_rules
(cs, parameters, sum_rules=None, alpha=1e-06, **ridge_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’
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]
ridge_kwargs (dict) – kwargs to be passed to sklearn Ridge
- Returns
constrained parameters
- Return type
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)