Chapter 5: Phonon Spectroscopy

Phonon spectroscopy encompasses a powerful suite of experimental techniques for measuring vibrational excitations in materials. From inelastic neutron scattering that directly probes phonon dispersion relations to Raman and infrared spectroscopies that reveal zone-center modes, these methods provide essential insights into lattice dynamics. This chapter explores the physical principles, selection rules, and practical applications of major phonon spectroscopy techniques, comparing their strengths and complementary capabilities for studying phonons across different energy and momentum scales.

5.1 Inelastic Neutron Scattering (INS)

Fundamental Principles

Inelastic neutron scattering (INS) is the gold standard technique for measuring phonon dispersion relations throughout the Brillouin zone. Neutrons interact with atomic nuclei via the strong force, and their wavelength and energy match typical phonon wavelengths and energies.

Why neutrons are ideal for phonon studies:
• Thermal neutron wavelength: λ ≈ 1-4 Å (comparable to interatomic spacing)
• Thermal neutron energy: E ≈ 5-100 meV (comparable to phonon energies)
• Direct momentum and energy conservation measurement
• Penetrating (bulk probe, not surface-sensitive)
• Isotope-specific via nuclear spin

Scattering Theory and Cross-Section

In an inelastic scattering event, a neutron transfers momentum $\hbar\mathbf{Q}$ and energy $\hbar\omega$ to the crystal:

$$\mathbf{Q} = \mathbf{k}_i - \mathbf{k}_f$$

$$\hbar\omega = E_i - E_f = \frac{\hbar^2}{2m_n}(k_i^2 - k_f^2)$$

where $\mathbf{k}_i, \mathbf{k}_f$ are initial and final neutron wavevectors, and $m_n$ is the neutron mass.

For phonon creation (energy loss by neutron):

$$\mathbf{Q} = \mathbf{q} + \mathbf{G}$$

$$\hbar\omega = \hbar\omega_{\mathbf{q},s}$$

where $\mathbf{q}$ is the phonon wavevector, $\mathbf{G}$ is a reciprocal lattice vector, and $s$ is the phonon branch index.

Differential scattering cross-section: The one-phonon inelastic scattering cross-section (coherent) is given by:

$$\frac{d^2\sigma}{d\Omega dE_f} = \frac{k_f}{k_i} \frac{(2\pi)^3}{v_0} \sum_{\mathbf{G}} \sum_{s} |F(\mathbf{Q})|^2 \frac{(\hbar Q)^2}{2M\omega_{\mathbf{q},s}} (n_{\mathbf{q},s} + \frac{1}{2} \pm \frac{1}{2}) \delta(\omega - \omega_{\mathbf{q},s})$$

where:

• $F(\mathbf{Q})$ is the structure factor (weighted sum of scattering lengths)
• $M$ is the average atomic mass
• $n_{\mathbf{q},s}$ is the Bose-Einstein occupation number
• ± corresponds to phonon creation (+) or annihilation (−)
• $v_0$ is the unit cell volume

Temperature dependence: The Bose-Einstein factor gives:

$$n_{\mathbf{q},s} = \frac{1}{e^{\hbar\omega_{\mathbf{q},s}/k_BT} - 1}$$

At higher temperatures, more phonons are thermally populated, enhancing anti-Stokes (phonon annihilation) scattering.

Coherent vs Incoherent Scattering

Neutron scattering has two contributions:

Coherent scattering length $b_{coh}$: Average scattering amplitude

$$b_{coh} = \langle b \rangle$$

Incoherent scattering length $b_{inc}$: Root-mean-square deviation

$$b_{inc}^2 = \langle b^2 \rangle - \langle b \rangle^2$$

Example: For hydrogen, $b_{coh} = -3.74$ fm, $b_{inc} = 25.27$ fm. The huge incoherent cross-section makes H-containing materials excellent for measuring phonon density of states.

Instrumentation: Triple-Axis Spectrometer (TAS)

The triple-axis spectrometer allows independent selection of incident and scattered neutron energies and scattering angle, enabling point-by-point mapping of $S(\mathbf{Q}, \omega)$.

Operating principle:

1. Monochromator: Selects neutron wavelength $\lambda_i$ via Bragg reflection
2. Sample: Oriented to set scattering vector $\mathbf{Q}$
3. Analyzer: Energy-analyzes scattered neutrons to select $\lambda_f$
4. Detector: Counts neutrons satisfying energy-momentum conservation

Advantages: High energy resolution (~0.1 meV), flexible $(\mathbf{Q}, \omega)$ space access

Disadvantages: Low count rate (scans one point at a time), requires large single crystals

Time-of-Flight (TOF) Spectrometer

Time-of-flight spectrometers measure neutron energy by travel time over a known distance, enabling simultaneous measurement over a range of energies.

Energy-time relation:

$$E = \frac{1}{2}m_n v^2 = \frac{1}{2}m_n \left(\frac{L}{t}\right)^2 = \frac{m_n L^2}{2t^2}$$

where $L$ is flight path and $t$ is time-of-flight.

Direct geometry TOF: Monochromatic incident beam, broad energy range detected simultaneously

Indirect geometry TOF: Broad incident spectrum, monochromatic final energy (higher resolution for low-energy phonons)

Advantages: High count rate, simultaneous broad energy range, efficient for powder samples

Disadvantages: Lower energy resolution than TAS, complex data reduction

Python Example: Simulating Neutron Scattering Intensity

import numpy as np
import matplotlib.pyplot as plt

# Simulate 1D phonon dispersion for acoustic and optical modes
def phonon_dispersion(q, omega_max_ac=30, omega_optical=80):
    """
    Calculate phonon frequencies (meV) for acoustic and optical branches.
    q: wavevector in units of pi/a (0 to 1 for Gamma to zone boundary)
    """
    omega_acoustic = omega_max_ac * np.sin(np.pi * q / 2)
    omega_optical = omega_optical * np.ones_like(q)
    return omega_acoustic, omega_optical

# Bose-Einstein distribution
def bose_einstein(omega, T=300):
    """
    Thermal occupation number at temperature T (K) and energy omega (meV).
    """
    k_B = 0.0862  # meV/K
    return 1.0 / (np.exp(omega / (k_B * T)) - 1.0 + 1e-10)

# Calculate scattering intensity (simplified one-phonon)
def neutron_intensity(q, omega, omega_phonon, T=300):
    """
    Simplified one-phonon scattering intensity.
    """
    n_BE = bose_einstein(omega_phonon, T)
    # Stokes (phonon creation): n+1
    # Anti-Stokes (phonon annihilation): n
    gamma = 2.0  # meV, phonon lifetime broadening

    # Lorentzian lineshape for phonon
    stokes = (n_BE + 1) * (gamma / np.pi) / ((omega - omega_phonon)**2 + gamma**2)
    antistokes = n_BE * (gamma / np.pi) / ((omega + omega_phonon)**2 + gamma**2)

    return stokes + antistokes

# Create scattering map
q_points = np.linspace(0, 1, 100)
omega_range = np.linspace(-100, 100, 200)
Q_grid, Omega_grid = np.meshgrid(q_points, omega_range)

# Calculate dispersion
omega_ac, omega_op = phonon_dispersion(q_points)

# Calculate scattering intensity
intensity_map = np.zeros_like(Q_grid)
for i, q in enumerate(q_points):
    omega_ac_q = omega_ac[i]
    omega_op_q = omega_op[i]

    for j, omega in enumerate(omega_range):
        intensity_map[j, i] = neutron_intensity(q, omega, omega_ac_q, T=300)
        intensity_map[j, i] += neutron_intensity(q, omega, omega_op_q, T=300)

# Plotting
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))

# Scattering intensity map
im = ax1.contourf(Q_grid, Omega_grid, np.log10(intensity_map + 1e-3),
                   levels=20, cmap='hot')
ax1.plot(q_points, omega_ac, 'c--', linewidth=2, label='Acoustic (theory)')
ax1.plot(q_points, -omega_ac, 'c--', linewidth=2)
ax1.plot(q_points, omega_op, 'b--', linewidth=2, label='Optical (theory)')
ax1.plot(q_points, -omega_op, 'b--', linewidth=2)
ax1.set_xlabel('Reduced wavevector q (π/a)', fontsize=12)
ax1.set_ylabel('Energy transfer ω (meV)', fontsize=12)
ax1.set_title('INS Scattering Intensity Map S(Q,ω)', fontsize=13, weight='bold')
ax1.legend()
ax1.grid(alpha=0.3)
cbar = plt.colorbar(im, ax=ax1)
cbar.set_label('log₁₀(Intensity)', fontsize=11)

# Constant-Q scan at q=0.5
q_fixed = 0.5
idx_q = np.argmin(np.abs(q_points - q_fixed))
intensity_slice = intensity_map[:, idx_q]

ax2.plot(omega_range, intensity_slice, 'r-', linewidth=2)
ax2.axvline(omega_ac[idx_q], color='c', linestyle='--', label=f'Acoustic: {omega_ac[idx_q]:.1f} meV')
ax2.axvline(-omega_ac[idx_q], color='c', linestyle='--')
ax2.axvline(omega_op[idx_q], color='b', linestyle='--', label=f'Optical: {omega_op[idx_q]:.1f} meV')
ax2.axvline(-omega_op[idx_q], color='b', linestyle='--')
ax2.set_xlabel('Energy transfer ω (meV)', fontsize=12)
ax2.set_ylabel('Intensity (arb. units)', fontsize=12)
ax2.set_title(f'Constant-Q Scan at q = {q_fixed}π/a', fontsize=13, weight='bold')
ax2.legend()
ax2.grid(alpha=0.3)

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

print("✓ Neutron scattering intensity map created")
print(f"✓ Acoustic mode at q=0.5: {omega_ac[idx_q]:.2f} meV")
print(f"✓ Optical mode at q=0.5: {omega_op[idx_q]:.2f} meV")

Key observations:

• Both Stokes (+ω) and anti-Stokes (−ω) scattering visible
• Intensity asymmetry due to Bose-Einstein factor (more Stokes at room temperature)
• Acoustic branch disperses from Γ point, optical branch nearly flat

5.2 Raman Scattering Spectroscopy

Quantum Theory of Raman Scattering

Raman scattering involves inelastic scattering of photons by phonons. A photon of energy $\hbar\omega_i$ creates (Stokes) or annihilates (anti-Stokes) a phonon of energy $\hbar\omega_{\mathbf{q}}$:

$$\hbar\omega_s = \hbar\omega_i \mp \hbar\omega_{\mathbf{q}}$$

where $\omega_s$ is the scattered photon frequency.

Raman scattering mechanism: Light induces polarization oscillations in the crystal. If the polarizability $\alpha$ depends on atomic positions (phonon displacement), modulation of polarizability scatters light at shifted frequencies.

The Raman scattering intensity is proportional to the square of the Raman tensor:

$$I_{Raman} \propto \left|\mathbf{e}_s \cdot \frac{\partial \alpha}{\partial u} \cdot \mathbf{e}_i\right|^2 (n_{\mathbf{q}} + 1)$$

for Stokes scattering, where $\mathbf{e}_i, \mathbf{e}_s$ are incident and scattered polarization vectors, and $\partial\alpha/\partial u$ is the Raman tensor (polarizability derivative with respect to phonon displacement).

Selection Rules for Raman Activity

Raman selection rule: A phonon mode is Raman-active if it modulates the polarizability tensor:

$$\frac{\partial \alpha_{ij}}{\partial Q_{\mathbf{q},s}} \neq 0$$

where $Q_{\mathbf{q},s}$ is the phonon normal coordinate.

Symmetry considerations:

• Centrosymmetric crystals: Raman and infrared activities are mutually exclusive (mutual exclusion rule)
• Non-centrosymmetric crystals: Modes can be both Raman and IR active
• Wavevector restriction: $\mathbf{q} \approx 0$ (photon momentum negligible compared to Brillouin zone), so only zone-center modes are probed

Example: Diamond (Oh symmetry):

• Optical phonon at Γ: triply degenerate T_2g mode → Raman-active
• Acoustic phonons: not Raman-active (no polarizability change for uniform translation)

First-Order vs Second-Order Raman Processes

Second-order processes:

• Overtones: $\mathbf{q}_1 = \mathbf{q}_2$ (2 phonons from same branch)
• Combinations: $\mathbf{q}_1 + \mathbf{q}_2 = 0$ (phonons from different branches)
• Much weaker intensity (~1% of first-order)
• Allows probing of non-Γ phonons (broader q-range accessible)

Example: Silicon Raman spectrum

• First-order: Sharp peak at 520 cm⁻¹ (zone-center optical phonon)
• Second-order: Broad features at ~300 cm⁻¹ and ~950 cm⁻¹ (two-phonon density of states)

Experimental Configuration

Key components:

• Laser: Monochromatic excitation (visible/NIR range)
• Notch filter: Rejects elastic Rayleigh scattering (10⁶× stronger than Raman)
• Spectrometer: High-resolution grating (typical resolution ~1 cm⁻¹)
• Raman shift: Reported in wavenumbers $\tilde{\nu} = 1/\lambda$ (cm⁻¹)

Advantages:

• Non-destructive, minimal sample preparation
• Spatial resolution ~1 μm (confocal micro-Raman)
• Polarization analysis reveals phonon symmetry
• Applicable to crystals, powders, thin films, liquids

Disadvantages:

• Only q ≈ 0 modes accessible (zone-center only)
• Weak signal (requires high laser power, can cause heating)
• Fluorescence can overwhelm Raman signal

Python Example: Simulating Raman Spectrum

import numpy as np
import matplotlib.pyplot as plt

# Define phonon modes for a hypothetical crystal
modes = [
    {'freq': 150, 'intensity': 30, 'width': 5, 'label': 'E₂g'},
    {'freq': 280, 'intensity': 15, 'width': 8, 'label': 'A₁g'},
    {'freq': 520, 'intensity': 100, 'width': 4, 'label': 'F₂g'},
    {'freq': 640, 'intensity': 25, 'width': 10, 'label': 'E₂g'},
]

# Second-order broad features
second_order = [
    {'freq': 300, 'intensity': 5, 'width': 50, 'label': '2TA'},
    {'freq': 950, 'intensity': 8, 'width': 80, 'label': '2TO'},
]

def lorentzian(x, x0, intensity, gamma):
    """Lorentzian lineshape for phonon peaks."""
    return intensity * (gamma**2) / ((x - x0)**2 + gamma**2)

# Raman shift range
raman_shift = np.linspace(50, 1100, 2000)  # cm^-1

# Calculate spectrum
spectrum = np.zeros_like(raman_shift)

# First-order peaks
for mode in modes:
    spectrum += lorentzian(raman_shift, mode['freq'], mode['intensity'], mode['width'])

# Second-order features
for mode in second_order:
    spectrum += lorentzian(raman_shift, mode['freq'], mode['intensity'], mode['width'])

# Add noise
noise = np.random.normal(0, 0.5, len(spectrum))
spectrum_noisy = spectrum + noise

# Temperature effect: calculate anti-Stokes intensity
T = 300  # K
k_B = 0.695  # cm^-1/K
anti_stokes_shift = -raman_shift  # Negative shift
anti_stokes_spectrum = np.zeros_like(raman_shift)

for mode in modes:
    n_BE = 1.0 / (np.exp(mode['freq'] / (k_B * T)) - 1)
    anti_stokes_spectrum += lorentzian(raman_shift, -mode['freq'],
                                        mode['intensity'] * n_BE, mode['width'])

# Plotting
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 10))

# Stokes spectrum with peak labels
ax1.plot(raman_shift, spectrum_noisy, 'b-', linewidth=1.5, alpha=0.7, label='Measured')
ax1.plot(raman_shift, spectrum, 'r-', linewidth=2, label='Fitted')

for mode in modes:
    ax1.annotate(mode['label'], xy=(mode['freq'], mode['intensity']),
                 xytext=(mode['freq'], mode['intensity'] + 10),
                 ha='center', fontsize=10, weight='bold',
                 arrowprops=dict(arrowstyle='->', color='red', lw=1.5))

ax1.fill_between(raman_shift, 0, spectrum, alpha=0.2, color='blue')
ax1.set_xlabel('Raman Shift (cm⁻¹)', fontsize=12)
ax1.set_ylabel('Intensity (arb. units)', fontsize=12)
ax1.set_title('First-Order Raman Spectrum (Stokes)', fontsize=14, weight='bold')
ax1.legend(fontsize=11)
ax1.grid(alpha=0.3)
ax1.set_xlim(50, 1100)

# Stokes vs anti-Stokes comparison
ax2.plot(raman_shift, spectrum, 'r-', linewidth=2, label='Stokes (phonon creation)')
ax2.plot(-raman_shift, anti_stokes_spectrum, 'b-', linewidth=2, label='Anti-Stokes (phonon annihilation)')
ax2.axvline(0, color='gray', linestyle='--', linewidth=1, label='Rayleigh (elastic)')
ax2.set_xlabel('Raman Shift (cm⁻¹)', fontsize=12)
ax2.set_ylabel('Intensity (arb. units)', fontsize=12)
ax2.set_title(f'Stokes vs Anti-Stokes at T = {T} K', fontsize=14, weight='bold')
ax2.legend(fontsize=11)
ax2.grid(alpha=0.3)
ax2.set_xlim(-700, 700)

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

# Calculate intensity ratio
mode_520 = modes[2]
n_BE_520 = 1.0 / (np.exp(mode_520['freq'] / (k_B * T)) - 1)
ratio = n_BE_520 / (n_BE_520 + 1)

print("✓ Raman spectrum simulation complete")
print(f"✓ Anti-Stokes/Stokes intensity ratio at {mode_520['freq']} cm⁻¹: {ratio:.3f}")
print(f"✓ Bose-Einstein occupation at {T} K: n = {n_BE_520:.3f}")

Key observations:

• Stokes scattering (positive shift) is much stronger than anti-Stokes at room temperature
• Anti-Stokes/Stokes ratio gives temperature via Bose-Einstein statistics
• Second-order features are broader and weaker than first-order peaks

5.3 Infrared Absorption Spectroscopy

Dipole Selection Rules

Infrared (IR) absorption occurs when a phonon mode modulates the electric dipole moment of the crystal. The electromagnetic wave couples to the oscillating dipole.

IR selection rule: A phonon mode is IR-active if it creates a time-varying electric dipole moment:

$$\frac{\partial \mathbf{P}}{\partial Q_{\mathbf{q},s}} \neq 0$$

where $\mathbf{P}$ is the electric polarization.

Key points:

• Only polar modes (modes with ionic displacement creating net dipole) are IR-active
• Acoustic phonons: Not IR-active (uniform translation doesn't change dipole)
• Wavevector restriction: $\mathbf{q} \approx 0$ (IR photon wavelength ~10 μm >> lattice constant)

Example: NaCl (rocksalt structure, cubic):

• Optical phonon at Γ: Na⁺ and Cl⁻ oscillate out-of-phase → dipole moment oscillates → IR-active
• This mode is NOT Raman-active (centrosymmetric crystal, mutual exclusion)

Polaritons in Polar Materials

In polar crystals, IR-active phonons couple to electromagnetic waves, forming phonon-polaritons (mixed light-phonon quasiparticles).

Polariton dispersion relation: For a simple polar crystal with one IR-active mode at frequency $\omega_{TO}$:

$$\omega^2(\mathbf{q}) = \frac{c^2 q^2}{2\epsilon_\infty} \left[1 + \epsilon_\infty \pm \sqrt{1 + 2\epsilon_\infty + \epsilon_\infty^2 + 4\epsilon_\infty \frac{\omega_{TO}^2}{c^2 q^2}(1 + \frac{\omega_{LO}^2 - \omega_{TO}^2}{\omega_{TO}^2})}\right]$$

This simplifies in two limits:

• Low q (photon-like): $\omega \approx cq/\sqrt{\epsilon_\infty}$ (dispersion of light in medium)
• High q (phonon-like): $\omega \approx \omega_{TO}$ or $\omega_{LO}$ (longitudinal optical phonon)

Lyddane-Sachs-Teller (LST) relation:

$$\frac{\omega_{LO}^2}{\omega_{TO}^2} = \frac{\epsilon_0}{\epsilon_\infty}$$

where $\epsilon_0$ is static dielectric constant and $\epsilon_\infty$ is high-frequency dielectric constant.

Physical interpretation:

• Polariton upper branch: photon-like at low q, asymptotes to $\omega_{LO}$ at high q
• Polariton lower branch: starts from $\omega_{TO}$, becomes photon-like at high q
• Stop band: Between $\omega_{TO}$ and $\omega_{LO}$, no propagating modes (total reflection)

Reflectivity and Kramers-Kronig Analysis

IR spectroscopy often measures reflectivity rather than transmission (especially for opaque materials). The complex dielectric function $\epsilon(\omega) = \epsilon_1(\omega) + i\epsilon_2(\omega)$ determines reflectivity:

$$R(\omega) = \left|\frac{\sqrt{\epsilon(\omega)} - 1}{\sqrt{\epsilon(\omega)} + 1}\right|^2$$

Kramers-Kronig relations allow extraction of phonon frequencies from reflectivity data (connects real and imaginary parts of $\epsilon(\omega)$ via causality).

Python Example: Polariton Dispersion

import numpy as np
import matplotlib.pyplot as plt

# Material parameters (example: GaP)
omega_TO = 367  # cm^-1, transverse optical phonon
omega_LO = 403  # cm^-1, longitudinal optical phonon
epsilon_inf = 9.1  # high-frequency dielectric constant
epsilon_0 = 11.1  # static dielectric constant

# Verify LST relation
LST_ratio = (omega_LO / omega_TO)**2
LST_expected = epsilon_0 / epsilon_inf
print(f"LST relation check: (ω_LO/ω_TO)² = {LST_ratio:.3f}, ε₀/ε_∞ = {LST_expected:.3f}")

# Wavevector range
c = 3e10  # cm/s, speed of light
q_max = 0.05  # cm^-1 (small compared to Brillouin zone ~10^8 cm^-1)
q = np.linspace(1e-6, q_max, 500)

# Polariton dispersion (solving coupled equations)
def polariton_dispersion(q, omega_TO, omega_LO, epsilon_inf, c):
    """
    Calculate upper and lower polariton branch frequencies.
    """
    A = c**2 * q**2 / epsilon_inf
    B = omega_TO**2 + omega_LO**2
    C = omega_TO**2 * omega_LO**2

    # Discriminant
    discriminant = B**2 - 4*C + 4*A*(B - omega_TO**2)

    omega_upper = np.sqrt((B + A + np.sqrt(discriminant)) / 2)
    omega_lower = np.sqrt((B + A - np.sqrt(discriminant)) / 2)

    return omega_lower, omega_upper

omega_lower, omega_upper = polariton_dispersion(q, omega_TO, omega_LO, epsilon_inf, c)

# Light line in medium
omega_light = c * q / np.sqrt(epsilon_inf)

# Plotting
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 6))

# Polariton dispersion
ax1.plot(q * 1e3, omega_upper, 'r-', linewidth=2.5, label='Upper polariton')
ax1.plot(q * 1e3, omega_lower, 'b-', linewidth=2.5, label='Lower polariton')
ax1.plot(q * 1e3, omega_light, 'k--', linewidth=1.5, label='Light line (c/√ε_∞)')
ax1.axhline(omega_TO, color='green', linestyle=':', linewidth=2, label=f'ω_TO = {omega_TO} cm⁻¹')
ax1.axhline(omega_LO, color='purple', linestyle=':', linewidth=2, label=f'ω_LO = {omega_LO} cm⁻¹')
ax1.fill_between(q * 1e3, omega_TO, omega_LO, alpha=0.2, color='gray', label='Stop band')

ax1.set_xlabel('Wavevector q (×10⁻³ cm⁻¹)', fontsize=12)
ax1.set_ylabel('Frequency ω (cm⁻¹)', fontsize=12)
ax1.set_title('Phonon-Polariton Dispersion in GaP', fontsize=14, weight='bold')
ax1.legend(fontsize=10, loc='upper left')
ax1.grid(alpha=0.3)
ax1.set_ylim(300, 450)

# Dielectric function and reflectivity
omega_range = np.linspace(300, 450, 500)
gamma = 2.0  # damping (cm^-1)

def dielectric_function(omega, omega_TO, omega_LO, epsilon_inf, gamma):
    """Complex dielectric function for polar phonon."""
    epsilon = epsilon_inf * (omega**2 - omega_LO**2 + 1j*gamma*omega) / \
                           (omega**2 - omega_TO**2 + 1j*gamma*omega)
    return epsilon

epsilon = dielectric_function(omega_range, omega_TO, omega_LO, epsilon_inf, gamma)
epsilon_real = np.real(epsilon)
epsilon_imag = np.imag(epsilon)

# Reflectivity
sqrt_epsilon = np.sqrt(epsilon)
reflectivity = np.abs((sqrt_epsilon - 1) / (sqrt_epsilon + 1))**2

ax2_twin = ax2.twinx()
ax2.plot(omega_range, epsilon_real, 'b-', linewidth=2, label='ε₁ (real)')
ax2.plot(omega_range, epsilon_imag, 'r-', linewidth=2, label='ε₂ (imaginary)')
ax2_twin.plot(omega_range, reflectivity, 'g-', linewidth=2.5, label='Reflectivity R')

ax2.axvline(omega_TO, color='green', linestyle=':', linewidth=1.5)
ax2.axvline(omega_LO, color='purple', linestyle=':', linewidth=1.5)
ax2.axhline(0, color='black', linestyle='-', linewidth=0.5)

ax2.set_xlabel('Frequency ω (cm⁻¹)', fontsize=12)
ax2.set_ylabel('Dielectric Function ε', fontsize=12, color='b')
ax2_twin.set_ylabel('Reflectivity R', fontsize=12, color='g')
ax2.set_title('Dielectric Function and IR Reflectivity', fontsize=14, weight='bold')
ax2.legend(loc='upper left', fontsize=10)
ax2_twin.legend(loc='upper right', fontsize=10)
ax2.grid(alpha=0.3)
ax2.tick_params(axis='y', labelcolor='b')
ax2_twin.tick_params(axis='y', labelcolor='g')

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

print(f"\n✓ Polariton dispersion calculated")
print(f"✓ Stop band (no propagation): {omega_TO} - {omega_LO} cm⁻¹")
print(f"✓ Maximum reflectivity in stop band: {reflectivity[np.argmax(reflectivity)]:.3f}")

Key observations:

• Upper and lower polariton branches avoid crossing (anti-crossing behavior)
• Stop band between $\omega_{TO}$ and $\omega_{LO}$ shows high reflectivity
• Dielectric function becomes negative in stop band (evanescent waves)

5.4 Other Phonon Spectroscopy Techniques

Electron Energy Loss Spectroscopy (EELS)

Principle: High-energy electrons (50-300 keV in TEM) lose energy by exciting phonons. Energy loss spectrum reveals phonon energies.

Energy-momentum relation:

$$\Delta E = \frac{\hbar^2}{2m_e}(k_i^2 - k_f^2) \approx \hbar\omega_{phonon}$$

$$\mathbf{q} = \mathbf{k}_i - \mathbf{k}_f$$

Advantages:

• Spatial resolution down to atomic scale (STEM-EELS)
• Can probe high-energy phonons (up to ~500 meV)
• Momentum-resolved measurements possible

Disadvantages:

• Requires ultrathin samples (~10-100 nm for TEM)
• Energy resolution limited by electron gun (~10-100 meV, worse than optical methods)
• Electron beam can damage sensitive materials

Typical applications:

• Mapping phonon modes in nanostructures
• Interface and surface phonons
• Local vibrational spectroscopy in defects

Inelastic X-ray Scattering (IXS)

Principle: Similar to neutron scattering, but uses high-energy X-rays (~10-20 keV). X-rays interact with electron density, which is modulated by phonon displacements.

Scattering cross-section:

$$\frac{d^2\sigma}{d\Omega dE_f} \propto |F(\mathbf{Q})|^2 S(\mathbf{Q}, \omega)$$

where $F(\mathbf{Q})$ is the X-ray structure factor and $S(\mathbf{Q}, \omega)$ is the dynamic structure factor.

Advantages:

• High momentum resolution (can access entire Brillouin zone)
• No neutron absorption issues (works for elements like B, Li that absorb neutrons)
• Small samples sufficient (~0.1 mm³)
• Combines with high-pressure diamond anvil cells

Disadvantages:

• Weaker scattering cross-section than neutrons (requires synchrotron)
• Energy resolution worse than neutron TAS (~1.5 meV typical)
• Expensive and limited beamtime availability

Typical applications:

• High-pressure phonon studies
• Light-element compounds (Li-ion battery materials, hydrogen storage)
• Acoustic phonons in liquids and amorphous materials

Time-Resolved Phonon Spectroscopy (Pump-Probe)

Principle: Ultrafast laser pulses (femtosecond) excite (pump) and probe phonon dynamics. Measures non-equilibrium phonon populations and lifetimes.

Experimental setup:

1. Pump pulse: Excites electrons, which relax and emit phonons (phonon generation)
2. Probe pulse: Delayed by time $\Delta t$, measures transient changes via reflectivity/transmission
3. Time-domain oscillations: Coherent phonon oscillations appear as periodic modulation

Key observables:

• Coherent phonons: Oscillations at phonon frequency (Fourier transform gives mode frequencies)
• Phonon lifetime: Exponential decay of oscillations gives anharmonic damping rate
• Phonon-phonon scattering: Temperature-dependent lifetimes reveal scattering channels

Advantages:

• Measures phonon dynamics in time domain (ultrafast, ~10 fs resolution)
• Access to phonon lifetimes and decay channels
• Can study phonon-electron coupling, phase transitions

Disadvantages:

• Limited to coherent, low-energy phonons (~THz range)
• Requires high-quality surfaces and intense laser pulses
• Complex data interpretation (multiple overlapping processes)

Typical applications:

• Phonon-mediated superconductivity dynamics
• Ultrafast phase transitions (photoinduced)
• Coherent control of phonons in quantum materials

5.5 Comparison of Phonon Spectroscopy Techniques

Complementary Use of Techniques

Strategy for complete phonon characterization:

1. INS (TAS or TOF): Measure complete phonon dispersion $\omega_{\mathbf{q},s}$ across Brillouin zone
2. Raman: Identify zone-center optical modes, determine symmetry via polarization analysis
3. IR: Measure IR-active polar modes, determine dielectric constants via LST relation
4. IXS: Complement neutron measurements for light-element or high-pressure systems
5. Pump-probe: Measure phonon lifetimes and anharmonic decay channels
6. EELS: Map local phonon modes in heterostructures and defects

Python Example: Technique Comparison Visualization

import numpy as np
import matplotlib.pyplot as plt

# Define accessible energy and q-ranges for different techniques
techniques = {
    'INS (TAS)': {'E_min': 0.1, 'E_max': 100, 'q_min': 0, 'q_max': 1.0, 'color': 'blue'},
    'INS (TOF)': {'E_min': 1, 'E_max': 500, 'q_min': 0.1, 'q_max': 1.0, 'color': 'cyan'},
    'Raman': {'E_min': 0.1, 'E_max': 300, 'q_min': 0, 'q_max': 0.01, 'color': 'red'},
    'Infrared': {'E_min': 1, 'E_max': 100, 'q_min': 0, 'q_max': 0.01, 'color': 'orange'},
    'IXS': {'E_min': 1, 'E_max': 200, 'q_min': 0, 'q_max': 1.0, 'color': 'green'},
    'EELS': {'E_min': 10, 'E_max': 500, 'q_min': 0.05, 'q_max': 0.5, 'color': 'purple'},
    'Pump-Probe': {'E_min': 0.01, 'E_max': 10, 'q_min': 0, 'q_max': 0.005, 'color': 'magenta'},
}

# Plot accessible regions
fig, ax = plt.subplots(figsize=(12, 8))

for i, (name, params) in enumerate(techniques.items()):
    # Rectangle showing accessible (E, q) range
    width = params['q_max'] - params['q_min']
    height = params['E_max'] - params['E_min']
    rect = plt.Rectangle((params['q_min'], params['E_min']), width, height,
                          alpha=0.3, facecolor=params['color'],
                          edgecolor=params['color'], linewidth=2,
                          label=name)
    ax.add_patch(rect)

# Typical phonon dispersion overlay
q_phonon = np.linspace(0, 1.0, 100)
omega_acoustic = 30 * np.sin(np.pi * q_phonon / 2)
omega_optical = 60 * np.ones_like(q_phonon)

ax.plot(q_phonon, omega_acoustic, 'k--', linewidth=2, label='Acoustic phonon')
ax.plot(q_phonon, omega_optical, 'k:', linewidth=2, label='Optical phonon')

ax.set_xlabel('Reduced wavevector q (π/a)', fontsize=13)
ax.set_ylabel('Phonon energy (meV)', fontsize=13)
ax.set_title('Accessible (Energy, Momentum) Ranges for Phonon Spectroscopy Techniques',
             fontsize=14, weight='bold')
ax.set_xlim(-0.05, 1.05)
ax.set_ylim(0.01, 600)
ax.set_yscale('log')
ax.legend(loc='upper left', fontsize=10, ncol=2)
ax.grid(alpha=0.3, which='both')

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

print("✓ Technique comparison plot created")
print("\nKey insights:")
print("  • INS provides full Brillouin zone coverage (best for dispersion)")
print("  • Raman/IR limited to q≈0 but high energy resolution")
print("  • IXS complements INS for light elements")
print("  • EELS offers spatial resolution with moderate energy resolution")
print("  • Pump-probe accesses ultrafast dynamics at low energies")

Summary

In this chapter, we explored the major experimental techniques for phonon spectroscopy:

• Inelastic Neutron Scattering (INS): The gold standard for measuring complete phonon dispersion relations. We learned the scattering cross-section formalism, difference between coherent and incoherent scattering, and operational principles of triple-axis and time-of-flight spectrometers.

• Raman Spectroscopy: Probes zone-center phonons via inelastic light scattering. We derived quantum selection rules, distinguished first-order and second-order processes, and understood how polarization analysis reveals phonon symmetry.

• Infrared Absorption: Measures IR-active polar phonon modes. We examined dipole selection rules, phonon-polariton coupling in polar crystals, and the Lyddane-Sachs-Teller relation connecting dielectric properties to phonon frequencies.

• Complementary Techniques: Electron energy loss spectroscopy (EELS) provides atomic-scale spatial resolution, inelastic X-ray scattering (IXS) enables high-pressure and light-element studies, and time-resolved pump-probe methods access phonon dynamics and lifetimes.

• Comparative Analysis: Each technique has unique strengths in energy resolution, momentum range, and sample requirements. Combining multiple techniques provides comprehensive phonon characterization from zone-center to full Brillouin zone, from static properties to ultrafast dynamics.

Understanding these experimental methods is essential for interpreting phonon measurements and designing experiments to probe specific aspects of lattice dynamics in materials.

Exercises