Phonon Interactions

1. Introduction: The Necessity of Anharmonicity

In Chapter 1, we developed the harmonic approximation where the lattice potential was expanded to second order in atomic displacements. While this approximation successfully explains many phonon properties (dispersion relations, density of states), it has fundamental limitations:

Real materials exhibit all of these phenomena, demonstrating that the harmonic approximation, while useful, is incomplete. The solution lies in including anharmonic terms in the lattice potential expansion.

1.1 Anharmonic Potential Energy

The total potential energy of a crystal can be expanded in powers of atomic displacements \(\mathbf{u}\) from equilibrium positions:

Where:

• \(V_0\) is the potential energy at equilibrium
• The linear term vanishes because we expand around the equilibrium (minimum energy) position
• The second-order term defines the harmonic approximation
• Higher-order terms introduce anharmonic effects

1.2 Physical Interpretation of Anharmonic Terms

Each order of anharmonicity has a distinct physical meaning:

2. Fermi's Golden Rule for Phonon Scattering

To calculate phonon scattering rates due to anharmonicity, we use time-dependent perturbation theory. The key tool is Fermi's golden rule, which gives the transition rate between quantum states.

2.1 Statement of Fermi's Golden Rule

For a system initially in state \(|i\rangle\) transitioning to final state \(|f\rangle\) due to a perturbation \(\hat{H}'\), the transition rate is:

Where:

• \(W_{i \to f}\) is the transition rate (probability per unit time)
• \(\langle f | \hat{H}' | i \rangle\) is the matrix element of the perturbation
• \(\delta(E_f - E_i)\) is the Dirac delta function enforcing energy conservation

2.2 Application to Phonon-Phonon Scattering

For phonons, we need to express the anharmonic potential in terms of phonon creation and annihilation operators. Recall from quantum mechanics that atomic displacements can be written as:

Where \(a_{\mathbf{q}s}\) and \(a_{\mathbf{q}s}^\dagger\) are phonon annihilation and creation operators. Substituting this into the third-order anharmonic term yields:

Where \(V_3\) is the Fourier-transformed third-order force constant. Expanding this product gives eight terms corresponding to different phonon processes (creation/annihilation combinations).

3. Three-Phonon Processes

The third-order anharmonic Hamiltonian allows for two fundamental types of three-phonon processes:

3.1 Decay Processes

In a decay process, one phonon spontaneously splits into two phonons:

This process corresponds to the term \(a_{\mathbf{q}_1}^\dagger a_{\mathbf{q}_2} a_{\mathbf{q}_3}\) in the Hamiltonian, which annihilates two phonons and creates one.

Conservation Laws:

• Energy conservation: \(\hbar\omega_{\mathbf{q}_1 s_1} = \hbar\omega_{\mathbf{q}_2 s_2} + \hbar\omega_{\mathbf{q}_3 s_3}\)
• Crystal momentum conservation: \(\mathbf{q}_1 = \mathbf{q}_2 + \mathbf{q}_3 + \mathbf{G}\)

The reciprocal lattice vector \(\mathbf{G}\) appears because crystal momentum is only defined modulo a reciprocal lattice vector. If \(\mathbf{G} = 0\), the process is called a normal process (N-process). If \(\mathbf{G} \neq 0\), it is an umklapp process (U-process).

3.2 Coalescence Processes

In a coalescence process, two phonons merge to form a single phonon:

This corresponds to the term \(a_{\mathbf{q}_1} a_{\mathbf{q}_2} a_{\mathbf{q}_3}^\dagger\) in the Hamiltonian.

Conservation Laws:

• Energy conservation: \(\hbar\omega_{\mathbf{q}_1 s_1} + \hbar\omega_{\mathbf{q}_2 s_2} = \hbar\omega_{\mathbf{q}_3 s_3}\)
• Crystal momentum conservation: \(\mathbf{q}_1 + \mathbf{q}_2 = \mathbf{q}_3 + \mathbf{G}\)

3.3 Scattering Rate Calculation

Using Fermi's golden rule, the scattering rate for a decay process is:

Where \(n_{\mathbf{q}s}\) is the Bose-Einstein occupation number:

The factors \((n + 1)\) come from the quantum mechanical matrix elements for phonon creation (stimulated emission plus spontaneous emission).

Similarly, for coalescence:

The total scattering rate is:

4. Normal vs Umklapp Processes

The distinction between normal (N) and umklapp (U) processes is crucial for understanding thermal conductivity.

4.1 Normal Processes (N-processes)

Normal processes conserve crystal momentum exactly (\(\mathbf{G} = 0\)):

4.2 Umklapp Processes (U-processes)

Umklapp processes (from German "umklappen" = to flip over) involve a reciprocal lattice vector \(\mathbf{G} \neq 0\):

4.3 Geometric Interpretation

The following diagram illustrates N and U processes in wavevector space:

Consider a decay process where \(\mathbf{q}_1\) splits into \(\mathbf{q}_2\) and \(\mathbf{q}_3\):

• N-process: If \(\mathbf{q}_2 + \mathbf{q}_3\) falls within the first Brillouin zone, we simply have \(\mathbf{q}_1 = \mathbf{q}_2 + \mathbf{q}_3\).
• U-process: If \(\mathbf{q}_2 + \mathbf{q}_3\) falls outside the first Brillouin zone, we must subtract a reciprocal lattice vector \(\mathbf{G}\) to bring it back: \(\mathbf{q}_1 = \mathbf{q}_2 + \mathbf{q}_3 - \mathbf{G}\).

4.4 Temperature Dependence and Thermal Conductivity

The rate of U-processes has strong temperature dependence:

Where \(\Theta_D\) is the Debye temperature and \(\beta\) is a numerical factor (typically \(\beta \approx 2\) to 3).

This leads to the characteristic behavior of thermal conductivity \(\kappa(T)\) in pure crystals:

5. Temperature Dependence of Phonon Scattering

The temperature dependence of phonon-phonon scattering arises from the Bose-Einstein occupation factors.

5.1 Occupation Number Temperature Dependence

The Bose-Einstein distribution is:

Limiting behaviors:

• Low temperature (\(k_B T \ll \hbar\omega\)): \(n \approx e^{-\hbar\omega/k_B T}\) (exponentially small)
• High temperature (\(k_B T \gg \hbar\omega\)): \(n \approx k_B T/\hbar\omega\) (classical limit)

5.2 Scattering Rate Temperature Scaling

The scattering rate for decay processes scales as:

At high temperatures where \(n \gg 1\):

For coalescence:

Therefore, in the high-temperature regime:

This leads to the classical thermal conductivity behavior \(\kappa \propto 1/T\) when U-processes dominate.

5.3 Low-Temperature Behavior

At low temperatures, only low-energy acoustic phonons are populated. The probability of U-processes is suppressed by:

This exponential suppression means that at very low \(T\), phonon-phonon scattering becomes negligible, and other mechanisms (boundary scattering, impurity scattering) dominate.

6. Fourth-Order Anharmonicity

While three-phonon processes dominate in most materials, fourth-order anharmonicity can be important in certain contexts.

6.1 Four-Phonon Scattering Processes

The fourth-order anharmonic Hamiltonian leads to processes like:

• \((\mathbf{q}_1, s_1) + (\mathbf{q}_2, s_2) \to (\mathbf{q}_3, s_3) + (\mathbf{q}_4, s_4)\) (two phonons scatter into two different phonons)
• \((\mathbf{q}_1, s_1) \to (\mathbf{q}_2, s_2) + (\mathbf{q}_3, s_3) + (\mathbf{q}_4, s_4)\) (one phonon decays into three)

These are typically weaker than three-phonon processes by a factor of \(\sim u/a\), but can become important:

• At very high temperatures
• When three-phonon processes are forbidden by selection rules
• In materials with high symmetry (where \(\Phi^{(3)}\) may vanish)
• For optical phonons in some materials

6.2 Recent Research: Four-Phonon Scattering in Thermal Conductivity

7. Connection to Thermal Resistance

Understanding phonon scattering is essential for predicting and engineering thermal conductivity.

7.1 Thermal Conductivity Formula

From the Boltzmann transport equation (covered in detail in Chapter 3), the lattice thermal conductivity can be expressed as:

Where:

• \(C_{\mathbf{q}s} = \hbar\omega_{\mathbf{q}s} \frac{\partial n}{\partial T}\) is the mode-specific heat capacity
• \(v_{\mathbf{q}s} = |\nabla_{\mathbf{q}} \omega_{\mathbf{q}s}|\) is the phonon group velocity
• \(\tau_{\mathbf{q}s}\) is the phonon lifetime from all scattering mechanisms

7.2 Matthiessen's Rule for Combined Scattering

When multiple scattering mechanisms are present (phonon-phonon, phonon-impurity, phonon-boundary), the total scattering rate is the sum:

This is Matthiessen's rule, valid when scattering events are uncorrelated.

7.3 Thermal Resistance Mechanisms Summary

8. Experimental Signatures of Anharmonicity

Anharmonic effects can be directly observed in various experimental techniques.

8.1 Spectroscopic Linewidths

In Raman, infrared, or neutron scattering spectroscopy, phonon peaks have finite width due to the finite lifetime:

This is a Lorentzian lineshape with FWHM = \(2\Gamma = 2\hbar/\tau\).

8.2 Phonon Frequency Shifts

Anharmonicity also causes phonon frequencies to shift with temperature:

The shift \(\Delta\omega(T)\) arises from:

• Volume expansion: Change in harmonic force constants due to lattice expansion
• Pure anharmonicity: Direct effect of anharmonic terms (phonon-phonon interactions)

Typically \(\Delta\omega < 0\) (phonons soften with increasing temperature).

8.3 Thermal Expansion

One of the most direct consequences of anharmonicity is thermal expansion. In the harmonic approximation, the lattice constant would be temperature-independent. The thermal expansion coefficient:

is directly proportional to the Grüneisen parameter \(\gamma\), which quantifies anharmonicity:

8.4 Inelastic Neutron Scattering

Neutron scattering can directly measure phonon dispersion curves. Anharmonic effects appear as:

• Asymmetric peak shapes (not perfect Lorentzians)
• Temperature-dependent peak positions
• Satellite peaks from higher-order processes
• Phonon decay pathways can sometimes be resolved

9. Computational Example: Phonon Scattering Rates

Let's implement a simplified calculation of three-phonon scattering rates using Python.

9.1 Simplified Model

We'll use a simplified isotropic model with Debye dispersion to demonstrate the concepts. For realistic calculations, one would use first-principles force constants.

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad
from scipy.constants import hbar, k as kB

class SimplifiedPhononScattering:
    """
    Simplified model for three-phonon scattering in an isotropic crystal.
    Uses Debye approximation and simplified anharmonic coupling.
    """

    def __init__(self, omega_D, v_sound, V0, gamma):
        """
        Parameters:
        -----------
        omega_D : float
            Debye frequency (rad/s)
        v_sound : float
            Sound velocity (m/s)
        V0 : float
            Cubic anharmonic coupling constant (J)
        gamma : float
            Grüneisen parameter (dimensionless)
        """
        self.omega_D = omega_D
        self.v_sound = v_sound
        self.V0 = V0
        self.gamma = gamma
        self.q_D = omega_D / v_sound  # Debye wavevector

    def omega(self, q):
        """Debye dispersion relation"""
        return self.v_sound * q

    def bose_einstein(self, omega, T):
        """Bose-Einstein occupation number"""
        x = hbar * omega / (kB * T)
        if x > 100:  # Avoid overflow
            return 0.0
        return 1.0 / (np.exp(x) - 1.0)

    def check_energy_conservation(self, omega1, omega2, omega3, tolerance=1e-6):
        """Check if energy is conserved (for decay: omega1 = omega2 + omega3)"""
        return np.abs(omega1 - omega2 - omega3) < tolerance

    def check_momentum_conservation(self, q1, q2, q3):
        """
        Check momentum conservation for decay: q1 = q2 + q3 (N-process)
        or q1 = q2 + q3 - G (U-process)
        Simplified 1D model: U-process if |q2 + q3| > q_D
        """
        q_sum = q2 + q3
        is_normal = np.abs(q1 - q_sum) < 1e-6
        is_umklapp = np.abs(q1 - (q_sum - 2*self.q_D)) < 1e-6
        return is_normal or is_umklapp, is_umklapp

    def scattering_rate_decay(self, q1, T):
        """
        Calculate decay scattering rate for phonon (q1) -> (q2) + (q3)
        Simplified analytical estimate using Debye approximation.
        """
        omega1 = self.omega(q1)

        # In full calculation, would integrate over all allowed (q2, q3) pairs
        # Here we use a simplified estimate proportional to phase space

        n_avg = self.bose_einstein(omega1/2, T)  # Average occupation

        # Simplified rate: Gamma ~ V0^2 * (n+1)^2 * DOS
        # With V0 ~ gamma * omega for cubic anharmonicity
        rate = (self.gamma * omega1 / hbar)**2 * (n_avg + 1)**2 * omega1

        # Add U-process suppression at low T
        # U-processes require high-energy phonons
        theta_D = hbar * self.omega_D / kB
        if T < theta_D:
            u_suppression = np.exp(-theta_D / (2 * T))
        else:
            u_suppression = 1.0

        return rate * u_suppression

    def scattering_rate_coalescence(self, q1, T):
        """
        Calculate coalescence scattering rate for (q2) + (q3) -> (q1)
        """
        omega1 = self.omega(q1)

        n_avg = self.bose_einstein(omega1/2, T)

        # Rate ~ V0^2 * n^2 * DOS
        rate = (self.gamma * omega1 / hbar)**2 * n_avg**2 * omega1

        return rate

    def total_scattering_rate(self, q, T):
        """Total three-phonon scattering rate"""
        return self.scattering_rate_decay(q, T) + self.scattering_rate_coalescence(q, T)

    def phonon_lifetime(self, q, T):
        """Phonon lifetime tau = 1/Gamma"""
        rate = self.total_scattering_rate(q, T)
        return 1.0 / rate if rate > 0 else np.inf

    def linewidth_fwhm(self, q, T):
        """Spectroscopic linewidth (FWHM) in energy units"""
        return 2 * hbar * self.total_scattering_rate(q, T)


# Example: Silicon-like parameters
omega_D_Si = 2 * np.pi * 15e12  # 15 THz Debye frequency
v_sound_Si = 6400  # m/s average sound velocity
V0_Si = 1e-20  # J (typical order of magnitude)
gamma_Si = 0.5  # Grüneisen parameter for Si

model = SimplifiedPhononScattering(omega_D_Si, v_sound_Si, V0_Si, gamma_Si)

# Temperature array
temperatures = np.linspace(50, 600, 50)

# Choose a phonon wavevector (middle of Brillouin zone)
q_test = 0.5 * model.q_D

# Calculate lifetime vs temperature
lifetimes = [model.phonon_lifetime(q_test, T) for T in temperatures]
linewidths = [model.linewidth_fwhm(q_test, T) / (1.602e-19 * 1e3) for T in temperatures]  # meV

# Plot results
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4))

ax1.plot(temperatures, np.array(lifetimes) * 1e12, 'b-', linewidth=2)
ax1.set_xlabel('Temperature (K)', fontsize=12)
ax1.set_ylabel('Phonon Lifetime (ps)', fontsize=12)
ax1.set_title('Phonon Lifetime vs Temperature', fontsize=13)
ax1.grid(True, alpha=0.3)
ax1.set_yscale('log')

ax2.plot(temperatures, linewidths, 'r-', linewidth=2)
ax2.set_xlabel('Temperature (K)', fontsize=12)
ax2.set_ylabel('Linewidth FWHM (meV)', fontsize=12)
ax2.set_title('Spectroscopic Linewidth vs Temperature', fontsize=13)
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('phonon_scattering_vs_temperature.png', dpi=150, bbox_inches='tight')
plt.show()

print(f"At 300 K:")
print(f"  Phonon lifetime: {model.phonon_lifetime(q_test, 300)*1e12:.2f} ps")
print(f"  Linewidth: {model.linewidth_fwhm(q_test, 300)/(1.602e-19*1e3):.2f} meV")

9.2 Expected Output

The code produces plots showing:

1. Phonon Lifetime vs Temperature: Lifetime decreases with increasing temperature, approximately as \(\tau \propto 1/T^2\) at high \(T\).
2. Linewidth vs Temperature: Linewidth increases with temperature, showing the anharmonic broadening effect observable in spectroscopy.

Typical values at 300 K for silicon:

• Phonon lifetime: ~1-10 ps
• Linewidth: ~0.1-1 meV

9.3 Extensions and Realistic Calculations

9.4 Exercise: Modify the Code

10. Summary

11. Exercises