ForceConstants

The ForceConstants object contains the force constants, supercell, and crystal structure information required to calculate phonon frequencies and eigenvectors at any arbitrary q via Fourier interpolation.

Reading From CASTEP

The force constants matrix and other required information can be read from a .castep_bin or .check file with ForceConstants.from_castep:

from euphonic import ForceConstants

filename = 'quartz/quartz.castep_bin'
fc = ForceConstants.from_castep(filename)

By default CASTEP may not write the force constants, if you receive an error saying the force constants could not be read, in the .param file ensure a PHONON_FINE_METHOD has been chosen e.g. PHONON_FINE_METHOD: interpolate, and set PHONON_WRITE_FORCE_CONSTANTS: true, then rerun CASTEP to trigger the force constants to be written.

Reading From Phonopy

When using Phonopy with Euphonic, it is recommended that all the required data (force constants, crystal structure, born charges if applicable) be collected in a single phonopy.yaml file. This can be done by running Phonopy with the --include-all flag or with INCLUDE_ALL = .TRUE. (phonopy >= 2.5.0 only).

Required information is read from Phonopy output files using ForceConstants.from_phonopy. A path keyword argument can be supplied (if the files are in another directory), and by default phonopy.yaml is read, but the filename can be changed with the summary_name keyword argument:

from euphonic import ForceConstants

fc = ForceConstants.from_phonopy(path='NaCl',
                                 summary_name='phonopy_fc.yaml')

If you are using an older version of Phonopy, the force constants and born charges can also be read from Phonopy plaintext or hdf5 files by specifying the fc_name and born_name keyword arguments:

from euphonic import ForceConstants

fc = ForceConstants.from_phonopy(path='NaCl',
                                 fc_name='force_constants.hdf5',
                                 born_name='BORN')

Calculating phonon frequencies/eigenvectors

Phonon frequencies and eigenvectors are calculated using ForceConstants.calculate_qpoint_phonon_modes (see the docstring for algorithm details). A Numpy array of q-points of shape (n_qpts, 3) must be provided, and a QpointPhononModes object is returned. A recommended q-point path for plotting bandstructures can be generated using seekpath:

import seekpath
import numpy as np
from euphonic import ForceConstants

# Read quartz data from quartz.castep_bin
filename = 'quartz/quartz.castep_bin'
fc = ForceConstants.from_castep(filename)

# Generate a recommended q-point path using seekpath
_, unique_atoms = np.unique(fc.crystal.atom_type, return_inverse=True)
structure = (fc.crystal.cell_vectors.magnitude,
             fc.crystal.atom_r, unique_atoms)
qpts = seekpath.get_explicit_k_path(structure)["explicit_kpoints_rel"]

# Calculate frequencies/eigenvectors
phonons = fc.calculate_qpoint_phonon_modes(qpts, asr='reciprocal')

Docstring

class ForceConstants(crystal, force_constants, sc_matrix, cell_origins, born=None, dielectric=None)

A class to read and store the data required for a phonon interpolation calculation from model (e.g. CASTEP) output, and calculate phonon frequencies/eigenvectors at arbitrary q-points via Fourier interpolation

Variables:
  • crystal (Crystal) – Lattice and atom information
  • force_constants ((n_cells_in_sc, 3*n_atoms, 3*n_atoms) float Quantity) – Force constants matrix
  • sc_matrix ((3, 3) int ndarray) – The supercell matrix
  • n_cells_in_sc (int) – Number of cells in the supercell
  • cell_origins ((n_cells_in_sc, 3) int ndarray) – The locations of the unit cells within the supercell
  • born ((n_atoms, 3, 3) float Quantity or None) – The Born charges for each atom
  • dielectric ((3, 3) float Quantity or None) – The dielectric permittivity tensor
__init__(crystal, force_constants, sc_matrix, cell_origins, born=None, dielectric=None)
Parameters:
  • crystal (Crystal) – Lattice and atom information
  • force_constants ((n_cells_in_sc, 3*n_atoms, 3*n_atoms) float Quantity) – Force constants matrix
  • sc_matrix ((3, 3) int ndarray) – The supercell matrix
  • cell_origins ((n_cells_in_sc, 3) int ndarray) – The locations of the unit cells within the supercell
  • born ((n_atoms, 3, 3) float Quantity, optional) – The Born charges for each atom
  • dielectric ((3, 3) float Quantity, optional) – The dielectric permittivity tensor
calculate_qpoint_phonon_modes(qpts, asr=None, dipole=True, eta_scale=1.0, splitting=True, insert_gamma=False, reduce_qpts=True, use_c=False, n_threads=1, fall_back_on_python=True)

Calculate phonon frequencies and eigenvectors at specified q-points from a force constants matrix via Fourier interpolation

Parameters:
  • qpts ((n_qpts, 3) float ndarray) – The q-points to interpolate onto
  • asr ({'realspace', 'reciprocal'}, optional) – Which acoustic sum rule correction to apply. ‘realspace’ applies the correction to the force constant matrix in real space. ‘reciprocal’ applies the correction to the dynamical matrix at every q-point
  • dipole (boolean, optional) – Calculates the dipole tail correction to the dynamical matrix at each q-point using the Ewald sum, if the Born charges and dielectric permitivitty tensor are present.
  • eta_scale (float, optional) – Changes the cutoff in real/reciprocal space for the dipole Ewald sum. A higher value uses more reciprocal terms. If tuned correctly this can result in performance improvements. See scripts/optimise_eta.py for help on choosing a good eta_scale.
  • splitting (boolean, optional) – Whether to calculate the LO-TO splitting at the gamma points. Only applied if dipole is True and the Born charges and dielectric permitivitty tensor are present.
  • insert_gamma (boolean, optional) – If splitting is True, this will insert gamma points into qpts to store the extra split frequencies. Note this means that the length of qpts in the output QpointPhononModes object will not necessarily be the same as the input qpts. If qpts already contains double gamma points where you want split frequencies, leave this as False.
  • reduce_qpts (boolean, optional) – Whether to use periodicity to reduce all q-points and only calculate for unique q-points within the 1st BZ. This won’t change the output but could increase performance.
  • use_c (boolean, optional) – Whether to use C instead of Python to calculate and diagonalise the dynamical matrix
  • n_threads (int, optional) – The number of OpenMP threads to use when looping over q-points in C. Only applicable if use_c=True
  • fall_back_on_python (boolean, optional) – If we cannot use the C extension, fall back on using python if this is true, else raise an ImportCError.
Returns:

A QpointPhononModes object containing the interpolated frequencies and eigenvectors at each q-point. Note that if there is LO-TO splitting, and insert_gamma=True, the number of input q-points may not be the same as in the output object
Return type:

QpointPhononModes
Raises:

ImportCError – If we have selected not to fall back on Python and cannot use the C extension

Notes

Phonon frequencies/eigenvectors are calculated at any q-point by Fourier interpolation of a force constants matrix. The force constants matrix is defined as [1]:

\[\phi_{\alpha, {\alpha}'}^{\kappa, {\kappa}'} = \frac{\delta^{2}E}{{\delta}u_{\kappa,\alpha}{\delta}u_{{\kappa}',{\alpha}'}}\]

Which gives the Dynamical matrix at q:

\[D_{\alpha, {\alpha}'}^{\kappa, {\kappa}'}(q) = \frac{1}{\sqrt{M_\kappa M_{\kappa '}}} \sum_{a}\phi_{\alpha, \alpha '}^{\kappa, \kappa '}e^{-iq\cdot r_a}\]

The eigenvalue equation for the dynamical matrix is then:

\[D_{\alpha, {\alpha}'}^{\kappa, {\kappa}'}(q) \epsilon_{q\nu\kappa\alpha} = \omega_{q\nu}^{2} \epsilon_{q\nu\kappa\alpha}\]

Where \(\nu\) runs over phonon modes, \(\kappa\) runs over atoms, \(\alpha\) runs over the Cartesian directions, \(a\) runs over unit cells in the supercell, \(u_{\kappa, \alpha}\) is the displacement of atom \(\kappa\) in direction \(\alpha\), \(M_{\kappa}\) is the mass of atom \(\kappa\), \(r_{a}\) is the vector to the origin of cell \(a\) in the supercell, \(\epsilon_{q\nu\kappa\alpha}\) are the eigevectors, and \(\omega_{q\nu}^{2}\) are the frequencies squared.

In polar materials, there is an additional long-ranged correction to the force constants matrix (applied if dipole=True) and a non-analytical correction at the gamma point [2] (applied if splitting=True).

[1]M.T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, Cambridge, 1993, 83-87
[2]
  1. Gonze, K. C. Charlier, D. C. Allan, M. P. Teter, Phys. Rev. B, 1994, 50, 13035-13038
to_dict()

Convert to a dictionary. See ForceConstants.from_dict for details on keys/values

Returns:
Return type:dict
to_json_file(filename)

Write to a JSON file. JSON fields are equivalent to ForceConstants.from_dict keys

Parameters:filename (str) – Name of the JSON file to write to
classmethod from_dict(d)

Convert a dictionary to a ForceConstants object

Parameters:d (dict) –

A dictionary with the following keys/values:

  • ’crystal’: dict, see Crystal.from_dict
  • ’force_constants’: (n_cells_in_sc, 3*crystal.n_atoms, 3*crystal.n_atoms) float ndarray
  • ’force_constants_unit’: str
  • ’sc_matrix’: (3,3) int ndarray
  • ’cell_origins’: (n_cells_in_sc, 3) int ndarray

There are also the following optional keys:

  • ’born’: (3*crystal.n_atoms, 3, 3) float ndarray
  • ’born_unit’: str
  • ’dielectric’: (3, 3) float ndarray
  • ’dielectric_unit’: str
Returns:
Return type:ForceConstants
classmethod from_json_file(filename)

Read from a JSON file. See ForceConstants.from_dict for required fields

Parameters:filename (str) – The file to read from
Returns:
Return type:ForceConstants
classmethod from_castep(filename)

Reads from a .castep_bin or .check file

Parameters:filename (str) – The path and name of the file to read
Returns:
Return type:ForceConstants
classmethod from_phonopy(path='.', summary_name='phonopy.yaml', born_name=None, fc_name='FORCE_CONSTANTS', fc_format=None)

Reads data from the phonopy summary file (default phonopy.yaml) and optionally born and force constants files. Only attempts to read from born or force constants files if these can’t be found in the summary file.

Parameters:
  • path (str, optional) – Path to directory containing the file(s)
  • summary_name (str, optional) – Filename of phonopy summary file, default phonopy.yaml. By default any information (e.g. force constants) read from this file takes priority
  • born_name (str, optional) – Name of the Phonopy file containing born charges and dielectric tensor (by convention in Phonopy this would be called BORN). Is only read if Born charges can’t be found in the summary_name file
  • fc_name (str, optional) – Name of file containing force constants. Is only read if force constants can’t be found in summary_name
  • fc_format ({'phonopy', 'hdf5'} str, optional) – Format of file containing force constants data. FORCE_CONSTANTS is type ‘phonopy’
Returns:
Return type:

ForceConstants