InfraredSpectrum

class InfraredSpectrum(md_trajectory, start_time=None, end_time=None, time_resolution=None, wavenumber_limit=None, temperature=None, path_length=None, smoothing_factor=None, info_panel=None)

Class for calculating the infrared vibrational spectrum from an MD simulation.

Parameters:

  • md_trajectory (MDTrajectory) – The MDtrajectory from which the infrared is calculated. If present the dipole measurement is used, otherwise information is taken from each image.

  • start_time (PhysicalQuantity of type time) – The start time.
    Default: 0.0 * fs.

  • end_time (PhysicalQuantity of type time) – The end time.
    Default: The last time frame.

  • time_resolution (PhysicalQuantity of type time) – The time interval between snapshots or measurements in the MD trajectory that are included in the analysis.

  • wavenumber_limit (PhysicalQuantity of type inverse length) – The maximum wavenumber displayed in the spectrum.
    Default: 5000 * cm^-1.

  • temperature (PhysicalQuantity of type temperature) – The temperature used in calculating absorption intensity.
    Default: 300 * Kelvin.

  • path_length (PhysicalQuantity of type length) – The length of sample the light passes through. Longer path lengths decreases the amount of light transmitted
    Default: 0.001 * cm.

  • smoothing_factor (int) – Apply smoothing to the spectrum by averaging n points.
    Default: 1.

  • info_panel (InfoPanel (Plot2D)) – Info panel to show the calculation progress.
    Default: No info panel.

absorption()
Returns:

The infrared absorption for each wavenumber in the spectrum.

Return type:

ndarray

absorptionCoefficient()
Returns:

The absorption coefficient for each wavenumber in the infrared spectrum.

Return type:

PhysicalQuantity of type inverse length

frequencies()
Returns:

The frequency of each datapoint in the infrared spectrum.

Return type:

PhysicalQuantity of type inverse time

transmission()
Returns:

The infrared transmission for each wavenumber in the spectrum.

Return type:

ndarray

wavenumbers()
Returns:

The wavenumber of each datapoint in the infrared spectrum.

Return type:

PhysicalQuantity of type inverse length

Usage Examples

Load an MDTrajectory containing a simulation of a liquid electrolyte, and calculate its infrared (IR) vibrational spectrum. The IR transmission through a sample with a path length of 100 micron is then plotted.


md_trajectory = nlread('electrolyte_trajectory.hdf5')[-1]

# Create the IR Spectrum
ir_spectrum = InfraredSpectrum(
    md_trajectory,
    wavenumber_limit=4000*cm**-1,
    temperature=300*Kelvin,
    path_length=0.01*cm,
    smoothing_factor=1,
)

# Create a model.
model = Plot.PlotModel(x_unit=cm**-1)
model.framing().setTitle('Electrolyte IR Spectra')
model.xAxis().setLabel('Wavenumber')
model.yAxis().setLabel('Transmission')
model.legend().setVisible(False)

# Create and add line.
line = Plot.Line(ir_spectrum.wavenumbers(), ir_spectrum.transmission())
line.setLabel('Spectrum')
line.setColor('blue')
model.addItem(line)

# Set x-axis limits.
model.setLimits('x', 4000*cm**-1, 500*cm**-1)

# Show the plot for interactive editing.
Plot.show(model)

infrared_spectrum.py

Notes

From a molecular dynamics trajectory, the infrared vibrational spectrum can be calculated by obtaining the Fourier transform of the system dipole autocorrelation function [1][2]. The adsorption coefficient \(\alpha(\omega)\) is given as:

\[\alpha(\omega) = \frac{\omega c \tanh\left(0.5 \beta \hbar \omega\right)} {3 \epsilon_0 \hbar n(\omega) V} \int_0^\infty \langle \mathbf{M}(0) \cdot \mathbf{M}(t) \rangle \cos(\omega t) dt\]

Here \(\omega\) is the vibrational frequency, \(\beta\) is the Boltzmann factor \(1/kT\), \(n(\omega)\) is the frequency dependent refractive index and \(V\) is the configuration volume. Note that in calculating the infrared spectrum the refractive index is assumed to be constant in the infrared frequency range. The total system dipole \(\mathbf{M}\) is also given as:

\[\mathbf{M} = \sum_{i=0}^N q_i\mathbf{r}_i\]

Here \(q_i\) and \(\mathbf{r}_i\) are the atomic charge and position respectively, with the sum being taken over all atoms in the system. Finally the transmission spectrum \(T(\omega)\) and absorption spectrum \(A(\omega)\) can be estimated with the relationship

\[T(\omega) = 1 - A(\omega) = \exp(-\alpha(\omega) l)\]

where \(l\) is a path length that represents the thickness of the sample.

In order to estimate the change in dipole moment during the simulation required to calculate the adsorption spectrum, partial charges must be set on the input configuration. Partial atomic charges can be added to a configuration using the setPartialCharges method. These charges are used to calculate the dipole in each image. Obtaining reliable vibrational spectra also requires saving the dipole moment at high frequency. As only the total system dipole is required, this is most efficiently done by saving the dipole as a measurement using the dipole_vector molecular dynamics measurement. Where appropriate dipole measurement data is contained within the input trajectory, this will be used in preference to calculating the dipole from each image. This allows using high-frequency dipole data with a minimal memory and computational overhead.

It should also be noted that since the analysis is based on the assigned partial charges, the relative size of different peaks may be dependent on the charge assignment. Of particular importance is the magnitude of the dipole moment change during each vibrational mode. Vibrations that do not result in a change in dipole moment do not result in peaks in the infrared spectrum.

The infrared spectrum analysis contains some optional arguments. The argument wavenumber_limit sets the largest wavenumber considered in the vibrational spectrum analysis. Similarly the argument temperature sets the temperature used in calculating the absorption coefficient, and is commonly set to the approximate temperature of the molecular dynamics simulation. The amount of infrared radiation absorbed or transmitted is controlled with the path_length variable. This can be used to enhance or reduce the size of peaks in the absorption spectrum to highlight important features. Finally the argument smoothing_factor controls the resolution of the final spectral data. Increasing the smoothing factor down-samples the resulting spectrum, averaging each value over the given number of points. This can help to reduce noise in the calculated spectrum.