Infrared intensities

Infrared is an extension of Vibrations, in addition to the vibrational modes, also the infrared intensities of the modes are calculated for an Atoms object.

class ase.vibrations.Infrared(atoms, indices=None, name='ir', delta=0.01, nfree=2, directions=None)[source]

Class for calculating vibrational modes and infrared intensities using finite difference.

The vibrational modes are calculated from a finite difference approximation of the Dynamical matrix and the IR intensities from a finite difference approximation of the gradient of the dipole moment. The method is described in:

D. Porezag, M. R. Pederson: “Infrared intensities and Raman-scattering activities within density-functional theory”, Phys. Rev. B 54, 7830 (1996)

The calculator object (calc) linked to the Atoms object (atoms) must have the attribute:

>>> calc.get_dipole_moment(atoms)

In addition to the methods included in the Vibrations class the Infrared class introduces two new methods; get_spectrum() and write_spectra(). The summary(), get_energies(), get_frequencies(), get_spectrum() and write_spectra() methods all take an optional method keyword. Use method=’Frederiksen’ to use the method described in:

T. Frederiksen, M. Paulsson, M. Brandbyge, A. P. Jauho: “Inelastic transport theory from first-principles: methodology and applications for nanoscale devices”, Phys. Rev. B 75, 205413 (2007)

atoms: Atoms object

The atoms to work on.

indices: list of int

List of indices of atoms to vibrate. Default behavior is to vibrate all atoms.

name: str

Name to use for files.

delta: float

Magnitude of displacements.

nfree: int

Number of displacements per degree of freedom, 2 or 4 are supported. Default is 2 which will displace each atom +delta and -delta in each cartesian direction.

directions: list of int

Cartesian coordinates to calculate the gradient of the dipole moment in. For example directions = 2 only dipole moment in the z-direction will be considered, whereas for directions = [0, 1] only the dipole moment in the xy-plane will be considered. Default behavior is to use the dipole moment in all directions.


>>> from import read
>>> from ase.calculators.vasp import Vasp
>>> from ase.vibrations import Infrared
>>> water = read('water.traj')  # read pre-relaxed structure of water
>>> calc = Vasp(prec='Accurate',
...             ediff=1E-8,
...             isym=0,
...             idipol=4,       # calculate the total dipole moment
...             dipol=water.get_center_of_mass(scaled=True),
...             ldipol=True)
>>> water.calc = calc
>>> ir = Infrared(water)
>>> ir.summary()
Mode    Frequency        Intensity
#    meV     cm^-1   (D/Å)^2 amu^-1
0   16.9i    136.2i     1.6108
1   10.5i     84.9i     2.1682
2    5.1i     41.1i     1.7327
3    0.3i      2.2i     0.0080
4    2.4      19.0      0.1186
5   15.3     123.5      1.4956
6  195.5    1576.7      1.6437
7  458.9    3701.3      0.0284
8  473.0    3814.6      1.1812
Zero-point energy: 0.573 eV
Static dipole moment: 1.833 D
Maximum force on atom in `equilibrium`: 0.0026 eV/Å

This interface now also works for calculator ‘siesta’, (added get_dipole_moment for siesta).


>>> #!/usr/bin/env python3
>>> from import read
>>> from ase.calculators.siesta import Siesta
>>> from ase.vibrations import Infrared
>>> bud = read('')
>>> calc = Siesta(label='bud',
...       meshcutoff=250 * Ry,
...       basis='DZP',
...       kpts=[1, 1, 1])
>>> calc.set_fdf('DM.MixingWeight', 0.08)
>>> calc.set_fdf('DM.NumberPulay', 3)
>>> calc.set_fdf('DM.NumberKick', 20)
>>> calc.set_fdf('DM.KickMixingWeight', 0.15)
>>> calc.set_fdf('SolutionMethod',      'Diagon')
>>> calc.set_fdf('MaxSCFIterations', 500)
>>> calc.set_fdf('PAO.BasisType',  'split')
>>> #50 meV = 0.003674931 * Ry
>>> calc.set_fdf('PAO.EnergyShift', 0.003674931 * Ry )
>>> calc.set_fdf('LatticeConstant', 1.000000 * Ang)
>>> calc.set_fdf('WriteCoorXmol',       'T')
>>> bud.calc = calc
>>> ir = Infrared(bud)
>>> ir.summary()
get_spectrum(start=800, end=4000, npts=None, width=4, type='Gaussian', method='standard', direction='central', intensity_unit='(D/A)2/amu', normalize=False)[source]

Get infrared spectrum.

The method returns wavenumbers in cm^-1 with corresponding absolute infrared intensity. Start and end point, and width of the Gaussian/Lorentzian should be given in cm^-1. normalize=True ensures the integral over the peaks to give the intensity.

write_spectra(out='ir-spectra.dat', start=800, end=4000, npts=None, width=10, type='Gaussian', method='standard', direction='central', intensity_unit='(D/A)2/amu', normalize=False)[source]

Write out infrared spectrum to file.

First column is the wavenumber in cm^-1, the second column the absolute infrared intensities, and the third column the absorbance scaled so that data runs from 1 to 0. Start and end point, and width of the Gaussian/Lorentzian should be given in cm^-1.