Optical Properties from DFT

Learning Objectives

By the end of this section you should be able to:

1. Translate any feature of an optical absorption spectrum into the band-structure language built in §5.1: explain colour, transparency, and absorption edges as statements about $E_n(\mathbf{k})$.
2. Write down — and read — the complex dielectric function $\varepsilon(\omega) = \varepsilon_1(\omega) + i\,\varepsilon_2(\omega)$ and connect it to refractive index $n(\omega)$, extinction coefficient $k(\omega)$, absorption coefficient $\alpha(\omega)$, and reflectivity $R(\omega)$.
3. Derive $\varepsilon_2(\omega)$ from Fermi's golden rule, and identify the three ingredients that sculpt every peak: the joint density of states, the momentum matrix element, and the broadening.
4. Use the Kramers–Kronig relation to argue why $\varepsilon_1(\omega)$ is uniquely determined by $\varepsilon_2(\omega)$ — and why a wrong gap in DFT poisons both real and imaginary parts at every frequency.
5. Run an optical calculation in VASP: set LOPTICS=.TRUE., choose NBANDS, NEDOS, and CSHIFT sensibly, and extract $\varepsilon_1, \varepsilon_2, n, k, \alpha$ from vasprun.xml.

One-line preview: the optical spectrum of a crystal is a one-dimensional histogram of vertical (k-conserving) jumps from occupied to unoccupied bands, weighted by how dipole-friendly each jump is. Everything else — refractive index, colour, reflectivity — is downstream.

What Optics Asks of a Crystal

Shine a laser of frequency $\omega$ on a crystal and three things can happen at the interface: a fraction of the light is reflected, a fraction is transmitted, and a fraction is absorbed inside the slab. The relative weights of these three outcomes are what every optical experiment ultimately measures, and all three are encoded in a single complex-valued function of one scalar variable:

the frequency-dependent dielectric function. The real part $\varepsilon_1$ describes how polarisable the electron cloud is at frequency $\omega$ — it controls phase velocity, refractive index, and the lensing of light. The imaginary part $\varepsilon_2$ describes how dissipative the crystal is at that frequency — it controls absorption, the Beer–Lambert decay of intensity, and ultimately the colour you see.

The Photon as a Vertical Arrow

From §5.1 you already know the geometry: a photon at visible frequency carries crystal momentum $k_\text{photon} \sim 10^{-3}\,\text{Å}^{-1}$, which on the Brillouin-zone scale of $\sim 1\,\text{Å}^{-1}$ is essentially zero. Conservation of crystal momentum therefore forces every photon-driven transition to be vertical in the band structure: $\mathbf{k}_\text{final} = \mathbf{k}_\text{initial}$. On a band diagram this is a strictly upward arrow whose tail sits on an occupied band and whose head sits on an empty band.

The arrow has two readings. Length is the photon energy $\hbar\omega$. Position is a wavevector $\mathbf{k}$ in the Brillouin zone. A photon of energy $\hbar\omega$finds a home in the crystal at any $\mathbf{k}$ where the gap between an occupied band $v$ and an empty band $c$ matches:

The set of $\mathbf{k}$ where this equation holds is a 2-D iso-energy surface inside the BZ. Sweep $\hbar\omega$ upward and the surface sweeps through the zone, gathering area. The amount of area gathered, weighted by the dipole strength of each jump, is the absorption spectrum. That single sentence is what the rest of this section unpacks.

Interactive — Vertical Transitions in a 1-D Band Structure

Below is a 1-D toy crystal with two parabolic bands. Drag the photon energy slider $\hbar\omega$: amber arrows appear at every $k$ where $E_c(k) - E_v(k) = \hbar\omega$. Drag the gap slider to widen or narrow the bands. Toggle Indirect to push the conduction-band minimum away from $\Gamma$ — the smallest verticalgap then sits above the true (indirect) gap, and there is a forbidden window where the photon is transparent even though there are electronic states available.

1. With indirect off, sweep $\hbar\omega$ from below the gap upward. The first arrow appears exactly at $\hbar\omega = E_g$ at $k = 0$. Above that, two arrows appear, symmetric about $k = 0$. That is the direct-gap absorption edge.
2. Switch indirect on. The arrows shift, and the smallest one no longer sits at $k = 0$ — the true minimum gap (red dashed line) is now diagonal. First-order (vertical) absorption begins at a higher energy than the true gap. The gap below this onset is bridged only with phonon assistance and is therefore weak.
3. Push $E_g$ to 0.5 eV (infrared), to 1.5 eV (visible red — like Mn:CdSe), and to 3.0 eV (UV). The colour of the crystal is not a separate property; it is the position of the first arrow.

The Dielectric Function: What DFT Actually Computes

Maxwell's equations in matter relate the electric displacement $\mathbf{D}$ and the electric field $\mathbf{E}$ through $\mathbf{D}(\omega) = \varepsilon_0\,\varepsilon(\omega)\,\mathbf{E}(\omega)$. For an anisotropic crystal $\varepsilon$ is a 3×3 tensor; cubic crystals like zincblende CdSe have only one independent component (the diagonal), uniaxial wurtzite CdSe has two (in-plane and out-of-plane). VASP writes the full tensor for you in any case.

The complex refractive index $\tilde n(\omega) = n(\omega) + i\,k(\omega)$ is the square root:

Solving this 2×2 algebra for $n$ and $k$ in terms of $\varepsilon_1, \varepsilon_2$ gives the formulas you will see in the VASP code below:

Two derived quantities then drop out by algebra. The absorption coefficient $\alpha(\omega) = 2\omega k(\omega)/c$ — Beer's law. The normal-incidence reflectivity $R(\omega) = \big[(n-1)^2 + k^2\big]/\big[(n+1)^2 + k^2\big]$ — the fraction of intensity bounced off a flat surface.

ε₂(ω) from Fermi's Golden Rule

The starting point is first-order perturbation theory for an electron in a Bloch state coupled to the vector potential of a photon. Fermi's golden rule gives the rate at which an electron in band $v$ at wavevector $\mathbf{k}$ jumps to band $c$ at the same $\mathbf{k}$ when bathed in a monochromatic field. Sum that rate over all initial states, divide by the photon flux, and you arrive at the textbook result:

Three pieces of physics live in that formula. Each one is a knob you can turn in a real DFT calculation.

The Joint Density of States — Where Peaks Come From

Set the matrix element to a constant and the prefactors aside. What survives in the formula is the joint density of states (JDOS):

Geometrically, $J_{cv}(\omega)$ measures the 2-D area in the Brillouin zone on which the vertical-transition energy equals $\hbar\omega$, divided by the BZ volume. It is to optics what the regular DOS $g(E)$ of §5.2 is to thermodynamics — the geometric factor that decides where the spectrum has structure.

Two universal features follow:

• Square-root edge. Near a direct, parabolic minimum gap, expand $E_c - E_v = E_g + \hbar^2 k^2 / 2\mu$ with reduced mass $\mu = m_e^* m_h^* / (m_e^* + m_h^*)$. The JDOS in 3-D is $J(\omega) \propto \sqrt{\hbar\omega - E_g}$ for $\hbar\omega \geq E_g$, and zero below. Every direct-gap absorption edge has this shape.
• Van Hove singularities. At wavevectors where $\nabla_{\mathbf{k}} (E_c - E_v) = 0$, the JDOS develops integrable kinks: jumps in 1-D, logarithmic divergences in 2-D, square-root cusps in 3-D. The named peaks of a real spectrum ($E_0, E_1, E_2$ in semiconductors) are van Hove singularities of the conduction– valence energy difference.

Interactive — Building ε(ω) from a JDOS

The widget below builds $\varepsilon_2(\omega)$ from a parabolic-JDOS model with adjustable gap $E_g$ and Lorentzian broadening $\Gamma$. It then computes $\varepsilon_1(\omega)$ numerically from $\varepsilon_2$ via the Kramers–Kronig integral — the two curves are not independent. Switch the matrix element between "direct-allowed" (the canonical $\sqrt{\hbar\omega - E_g}$ edge) and "forbidden" (a softer $(\hbar\omega - E_g)^{3/2}$ onset) to see how the spectrum reshapes.

Pay attention to four things as you slide:

1. The threshold. $\varepsilon_2(\omega)$ is exactly zero for $\hbar\omega < E_g$. The crystal is transparent below the gap. This is why pure silicon ($E_g = 1.1$ eV) is opaque to visible light but transparent in the mid-infrared.
2. The shape near the threshold. Allowed transitions give a vertical-tangent square-root rise; forbidden ones give a horizontal-tangent $(E - E_g)^{3/2}$ rise. That is a fingerprint you can read off an experimental absorption spectrum.
3. The Kramers–Kronig coupling. $\varepsilon_1$ develops a peak just below the absorption onset and dips to a minimum just above — this is the universal dispersion shape every refractive-index curve shows near an absorption band.
4. Broadening. Increasing $\Gamma$ smooths peaks but conserves the area under $\varepsilon_2(\omega)$ (the f-sum rule). In a real DFT calculation, $\Gamma$ is set by the CSHIFT tag.

Kramers–Kronig: Why ε₁ Is Not Independent of ε₂

Causality — the crystal cannot respond before the field arrives — forces $\varepsilon(\omega)$ to be analytic in the upper half complex-frequency plane. A standard contour argument then gives the Kramers–Kronig relation:

Two consequences worth burning in:

• You only need ε₂. If you have computed $\varepsilon_2$ well, you do not need a separate calculation for $\varepsilon_1$. VASP applies the KK transform internally to give you both.
• Errors propagate everywhere. A wrong absorption edge — say, the notorious DFT band-gap underestimation that we will fix with hybrids in §6.4 — is not a localised error. It poisons $\varepsilon_1$ at everyfrequency. Static dielectric constants, refractive indices, Hamaker coefficients of dispersion forces — all of them inherit the gap error through KK.

Selection Rules and the Momentum Matrix Element

The matrix element $\langle \psi_{c\mathbf{k}} | \mathbf{p} | \psi_{v\mathbf{k}} \rangle$ is the genetic sequence of the spectrum. Two of its features matter most:

• Symmetry zeros. If the valence and conduction states at $\mathbf{k}$ belong to irreducible representations whose direct product does not contain the representation of $\mathbf{p}$, the matrix element vanishes by symmetry. This is the same selection rule you met in the H-atom — no $s \to s$ dipole transition — applied to crystal point groups and §2.7 direct products.
• Orbital character. The matrix element peaks when the valence orbital and conduction orbital share spatial overlap and have opposite parity (s ↔ p, p ↔ d). In CdSe the VBM is mostly Se-4p and the CBM is mostly Cd-5s; this s-p product islarge, which is one reason CdSe absorbs strongly.

From ε(ω) to What an Experimentalist Sees

A spectroscopist does not measure $\varepsilon$ directly; they measure intensity ratios. Three derived quantities bridge the theorist's$\varepsilon$ to those measurements:

Notice the round trip. DFT computes $\varepsilon$; algebra produces $n, k, \alpha, R$; the experimentalist measures one of those. To compare with experiment you simply look up the right algebraic combination of $\varepsilon_1$ and $\varepsilon_2$. That is why optical-output figures in the VASP literature usually plot the same underlying spectrum in three or four guises on the same page.

VASP — LOPTICS, NBANDS, and the Optical Tensor

Optical calculations in VASP run on top of a converged self- consistent ground state, exactly like the band-structure recipe of §5.1. The two extras you need are an LOPTICS flag and a much larger NBANDS than usual, because the sum over conduction states must be carried high enough into the empty manifold for the spectrum to be saturated up to your largest photon energy of interest.

Step 1 — converged SCF (as before)

1# INCAR  (step 1: self-consistent ground state)
2SYSTEM = CdSe zincblende — SCF
3PREC   = Accurate
4ENCUT  = 400
5EDIFF  = 1E-7        # tighter than usual — KK is sensitive to noise
6ISMEAR = 0
7SIGMA  = 0.05
8LCHARG = .TRUE.      # write CHGCAR for step 2
9LWAVE  = .TRUE.      # write WAVECAR for step 2

Step 2 — non-SCF + LOPTICS

1# INCAR  (step 2: optical / dielectric tensor)
2SYSTEM   = CdSe zincblende — optics
3PREC     = Accurate
4ENCUT    = 400
5ICHARG   = 11        # read CHGCAR; non-self-consistent
6ISMEAR   = 0
7SIGMA    = 0.05
8NBANDS   = 96        # ≈ 4–6× number of valence bands; converge!
9LOPTICS  = .TRUE.    # compute frequency-dependent dielectric tensor
10NEDOS    = 2000      # frequency-grid points for ε(ω)
11CSHIFT   = 0.10      # Lorentzian broadening Γ in eV
12LREAL    = .FALSE.   # full reciprocal-space projectors — accuracy first

1# KPOINTS  (step 2: dense uniform grid; bigger than the SCF grid)
2Γ-centred dense grid for optical convergence
30
4Gamma
512 12 12
60  0  0

Reading the output: vasprun.xml → pymatgen

Everything ends up in vasprun.xml. Pymatgen exposes it as a 3-tuple $(\omega_i, \varepsilon_1^{ab}(\omega_i), \varepsilon_2^{ab}(\omega_i))$ with six independent tensor components per frequency. The script below reads that tuple, derives $n, k, \alpha$ from the formulas of the previous section, and saves a publication- ready figure in eight executable lines.

1from pymatgen.io.vasp import Vasprun
2import numpy as np
3import matplotlib.pyplot as plt
4
5vr   = Vasprun("vasprun.xml", parse_potcar_file=False)
6diel = vr.dielectric                              # (energies, eps1_tensor, eps2_tensor)
7
8energies = np.array(diel[0])                      # photon energies ω in eV
9eps1_xx  = np.array([t[0] for t in diel[1]])      # real part — xx component
10eps2_xx  = np.array([t[0] for t in diel[2]])      # imaginary part — xx component
11
12mod    = np.sqrt(eps1_xx**2 + eps2_xx**2)
13n_xx   = np.sqrt((mod + eps1_xx) / 2.0)           # refractive index
14k_xx   = np.sqrt((mod - eps1_xx) / 2.0)           # extinction coefficient
15alpha  = (2 * energies * k_xx) / 1.973e-5         # absorption coefficient (1/cm)
16
17plt.plot(energies, eps2_xx, label=r"$\varepsilon_2(\omega)$")
18plt.plot(energies, alpha,   label=r"$\alpha(\omega)$")
19plt.xlabel("photon energy ℏω (eV)"); plt.legend(); plt.savefig("optics.png", dpi=200)

What Standard DFT Optics Misses — Excitons and the Gap

Everything we have computed sits inside the independent-particle approximation: an electron and the hole it leaves behind do not interact. Real life is not so forgiving. The ejected electron sees a Coulomb attraction to its own hole, and the bound electron–hole pair — an exciton — has internal binding energy of order meV (Wannier excitons in GaAs) up to hundreds of meV (Frenkel excitons in molecular crystals). In an absorption spectrum, excitons show up as sharp below-gap peaks and a redistribution of oscillator strength piling up near the gap. Standard DFT-LOPTICS misses both effects.

Take the result of this section as the floor of optical accuracy: cheap, fast, and good enough to confirm directness, find peak positions to ~10–20%, and compute static dielectric constants within a factor of two. For Mn-doped CdSe quantum dots — Chapter 7 — independent-particle DFT will tell us where the absorption begins and confirm the directness of the gap; Chapter 6's hybrid functionals will then put the band gap in the right place; a Bethe–Salpeter calculation (briefly previewed in §5.10) is what you run if you need exciton-binding energies to within meV.

Summary

• A photon in a crystal is a vertical arrow in the band structure: head on a conduction band, tail on a valence band, at the same $\mathbf{k}$, with length $\hbar\omega$. The set of allowed $\mathbf{k}$ at each $\omega$ defines an iso-energy surface in the BZ.
• The complex dielectric function $\varepsilon(\omega) = \varepsilon_1 + i\varepsilon_2$ encodes everything: refractive index, absorption, reflectivity, colour, plasmons. $\varepsilon_2$ is the dissipative channel; $\varepsilon_1$ follows from it via Kramers–Kronig.
• Fermi's golden rule gives $\varepsilon_2 \propto (1/\omega^2)\sum_{c,v}\int |\langle p \rangle|^2 \delta(E_c - E_v - \hbar\omega)$. Three ingredients: the JDOS (geometry), the matrix element (chemistry), the energy-conservation δ (band structure).
• Direct-allowed parabolic edges have a universal $\sqrt{\hbar\omega - E_g}$ onset; van Hove singularities of the JDOS are where the named optical peaks ($E_0, E_1, E_2$) of every textbook spectrum live.
• In VASP: LOPTICS=.TRUE., large NBANDS, dense k-grid, sensible CSHIFT. vasprun.xml hands you the full tensor; pymatgen turns it into $n, k, \alpha, R$ in eight lines.
• Standard DFT-LOPTICS misses excitons and underestimates the gap. Use it as a fast, qualitative spectrum; reach for hybrids (§6.4) for the gap and BSE (§5.10) for excitons when accuracy demands.

Coming next: §5.6 — Magnetic Ordering and Spin — where we duplicate the band structure into spin-up and spin-down channels, and learn how the same vertical-arrow logic tells us which materials are ferromagnetic, antiferromagnetic, or non-magnetic.