QpointPhononModes¶

The QpointPhononModes object contains precalculated phonon frequencies and eigenvectors at certain q-points.

Reading From CASTEP
Reading From Phonopy
From Force Constants
Reordering frequencies
Calculating The Coherent Neutron Structure Factor
Calculating The Debye-Waller Exponent
Calculating Density of States
Docstring

Reading From CASTEP ¶

Phonon frequencies and eigenvectors can be read from a .phonon file using QpointPhononModes.from_castep.

from euphonic import QpointPhononModes

filename = 'quartz.phonon'
phonons = QpointPhononModes.from_castep(filename)

Reading From Phonopy ¶

Phonopy should be run with the --eigvecs flag, or EIGENVECTORS = .TRUE. for use with Euphonic.

Using QpointPhononModes.from_phonopy Euphonic can read frequencies and eigenvectors from Phonopy files with the following default names:

mesh.yaml/mesh.hdf5
qpoints.yaml/qpoints.hdf5
bands.yaml/bands.hdf5

The file to be read can be specified with the phonon_name argument. Some of these files do not include the crystal information, so it must be read from a phonopy.yaml file, which can be specified with the summary_name argument. A path can also be specified.

from euphonic import QpointPhononModes

phonons = QpointPhononModes.from_phonopy(path='NaCl', phonon_name='mesh.hdf5')

From Force Constants ¶

See Force Constants

Reordering frequencies ¶

The stored frequencies can be reordered by comparing eigenvectors using QpointPhononModes.reorder_frequencies. This reordering can be seen the plotting dispersion (see Plotting)

from euphonic import QpointPhononModes

phonons = QpointPhononModes.from_castep('quartz.phonon')
phonons.reorder_frequencies()

Calculating The Coherent Neutron Structure Factor ¶

The neutron structure factor can be calculated for each branch and q-point using QpointPhononModes.calculate_structure_factor (see the docstring for algorithm details). A dictionary containing the coherent neutron scattering lengths for each atom must be provided. The units of each scattering length in the dictionary must be specified by providing them as a pint.Quantity, see Units for details. Inclusion of the Debye-Waller factor is optional, and can be provided in the dw keyword argument, see Calculating The Debye-Waller Exponent. From this function a StructureFactor object is returned.

The following example shows a full calculation from the force constants to the structure factor with Debye-Waller:

import seekpath
import numpy as np
from euphonic import ureg, QpointPhononModes, ForceConstants
from euphonic.util import mp_grid

# Read the force constants
fc = ForceConstants.from_castep('quartz.castep_bin')

# Generate a recommended q-point path to calculate the structure factor on
# 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 for the q-point path
phonons = fc.calculate_qpoint_phonon_modes(qpts, asr='reciprocal')

# For the Debye-Waller calculation, generate and calculate
# frequencies/eigenvectors on a grid (generate a Monkhorst-Pack grid of
# q-points using the mp-grid helper function)
q_grid = mp_grid([5,5,5])
phonons_grid = fc.calculate_qpoint_phonon_modes(q_grid, asr='reciprocal')
# Now calculate the Debye-Waller exponent
temperature = 5*ureg('K')
dw = phonons_grid.calculate_debye_waller(temperature)

# Calculate the structure factor for each q-point in phonons. A
# StructureFactor object is returned
fm = ureg('fm')
scattering_lengths = {'Si': 4.1491*fm, 'O': 5.803*fm}
sf = phonons.calculate_structure_factor(scattering_lengths, dw=dw)

Calculating The Debye-Waller Exponent ¶

The Debye-Waller factor is an optional part of the structure factor calculation. The exponent part of the Debye-Waller factor is independent of Q and should be precalculated using QpointPhononModes.calculate_debye_waller (see the docstring for algorithm details). This requires a QpointPhononModes object calculated on a grid of q-points and a temperature, and returns a DebyeWaller object. The Debye-Waller exponent can be calculated by:

from euphonic import ureg, QpointPhononModes

phonons = QpointPhononModes.from_castep('quartz-grid.phonon')
temperature = 5*ureg('K')
dw = phonons.calculate_debye_waller(temperature)

Calculating Density of States ¶

Density of states can be calculated using QpointPhononModes.calculate_dos. This requires an array of energy bin edges, with the units specified by wrapping it as a pint.Quantity (see Units for details). This function returns a generic Spectrum1D object. For example:

from euphonic import ureg, QpointPhononModes
import numpy as np

phonons = QpointPhononModes.from_castep('quartz.phonon')

# Create an array of energy bins 0 - 100 in meV
energy_bins = np.arange(0, 101, 1)*ureg('meV')

# Calculate dos
dos = phonons.calculate_dos(energy_bins)

Docstring ¶

class QpointPhononModes(crystal, qpts, frequencies, eigenvectors, weights=None)¶

A class to read and store vibrational data from model (e.g. CASTEP) output files

Variables:

crystal (Crystal) – Lattice and atom information
n_qpts (int) – Number of q-points in the object
qpts ((n_qpts, 3) float ndarray) – Q-point coordinates, in fractional coordinates of the reciprocal lattice
weights ((n_qpts,) float ndarray) – The weight for each q-point
frequencies ((n_qpts, 3*crystal.n_atoms) float Quantity) – Phonon frequencies per q-point and mode
eigenvectors ((n_qpts, 3*crystal.n_atoms, crystal.n_atoms, 3) complex ndarray) – Dynamical matrix eigenvectors

__init__(crystal, qpts, frequencies, eigenvectors, weights=None)¶

Parameters:

crystal (Crystal) – Lattice and atom information
qpts ((n_qpts, 3) float ndarray) – Q-point coordinates
frequencies ((n_qpts, 3*crystal.n_atoms) float Quantity) – Phonon frequencies, ordered according to increasing q-point number. Default units meV
eigenvectors ((n_qpts, 3*crystal.n_atoms, crystal.n_atoms, 3) complex ndarray) – Dynamical matrix eigenvectors
weights ((n_qpts,) float ndarray, optional) – The weight for each q-point. If None, equal weights are assumed

reorder_frequencies(reorder_gamma=True)¶

By doing a dot product of eigenvectors at adjacent q-points, determines which modes are most similar and reorders the frequencies at each q-point

Parameters:	reorder_gamma (bool) – Whether to reorder frequencies at gamma-equivalent points. If an analytical correction has been applied at the gamma points (i.e LO-TO splitting) mode assignments can be incorrect at adjacent q-points where the correction hasn’t been applied. So you might not want to reorder at gamma for some materials

calculate_structure_factor(scattering_lengths, dw=None)¶

Calculate the one phonon inelastic scattering for neutrons at each q-point

Parameters:	scattering_lengths (dictionary of float Quantity) – Dictionary of spin and isotope averaged coherent scattering length for each element in the structure in e.g. {‘O’: 5.803ureg(‘fm’), ‘Zn’: 5.680ureg(‘fm’)} dw (DebyeWaller) – A DebyeWaller exponent object
Returns:	sf – An object containing the structure factor for each q-point and phonon mode
Return type:	StructureFactor

Notes

The structure factor is defined as [1]:

\[F(Q, \nu) = \frac{b_\kappa}{M_{\kappa}^{1/2}\omega_{q\nu}^{1/2}} \ [Q\cdot\epsilon_{q\nu\kappa\alpha}]e^{Q{\cdot}r_\kappa}e^{-W}\]

Where \(\nu\) runs over phonon modes, \(\kappa\) runs over atoms, \(\alpha\) runs over the Cartesian directions, \(b_\kappa\) is the coherent neutron scattering length, \(M_{\kappa}\) is the atom mass, \(r_{\kappa}\) is the vector to atom \(\kappa\) in the unit cell, \(\epsilon_{q\nu\kappa\alpha}\) are the eigevectors, \(\omega_{q\nu}\) are the frequencies and \(e^{-W}\) is the Debye-Waller factor. Note that a factor N for the number of unit cells in the sample hasn’t been included, so the returned structure factor is per unit cell.

[1]	M.T. Dove, Structure and Dynamics, Oxford University Press, Oxford, 2003, 225-226

calculate_debye_waller(temperature)¶

Calculate the 3 x 3 Debye-Waller exponent for each atom over the q-points contained in this object

Parameters:	temperature (float Quantity) – Temperature
Returns:	dw – An object containing the 3x3 Debye-Waller exponent for each atom
Return type:	DebyeWaller

Notes

As part of the structure factor calculation, the anisotropic Debye-Waller factor is defined as:

\[e^{-W} = e^{-\sum_{\alpha\beta}\frac{W^{\kappa}_{\alpha\beta}Q_{\alpha}Q_{\beta}}{2}}\]

The Debye-Waller exponent is defined as \(W^{\kappa}_{\alpha\beta}\) and is independent of Q, so for efficiency can be precalculated to be used in the structure factor calculation. The Debye-Waller exponent is calculated by [2]

\[W^{\kappa}_{\alpha\beta} = \frac{1}{2N_{q}M_{\kappa}} \sum_{BZ}\frac{\epsilon_{q\nu\kappa\alpha}\epsilon^{*}_{q\nu\kappa\beta}} {\omega_{q\nu}} coth(\frac{\omega_{q\nu}}{2k_BT})\]

Where \(\nu\) runs over phonon modes, \(\kappa\) runs over atoms, \(\alpha,\beta\) run over the Cartesian directions, \(M_{\kappa}\) is the atom mass, \(\epsilon_{q\nu\kappa\alpha}\) are the eigenvectors, \(\omega_{q\nu}\) are the frequencies, \(\sum_{BZ}\) is a sum over the 1st Brillouin Zone, and \(N_q\) is the number of q-point samples in the BZ.

[2]	G.L. Squires, Introduction to the Theory of Thermal Neutron Scattering, Dover Publications, New York, 1996, 34-37

calculate_dos(dos_bins)¶

Calculates a density of states

Parameters:	dos_bins ((n_ebins + 1,) float Quantity) – The energy bin edges to use for calculating the DOS
Returns:	dos – A spectrum containing the energy bins on the x-axis and dos on the y-axis
Return type:	Spectrum1D

to_dict()¶

Convert to a dictionary. See QpointPhononModes.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 QpointPhononModes.from_dict keys

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

classmethod from_dict(d)¶

Convert a dictionary to a QpointPhononModes object

Parameters:

d (dict) –

A dictionary with the following keys/values:

’crystal’: dict, see Crystal.from_dict
’qpts’: (n_qpts, 3) float ndarray
’frequencies’: (n_qpts, 3*crystal.n_atoms) float ndarray
’frequencies_unit’: str
’eigenvectors’: (n_qpts, 3*crystal.n_atoms, crystal.n_atoms, 3) complex ndarray

There are also the following optional keys:

’weights’: (n_qpts,) float ndarray

Returns:

Return type: QpointPhononModes

classmethod from_json_file(filename)¶

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

Parameters:	filename (str) – The file to read from
Returns:
Return type:	QpointPhononModes

classmethod from_castep(filename)¶

Reads precalculated phonon mode data from a CASTEP .phonon file

Parameters:	filename (str) – The path and name of the .phonon file to read
Returns:
Return type:	QpointPhononModes

classmethod from_phonopy(path='.', phonon_name='band.yaml', phonon_format=None, summary_name='phonopy.yaml')¶

Reads precalculated phonon mode data from a Phonopy mesh/band/qpoints.yaml/hdf5 file. May also read from phonopy.yaml for structure information.

Parameters:	path (str, optional) – Path to directory containing the file(s) phonon_name (str, optional) – Name of Phonopy file including the frequencies and eigenvectors phonon_format ({'yaml', 'hdf5'} str, optional) – Format of the phonon_name file if it isn’t obvious from the phonon_name extension summary_name (str, optional) – Name of Phonopy summary file to read the crystal information from. Crystal information in the phonon_name file takes priority, but if it isn’t present, crystal information is read from summary_name instead
Returns:
Return type:	QpointPhononModes