Mode projection

The mode projection functionality in dynasor is mainly handled by two objects: the dynasor.ModeProjector class and dynasor.project_modes function. From the dynasor.ModeProjector there is access to data-objects representing a q-point dynasor.modes.qpoint.QPoint and from the q-point there is access to an object representing a particular band dynasor.modes.band.Band at that q-point. In addition, simple wrappers around the coordinates Q, P and F exists via dynasor.modes.complex_coordinate.ComplexCoordinate to easily set the amplitude and phase of a mode while preserving the \(Q(-q)=Q^*(q)\) symmetries. Internally dynasor wraps the primitive cell dynasor.modes.atoms.Prim and supercell dynasor.modes.atoms.Supercell. As a user only the dynasor.ModeProjector and dynasor.project_modes should need to be imported.

class dynasor.ModeProjector(primitive, supercell, force_constants)

The ModeProjector maps between real atomic displacements u and complex mode coordinates Q.

Some special python methods are implemented. The __str__ and __repr__ provides useful info. QPoint objects are representations of a single q-point and associated information and can be accessed either by call providing a reduced wavevector

>>> mp((1/2, 0, 0))  

or by index corresponding to reduced q-point accessible from ModeProjector.q_reduced

>>> mp[2]  

In addition to mode coordinates Q the class can also map the atomic velocities v to mode momenta P as well as atomic forces f to mode forces F. The class can also map back. This mapping is done using getters and setters. Internally only Q, P and F are stored

>>> Q = mp.get_Q()  
>>> mp.set_u(u)  

In addition, the forces corresponding to the harmonic forces can be accessed by ModeProjector.get_f_harmonic() and ModeProjector.get_F_harmonic(). For ASE Atoms objects the displacments etc. can be updated and applied by

>>> mp.update_from_atoms(atoms)  
>>> atoms = mp.get_atoms(harmonic_forces=False)  

The shapes for each property uses the follwoing varaibles

• N: Number of atoms in supercell

• Nu: unit cells (=N/Np)

• Np: primitive basis atoms (=N/Nu)

• Nb: bands (=Np*3)

• Nq: q-points (=Nu)

Please consult the documentation or the specific getters and setters to see the exact transformations used.

Units

The internal units in dynasor is Å, fs and eV. All frequencies are angular (a.k.a. “physicist’s” convention with 2π included). These are the units dynasor will expect and return. In e.g. print functions conventional units like fs, Å, THz, Da, meV will often be used.

Mass: The internal unit choice means that the mass unit is not Dalton (a.k.a. amu or u) but instead 0.009648533290731906 Da. We denote this with dmu (dynasor mass unit). In other words one Dalton is equal to 103.64269572045423 dmu units. This is typically not a problem since the masses provided via e.g. ASE Atoms objects are converted internally.

Waves: dynasor will always communicate and expect spatial (angular) frequencies in rad/Å and temporal (angular) frequencies in rad/fs. This follows the “physicist’s” convention as the factor of 2π is included in the wave vectors. For instance the wavelength is given by λ=2π/q.

Mode amplitudes: Usually the mode amplitudes is communicated in Å√dmu a.k.a. fs√eV.

Velocities: For modes the momenta (Å√dmu/fs or just √eV) is used but for atoms the velocities (Å/fs) are used.

Mode forces: The force is just defined as the rate of change of the canonical momenta so the unit would be Å√dmu/fs² (or √eV/fs)

Internal arrays

For the curious, the internal data arrays are

• primitive, supercell, force_constants (input)

• _q, q_minus (reduced q-points and which q-points are related by inversion)

• _D, _w2, _W (dynamical matrices, frequencies (ev/Ų/Da), polarization vectors)

• _X (eigenmodes which are mass weighted “polarization vectors” in the supercell)

property eigenmodes

The eigenmodes in the supercell as (Nq, Nb, N, 3)-array

The eigenmodes with masses included so that

Q = X u

where u is the supercell displacements

get_F()

The mode force amplitudes in eV/Å√dmu

get_F_harmonic()

The harmonic mode forces, computed as -omega^2 * Q, in eV/Å√dmu

get_P()

The mode momentum amplitudes in √eV

get_Q()

The mode coordinate amplitudes in Å√dmu

get_atoms(harmonic_forces=False)

Returns ase Atoms object with displacement, velocities, forces and harmonic energies set

Parameters

harmonic_forces – Whether the forces should be taken from the internal F or via -w2 Q

Return type

Atoms

get_f()

The atomic forces in eV/Å

get_f_harmonic()

The harmonic atomic forces for the current displacements

get_u()

The atomic displacements in Å

get_v()

The atomic velocities in Å/fs

property kinetic_energies

Kinetic energy per mode as (Nq, Nb)-array

The kinetic energies are defined as \(1/2PP^*\). Should equal \(1/2k_BT\) in equilibrium.`

property omegas

The frequencies of each mode in (rad/fs).

Negative frequencies correspond to imaginary frequencies

property polarizations

The polarization vectors for each mode (Nq, Nb, Np, 3)

property potential_energies

Potential energy per mode as (Nq, Nb)-array

The potential energies are defined as \(1/2QQ^*\) and should equal \(1/2k_BT\) in equilibrium for a harmonic system.`

property q_cartesian

The q-points in cartesian coordinates with unit of rad/Å (2π included)

property q_minus

The index of the corresponding counter-propagating mode (-q)

property q_reduced

The q-points in reduced coordinates.

A zone boundary mode would be e.g. (1/2, 0, 0)

classmethod read(file_name)

Return mp instance from pickle file that was saved using mp.write

Return type

ModeProjector

set_F(F)

Sets the internal mode forces F.

The functions ensures \(F(-q)=F^*(q)\)

set_P(P)

Sets the internal mode momenta P.

The functions ensures \(P(-q)=P^*(q)\)

set_Q(Q)

Sets the internal mode coordinates Q

The functions ensures Q(-q)=Q*(q)

set_f(f)

Sets the internal mode forces F given the forces f

F = X.conj() * f / m

Parameters

f – The atomic forces in eV/Å

set_u(u)

Sets the internal mode coordinates Q given the atomic displacements u

Q = X * u

Parameters

u – The atomic displacements in Å

set_v(v)

Sets the internal mode momenta P given the velocities v

P = X.conj() * v

Parameters

v – The atomic velocities in fs/Å

update_from_atoms(atoms)

Updates the ModeProjector with displacments, velocities and forces from an ASE Atoms object

Checks for attached calculator in first place and next for forces array.

If no data sets corresponding array to zeros.

The masses and velocities are converted to dynasor units internally.

property virial_energies

The virial energies per mode as (Nq, Nb)-array

The virial energies are defined here as \(-1/2QF\) which should have an expectation value of \(1/2k_BT\) per mode in equilibrium. For a harmonic system this is simply equal to the potential energy. This means that that the virial energy can be used to monitor the anharmonicity or define a measure of the potential energy.

write(file_name)

Uses pickle to write mp to file

class dynasor.modes.qpoint.QPoint(q, mp)

Representation of a single q-point and properties

The bands can be accessed by e.g. qp[2] to get band index 2 in the form of a dynasor.modes.band.Band object

property eigenmodes

Slice, see dynasor.ModeProjector

property index

q-point index corresponding to dynasor.ModeProjector.q_reduced

property is_real

If the q-point has purely real mode coordinates, ‘q=-q’

property kinetic_energies

Slice, see dynasor.ModeProjector

property omegas

Slice, see dynasor.ModeProjector

property polarizations

Slice, see dynasor.ModeProjector

property potential_energies

Slice, see dynasor.ModeProjector

property q_cartesian

Slice, see dynasor.ModeProjector

property q_minus: QPoint

The corresponding counter-propagating mode

property q_reduced

Slice, see dynasor.ModeProjector

property virial_energies

Slice, see dynasor.ModeProjector

property wavelength

Wavelength of mode in Å

property wavenumber

Wavenumber of mode in rad/Å

class dynasor.modes.band.Band(q_index, band_index, mp)

Represents properties of a single band belonging to a q-point

To modify the coordinates Q, P and F use band.Q and see doc for dynasor.modes.complex_coordinate.ComplexCoordinate

property F: F

The mode force

property P: P

The mode momentum

property Q: Q

The mode coordinate

property eigenmode

Slice, see doc for ModeProjector

property index

The index of the band at the q-point

property kinetic_energy

Slice, see doc for ModeProjector

property omega

Slice, see doc for ModeProjector

property polarization

Slice, see doc for ModeProjector

property potential_energy

Slice, see doc for ModeProjector

property qpoint: QPoint

The q-point to which the band belongs

property virial_energy

Slice, see doc for ModeProjector

class dynasor.modes.complex_coordinate.ComplexCoordinate(q, s, mp)

Class to work with the in general complex coordinates of lattice dynamics.

Can be cast to complex by complex(cc)

Example

>>> cc.complex = 1.7j  
>>> cc.amplitude = 3.8  
>>> cc.phase = 2*np.pi  

property amplitude

The amplitude of the complex coordinate

property band: Band

The band to which this coordinate belongs

property complex

The complex coordinate as a complex number

property phase

The phase of the complex coordinate

dynasor.project_modes(traj, modes, ideal_supercell, check_mic=True)

Projects an atomic trajectory onto set of phonon modes

Parameters

• traj (Trajectory) – Input trajectory

• modes (ndarray[Any, dtype[float]]) – Modes to project on (1, Nlambda, N, 3) array where Nm is the number of modes and N is the number of atoms in the supercell

• ideal_supercell (Atoms) – Used to find atomic displacements and should correspond to the ideal structure. Be careful not to mess up the permutation

• check_mic – Whether to wrap the displacements or not, faster if no wrap.

Return type

Tuple[ndarray[Any, dtype[float]], ndarray[Any, dtype[float]]]

Returns

• Q – mode coordinates, complex ndarray (length of traj, number of modes)

• P – mode momenta, complex ndarray (length of traj, number of modes)

class dynasor.modes.atoms.Prim(atoms)

class dynasor.modes.atoms.Supercell(supercell, prim)

The supercell takes care of some mappings between the primitive and repeated structure

In particular the P-matrix connecting the cells as well as the offset-index of each atom is calculated.

Note that the positions can not be revovered as offset x cell + basis since the atoms gets wrapped.

Parameters

• supercell – some ideal repetition of the prim cell and possible wrapping as either ASEAtoms

• prim – primitive cell as either ASEAtoms

property P

P-matrix is defined as dot(P, prim.cell) = supercell.cell

property P_inv

Inverse of P

property indices

The basis index of each atom

property n_cells

Number of unit cells

property offsets

The offset of each atom