PhononTransmissionSpectrum

class PhononTransmissionSpectrum(dynamical_matrix, configuration=None, energies=None, qpoints=None, qpoints_weights=None, self_energy_calculator=None, infinitesimal=None)

Class for representing the phonon transmission spectrum for a given device configuration.

Parameters:
  • dynamical_matrix (DynamicalMatrix) – The dynamical matrix to use for the calculation.

  • configuration (DeviceConfiguration) – The device configuration for which the vibrational modes should be calculated. If no calculator is attached to the configuration, the calculator from dynamical_matrix is used. Currently, configuration must be the same as the configuration from dynamical_matrix and is therefore not needed.
    Default: The configuration from dynamical_matrix

  • energies (PhysicalQuantity of type energy.) – The energies for which the transmission spectrum should be calculated.
    Default: numpy.linspace(0, 0.5, 101)*eV

  • qpoints – The q-points for which to calculate the transmission spectrum. This can either be given as a q-point grid or a list of fractional q-points, e.g., [[0.0, 0.0, 0.0], [0.0, 0.0, 0.1], ...]. Note that the q-points must be in the xy-plane.
    Default: MonkhorstPackGrid(na, nb) where (na, nb) is the sampling used for the self-consistent calculation.

  • qpoints_weights – The weight of each q-point for which the transmission spectrum should be calculated.
    Default: The weights corresponding to the MonkhorstPackGrid, or a list of [1.0, 1.0, …] if q-points are specified as floats.

  • self_energy_calculator – The SelfEnergyCalculator to be used for the transmission spectrum.
    Default: RecursionSelfEnergy(storage_strategy=NoStorage()) :type self_energy_calculator: DirectSelfEnergy | RecursionSelfEnergy | SparseRecursionSelfEnergy | KrylovSelfEnergy

  • infinitesimal (PhysicalQuantity of type energy.) – Small positive energy, used to move the transmission calculation away from the real axis. This is only relevant for recursion-style self-energy calculators.
    Default: 1.0e-6*eV

adaptive()
Returns:

True if adaptive grids are used, otherwise False.

Return type:

bool

energies()
Returns:

The energies used in this transmission spectrum

Return type:

PhysicalQuantity of type energy

evaluate(qpoints=None)

Evaluates the q-point averaged transmission spectrum.

Parameters:

qpoints (MonkhorstPackGrid | RegularKpointGrid | list(n) of list(3) of float.) – The set of q-points to average over. The q-points must be compatible with the q-points used in the constructor.
Default: Average over all q-points.

Returns:

The q-point averaged transmission spectrum.

Return type:

array(n_energies)

infinitesimal()
Returns:

The infinitesimal used for calculating the transmission.

Return type:

PhysicalQuantity of type energy

metatext()
Returns:

The metatext of the object or None if no metatext is present.

Return type:

str | None

nlprint(stream=None)

Print a string containing an ASCII table useful for plotting the AnalysisSpin object.

Parameters:

stream (python stream) – The stream the table should be written to.
Default: NLPrintLogger()

qpoints()
Returns:

The q-points used in this transmission spectrum. The shape of the array is (number_of_qpoints, 3)

Return type:

numpy.ndarray

qpointsWeights()
Returns:

The weights of each q-point used in this transmission spectrum.

Return type:

list.

setMetatext(metatext)

Set a given metatext string on the object.

Parameters:

metatext (str | None) – The metatext string that should be set. A value of “None” can be given to remove the current metatext.

thermalConductance(phonon_temperature=None)

Calculate the thermal conductance of the phonon spectrum.

Parameters:

phonon_temperature (Positive PhysicalQuantity of type temperature.) – The temperature of the phonon system.
Default: 300*Kelvin

Returns:

The thermal conductance of the phonon spectrum.

Return type:

PhysicalQuantity with unit Joule/Second/Kelvin

transmission()

Return the transmission coefficients of the phonon transmission spectrum.

Returns:

The transmission coefficients for each energy and q-point. The shape of the array is (number_of_energies, number_of_qpoints)

Return type:

numpy.ndarray

uniqueString()

Return a unique string representing the state of the object.

Usage Examples

Calculate the PhononTransmissionSpectrum for a generic bulk silicon device:

# -------------------------------------------------------------
# Left electrode
# -------------------------------------------------------------

# Set up lattice
vector_a = [5.4306, 0.0, 0.0]*Angstrom
vector_b = [0.0, 5.4306, 0.0]*Angstrom
vector_c = [0.0, 0.0, 5.4306]*Angstrom
left_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
left_electrode_elements = [Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                           Silicon]

# Define coordinates
left_electrode_coordinates = [[ 0.      ,  0.      ,  0.678825],
                              [ 2.7153  ,  2.7153  ,  0.678825],
                              [ 1.35765 ,  1.35765 ,  2.036475],
                              [ 4.07295 ,  4.07295 ,  2.036475],
                              [ 2.7153  ,  0.      ,  3.394125],
                              [ 0.      ,  2.7153  ,  3.394125],
                              [ 4.07295 ,  1.35765 ,  4.751775],
                              [ 1.35765 ,  4.07295 ,  4.751775]]*Angstrom

# Set up configuration
left_electrode = BulkConfiguration(
    bravais_lattice=left_electrode_lattice,
    elements=left_electrode_elements,
    cartesian_coordinates=left_electrode_coordinates
    )

# -------------------------------------------------------------
# Right electrode
# -------------------------------------------------------------

# Set up lattice
vector_a = [5.4306, 0.0, 0.0]*Angstrom
vector_b = [0.0, 5.4306, 0.0]*Angstrom
vector_c = [0.0, 0.0, 5.4306]*Angstrom
right_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
right_electrode_elements = [Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                            Silicon]

# Define coordinates
right_electrode_coordinates = [[ 0.      ,  0.      ,  0.678825],
                               [ 2.7153  ,  2.7153  ,  0.678825],
                               [ 1.35765 ,  1.35765 ,  2.036475],
                               [ 4.07295 ,  4.07295 ,  2.036475],
                               [ 2.7153  ,  0.      ,  3.394125],
                               [ 0.      ,  2.7153  ,  3.394125],
                               [ 4.07295 ,  1.35765 ,  4.751775],
                               [ 1.35765 ,  4.07295 ,  4.751775]]*Angstrom

# Set up configuration
right_electrode = BulkConfiguration(
    bravais_lattice=right_electrode_lattice,
    elements=right_electrode_elements,
    cartesian_coordinates=right_electrode_coordinates
    )

# -------------------------------------------------------------
# Central region
# -------------------------------------------------------------

# Set up lattice
vector_a = [5.4306, 0.0, 0.0]*Angstrom
vector_b = [0.0, 5.4306, 0.0]*Angstrom
vector_c = [0.0, 0.0, 21.7224]*Angstrom
central_region_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
central_region_elements = [Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                           Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                           Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                           Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
                           Silicon, Silicon, Silicon, Silicon]

# Define coordinates
central_region_coordinates = [[  0.      ,   0.      ,   0.678825],
                              [  2.7153  ,   2.7153  ,   0.678825],
                              [  1.35765 ,   1.35765 ,   2.036475],
                              [  4.07295 ,   4.07295 ,   2.036475],
                              [  2.7153  ,   0.      ,   3.394125],
                              [  0.      ,   2.7153  ,   3.394125],
                              [  4.07295 ,   1.35765 ,   4.751775],
                              [  1.35765 ,   4.07295 ,   4.751775],
                              [  0.      ,   0.      ,   6.109425],
                              [  2.7153  ,   2.7153  ,   6.109425],
                              [  1.35765 ,   1.35765 ,   7.467075],
                              [  4.07295 ,   4.07295 ,   7.467075],
                              [  2.7153  ,   0.      ,   8.824725],
                              [  0.      ,   2.7153  ,   8.824725],
                              [  4.07295 ,   1.35765 ,  10.182375],
                              [  1.35765 ,   4.07295 ,  10.182375],
                              [  0.      ,   0.      ,  11.540025],
                              [  2.7153  ,   2.7153  ,  11.540025],
                              [  1.35765 ,   1.35765 ,  12.897675],
                              [  4.07295 ,   4.07295 ,  12.897675],
                              [  2.7153  ,   0.      ,  14.255325],
                              [  0.      ,   2.7153  ,  14.255325],
                              [  4.07295 ,   1.35765 ,  15.612975],
                              [  1.35765 ,   4.07295 ,  15.612975],
                              [  0.      ,   0.      ,  16.970625],
                              [  2.7153  ,   2.7153  ,  16.970625],
                              [  1.35765 ,   1.35765 ,  18.328275],
                              [  4.07295 ,   4.07295 ,  18.328275],
                              [  2.7153  ,   0.      ,  19.685925],
                              [  0.      ,   2.7153  ,  19.685925],
                              [  4.07295 ,   1.35765 ,  21.043575],
                              [  1.35765 ,   4.07295 ,  21.043575]]*Angstrom

# Set up configuration
central_region = BulkConfiguration(
    bravais_lattice=central_region_lattice,
    elements=central_region_elements,
    cartesian_coordinates=central_region_coordinates
    )

device_configuration = DeviceConfiguration(
    central_region,
    [left_electrode, right_electrode]
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------

potentialSet = StillingerWeber_Si_1985()
calculator = TremoloXCalculator(parameters=potentialSet)

device_configuration.setCalculator(calculator)
device_configuration.update()

# -------------------------------------------------------------
# Phonon transmission spectrum
# -------------------------------------------------------------
phonon_transmission_spectrum = PhononTransmissionSpectrum(
    configuration=device_configuration,
    energies=numpy.linspace(0,0.5,101)*eV,
    qpoints=MonkhorstPackGrid(5,5),
    infinitesimal=1e-06*eV,
    self_energy_calculator=RecursionSelfEnergy(),
    )

nlsave('si_phonon_transmission.nc', phonon_transmission_spectrum)

phonon_transmission.py

Calculate the thermal conductance for different phonon temperatures:

temp = numpy.linspace(100, 1000, 19)
conductances = []

for t in temp:
    conductance = phonon_transmission_spectrum.thermalConductance(t*Kelvin)
    conductances.append(conductance)

print("Phonon temperature   Thermal Conductance")
output_format = "      %6.0f K              %s"
for i in range(len(temp)):
    print(output_format % (temp[i], conductances[i]))

thermal_conductance.py

Notes

The calculation of the phonon transmission essentially follows the formalism employed to calculate the electronic TransmissionSpectrum, by using the substitutes:

\[E \mathbf{S}_C \rightarrow \omega^2\mathbf{M} \, ,\]

and

\[\mathbf{H}_C \rightarrow \mathbf{K}_C \, ,\]

where \(\mathbf{M}\) is the diagonal matrix with the element masses and \(\mathbf{K}_C\) is the dynamical matrix of the central region. Similar replacements are made for the electrodes. For more details see, e.g. Ref. [1].

Note

The phonon transmission calculated via the non-equilibrium Green’s function formalism gives the ideal, ballistic contributions to the thermal transport through the device. It does not account for anharmonic effects, such as phonon-phonon-scattering. To take such effects into account as well, you should consider other approaches, e.g. based on molecular dynamics simulations.

To export the data of a transmission spectrum, use the method nlsave.

The details and parameters of the force constant calculation to obtain the dynamical matrix, can be found under DynamicalMatrix.

Unit of transmission

A fully transmitting channel contributes \(1\). See also, TransmissionEigenvalues.

For a PhononTransmissionSpectrum, the first index (i.e. the spin index) in the list, that is returned by transmission(), must always be 0, as phonons do not have a spin.

Current

The heat current can be calculated from the transmission coefficients,

\[J(T_L, T_R) = \frac{\hbar}{2 \pi} \int_0^\infty \omega T(\omega) \left[ n_B(T_L) - n_B (T_R) \right] d\omega \, ,\]

where \(n_B=(e^{\hbar \omega/k_B T}-1)^{-1}\) is the Bose-Einstein distribution of the phonons, \(T_{L/R}\) are the phonon temperatures of the left/right electrode, and \(T(\omega)\) is the transmission coefficient at frequency \(\omega\) .

Thermal conductance

In the limit of small temperature differences, \(\Delta T\), between the left and right electrodes, with \(T_L = T+\Delta T/2\) and \(T_R = T-\Delta T/2\), the thermal conductance at average temperature, \(T\), is calculated from the phonon transmission as follows:

\[\kappa_{ph} (T) = \frac{\hbar}{2 \pi k_B T^2} \int_0^\infty \omega^2 T (\omega)\frac{e^{\hbar\omega/k_B T}}{(e^{\hbar\omega/k_B T}-1)^2} d\omega.\]