SpecificHeatCapacity

class SpecificHeatCapacity(md, start_time=None, end_time=None, atom_selection=None, time_resolution=None, min_temperature=PhysicalQuantity(0.0, K), max_temperature=PhysicalQuantity(1000.0, K), temperature_resolution=None, info_panel=None)

Class for calculating the specific heat capacity of a system from a MD simulation.

Parameters:

  • md_trajectory (MDTrajectory) – The MDtrajectory to calculate the specific heat capacity for.

  • 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

  • atom_selection (PeriodicTableElement | str | list of ints) – Only include contributions from this selection. The atoms can be selected by element i.e. PeriodicTableElement, tag or a list of atomic indices.
    Default: All atoms.

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

  • min_temperature (PhysicalQuantity of type temperature) – The minimum temperature.
    Default: 0*Kelvin

  • max_temperature (PhysicalQuantity of type temperature) – The maximum temperature.
    Default: 1000*Kelvin

  • temperature_resolution (PhysicalQuantity of type temperature) – The difference between individual temperature points.
    Default: 100 points

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

frequencies()
Returns:

The frequencies associated with the vibrational density of states.

Return type:

PhysicalQuantity of type 1/time

specificHeatCapacity()
Returns:

The specific heat capacity calculated over the requested temperature range.

Return type:

PhysicalQuantity of type energy per mass per temperature

temperatures()
Returns:

The array of temperature values that the specific heat capacity has been evaluated at.

Return type:

PhysicalQuantity of type temperature

vibrationalDensityOfStates()
Returns:

The vibrational density-of-states, calculated from the Fourier transformation of the velocities.

Return type:

PhysicalQuantity of type time

Usage Examples

This example calculates the specific heat capacity of bulk silicon using classical molecular dynamics:

md_trajectory = MolecularDynamics(
    bulk_configuration,
    constraints=constraints,
    trajectory_filename='si_bulk_specific_heat.hdf5',
    steps=5000,
    log_interval=1,
    method=method
)

bulk_configuration = md_trajectory.lastImage()

# Calculate specific heat capacity.
specific_heat_capacity = SpecificHeatCapacity(
    md_trajectory,
    min_temperature=0.0*Kelvin,
    max_temperature=1000.0*Kelvin,
)

# Plot specific heat capacity.
import pylab
temperatures = specific_heat_capacity.temperatures().inUnitsOf(Kelvin)
specific_heats = specific_heat_capacity.specificHeatCapacity().inUnitsOf(Joule/gram/Kelvin)
pylab.plot(temperatures, specific_heats)
pylab.xlabel('Temperature (K)')
pylab.ylabel('Specific Heat Capacity (J/g/K)')
pylab.savefig('si_bulk_specific_heat.png')

si_bulk_specific_heat.py

The resulting plot of heat capacity versus temperature shows the correct limiting behavior and low and high temperature. This is achieved, when calculating the heat capacity, by modeling the system as a collection of quantum harmonic oscillators, whose frequencies are determined from the vibrational density of states.

../../../_images/si_bulk_specific_heat.png

Notes

The heat capacity, is calculated by modeling the system as a collection of quantum harmonic oscillators, whose frequencies are determined using the vibrational density of states. [1] The specific heat, as a function of temperature, is

\[C_{\mathrm{v}}(T) = \frac{h^2}{m k_{\mathrm{B}} T^2} \int_0^\infty \frac{\nu^2 \exp(h \nu / k_{\mathrm{B}} T)}{(\exp(h \nu/k_{\mathrm{B}} T) -1)^2} S(\nu) \mathrm{d}\nu\]

where, \(T\) is the temperature, \(\nu\) is the vibrational frequency, \(h\) is Planck’s constant, \(m\) is the total mass of the atoms in the unit cell, \(k_{\mathrm{B}}\) is Boltzmann’s constant, and \(S(\nu)\) is the vibrational density of states (VDOS). See the VibrationalDensityOfStates notes for details on how the VDOS is calculated.