From 7d38601f2380bf7b8ad57816ab3695269794a4a5 Mon Sep 17 00:00:00 2001
From: Anton Akhmerov <anton.akhmerov@gmail.com>
Date: Mon, 14 Jan 2019 23:52:02 +0100
Subject: [PATCH] first attempt at splitting

---
 mkdocs.yml                                    |  28 ++-
 src/{lecture_5.md => 10_xray.md}              | 104 --------
 src/11_nearly_free_electron_model.md          | 146 ++++++++++++
 ...2_band_structures_in_higher_dimensions.md} | 137 +----------
 src/13_semiconductors.md                      | 222 +++++++++++++++++
 ...{lecture_7.md => 14_doping_and_devices.md} | 205 +---------------
 src/1_einstein_model.md                       | 119 ++++++++++
 src/{lecture_1.md => 2_debye_model.md}        | 115 +--------
 src/3_drude_model.md                          |  86 +++++++
 src/{lecture_2.md => 4_sommerfeld_model.md}   |  79 +-----
 src/{lecture_3.md => 5_atoms_and_lcao.md}     | 197 ---------------
 src/6_bonds_and_spectra.md                    | 208 ++++++++++++++++
 src/{lecture_4.md => 7_tight_binding.md}      | 212 -----------------
 src/8_many_atoms.md                           | 224 ++++++++++++++++++
 src/9_crystal_structure.md                    | 103 ++++++++
 15 files changed, 1139 insertions(+), 1046 deletions(-)
 rename src/{lecture_5.md => 10_xray.md} (56%)
 create mode 100644 src/11_nearly_free_electron_model.md
 rename src/{lecture_6.md => 12_band_structures_in_higher_dimensions.md} (53%)
 create mode 100644 src/13_semiconductors.md
 rename src/{lecture_7.md => 14_doping_and_devices.md} (53%)
 create mode 100644 src/1_einstein_model.md
 rename src/{lecture_1.md => 2_debye_model.md} (60%)
 create mode 100644 src/3_drude_model.md
 rename src/{lecture_2.md => 4_sommerfeld_model.md} (68%)
 rename src/{lecture_3.md => 5_atoms_and_lcao.md} (53%)
 create mode 100644 src/6_bonds_and_spectra.md
 rename src/{lecture_4.md => 7_tight_binding.md} (50%)
 create mode 100644 src/8_many_atoms.md
 create mode 100644 src/9_crystal_structure.md

diff --git a/mkdocs.yml b/mkdocs.yml
index 1bc33959..12ce6b19 100644
--- a/mkdocs.yml
+++ b/mkdocs.yml
@@ -6,13 +6,27 @@ site_description: |
   Lecture notes for the TU Delft course TN2844 Solid state physics
 nav:
   - Intro: 'index.md'
-  - Phonons and specific Heat: 'lecture_1.md'
-  - Free electron model: 'lecture_2.md'
-  - Atoms and bonds: 'lecture_3.md'
-  - Electrons and phonons in 1D: 'lecture_4.md'
-  - Crystal structure and diffraction: 'lecture_5.md'
-  - Tight binding and nearly free electrons: 'lecture_6.md'
-  - Semiconductors: 'lecture_7.md'
+  - Phonons and specific Heat:
+      - Einstein model: '1_einstein_model.md'
+      - Debye model: '2_debye_model.md'
+  - Free electron model:
+      - Drude model: '3_drude_model.md'
+      - Sommerfeld model: '4_sommerfeld_model.md'
+  - Atoms and bonds:
+      - LCAO model: '5_atoms_and_lcao.md'
+      - Bonds and spectra: '6_bonds_and_spectra.md'
+  - Electrons and phonons in 1D:
+    - Tight-binding model: '7_tight_binding.md'
+    - Many atoms per unit cell: '8_many_atoms.md'
+  - Crystal structure and diffraction:
+    - Crystal structure: '9_crystal_structure.md'
+    - X-ray diffraction: '10_xray.md'
+  - Tight binding and nearly free electrons:
+    - Nearly free electron model: '11_nearly_free_electron_model.md'
+    - Band structures in 2D: '12_band_structures_in_higher_dimensions.md'
+  - Semiconductors:
+    - Basic principles: '13_semiconductors.md'
+    - Doping and devices: '14_doping_and_devices.md'
   - Attic:
     - Magnetism: 'lecture_8.md'
     - Fermi surface periodic table: 'fermi_surfaces.md'
diff --git a/src/lecture_5.md b/src/10_xray.md
similarity index 56%
rename from src/lecture_5.md
rename to src/10_xray.md
index b1c6920b..5b9f7edc 100644
--- a/src/lecture_5.md
+++ b/src/10_xray.md
@@ -1,107 +1,3 @@
-# Lecture 5 – Crystal structure and diffraction
-
-_based on chapters 12–14, (up to and including 14.2) of the book_  
-Exercises 12.3, 12.4, 13.3, 13.4, 14.2
-
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - Describe any crystal using crystallographic terminology, and interpret this terminology
-    - Compute the volume filling fraction given a crystal structure
-    - Determine the primitive, conventional, and Wigner-Seitz unit cells of a given lattice
-    - Determine the Miller planes of a given lattice
-
-### Crystal classification
-
-- **_Lattice_**
-    + periodic pattern of *lattice points*, which all have an identical view
-    + lattice points are not necessarily the same as atom positions
-    + there can be multiple atoms per lattice point
-    + freedom of translation
-    + multiple lattices with different point densities possible
-- **_Lattice vectors_**
-    + from lattice point to lattice point
-    + $N$ vectors for $N$ dimensions
-    + multiple combinations possible
-    + not all combinations provide full coverage
-- **_Unit cell_**
-    + spanned by lattice vectors
-    + has 4 corners in 2D, 8 corners in 3D
-    + check if copying unit cell along lattice vectors gives full lattice
-- **_Primitive unit cell_**
-    + smallest possible $\rightarrow$ no identical points skipped
-    + not always most practical choice
-- **_Basis_**
-    + only now we care about the contents (i.e. atoms)
-    + gives element and position of atoms
-    + properly count partial atoms $\rightarrow$ choose which belongs to unit cell
-    + positions in terms of lattice vectors, *not* Cartesian coordinates!
-
-### Example: graphite
-
-![](figures/graphite_mod.svg)
-
-1. Choose origin (can be atom, not necessary)
-2. Find other lattice points that are identical
-3. Choose lattice vectors, either primitive (red) or not primitive (blue)
-    - lengths of lattice vectors and angle(s) between them fully define the crystal lattice
-    - for graphite: $|{\bf a}_1|=|{\bf a}_2|$ = 0.246 nm = 2.46 Ã…, $\gamma$ = 60$^{\circ}$
-4. Specify basis
-    - using ${\bf a}_1$ and ${\bf a}_2$: C$(0,0)$, C$\left(\frac{2}{3},\frac{2}{3}\right)$
-    - using ${\bf a}_1$ and ${\bf a}_{2}'$: C$(0,0)$, C$\left(0,\frac{1}{3}\right)$, C$\left(\frac{1}{2},\frac{1}{2}\right)$, C$\left(\frac{1}{2},\frac{5}{6}\right)$
-
-An alternative type of unit cell is the _Wigner-Seitz cell_: the collection of all points that are closer to one specific lattice point than to any other lattice point. You form this cell by taking all the perpendicular bisectrices or lines connecting a lattice point to its neighboring lattice points.
-
-### Stacking of atoms
-What determines what crystal structure a material adopts? $\rightarrow$ To good approximation, atoms are solid incrompressible spheres that attract each other. How will these organize?
-
-We start with the densest possible packing in 2D:
-![](figures/packing.svg)
-Will the second layer go on sites A, B or C?
-
-ABCABC stacking $\rightarrow$ _cubic close packed_, also known as _face centered cubic_ (fcc):
-
-![](figures/stacking.svg)
-
-- One atom on the center of each side-plane: 'a dice that always throws 1'
-- Conventional unit cell $\neq$ primitive unit cell
-- Cyclic ABC $\rightarrow$ all atoms identical $\rightarrow$ 1 atom per primitive unit cell
-- Conventional cell: $8\times\frac{1}{8}+6\times\frac{1}{2}=1+3=4$ atoms
-
-Examples of fcc crystals: Al, Cu, Ni, Pb, Au and Ag.
-
-### Filling factor
-Filling factor = # of atoms per cell $\times$ volume of 1 atom / volume of cell
-
-$=\frac{4\times\frac{4}{3}\pi R^3}{a^3}=\frac{1}{6}\sqrt{2}\pi\approx 0.74$, where we have used that $\sqrt{2}a=4R$.
-
-Compare this to _body centered cubic_ (bcc), which consists of a cube with atoms on the corners and one atom in the center (Fig. 1.12): filling factor = 0.68. Examples of bcc crystals: Fe, Na, W, Nb.
-
-Question: is 74% the largest possible filling factor? $\rightarrow$ Kepler conjecture (1571 – 1630). Positive proof by Hales _et al._ in 2015!
-
-Crystal structures that are related to fcc:
-
-1. ionic crystals (Fig. 1.13), e.g. NaCl
-2. zincblende (Fig. 1.15), e.g. diamond
-
-ABABAB stacking $\rightarrow$ _hexagonally close-packed_ (hcp), e.g. Co, Zn. In this case there is no cubic symmetry (Fig. 1.11).
-
-### Miller planes
-We start with a simple cubic lattice:
-
-![](figures/cubic_mod.svg)
-
-$|{\bf a}_1|=|{\bf a}_2|=|{\bf a}_3|\equiv a$ (_lattice constant_)
-
-The plane designated by Miller indices $(u,v,w)$ intersects lattice vector ${\bf a}_1$ at $\frac{|{\bf a}_1|}{u}$, ${\bf a}_2$ at $\frac{|{\bf a}_2|}{v}$ and ${\bf a}_3$ at $\frac{|{\bf a}_3|}{w}$.
-
-![](figures/miller.svg)
-
-Miller index 0 means that the plane is parallel to that axis (intersection at "$\frac{|{\bf a}_3|}{0}=\infty$"). A bar above a Miller index means intersection at a negative coordinate.
-
-If a crystal is symmetric under $90^\circ$ rotations, then $(100)$, $(010)$ and $(001)$ are physically indistinguishable. This is indicated with $\{100\}$. $[100]$ is a vector. In a cubic crystal, $[100]$ is perpendicular to $(100)$ $\rightarrow$ proof in problem set.
-
 !!! summary "Learning goals"
 
     After this lecture you will be able to:
diff --git a/src/11_nearly_free_electron_model.md b/src/11_nearly_free_electron_model.md
new file mode 100644
index 00000000..6ea31457
--- /dev/null
+++ b/src/11_nearly_free_electron_model.md
@@ -0,0 +1,146 @@
+```python
+import numpy as np
+import plotly.offline as py
+import plotly.graph_objs as go
+
+pi = np.pi
+```
+
+# Tight binding and nearly free electrons
+_(based on chapters 15–16 of the book)_  
+Exercises 15.1, 15.3, 15.4, 16.1, 16.2
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - formulate a general way of computing the electron band structure - the **Bloch theorem**.
+    - recall that in a periodic potential all electron states are Bloch waves.
+    - derive the electron band structure when the interaction with the lattice is weak using the **Nearly free electron model**.
+
+Let's summarize what we learned about electrons so far:
+
+* Free electrons form a Fermi sea ([lecture 4](4_sommerfeld_model.md))
+* Isolated atoms have discrete orbitals ([lecture 5](5_atoms_and_lcao.md))
+* When orbitals hybridize we get *LCAO* or *tight-binding* band structures ([lecture 7](7_tight_binding.md))
+
+The nearly free electron model (the topic of this lecture) helps to understand the relation between tight-binding and free electron models. It describes the properties of metals.
+
+These different models can be organized as a function of the strength of the lattice potential $V(x)$:
+
+![](figures/models.svg)
+
+## Bloch theorem
+
+> All Hamiltonian eigenstates in a crystal have the form
+> $$ \psi_n(\mathbf{r}) = u_n(\mathbf{r})e^{i\mathbf{kr}} $$
+> with $u_n(\mathbf{r})$ having the same periodicity as the lattice potential $V(\mathbf{r})$, and index $n$ labeling electron bands with energies $E_n(\mathbf{k})$.
+
+In other words: any electron wave function in a crystal is a product of a periodic part that describes electron motion within a unit cell and a plane wave.
+
+### Extra remarks
+
+The wave function $u_n(\mathbf{r})e^{i\mathbf{kr}}$ is called a **Bloch wave**.
+
+The $u_n(\mathbf{r})$ part is some unknown function. To calculate it we need to solve the Schrödinger equation. It is hard in general, but there are two limits when $U$ is "weak" and $U$ is "large" that provide us with most intuition.
+
+If we change $\mathbf{k}$ by a reciprocal lattice vector $\mathbf{k} \rightarrow \mathbf{k} + h\mathbf{b}_1 + k\mathbf{b}_2 + l\mathbf{b}_3$, and we change $u_n(\mathbf{r}) \rightarrow u_n(\mathbf{r})\exp\left[-h\mathbf{b}_1 - k\mathbf{b}_2 - l\mathbf{b}_3\right]$ (also periodic!), we obtain the same wave function. Therefore energies of all bands $E_n(\mathbf{k})$ are periodic in reciprocal space with the periodicity of the reciprocal lattice.
+
+Bloch theorem is extremely similar to the ansatz we used in [1D](7_tight_binding.md), and to the description of the [X-ray scattering](10_xray.md).
+
+## Nearly free electron model
+
+In free electron model $E = \hbar^2 k^2/2m$.
+
+* There is only one band
+* The band structure is not periodic in $k$-space
+* In other words the Brillouin zone is infinite in momentum space
+
+Within the **nearly free electron model** we want to start from the dispersion relation of the free electrons and consider the effect of introducing a weak lattice potential. The logic is very similar to getting optical and acoustic phonon branches by changing atom masses (and therefore reducing the size of the Brillouin zone).
+
+In our situation:
+
+![](figures/nearly_free_electron_bands.svg)
+
+The band gaps open where two copies of the free electron dispersion cross.
+
+### Avoided level crossing
+
+*Remark: this is an important concept in quantum mechanics, based on the perturbation theory. You will only learn it later in QMIII, so we will need to postulate some important facts.*
+
+Let's focus on the first crossing. The momentum near it is $k = \pi/a + \delta k$ and we have two copies of the original band structure coming together. One with $\psi_+ \propto e^{i\pi x/a}$, another with $\psi_- \propto e^{-i\pi x/a}$. Near the crossing the wave function is the linear superposition of $\psi_+$ and $\psi_-$: $\psi = \alpha \psi_+ + \beta \psi_-$. We actually used almost the same form of the wave function in LCAO, except instead of $\psi_\pm$ we used the orbitals $\phi_1$ and $\phi_2$ there.
+
+Without the lattice potential we can approximate the Hamiltonian of these two states as follows:
+$$H\begin{pmatrix}\alpha \\ \beta \end{pmatrix} =
+\begin{pmatrix} E_0 + v \hbar \delta k & 0 \\ 0 & E_0 - v \hbar \delta k\end{pmatrix}
+\begin{pmatrix}\alpha \\ \beta \end{pmatrix}.
+$$
+
+Here we used $\delta p = \hbar \delta k$, and we expanded the quadratic function into a linear plus a small correction.  
+
+??? question "calculate $E_0$ and the velocity $v$"
+    The edge of the Brilloin zone has $k = \pi/a$. Substituting this in the free electron dispersion $E = \hbar^2 k^2/2m$ we get $E_0 = \hbar^2 \pi^2/2 m a^2$, and $v=\hbar p/m=\hbar \pi/ma$.
+
+Without $V(x)$ the two wave functions $\psi_+$ and $\psi_-$ are independent since they have a different momentum. When $V(x)$ is present, it may couple these two states.
+
+So in presence of $V(x)$ the Hamiltonian becomes
+
+$$
+H\begin{pmatrix}\alpha \\ \beta \end{pmatrix} =
+\begin{pmatrix} E_0 + v \hbar \delta k & W \\ W^* & E_0 - v \hbar \delta k\end{pmatrix}
+\begin{pmatrix}\alpha \\ \beta \end{pmatrix},
+$$
+
+Here the coupling strength $W = \langle \psi_+ | V(x) | \psi_- \rangle$ is the matrix element of the potential between two states. *(This where we need to apply the perturbation theory, and this is very similar to the LCAO Hamiltonian)*.
+
+??? question "how does our solution satisfy the Bloch theorem? What is $u(x)$ in this case?"
+    The wave function has a form $\psi(x) = \alpha \exp[ikx] + \beta \exp[i(k - 2\pi/a)x]$
+    (here $k = \pi/a + \delta k$). Choosing $u(x) = \alpha + \beta \exp(2\pi i x/a)$ we see
+    that $\psi(x) = u(x) \exp(ikx)$.
+
+#### Dispersion relation near the avoided level crossing
+
+We need to diagonalize a 2x2 matrix Hamiltonian. The answer is
+$$ E(\delta k) = E_0 \pm \sqrt{v^2\hbar^2\delta k^2 + |W|^2}$$
+
+Check out section 15.1.1 of the book for the details of this calculation.
+
+#### Physical meaning of $W$
+
+Now we expand the definition of $W$:
+
+$$W = \langle \psi_+ | V(x) | \psi_- \rangle = \int_0^{a} dx \left[e^{i\pi x/a}\right]^* V(x) \left[e^{-i\pi x/a}\right] = \int_0^a e^{2\pi i x /a} V(x) = V_1$$
+
+Here $V_1$ is the first Fourier component of $V(x)$ (using a complex Fourier transform).
+
+$$ V(x) = \sum_{n=-\infty}^{\infty} V_n e^{2\pi i n x/a}$$
+
+#### Crossings between the higher bands
+
+Everything we did applies to the crossings at higher energies, only there we would get higher Fourier components of $V(x)$: $V_2$ for the crossing between the second and third band, $V_3$ for the crossing between third and 4th, etc.
+
+### Repeated vs reduced vs extended Brillouin zone
+
+All different ways to **plot** the same dispersion relation (no difference in physical information).
+
+Repeated BZ (all possible Bloch bands):
+
+![](figures/nearly_free_electron_bands.svg)
+
+* Contains redundant information
+* May be easier to count/follow the bands
+
+Reduced BZ (all bands within 1st BZ):
+
+![](figures/reduced_nearly_free_electron_bands.svg)
+
+* No redundant information
+* Hard to relate to original dispersion
+
+Extended BZ (n-th band within n-th BZ):
+
+![](figures/extended_nearly_free_electron_bands.svg)
+
+* No redundant information
+* Easy to relate to free electron model
+* Contains discontinuities
diff --git a/src/lecture_6.md b/src/12_band_structures_in_higher_dimensions.md
similarity index 53%
rename from src/lecture_6.md
rename to src/12_band_structures_in_higher_dimensions.md
index 430a3439..476eae33 100644
--- a/src/lecture_6.md
+++ b/src/12_band_structures_in_higher_dimensions.md
@@ -10,148 +10,13 @@ pi = np.pi
 _(based on chapters 15–16 of the book)_  
 Exercises 15.1, 15.3, 15.4, 16.1, 16.2
 
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - formulate a general way of computing the electron band structure - the **Bloch theorem**.
-    - recall that in a periodic potential all electron states are Bloch waves.
-    - derive the electron band structure when the interaction with the lattice is weak using the **Nearly free electron model**.
-
-Let's summarize what we learned about electrons so far:
-
-* Free electrons form a Fermi sea ([lecture 2](lecture_1.md))
-* Isolated atoms have discrete orbitals ([lecture 3](lecture_3.md))
-* When orbitals hybridize we get *LCAO* or *tight-binding* band structures ([lecture 4](lecture_4.md))
-
-The nearly free electron model (the topic of this lecture) helps to understand the relation between tight-binding and free electron models. It describes the properties of metals.
-
-These different models can be organized as a function of the strength of the lattice potential $V(x)$:
-
-![](figures/models.svg)
-
-## Bloch theorem
-
-> All Hamiltonian eigenstates in a crystal have the form
-> $$ \psi_n(\mathbf{r}) = u_n(\mathbf{r})e^{i\mathbf{kr}} $$
-> with $u_n(\mathbf{r})$ having the same periodicity as the lattice potential $V(\mathbf{r})$, and index $n$ labeling electron bands with energies $E_n(\mathbf{k})$.
-
-In other words: any electron wave function in a crystal is a product of a periodic part that describes electron motion within a unit cell and a plane wave.
-
-### Extra remarks
-
-The wave function $u_n(\mathbf{r})e^{i\mathbf{kr}}$ is called a **Bloch wave**.
-
-The $u_n(\mathbf{r})$ part is some unknown function. To calculate it we need to solve the Schrödinger equation. It is hard in general, but there are two limits when $U$ is "weak" and $U$ is "large" that provide us with most intuition.
-
-If we change $\mathbf{k}$ by a reciprocal lattice vector $\mathbf{k} \rightarrow \mathbf{k} + h\mathbf{b}_1 + k\mathbf{b}_2 + l\mathbf{b}_3$, and we change $u_n(\mathbf{r}) \rightarrow u_n(\mathbf{r})\exp\left[-h\mathbf{b}_1 - k\mathbf{b}_2 - l\mathbf{b}_3\right]$ (also periodic!), we obtain the same wave function. Therefore energies of all bands $E_n(\mathbf{k})$ are periodic in reciprocal space with the periodicity of the reciprocal lattice.
-
-Bloch theorem is extremely similar to the ansatz we used in [1D](lecture_4.md), and to the description of the [X-ray scattering](lecture_5.md).
-
-## Nearly free electron model
-
-In free electron model $E = \hbar^2 k^2/2m$.
-
-* There is only one band
-* The band structure is not periodic in $k$-space
-* In other words the Brillouin zone is infinite in momentum space
-
-Within the **nearly free electron model** we want to start from the dispersion relation of the free electrons and consider the effect of introducing a weak lattice potential. The logic is very similar to getting optical and acoustic phonon branches by changing atom masses (and therefore reducing the size of the Brillouin zone).
-
-In our situation:
-
-![](figures/nearly_free_electron_bands.svg)
-
-The band gaps open where two copies of the free electron dispersion cross.
-
-### Avoided level crossing
-
-*Remark: this is an important concept in quantum mechanics, based on the perturbation theory. You will only learn it later in QMIII, so we will need to postulate some important facts.*
-
-Let's focus on the first crossing. The momentum near it is $k = \pi/a + \delta k$ and we have two copies of the original band structure coming together. One with $\psi_+ \propto e^{i\pi x/a}$, another with $\psi_- \propto e^{-i\pi x/a}$. Near the crossing the wave function is the linear superposition of $\psi_+$ and $\psi_-$: $\psi = \alpha \psi_+ + \beta \psi_-$. We actually used almost the same form of the wave function in LCAO, except instead of $\psi_\pm$ we used the orbitals $\phi_1$ and $\phi_2$ there.
-
-Without the lattice potential we can approximate the Hamiltonian of these two states as follows:
-$$H\begin{pmatrix}\alpha \\ \beta \end{pmatrix} =
-\begin{pmatrix} E_0 + v \hbar \delta k & 0 \\ 0 & E_0 - v \hbar \delta k\end{pmatrix}
-\begin{pmatrix}\alpha \\ \beta \end{pmatrix}.
-$$
-
-Here we used $\delta p = \hbar \delta k$, and we expanded the quadratic function into a linear plus a small correction.  
-
-??? question "calculate $E_0$ and the velocity $v$"
-    The edge of the Brilloin zone has $k = \pi/a$. Substituting this in the free electron dispersion $E = \hbar^2 k^2/2m$ we get $E_0 = \hbar^2 \pi^2/2 m a^2$, and $v=\hbar p/m=\hbar \pi/ma$.
-
-Without $V(x)$ the two wave functions $\psi_+$ and $\psi_-$ are independent since they have a different momentum. When $V(x)$ is present, it may couple these two states.
-
-So in presence of $V(x)$ the Hamiltonian becomes
-
-$$
-H\begin{pmatrix}\alpha \\ \beta \end{pmatrix} =
-\begin{pmatrix} E_0 + v \hbar \delta k & W \\ W^* & E_0 - v \hbar \delta k\end{pmatrix}
-\begin{pmatrix}\alpha \\ \beta \end{pmatrix},
-$$
-
-Here the coupling strength $W = \langle \psi_+ | V(x) | \psi_- \rangle$ is the matrix element of the potential between two states. *(This where we need to apply the perturbation theory, and this is very similar to the LCAO Hamiltonian)*.
-
-??? question "how does our solution satisfy the Bloch theorem? What is $u(x)$ in this case?"
-    The wave function has a form $\psi(x) = \alpha \exp[ikx] + \beta \exp[i(k - 2\pi/a)x]$
-    (here $k = \pi/a + \delta k$). Choosing $u(x) = \alpha + \beta \exp(2\pi i x/a)$ we see
-    that $\psi(x) = u(x) \exp(ikx)$.
-
-#### Dispersion relation near the avoided level crossing
-
-We need to diagonalize a 2x2 matrix Hamiltonian. The answer is
-$$ E(\delta k) = E_0 \pm \sqrt{v^2\hbar^2\delta k^2 + |W|^2}$$
-
-Check out section 15.1.1 of the book for the details of this calculation.
-
-#### Physical meaning of $W$
-
-Now we expand the definition of $W$:
-
-$$W = \langle \psi_+ | V(x) | \psi_- \rangle = \int_0^{a} dx \left[e^{i\pi x/a}\right]^* V(x) \left[e^{-i\pi x/a}\right] = \int_0^a e^{2\pi i x /a} V(x) = V_1$$
-
-Here $V_1$ is the first Fourier component of $V(x)$ (using a complex Fourier transform).
-
-$$ V(x) = \sum_{n=-\infty}^{\infty} V_n e^{2\pi i n x/a}$$
-
-#### Crossings between the higher bands
-
-Everything we did applies to the crossings at higher energies, only there we would get higher Fourier components of $V(x)$: $V_2$ for the crossing between the second and third band, $V_3$ for the crossing between third and 4th, etc.
-
-### Repeated vs reduced vs extended Brillouin zone
-
-All different ways to **plot** the same dispersion relation (no difference in physical information).
-
-Repeated BZ (all possible Bloch bands):
-
-![](figures/nearly_free_electron_bands.svg)
-
-* Contains redundant information
-* May be easier to count/follow the bands
-
-Reduced BZ (all bands within 1st BZ):
-
-![](figures/reduced_nearly_free_electron_bands.svg)
-
-* No redundant information
-* Hard to relate to original dispersion
-
-Extended BZ (n-th band within n-th BZ):
-
-![](figures/extended_nearly_free_electron_bands.svg)
-
-* No redundant information
-* Easy to relate to free electron model
-* Contains discontinuities
-
 !!! summary "Learning goals"
 
     After this lecture you will be able to:
 
     - examine 1D and 2D band structures and argue if you expect the corresponding material to be an insulator/semiconductor or a conductor.
     - describe how the light absorption spectrum of a material relates to its band structure.
-    
+
 ## Band structure
 
 How are material properties related to the band structure?
diff --git a/src/13_semiconductors.md b/src/13_semiconductors.md
new file mode 100644
index 00000000..feb4c260
--- /dev/null
+++ b/src/13_semiconductors.md
@@ -0,0 +1,222 @@
+```python
+from matplotlib import pyplot
+
+import numpy as np
+from scipy.optimize import curve_fit
+from scipy.integrate import quad
+
+from common import draw_classic_axes, configure_plotting
+
+configure_plotting()
+
+def sqrt_plus(x):
+    return np.sqrt(x * (x >= 0))
+
+# Band structure parameters.
+E_V, E_C, E_F = -1.2, 1.8, .4
+E_D, E_A = E_C - .7, E_V + .5
+m_h, m_e = 1, .5
+```
+
+# Semiconductor physics
+_(based on chapters 17–18 of the book)_  
+Exercises 17.1–17.4, 18.1, 18.2
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - Describe the concept of holes, and apply this concept to describe the properties of semiconductors
+    - Compute the density of states of electrons in semiconductors
+    - Compute the density of charge carriers and chemical potential as a function of temperature
+
+## Review of band structure properties
+
+* Group velocity $v=\hbar^{-1}\partial E(k)/\partial k$
+* Effective mass $m_{eff} = \hbar^2\left(d^2 E(k)/dk^2\right)^{-1}$
+* Density of states $g(E) = \sum_{\textrm{FS}}dk \times (dn/dk) \times (dk/dE)$ the amount of states per infinitesimal interval of energy at each energy.
+
+A simple check that everything is correct is the free electron model:
+
+$$H = \hbar^2 k^2/2m$$
+
+The velocity is $\hbar^{-1}\partial E(k)/\partial k = \hbar k / m \equiv p/m$.  
+The effective mass is $m_{eff} = \hbar^2\left(d^2 E(k)/dk^2\right)^{-1} = m$.
+
+So in the simplest case the definitions match the usual expressions.
+
+#### Why are these properties important?
+
+* Group velocity descibes how quickly electrons interacting with the lattice move.
+* Effective mass tells us how hard it is to *accelerate* the particles and it enters $g(E)$ for a parabolic band
+* Density of states is required to determine the Fermi level, heat capacity, etc.
+
+## Filled vs empty bands
+
+A completely filled band is very similar to a completely empty band.
+
+In a filled band $n(E)=1$ because $|E - E_F| \gg kT$. In an empty band $n(E)=0$.
+
+Heat capacity $C_v = \tfrac{d}{dT}\int_{-\infty}^\infty E\times g(E) \times dE\times n_F(E, T) = 0$.
+
+A completely filled band carries no electric current:
+
+$$
+j = 2e \int_{-\pi/a}^{\pi/a} v(k) dk = 2e \int_{-\pi/a}^{\pi/a} (dE/dk)\times dk = 2e [E(-\pi/a) - E(\pi/a)] = 0
+$$
+
+## From electrons to holes
+
+A completely filled band $\approx$ completely empty band. The idea of introducing **holes** is to transform an *almost* completely filled band $\Rightarrow$ almost empty one. Instead of describing a lot of electrons that are present, we focus on those that are absent.
+
+Definition: a **hole** is a state of a completely filled band with one particle missing.
+
+![](figures/holes.svg)
+
+In this schematic we can either say that 8×2 electron states are occupied (the system has 8×2 electrons counting spin), or 10×2 hole states are occupied. A useful analogy to remember: glass half-full or glass half-empty.
+
+## Properties of holes
+
+The probability for an electron state to be occupied in equilibrium is $f(E)$:
+
+$$f(E) = \frac{1}{e^{(E-E_F)/kT} + 1}$$
+
+The probability for a hole state to be occupied:
+
+$$1 - f(E) = 1 - \frac{1}{e^{(E-E_F)/kT} + 1} = \frac{1}{e^{(-E+E_F)/kT} + 1}$$
+
+therefore for holes both energy $E_h= -E$ and $E_{F,h} = -E_F$.
+
+The **momentum** $p_h$ of a hole should give the correct total momentum of a partially filled band if one sums momenta of all holes. Therefore $p_h = -\hbar k$, where $k$ is the wave vector of the electron.
+
+Similarly, the total **charge** should be the same regardless of whether we count electrons or holes, so holes have a positive charge $+e$ (electrons have $-e$).
+
+On the other hand, the velocity of a hole is **the same**:
+$$\frac{dE_h}{dp_h} = \frac{-dE}{-d\hbar k} = \frac{dE}{dp}$$
+
+Finally, we derive the hole effective mass from the equations of motion:
+
+$$m_h \frac{d v}{d t} = +e (E + v\times B)$$
+
+Comparing with
+
+$$m_e \frac{d v}{d t} = -e (E + v\times B)$$
+
+we get $m_h = -m_e$.
+
+## Semiconductors: materials with two bands.
+
+In semiconductors the Fermi level is between two bands. The unoccupied band is the **conduction band**, the occupied one is the **valence band**. In the conduction band the **charge carriers** (particles carrying electric current) are electrons, in the valence band they are holes.
+
+We can control the position of the Fermi level (or create additional excitations) making semiconductors conduct when needed.
+
+Only the bottom of the conduction band has electrons and the top of the valence band has holes because the temperature is always smaller than the width of the band.
+
+Therefore we can approximate the dispersion relation of both bands as parabolic.
+
+![](figures/semiconductor.svg)
+
+Or in other words
+
+$$E_e = E_G + \frac{\hbar^2k^2}{2m_e},$$
+$$E_h = \frac{\hbar^2k^2}{2m_h}.$$
+
+Here $E_G$ is the band gap, and the top of the valence band is at $E=0$.
+
+Observe that because we are describing particles in the valence band as holes, $m_h > 0$ and $E_h > 0$.
+
+??? question "a photon gives a single electron enough energy to move from the valence band to the conduction band. How many particles does this process create?"
+    Two: one electron and one hole.
+
+## Semiconductor density of states and Fermi level
+
+### Part 1: pristine semiconductor
+
+```python
+E = np.linspace(-3, 3, 1000)
+fig, ax = pyplot.subplots()
+
+n_F = 1/(np.exp(2*(E - E_F)) + 1)
+g_e = m_e * sqrt_plus(E - E_C)
+g_h = m_h * sqrt_plus(E_V - E)
+ax.plot(E, g_h, label="$g_e$")
+ax.plot(E, g_e, label="$g_h$")
+ax.plot(E, 10 * g_h * (1-n_F), label="$n_h$")
+ax.plot(E, 10 * g_e * n_F, label="$n_e$")
+ax.plot(E, n_F, label="$n_F$", linestyle='dashed')
+ax.set_ylim(top=1.5)
+
+ax.set_xlabel('$E$')
+ax.set_ylabel('$g$')
+ax.set_xticks([E_V, E_C, E_F])
+ax.set_xticklabels(['$E_V$', '$E_C$', '$E_F$'])
+ax.legend()
+draw_classic_axes(ax, xlabeloffset=.2)
+```
+
+Electrons
+$$ E = E_G + {p^2}/{2m_e}$$
+
+$$ g(E) = (2m_e)^{3/2}\sqrt{E-E_G}/2\pi^2\hbar^3$$
+
+Holes
+$$ E_h = {p^2}/{2m_h}$$
+
+$$ g(E_h) = (2m_h)^{3/2}\sqrt{E_h}/2\pi^2\hbar^3$$
+
+**The key algorithm of describing the state of a semiconductor:**
+
+1. Write down the density of states, assuming a certain position of the Fermi level
+2. Calculate the total amount of electrons and holes, equate the difference to the total amount of electrons $-$ holes available.
+3. Use physics intuition to simplify the equations (this is important!)
+4. Find $E_F$ and concentrations of electrons and holes
+
+Applying the algorithm:
+
+$$n_h = \int_0^\infty f_h(E+E_F)g_h(E)dE = \int_0^\infty\frac{(2m_h)^{3/2}}{2\pi^2\hbar^3} \sqrt{E}\frac{1}{e^{(E+E_F)/kT}+1}dE$$
+
+$$n_e = \int_{E_G}^\infty f_h(E-E_F)g_e(E)dE = \int_{E_G}^\infty\frac{(2m_e)^{3/2}}{2\pi^2\hbar^3} \sqrt{E-E_G}\frac{1}{e^{(E-E_F)/kT}+1}dE$$
+
+We need to solve $n_e = n_h$
+
+Simplification:
+Fermi level is far from both bands $E_F \gg kT$ and $E_G - E_F \gg kT$
+
+Therefore Fermi-Dirac distribution is approximately similar to Boltzmann distribution.
+
+$$f(E\pm E_F) \approx e^{-(E\pm E_F)/kT}$$
+
+Now we can calculate $n_e$ and $n_h$:
+
+$$n_h \approx \frac{V(2m_h)^{3/2}}{2\pi^2\hbar^3}e^{-E_F/kT} \int_0^\infty\sqrt{E}e^{-E/kT}dE =
+N_V e^{-E_F/kT},$$
+
+with
+
+$$N_V = 2\left(\frac{2\pi m_h kT}{h^2}\right)^{3/2}$$
+
+the number of holes with energy $E<T$ (compare with the rule above).
+
+??? question "how large is $N_V$ at room temperature? (hard question)"
+    If $kT \sim 1\textrm{eV}$ (the typical energy size of a band), then electrons in the whole band may be excited and $N_V \sim 1$. On the other hand, $N_V \sim T^{3/2}$ Therefore $N_V \sim (kT/1 \textrm{eV})^{3/2}\sim 1\%$.
+
+Similarly for electrons:
+
+$$n_e = N_C e^{-(E_G - E_F)/kT},\quad N_C = 2\left(\frac{2\pi m_e kT}{h^2}\right)^{3/2}$$
+
+
+Combining everything together:
+
+$$n_h \approx N_V e^{-E_F/kT} = N_C e^{-(E_G-E_F)/kT} \approx n_e$$
+
+Solving for $E_F$:
+
+$$E_F = \frac{E_G}{2} - \frac{3}{4}kT\log(m_e/m_h)$$
+
+An extra observation: regardless of where $E_F$ is located, $n_e n_h = N_C N_V e^{-E_G/kT} \equiv n_i^2$.
+
+$n_i$ is the **intrinsic carrier concentration**, and for a pristine semiconductor $n_e = n_h = n_i$.
+
+> The equation
+> $$n_e n_h = n_i^2$$
+> is the **law of mass action**. The name is borrowed from chemistry, and describes the equilibrium concentration of two reagents in a reaction $A+B \leftrightarrow AB$. Here electrons and hole constantly split and recombine.
diff --git a/src/lecture_7.md b/src/14_doping_and_devices.md
similarity index 53%
rename from src/lecture_7.md
rename to src/14_doping_and_devices.md
index 335b8a39..ebeda3db 100644
--- a/src/lecture_7.md
+++ b/src/14_doping_and_devices.md
@@ -18,209 +18,6 @@ E_D, E_A = E_C - .7, E_V + .5
 m_h, m_e = 1, .5
 ```
 
-# Semiconductor physics
-_(based on chapters 17–18 of the book)_  
-Exercises 17.1–17.4, 18.1, 18.2
-
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - Describe the concept of holes, and apply this concept to describe the properties of semiconductors
-    - Compute the density of states of electrons in semiconductors
-    - Compute the density of charge carriers and chemical potential as a function of temperature
-
-## Review of band structure properties
-
-* Group velocity $v=\hbar^{-1}\partial E(k)/\partial k$
-* Effective mass $m_{eff} = \hbar^2\left(d^2 E(k)/dk^2\right)^{-1}$
-* Density of states $g(E) = \sum_{\textrm{FS}}dk \times (dn/dk) \times (dk/dE)$ the amount of states per infinitesimal interval of energy at each energy.
-
-A simple check that everything is correct is the free electron model:
-
-$$H = \hbar^2 k^2/2m$$
-
-The velocity is $\hbar^{-1}\partial E(k)/\partial k = \hbar k / m \equiv p/m$.  
-The effective mass is $m_{eff} = \hbar^2\left(d^2 E(k)/dk^2\right)^{-1} = m$.
-
-So in the simplest case the definitions match the usual expressions.
-
-#### Why are these properties important?
-
-* Group velocity descibes how quickly electrons interacting with the lattice move.
-* Effective mass tells us how hard it is to *accelerate* the particles and it enters $g(E)$ for a parabolic band
-* Density of states is required to determine the Fermi level, heat capacity, etc.
-
-## Filled vs empty bands
-
-A completely filled band is very similar to a completely empty band.
-
-In a filled band $n(E)=1$ because $|E - E_F| \gg kT$. In an empty band $n(E)=0$.
-
-Heat capacity $C_v = \tfrac{d}{dT}\int_{-\infty}^\infty E\times g(E) \times dE\times n_F(E, T) = 0$.
-
-A completely filled band carries no electric current:
-
-$$
-j = 2e \int_{-\pi/a}^{\pi/a} v(k) dk = 2e \int_{-\pi/a}^{\pi/a} (dE/dk)\times dk = 2e [E(-\pi/a) - E(\pi/a)] = 0
-$$
-
-## From electrons to holes
-
-A completely filled band $\approx$ completely empty band. The idea of introducing **holes** is to transform an *almost* completely filled band $\Rightarrow$ almost empty one. Instead of describing a lot of electrons that are present, we focus on those that are absent.
-
-Definition: a **hole** is a state of a completely filled band with one particle missing.
-
-![](figures/holes.svg)
-
-In this schematic we can either say that 8×2 electron states are occupied (the system has 8×2 electrons counting spin), or 10×2 hole states are occupied. A useful analogy to remember: glass half-full or glass half-empty.
-
-## Properties of holes
-
-The probability for an electron state to be occupied in equilibrium is $f(E)$:
-
-$$f(E) = \frac{1}{e^{(E-E_F)/kT} + 1}$$
-
-The probability for a hole state to be occupied:
-
-$$1 - f(E) = 1 - \frac{1}{e^{(E-E_F)/kT} + 1} = \frac{1}{e^{(-E+E_F)/kT} + 1}$$
-
-therefore for holes both energy $E_h= -E$ and $E_{F,h} = -E_F$.
-
-The **momentum** $p_h$ of a hole should give the correct total momentum of a partially filled band if one sums momenta of all holes. Therefore $p_h = -\hbar k$, where $k$ is the wave vector of the electron.
-
-Similarly, the total **charge** should be the same regardless of whether we count electrons or holes, so holes have a positive charge $+e$ (electrons have $-e$).
-
-On the other hand, the velocity of a hole is **the same**:
-$$\frac{dE_h}{dp_h} = \frac{-dE}{-d\hbar k} = \frac{dE}{dp}$$
-
-Finally, we derive the hole effective mass from the equations of motion:
-
-$$m_h \frac{d v}{d t} = +e (E + v\times B)$$
-
-Comparing with
-
-$$m_e \frac{d v}{d t} = -e (E + v\times B)$$
-
-we get $m_h = -m_e$.
-
-## Semiconductors: materials with two bands.
-
-In semiconductors the Fermi level is between two bands. The unoccupied band is the **conduction band**, the occupied one is the **valence band**. In the conduction band the **charge carriers** (particles carrying electric current) are electrons, in the valence band they are holes.
-
-We can control the position of the Fermi level (or create additional excitations) making semiconductors conduct when needed.
-
-Only the bottom of the conduction band has electrons and the top of the valence band has holes because the temperature is always smaller than the width of the band.
-
-Therefore we can approximate the dispersion relation of both bands as parabolic.
-
-![](figures/semiconductor.svg)
-
-Or in other words
-
-$$E_e = E_G + \frac{\hbar^2k^2}{2m_e},$$
-$$E_h = \frac{\hbar^2k^2}{2m_h}.$$
-
-Here $E_G$ is the band gap, and the top of the valence band is at $E=0$.
-
-Observe that because we are describing particles in the valence band as holes, $m_h > 0$ and $E_h > 0$.
-
-??? question "a photon gives a single electron enough energy to move from the valence band to the conduction band. How many particles does this process create?"
-    Two: one electron and one hole.
-
-## Semiconductor density of states and Fermi level
-
-### Part 1: pristine semiconductor
-
-```python
-E = np.linspace(-3, 3, 1000)
-fig, ax = pyplot.subplots()
-
-n_F = 1/(np.exp(2*(E - E_F)) + 1)
-g_e = m_e * sqrt_plus(E - E_C)
-g_h = m_h * sqrt_plus(E_V - E)
-ax.plot(E, g_h, label="$g_e$")
-ax.plot(E, g_e, label="$g_h$")
-ax.plot(E, 10 * g_h * (1-n_F), label="$n_h$")
-ax.plot(E, 10 * g_e * n_F, label="$n_e$")
-ax.plot(E, n_F, label="$n_F$", linestyle='dashed')
-ax.set_ylim(top=1.5)
-
-ax.set_xlabel('$E$')
-ax.set_ylabel('$g$')
-ax.set_xticks([E_V, E_C, E_F])
-ax.set_xticklabels(['$E_V$', '$E_C$', '$E_F$'])
-ax.legend()
-draw_classic_axes(ax, xlabeloffset=.2)
-```
-
-Electrons
-$$ E = E_G + {p^2}/{2m_e}$$
-
-$$ g(E) = (2m_e)^{3/2}\sqrt{E-E_G}/2\pi^2\hbar^3$$
-
-Holes
-$$ E_h = {p^2}/{2m_h}$$
-
-$$ g(E_h) = (2m_h)^{3/2}\sqrt{E_h}/2\pi^2\hbar^3$$
-
-**The key algorithm of describing the state of a semiconductor:**
-
-1. Write down the density of states, assuming a certain position of the Fermi level
-2. Calculate the total amount of electrons and holes, equate the difference to the total amount of electrons $-$ holes available.
-3. Use physics intuition to simplify the equations (this is important!)
-4. Find $E_F$ and concentrations of electrons and holes
-
-Applying the algorithm:
-
-$$n_h = \int_0^\infty f_h(E+E_F)g_h(E)dE = \int_0^\infty\frac{(2m_h)^{3/2}}{2\pi^2\hbar^3} \sqrt{E}\frac{1}{e^{(E+E_F)/kT}+1}dE$$
-
-$$n_e = \int_{E_G}^\infty f_h(E-E_F)g_e(E)dE = \int_{E_G}^\infty\frac{(2m_e)^{3/2}}{2\pi^2\hbar^3} \sqrt{E-E_G}\frac{1}{e^{(E-E_F)/kT}+1}dE$$
-
-We need to solve $n_e = n_h$
-
-Simplification:
-Fermi level is far from both bands $E_F \gg kT$ and $E_G - E_F \gg kT$
-
-Therefore Fermi-Dirac distribution is approximately similar to Boltzmann distribution.
-
-$$f(E\pm E_F) \approx e^{-(E\pm E_F)/kT}$$
-
-Now we can calculate $n_e$ and $n_h$:
-
-$$n_h \approx \frac{V(2m_h)^{3/2}}{2\pi^2\hbar^3}e^{-E_F/kT} \int_0^\infty\sqrt{E}e^{-E/kT}dE =
-N_V e^{-E_F/kT},$$
-
-with
-
-$$N_V = 2\left(\frac{2\pi m_h kT}{h^2}\right)^{3/2}$$
-
-the number of holes with energy $E<T$ (compare with the rule above).
-
-??? question "how large is $N_V$ at room temperature? (hard question)"
-    If $kT \sim 1\textrm{eV}$ (the typical energy size of a band), then electrons in the whole band may be excited and $N_V \sim 1$. On the other hand, $N_V \sim T^{3/2}$ Therefore $N_V \sim (kT/1 \textrm{eV})^{3/2}\sim 1\%$.
-
-Similarly for electrons:
-
-$$n_e = N_C e^{-(E_G - E_F)/kT},\quad N_C = 2\left(\frac{2\pi m_e kT}{h^2}\right)^{3/2}$$
-
-
-Combining everything together:
-
-$$n_h \approx N_V e^{-E_F/kT} = N_C e^{-(E_G-E_F)/kT} \approx n_e$$
-
-Solving for $E_F$:
-
-$$E_F = \frac{E_G}{2} - \frac{3}{4}kT\log(m_e/m_h)$$
-
-An extra observation: regardless of where $E_F$ is located, $n_e n_h = N_C N_V e^{-E_G/kT} \equiv n_i^2$.
-
-$n_i$ is the **intrinsic carrier concentration**, and for a pristine semiconductor $n_e = n_h = n_i$.
-
-> The equation
-> $$n_e n_h = n_i^2$$
-> is the **law of mass action**. The name is borrowed from chemistry, and describes the equilibrium concentration of two reagents in a reaction $A+B \leftrightarrow AB$. Here electrons and hole constantly split and recombine.
-
 !!! summary "Learning goals"
 
     After this lecture you will be able to:
@@ -347,7 +144,7 @@ Additional information can be obtained using Hall effect. However Hall effect is
 
 ### Light adsorption
 
-See [previous lecture](lecture_6.md#light-adsorption)
+See [previous lecture](.md#light-adsorption)
 
 ## Combining semiconductors: $pn$-junction
 
diff --git a/src/1_einstein_model.md b/src/1_einstein_model.md
new file mode 100644
index 00000000..df03afd1
--- /dev/null
+++ b/src/1_einstein_model.md
@@ -0,0 +1,119 @@
+```python
+from matplotlib import pyplot
+
+import numpy as np
+from scipy.optimize import curve_fit
+from scipy.integrate import quad
+
+from common import draw_classic_axes, configure_plotting
+
+configure_plotting()
+```
+
+# Lecture 1 – Phonons and specific Heat
+_(based on chapter 2 of the book)_  
+Exercises: 2.3, 2.4, 2.5, 2.6, 2.8
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - Explain quantum mechanical effects on the heat capacity of solids (Einstein model)
+    - Compute the expected particle number, energy, and heat capacity of a quantum harmonic oscillator (a single boson)
+    - Write down the total thermal energy of a material
+
+### Einstein model
+Before solid state physics: heat capacity per atom $C=3k_{\rm B}$ (Dulong-Petit). Each atom is (classical) harmonic oscillator in three directions. Experiments showed that this law breaks down at low temperatures, where $C$ reduces to zero ($C\propto T^3$).
+
+This can be explained by considering a _quantum_ harmonic oscillator:
+![](figures/harmonic.svg)
+
+$$\varepsilon_n=\left(n+\frac{1}{2}\right)\hbar\omega$$
+
+Phonons are bosons $\Rightarrow$ they follow Bose-Einstein statistics.
+
+$$
+n(\omega,T)=\frac{1}{ {\rm e}^{\hbar\omega/k_{\rm B}T}-1}\Rightarrow\bar{\varepsilon}=\frac{1}{2}\hbar\omega+\frac{\hbar\omega}{ {\rm e}^{\hbar\omega/k_{\rm B}T}-1}
+$$
+
+```python
+fig, (ax, ax2) = pyplot.subplots(ncols=2, figsize=(10, 5))
+omega = np.linspace(0.1, 2)
+ax.plot(omega, 1/(np.exp(omega) - 1))
+ax.set_ylim(0, top=3)
+ax.set_xlim(left=0)
+ax.set_xlabel(r'$\hbar \omega$')
+ax.set_xticks([0, 1])
+ax.set_xticklabels(['$0$', '$k_B T$'])
+ax.set_ylabel('$n$')
+ax.set_yticks([1, 2])
+ax.set_yticklabels(['$1$', '$2$'])
+draw_classic_axes(ax, xlabeloffset=.2)
+T = np.linspace(0.01, 2)
+ax2.plot(T, 1/2 + 1/(np.exp(1/T)-1))
+ax2.set_ylim(bottom=0)
+ax2.set_xlabel('$k_B T$')
+ax2.set_xticks([0, 1])
+ax2.set_xticklabels(['$0$', r'$\hbar \omega$'])
+ax2.set_ylabel(r"$\bar\varepsilon$")
+ax2.set_yticks([1/2])
+ax2.set_yticklabels([r'$\hbar\omega/2$'])
+draw_classic_axes(ax2, xlabeloffset=.15)
+```
+
+The term $\frac{1}{2}\hbar\omega$ is the _zero point energy_, which follows from the uncertainty principle.
+
+In order to calculate the heat capacity per atom $C$, we need to differentiate $\bar{\varepsilon}$ to $T$.
+
+$$
+\begin{multline}
+C = \frac{\partial\bar{\varepsilon}}{\partial T}
+= -\frac{\hbar\omega}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}\frac{\partial}{\partial T}\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)\\
+= \frac{\hbar^2\omega^2}{k_{\rm B}T^2}\frac{ {\rm e}^{\hbar\omega/k_{\rm B}T}}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}
+=k_{\rm B}\left(\frac{\hbar\omega}{k_{\rm B}T}\right)^2\frac{ {\rm e}^{\hbar\omega/k_{\rm B}T}}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}
+\end{multline}
+$$
+
+```python
+def c_einstein(T, T_E=1):
+    x = T_E / T
+    return 3 * x**2 * np.exp(x) / (np.exp(x) - 1)**2
+
+T = np.linspace(0.01, 1.5, 500)
+fig, ax = pyplot.subplots()
+
+ax.plot(T, c_einstein(T)/3)
+ax.fill_between(T, c_einstein(T)/3, 1, alpha=0.5)
+
+ax.set_ylim(bottom=0, top=1.2)
+ax.set_xlabel('$T$')
+ax.set_ylabel(r'$\omega$')
+ax.set_xticks([1])
+ax.set_xticklabels([r'$\hbar \omega/k_B$'])
+ax.set_yticks([1])
+ax.set_yticklabels(['$k_B$'])
+pyplot.hlines([1], 0, 1.5, linestyles='dashed')
+draw_classic_axes(ax)
+```
+
+The dashed line signifies the classical value, $k_{\rm B}$. Shaded area $=\frac{1}{2}\hbar\omega$, the zero point energy that cannot be removed through cooling.
+
+This is for just one atom. In order to obtain the heat capacity of a full material, we would have to multiply $C$ (or $\bar{\varepsilon}$) by $3N$, _i.e._ the number of harmonic oscillators according to Einstein model.
+
+```python
+# Data from Einstein's paper
+T = [222.4, 262.4, 283.7, 306.4, 331.3, 358.5, 413.0, 479.2, 520.0, 879.7, 1079.7, 1258.0]
+c = [0.384, 0.578, 0.683, 0.798, 0.928, 1.069, 1.343, 1.656, 1.833, 2.671, 2.720, 2.781]
+
+fit = curve_fit(c_einstein, T, c, 1000)
+T_E = fit[0][0]
+delta_T_E = np.sqrt(fit[1][0, 0])
+
+fig, ax = pyplot.subplots()
+ax.scatter(T, c)
+temps = np.linspace(10, T[-1], 100)
+ax.plot(temps, c_einstein(temps, T_E));
+ax.set_xlabel('$T[K]$')
+ax.set_ylabel('$C/k_B$')
+ax.set_ylim((0, 3));
+```
diff --git a/src/lecture_1.md b/src/2_debye_model.md
similarity index 60%
rename from src/lecture_1.md
rename to src/2_debye_model.md
index 9c49766f..abd2b2a0 100644
--- a/src/lecture_1.md
+++ b/src/2_debye_model.md
@@ -10,116 +10,18 @@ from common import draw_classic_axes, configure_plotting
 configure_plotting()
 ```
 
-# Lecture 1 – Phonons and specific Heat
-_(based on chapter 2 of the book)_  
-Exercises: 2.3, 2.4, 2.5, 2.6, 2.8
-
 !!! summary "Learning goals"
 
     After this lecture you will be able to:
 
-    - Explain quantum mechanical effects on the heat capacity of solids (Einstein model)
-    - Compute the expected particle number, energy, and heat capacity of a quantum harmonic oscillator (a single boson)
-    - Write down the total thermal energy of a material
-
-### Einstein model
-Before solid state physics: heat capacity per atom $C=3k_{\rm B}$ (Dulong-Petit). Each atom is (classical) harmonic oscillator in three directions. Experiments showed that this law breaks down at low temperatures, where $C$ reduces to zero ($C\propto T^3$).
-
-This can be explained by considering a _quantum_ harmonic oscillator:
-![](figures/harmonic.svg)
-
-$$\varepsilon_n=\left(n+\frac{1}{2}\right)\hbar\omega$$
-
-Phonons are bosons $\Rightarrow$ they follow Bose-Einstein statistics.
-
-$$
-n(\omega,T)=\frac{1}{ {\rm e}^{\hbar\omega/k_{\rm B}T}-1}\Rightarrow\bar{\varepsilon}=\frac{1}{2}\hbar\omega+\frac{\hbar\omega}{ {\rm e}^{\hbar\omega/k_{\rm B}T}-1}
-$$
-
-```python
-fig, (ax, ax2) = pyplot.subplots(ncols=2, figsize=(10, 5))
-omega = np.linspace(0.1, 2)
-ax.plot(omega, 1/(np.exp(omega) - 1))
-ax.set_ylim(0, top=3)
-ax.set_xlim(left=0)
-ax.set_xlabel(r'$\hbar \omega$')
-ax.set_xticks([0, 1])
-ax.set_xticklabels(['$0$', '$k_B T$'])
-ax.set_ylabel('$n$')
-ax.set_yticks([1, 2])
-ax.set_yticklabels(['$1$', '$2$'])
-draw_classic_axes(ax, xlabeloffset=.2)
-T = np.linspace(0.01, 2)
-ax2.plot(T, 1/2 + 1/(np.exp(1/T)-1))
-ax2.set_ylim(bottom=0)
-ax2.set_xlabel('$k_B T$')
-ax2.set_xticks([0, 1])
-ax2.set_xticklabels(['$0$', r'$\hbar \omega$'])
-ax2.set_ylabel(r"$\bar\varepsilon$")
-ax2.set_yticks([1/2])
-ax2.set_yticklabels([r'$\hbar\omega/2$'])
-draw_classic_axes(ax2, xlabeloffset=.15)
-```
-
-The term $\frac{1}{2}\hbar\omega$ is the _zero point energy_, which follows from the uncertainty principle.
-
-In order to calculate the heat capacity per atom $C$, we need to differentiate $\bar{\varepsilon}$ to $T$.
-
-$$
-\begin{multline}
-C = \frac{\partial\bar{\varepsilon}}{\partial T}
-= -\frac{\hbar\omega}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}\frac{\partial}{\partial T}\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)\\
-= \frac{\hbar^2\omega^2}{k_{\rm B}T^2}\frac{ {\rm e}^{\hbar\omega/k_{\rm B}T}}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}
-=k_{\rm B}\left(\frac{\hbar\omega}{k_{\rm B}T}\right)^2\frac{ {\rm e}^{\hbar\omega/k_{\rm B}T}}{\left({\rm e}^{\hbar\omega/k_{\rm B}T}-1\right)^2}
-\end{multline}
-$$
-
-```python
-def c_einstein(T, T_E=1):
-    x = T_E / T
-    return 3 * x**2 * np.exp(x) / (np.exp(x) - 1)**2
-
-T = np.linspace(0.01, 1.5, 500)
-fig, ax = pyplot.subplots()
-
-ax.plot(T, c_einstein(T)/3)
-ax.fill_between(T, c_einstein(T)/3, 1, alpha=0.5)
-
-ax.set_ylim(bottom=0, top=1.2)
-ax.set_xlabel('$T$')
-ax.set_ylabel(r'$\omega$')
-ax.set_xticks([1])
-ax.set_xticklabels([r'$\hbar \omega/k_B$'])
-ax.set_yticks([1])
-ax.set_yticklabels(['$k_B$'])
-pyplot.hlines([1], 0, 1.5, linestyles='dashed')
-draw_classic_axes(ax)
-```
-
-The dashed line signifies the classical value, $k_{\rm B}$. Shaded area $=\frac{1}{2}\hbar\omega$, the zero point energy that cannot be removed through cooling.
-
-This is for just one atom. In order to obtain the heat capacity of a full material, we would have to multiply $C$ (or $\bar{\varepsilon}$) by $3N$, _i.e._ the number of harmonic oscillators according to Einstein model.
+    - Describe the concept of reciprocal space and allowed momenta
+    - Write down the total energy of phonons given the temperature and the dispersion relation
+    - Estimate heat capacity due to phonons in high temperature and low temperature regimes of the Debye model
 
 ### Debye model
 The Einstein model explained the experimental data quite well, but still slightly underestimated the observed values of $C$. Apparently the "each atom is an oscillator"-idea is too simplistic.
 
-```python
-# Data from Einstein's paper
-T = [222.4, 262.4, 283.7, 306.4, 331.3, 358.5, 413.0, 479.2, 520.0, 879.7, 1079.7, 1258.0]
-c = [0.384, 0.578, 0.683, 0.798, 0.928, 1.069, 1.343, 1.656, 1.833, 2.671, 2.720, 2.781]
 
-fit = curve_fit(c_einstein, T, c, 1000)
-T_E = fit[0][0]
-delta_T_E = np.sqrt(fit[1][0, 0])
-
-fig, ax = pyplot.subplots()
-ax.scatter(T, c)
-temps = np.linspace(10, T[-1], 100)
-ax.plot(temps, c_einstein(temps, T_E));
-ax.set_xlabel('$T[K]$')
-ax.set_ylabel('$C/k_B$')
-ax.set_ylim((0, 3));
-```
 
 Peter Debye (1884 – 1966) suggested to instead consider _normal modes_: sound waves that propagate through a solid with velocity $v=\omega/k$, where $k=2\pi/\lambda$ is the _wave number_. Instead of just multiplying $\bar{\varepsilon}$ by $3N$, we now use:
 
@@ -127,13 +29,6 @@ $$E=\int\limits_0^\infty\left(\frac{1}{2}\hbar\omega+\frac{\hbar\omega}{ {\rm e}
 
 $g(\omega)$ is the _density of states_: the number of normal modes found at each position along the $\omega$-axis. How do we calculate $g(\omega)$?
 
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - Describe the concept of reciprocal space and allowed momenta
-    - Write down the total energy of phonons given the temperature and the dispersion relation
-    - Estimate heat capacity due to phonons in high temperature and low temperature regimes of the Debye model
 
 #### Reciprocal space, periodic boundary conditions
 Each normal mode can be described by a _wave vector_ ${\bf k}$. A wave vector represents a point in _reciprocal space_ or _k-space_. We can find $g(\omega)$ by counting the number of normal modes in k-space and then converting those to $\omega$.
@@ -207,6 +102,10 @@ where $x=\frac{\hbar\omega}{k_{\rm B}T}$ and $\Theta_{\rm D}\equiv\frac{\hbar\om
 def integrand(y):
     return y**4 * np.exp(y) / (np.exp(y) - 1)**2
 
+def c_einstein(T, T_E=1):
+    x = T_E / T
+    return 3 * x**2 * np.exp(x) / (np.exp(x) - 1)**2
+
 @np.vectorize
 def c_debye(T, T_D=1):
     x = T / T_D
diff --git a/src/3_drude_model.md b/src/3_drude_model.md
new file mode 100644
index 00000000..aed19016
--- /dev/null
+++ b/src/3_drude_model.md
@@ -0,0 +1,86 @@
+```python
+from matplotlib import pyplot
+
+import numpy as np
+
+from common import draw_classic_axes, configure_plotting
+
+configure_plotting()
+```
+
+# Drude theory and the free electron model
+_(based on chapters 3–4 of the book)_  
+Exercises 3.1, 3.3, 4.2, 4.5, 4.6, 4.7
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - discuss the basics of 'Drude theory', which describes electron motion in metals.
+    - use Drude theory to analyze transport of electrons through conductors in electric and magnetic fields.
+    - describe central terms such as the mobility and the Hall resistance.
+
+### Drude theory
+Ohm's law states that $V=IR=I\rho\frac{l}{A}$. In this lecture we will investigate where this law comes from. We will use the theory developed by Paul Drude in 1900, which is based on three assumptions:
+
+- Electrons have an average scattering time $\tau$.
+- At each scattering event an electron returns to momentum ${\bf p}=0$.
+- In-between scattering events electrons respond to the Lorentz force ${\bf F}_{\rm L}=-e\left({\bf E}+{\bf v}\times{\bf B}\right)$.
+
+For now we will consider only an electric field (_i.e._ ${\bf B}=0$). What velocity do electrons acquire in-between collisions?
+
+$$
+{\bf v}=-\int_0^\tau\frac{e{\bf E}}{m_{\rm e}}{\rm d}t=-\frac{e\tau}{m_{\rm e}}{\bf E}=-\mu{\bf E}
+$$
+
+Here we have defined the quantity $\mu\equiv e\tau/m_{\rm e}$, which is the _mobility_. If we have a density $n$ of electrons in our solid, the current density ${\bf j}$ [A/m$^2$] then becomes:
+
+$$
+{\bf j}=-en{\bf v}=\frac{n e^2\tau}{m_{\rm e}}{\bf E}=\sigma{\bf E}\ ,\ \ \sigma=\frac{ne^2\tau}{m_{\rm e}}=ne\mu
+$$
+
+$\sigma$ is the conductivity, which is the inverse of resistivity: $\rho=\frac{1}{\sigma}$. If we now take $j=\frac{I}{A}$ and $E=\frac{V}{l}$, we retrieve Ohm's Law: $\frac{I}{A}=\frac{V}{\rho l}$.
+
+Scattering is caused by collisions with:
+
+- Phonons: $\tau_{\rm ph}(T)$ ($\tau_{\rm ph}\rightarrow\infty$ as $T\rightarrow 0$)
+- Impurities/vacancies: $\tau_0$
+
+Scattering rate $\frac{1}{\tau}$:
+
+$$
+\frac{1}{\tau}=\frac{1}{\tau_{\rm ph}(T)}+\frac{1}{\tau_0}\ \Rightarrow\ \rho=\frac{1}{\sigma}=\frac{m}{ne^2}\left( \frac{1}{\tau_{\rm ph}(T)}+\frac{1}{\tau_0} \right)\equiv \rho_{\rm ph}(T)+\rho_0
+$$
+
+![](figures/matthiessen.svg)
+
+_Matthiessen's Rule_ (1864). Solid (dashed) curve: $\rho(T)$ for a pure (impure) crystal.
+
+How fast do electrons travel through a copper wire? Let's take $E$ = 1 volt/m, $\tau$ ~ 25 fs (Cu, $T=$ 300 K).
+
+$\rightarrow v=\mu E=\frac{e\tau}{m_{\rm e}}E=\frac{10^{-19}\times 2.5\times 10^{-14}}{10^{-30}}=2.5\times10^{-3}=2.5$ mm/s ! (= 50 $\mu$m @ 50 Hz AC)
+
+### Hall effect
+Consider a conductive wire in a magnetic field ${\bf B} \rightarrow$ electrons are deflected in a direction perpendicular to ${\bf B}$ and ${\bf j}$.
+
+![](figures/hall_effect.svg)
+
+${\bf E}_{\rm H}$ = _Hall voltage_, caused by the Lorentz force.
+
+In equilibrium, assuming that the average velocity becomes zero after every collision: $\frac{mv_x}{\tau}=-eE$
+
+The $y$-component of the Lorentz force $-e{\bf v}_x\times{\bf B}$ is being compensated by the Hall voltage ${\bf E}_{\rm H}={\bf v}_x\times{\bf B}=\frac{1}{ne}{\bf j}\times{\bf B}$. The total electric field then becomes
+
+$$
+{\bf E}=\left(\frac{1}{ne}{\bf j}\times{\bf B}+\frac{m}{ne^2\tau}{\bf j}\right)
+$$
+
+We now introduce the _resistivity matrix_ $\tilde{\rho}$ as ${\bf E}=\tilde{\rho}{\bf j}$, where the diagonal elements are simply $\rho_{xx}=\rho_{yy}=\rho_{zz}=\frac{m}{ne^2\tau}$. The off-diagonal element $\rho_{xy}$ gives us:
+
+$$
+\rho_{xy}=\frac{B}{ne}\equiv -R_{\rm H}B
+$$
+
+where $R_{\rm H}=-\frac{1}{ne}$ is the _Hall resistance_. So by measuring the Hall resistance, we can obtain $n$, the density of free electrons in a material.
+
+While most materials have $R_{\rm H}>0$, interestingly some materials are found to have $R_{\rm H}<0$. This would imply that the charge carriers either have a positive charge, or a negative mass. We will see later (chapter 17) how to interpret this.
diff --git a/src/lecture_2.md b/src/4_sommerfeld_model.md
similarity index 68%
rename from src/lecture_2.md
rename to src/4_sommerfeld_model.md
index 056e6163..905c15a1 100644
--- a/src/lecture_2.md
+++ b/src/4_sommerfeld_model.md
@@ -8,89 +8,12 @@ from common import draw_classic_axes, configure_plotting
 configure_plotting()
 ```
 
-# Lectures 2A & 2B – Drude theory and the free electron model
-_(based on chapters 3–4 of the book)_  
-Exercises 3.1, 3.3, 4.2, 4.5, 4.6, 4.7
-
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - discuss the basics of 'Drude theory', which describes electron motion in metals.
-    - use Drude theory to analyze transport of electrons through conductors in electric and magnetic fields.
-    - describe central terms such as the mobility and the Hall resistance.
-
-### Drude theory
-Ohm's law states that $V=IR=I\rho\frac{l}{A}$. In this lecture we will investigate where this law comes from. We will use the theory developed by Paul Drude in 1900, which is based on three assumptions:
-
-- Electrons have an average scattering time $\tau$.
-- At each scattering event an electron returns to momentum ${\bf p}=0$.
-- In-between scattering events electrons respond to the Lorentz force ${\bf F}_{\rm L}=-e\left({\bf E}+{\bf v}\times{\bf B}\right)$.
-
-For now we will consider only an electric field (_i.e._ ${\bf B}=0$). What velocity do electrons acquire in-between collisions?
-
-$$
-{\bf v}=-\int_0^\tau\frac{e{\bf E}}{m_{\rm e}}{\rm d}t=-\frac{e\tau}{m_{\rm e}}{\bf E}=-\mu{\bf E}
-$$
-
-Here we have defined the quantity $\mu\equiv e\tau/m_{\rm e}$, which is the _mobility_. If we have a density $n$ of electrons in our solid, the current density ${\bf j}$ [A/m$^2$] then becomes:
-
-$$
-{\bf j}=-en{\bf v}=\frac{n e^2\tau}{m_{\rm e}}{\bf E}=\sigma{\bf E}\ ,\ \ \sigma=\frac{ne^2\tau}{m_{\rm e}}=ne\mu
-$$
-
-$\sigma$ is the conductivity, which is the inverse of resistivity: $\rho=\frac{1}{\sigma}$. If we now take $j=\frac{I}{A}$ and $E=\frac{V}{l}$, we retrieve Ohm's Law: $\frac{I}{A}=\frac{V}{\rho l}$.
-
-Scattering is caused by collisions with:
-
-- Phonons: $\tau_{\rm ph}(T)$ ($\tau_{\rm ph}\rightarrow\infty$ as $T\rightarrow 0$)
-- Impurities/vacancies: $\tau_0$
-
-Scattering rate $\frac{1}{\tau}$:
-
-$$
-\frac{1}{\tau}=\frac{1}{\tau_{\rm ph}(T)}+\frac{1}{\tau_0}\ \Rightarrow\ \rho=\frac{1}{\sigma}=\frac{m}{ne^2}\left( \frac{1}{\tau_{\rm ph}(T)}+\frac{1}{\tau_0} \right)\equiv \rho_{\rm ph}(T)+\rho_0
-$$
-
-![](figures/matthiessen.svg)
-
-_Matthiessen's Rule_ (1864). Solid (dashed) curve: $\rho(T)$ for a pure (impure) crystal.
-
-How fast do electrons travel through a copper wire? Let's take $E$ = 1 volt/m, $\tau$ ~ 25 fs (Cu, $T=$ 300 K).
-
-$\rightarrow v=\mu E=\frac{e\tau}{m_{\rm e}}E=\frac{10^{-19}\times 2.5\times 10^{-14}}{10^{-30}}=2.5\times10^{-3}=2.5$ mm/s ! (= 50 $\mu$m @ 50 Hz AC)
-
-### Hall effect
-Consider a conductive wire in a magnetic field ${\bf B} \rightarrow$ electrons are deflected in a direction perpendicular to ${\bf B}$ and ${\bf j}$.
-
-![](figures/hall_effect.svg)
-
-${\bf E}_{\rm H}$ = _Hall voltage_, caused by the Lorentz force.
-
-In equilibrium, assuming that the average velocity becomes zero after every collision: $\frac{mv_x}{\tau}=-eE$
-
-The $y$-component of the Lorentz force $-e{\bf v}_x\times{\bf B}$ is being compensated by the Hall voltage ${\bf E}_{\rm H}={\bf v}_x\times{\bf B}=\frac{1}{ne}{\bf j}\times{\bf B}$. The total electric field then becomes
-
-$$
-{\bf E}=\left(\frac{1}{ne}{\bf j}\times{\bf B}+\frac{m}{ne^2\tau}{\bf j}\right)
-$$
-
-We now introduce the _resistivity matrix_ $\tilde{\rho}$ as ${\bf E}=\tilde{\rho}{\bf j}$, where the diagonal elements are simply $\rho_{xx}=\rho_{yy}=\rho_{zz}=\frac{m}{ne^2\tau}$. The off-diagonal element $\rho_{xy}$ gives us:
-
-$$
-\rho_{xy}=\frac{B}{ne}\equiv -R_{\rm H}B
-$$
-
-where $R_{\rm H}=-\frac{1}{ne}$ is the _Hall resistance_. So by measuring the Hall resistance, we can obtain $n$, the density of free electrons in a material.
-
-While most materials have $R_{\rm H}>0$, interestingly some materials are found to have $R_{\rm H}<0$. This would imply that the charge carriers either have a positive charge, or a negative mass. We will see later (chapter 17) how to interpret this.
-
 ### Sommerfeld theory (free electron model)
 
 !!! summary "Learning goals"
 
     After this lecture you will be able to:
-    
+
     - calculate the electron density of states in 1D, 2D, and 3D using the Sommerfeld free-electron model.
     - express the number and energy of particles in a system in terms of integrals over k-space.
     - use the Fermi distribution to extend the previous learning goal to finite T.
diff --git a/src/lecture_3.md b/src/5_atoms_and_lcao.md
similarity index 53%
rename from src/lecture_3.md
rename to src/5_atoms_and_lcao.md
index b5a913c5..166e21a2 100644
--- a/src/lecture_3.md
+++ b/src/5_atoms_and_lcao.md
@@ -142,200 +142,3 @@ An occupied $\psi_-$ increases the molecule energy as the atoms move closer ⇒
 Therefore if each atom has a single electron in the outermost shell, these atoms attract, if there are 0 or 2 electrons, the net electron force cancels (but electrostatic repulsion remains).
 
 *[LCAO]: Linear Combination of Atomic Orbitals
-
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - Explain the origins of interatomic forces
-    - Compute vibrational spectra of small molecules in 1D
-    - Formulate Hamiltonians and equations of motion of bulk materials (but not yet solve them)
-
-## Adding repulsion
-
-If bringing two atoms closer would keep increasing energy, any covalent bond would collapse; eventually the two atoms must start repelling (at least when the nuclei get close, but really already much earlier).
-
-![](figures/bonding_with_repulsion.svg)
-
-### Van der Waals bond
-
-*Not our focus, but we will discuss it for completeness*
-
-Originates from attraction between dipole moments of two atoms:
-
-If one atom has a dipole moment $p_1$, it creates electric field
-
-$$ E = \frac{p_1}{4\pi \varepsilon_0 r^3} $$
-
-at the position of another atom. The other atom then develops a dipole moment $p_2 = \chi E$ with $\chi$ the *polarizability* of the atom.
-
-This results in an attractive potential $U(r) = -E p_2 \sim r^{-6}$.
-
-This attraction is much weaker than covalent bonds but drops slower with increasing distance. (Q: how fast does the interaction of covalent bonds drop?)
-
-Van der Waals bonds hold layers of covalently bonded carbon atoms together when forming graphite:
-
-![Graphite atomic layers](https://upload.wikimedia.org/wikipedia/commons/thumb/5/5f/Graphite-layers-side-3D-balls.png/320px-Graphite-layers-side-3D-balls.png)  
-(image source: Wikipedia)
-
-## First steps towards phonons
-
-The atomic interaction is minimized when the distance between the atoms is at the minimum of the potential: $U = U_0$ when $\delta r = a$.
-
-![](figures/interatomic_interaction.svg)
-
-Near the minimum, the potential is approximately parabolic:
-
-$$
-U = U_0 + \frac{\kappa}{2} (\delta r - a)^2 - \frac{\kappa_3}{3!} (\delta r - a)^3 + \ldots
-$$
-
-This is a harmonic oscillator potential with the higher order term $\kappa_3$ providing anharmonicity.
-
-### Rigidity
-
-If we stretch a material with length $L$ containing $N=L/a$ atoms by a lenght $\delta L$, the displacement of each atom is $\delta r - a = \delta L/N$.
-
-The returning force is
-
-$$
-F = - \frac{d U}{d(\delta r)} (\delta r - a) = \kappa a \delta L / L,
-$$
-
-which leads us to the material *compressibility*
-
-$$
-\beta \equiv -\frac{1}{L} \frac{\partial L}{\partial F} = \frac{1}{\kappa a}.
-$$
-
-Then using compressibility we can derive the wave equation and obtain the speed of sound in the material relating it to the bond strength.
-
-### Thermal expansion
-
-The nonlinearity $\kappa_3$ means that the potential grows slower when the atoms are further apart than when they are closer to each other.
-
-Therefore the larger the energy of the atoms, the more their average distance increases $\Rightarrow$ we have a model of thermal expansion.
-
-## Looking ahead: multiple atoms
-
-Our aim is to understand electrons and phonons in solids containing $N\to\infty$ atoms.
-
-We now have a microscopic model of what happens when there are two atoms; let's try to see what happens when we take several.
-
-### Phonons
-
-Our plan:
-
-* Consider only *harmonic potential* acting between atoms
-* Write down equations of motion
-* Compute normal modes
-
-For simplicity we consider 1D motion, and let's start with 3 atoms:
-
-![](figures/phonons2.svg)
-
-Newton's equations of motion:
-
-$$
-\begin{aligned}
-m \ddot{x}_1 &= - \kappa (x_1 - x_2) \\
-m \ddot{x}_2 &= -\kappa (x_2 - x_1) -\kappa (x_2 - x_3) \\
-m \ddot{x}_3 &= - \kappa (x_3 - x_2)
-\end{aligned},
-$$
-
-Or in matrix form:
-
-$$
-m \ddot{x} = -\kappa
-\begin{pmatrix}
-1 & -1 & 0\\
--1 & 2 & -1\\
-0 & -1 & 1
-\end{pmatrix}x
-$$
-
-We search for *normal modes*: patterns of motion that are periodic and have a fixed frequency: $x(t) = x_0 e^{i\omega t}$.
-
-Substituting into the equations of motion we get an eigenvalue problem:
-
-$$
-\omega^2 x_0 = \sqrt{\frac{\kappa}{m}}
-\begin{pmatrix}
-1 & -1 & 0\\
--1 & 2 & -1\\
-0 & -1 & 1
-\end{pmatrix}x_0.
-$$
-
-The solutions of this eigenvalue problem are the phonon modes that we occupy.
-
-### Electrons
-
-Same as a single covalent bond, only more atoms in a line. Considering 3 atoms:
-
-$$
-E \begin{pmatrix}
-c_1 \\ c_2 \\ c_3
-\end{pmatrix} =
-\begin{pmatrix}
-E_0 & t & 0 \\
-t & E_0 & t \\
-0 & t & E_0
-\end{pmatrix}
-\begin{pmatrix}
-c_1 \\ c_2 \\ c_3
-\end{pmatrix}
-$$
-
-### Numerical test
-
-Diagonalizing large matrices is unwieldy, but let's try and check it numerically to see if we notice a trend.
-
-Eigenfrequencies of 3 atoms: `[0.0 1.0 1.732050]`
-
-```python
-def DOS_finite_phonon_chain(n):
-    rhs = 2 * np.eye(n) - np.eye(n, k=1) - np.eye(n, k=-1)
-    rhs[0, 0] -= 1
-    rhs[-1, -1] -= 1
-    pyplot.figure()
-    pyplot.hist(np.sqrt(np.abs(np.linalg.eigvalsh(rhs))), bins=30)
-    pyplot.xlabel("$\omega$")
-    pyplot.ylabel("Number of levels")
-
-DOS_finite_phonon_chain(3)
-```
-
-Energies of 3 orbitals: `[-1.41421356 0.0  1.41421356]`
-
-```python
-def DOS_finite_electron_chain(n):
-    rhs = - np.eye(n, k=1) - np.eye(n, k=-1)
-    pyplot.figure()
-    pyplot.hist(np.linalg.eigvalsh(rhs), bins=30)
-    pyplot.xlabel("$E$")
-    pyplot.ylabel("Number of levels")
-
-DOS_finite_electron_chain(3)
-```
-
-### From 3 atoms to 300
-
-Frequencies of many phonons:
-
-```python
-DOS_finite_phonon_chain(300)
-```
-
-Energies of many electrons:
-
-```python
-DOS_finite_electron_chain(300)
-```
-## Summary
-
-* Electrons in atoms occupy "shells" in the energetically favorable order
-* When two atoms come close, electrons occupy molecular orbitals and bind atoms together.
-* Motion of electrons makes atoms attract
-* Oscillatory motion of atoms and hopping of electrons between atoms give rise to the macroscopic behavior of the materials (next week)
diff --git a/src/6_bonds_and_spectra.md b/src/6_bonds_and_spectra.md
new file mode 100644
index 00000000..304926ca
--- /dev/null
+++ b/src/6_bonds_and_spectra.md
@@ -0,0 +1,208 @@
+```python
+from matplotlib import pyplot
+
+import numpy as np
+
+from common import draw_classic_axes, configure_plotting
+
+configure_plotting()
+
+pi = np.pi
+```
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - Explain the origins of interatomic forces
+    - Compute vibrational spectra of small molecules in 1D
+    - Formulate Hamiltonians and equations of motion of bulk materials (but not yet solve them)
+
+## Adding repulsion
+
+If bringing two atoms closer would keep increasing energy, any covalent bond would collapse; eventually the two atoms must start repelling (at least when the nuclei get close, but really already much earlier).
+
+![](figures/bonding_with_repulsion.svg)
+
+### Van der Waals bond
+
+*Not our focus, but we will discuss it for completeness*
+
+Originates from attraction between dipole moments of two atoms:
+
+If one atom has a dipole moment $p_1$, it creates electric field
+
+$$ E = \frac{p_1}{4\pi \varepsilon_0 r^3} $$
+
+at the position of another atom. The other atom then develops a dipole moment $p_2 = \chi E$ with $\chi$ the *polarizability* of the atom.
+
+This results in an attractive potential $U(r) = -E p_2 \sim r^{-6}$.
+
+This attraction is much weaker than covalent bonds but drops slower with increasing distance. (Q: how fast does the interaction of covalent bonds drop?)
+
+Van der Waals bonds hold layers of covalently bonded carbon atoms together when forming graphite:
+
+![Graphite atomic layers](https://upload.wikimedia.org/wikipedia/commons/thumb/5/5f/Graphite-layers-side-3D-balls.png/320px-Graphite-layers-side-3D-balls.png)  
+(image source: Wikipedia)
+
+## First steps towards phonons
+
+The atomic interaction is minimized when the distance between the atoms is at the minimum of the potential: $U = U_0$ when $\delta r = a$.
+
+![](figures/interatomic_interaction.svg)
+
+Near the minimum, the potential is approximately parabolic:
+
+$$
+U = U_0 + \frac{\kappa}{2} (\delta r - a)^2 - \frac{\kappa_3}{3!} (\delta r - a)^3 + \ldots
+$$
+
+This is a harmonic oscillator potential with the higher order term $\kappa_3$ providing anharmonicity.
+
+### Rigidity
+
+If we stretch a material with length $L$ containing $N=L/a$ atoms by a lenght $\delta L$, the displacement of each atom is $\delta r - a = \delta L/N$.
+
+The returning force is
+
+$$
+F = - \frac{d U}{d(\delta r)} (\delta r - a) = \kappa a \delta L / L,
+$$
+
+which leads us to the material *compressibility*
+
+$$
+\beta \equiv -\frac{1}{L} \frac{\partial L}{\partial F} = \frac{1}{\kappa a}.
+$$
+
+Then using compressibility we can derive the wave equation and obtain the speed of sound in the material relating it to the bond strength.
+
+### Thermal expansion
+
+The nonlinearity $\kappa_3$ means that the potential grows slower when the atoms are further apart than when they are closer to each other.
+
+Therefore the larger the energy of the atoms, the more their average distance increases $\Rightarrow$ we have a model of thermal expansion.
+
+## Looking ahead: multiple atoms
+
+Our aim is to understand electrons and phonons in solids containing $N\to\infty$ atoms.
+
+We now have a microscopic model of what happens when there are two atoms; let's try to see what happens when we take several.
+
+### Phonons
+
+Our plan:
+
+* Consider only *harmonic potential* acting between atoms
+* Write down equations of motion
+* Compute normal modes
+
+For simplicity we consider 1D motion, and let's start with 3 atoms:
+
+![](figures/phonons2.svg)
+
+Newton's equations of motion:
+
+$$
+\begin{aligned}
+m \ddot{x}_1 &= - \kappa (x_1 - x_2) \\
+m \ddot{x}_2 &= -\kappa (x_2 - x_1) -\kappa (x_2 - x_3) \\
+m \ddot{x}_3 &= - \kappa (x_3 - x_2)
+\end{aligned},
+$$
+
+Or in matrix form:
+
+$$
+m \ddot{x} = -\kappa
+\begin{pmatrix}
+1 & -1 & 0\\
+-1 & 2 & -1\\
+0 & -1 & 1
+\end{pmatrix}x
+$$
+
+We search for *normal modes*: patterns of motion that are periodic and have a fixed frequency: $x(t) = x_0 e^{i\omega t}$.
+
+Substituting into the equations of motion we get an eigenvalue problem:
+
+$$
+\omega^2 x_0 = \sqrt{\frac{\kappa}{m}}
+\begin{pmatrix}
+1 & -1 & 0\\
+-1 & 2 & -1\\
+0 & -1 & 1
+\end{pmatrix}x_0.
+$$
+
+The solutions of this eigenvalue problem are the phonon modes that we occupy.
+
+### Electrons
+
+Same as a single covalent bond, only more atoms in a line. Considering 3 atoms:
+
+$$
+E \begin{pmatrix}
+c_1 \\ c_2 \\ c_3
+\end{pmatrix} =
+\begin{pmatrix}
+E_0 & t & 0 \\
+t & E_0 & t \\
+0 & t & E_0
+\end{pmatrix}
+\begin{pmatrix}
+c_1 \\ c_2 \\ c_3
+\end{pmatrix}
+$$
+
+### Numerical test
+
+Diagonalizing large matrices is unwieldy, but let's try and check it numerically to see if we notice a trend.
+
+Eigenfrequencies of 3 atoms: `[0.0 1.0 1.732050]`
+
+```python
+def DOS_finite_phonon_chain(n):
+    rhs = 2 * np.eye(n) - np.eye(n, k=1) - np.eye(n, k=-1)
+    rhs[0, 0] -= 1
+    rhs[-1, -1] -= 1
+    pyplot.figure()
+    pyplot.hist(np.sqrt(np.abs(np.linalg.eigvalsh(rhs))), bins=30)
+    pyplot.xlabel("$\omega$")
+    pyplot.ylabel("Number of levels")
+
+DOS_finite_phonon_chain(3)
+```
+
+Energies of 3 orbitals: `[-1.41421356 0.0  1.41421356]`
+
+```python
+def DOS_finite_electron_chain(n):
+    rhs = - np.eye(n, k=1) - np.eye(n, k=-1)
+    pyplot.figure()
+    pyplot.hist(np.linalg.eigvalsh(rhs), bins=30)
+    pyplot.xlabel("$E$")
+    pyplot.ylabel("Number of levels")
+
+DOS_finite_electron_chain(3)
+```
+
+### From 3 atoms to 300
+
+Frequencies of many phonons:
+
+```python
+DOS_finite_phonon_chain(300)
+```
+
+Energies of many electrons:
+
+```python
+DOS_finite_electron_chain(300)
+```
+## Summary
+
+* Electrons in atoms occupy "shells" in the energetically favorable order
+* When two atoms come close, electrons occupy molecular orbitals and bind atoms together.
+* Motion of electrons makes atoms attract
+* Oscillatory motion of atoms and hopping of electrons between atoms give rise to the macroscopic behavior of the materials (next week)
diff --git a/src/lecture_4.md b/src/7_tight_binding.md
similarity index 50%
rename from src/lecture_4.md
rename to src/7_tight_binding.md
index 9728744c..cd106855 100644
--- a/src/lecture_4.md
+++ b/src/7_tight_binding.md
@@ -250,215 +250,3 @@ $$
 You can get to this result immediately if you remember the derivative of arccosine. Otherwise you need to go the long way: compute $dE/dk$ as a function of $k$, express $k$ through $E$ as we did above, and take the inverse.
 
 Also a sanity check: when the energy is close to the bottom of the band, $E = E_0 + 2t + \delta E$ we get $g(E) \sim E^{-1/2}$, as we would expect in 1D.
-
-## More complex systems
-
-!!! summary "Learning goals"
-
-    After this lecture you will be able to:
-
-    - formulate equations of motion for electrons or phonons in 1D for systems with more than one degree of freedom per unit cell.
-    - solve these equations to arrive at the dispersion relation.
-    - derive the group velocity, effective mass, and density of states.
-    - explain what happens with the band structure when the periodicity of the lattice is increased or reduced.
-    
-### More hoppings in LCAO
-
-Consider electrons in a 1D atomic chain again
-
-![](figures/lattice_potential.svg)
-
-Let's relax the assumption that there is no hopping between atoms that are further than nearest neighbors:
-
-$$ \langle \phi_n | H | \phi_{n+2}\rangle \equiv t' ≠ 0.$$
-
-We now have more terms in the Schrödinger equation:
-
-$$ E c_n = E_0 c_n + t c_{n-1} + t c_{n+1} + t' c_{n-2} + t' c_{n+2},$$
-
-but the same Ansatz as before should work: $ c_n = c_0 \exp(ikna) $.
-
-Indeed, after substituting we get:
-
-$$ E = E_0 + 2t\cos(ka) + 2t'\cos(2ka).$$
-
-*(check that!)*
-
-Therefore including more hopping terms works out of the box using the same Ansatz.
-
-### More degrees of freedom per unit cell:
-
-![](figures/phonons5.svg)
-
-Before we used that all atoms are similar, but this doesn't work anymore. We need to generalize our anzatz. We label all degrees of freedom in each **unit cell** (a repeated element of the system).
-
-For two atoms in a unit cell we have displacements $u_{1,n}$ and $u_{2,n}$. Here the first index is the atom number within the unit cell, second is the unit cell number. The atom with mass $M$ is the atom 1, and the atom with mass $m$ is the atom 2.
-
-Having the degrees of freedom let's write down the equations of motion:
-
-$$
-\begin{aligned}
-M\ddot{u}_{1,n}=\kappa(u_{2,n}-2u_{1,n}+u_{2,n-1})\\
-m\ddot{u}_{2,n}=\kappa(u_{1, n} - 2u_{2,n}+u_{1,n+1}),
-\end{aligned}
-$$
-
-The new Ansatz is conceptually the same as before: all unit cells should be the same up to a phase factor:
-
-$$
-\begin{pmatrix}
-u_{1,n}\\
-u_{2,n}
-\end{pmatrix} =
-e^{i\omega t - ik na}
-\begin{pmatrix}
-u_{1}\\
-u_{2}
-\end{pmatrix}.
-$$
-
-Substituting this ansatz into the equations of motion we end up with an eigenvalue problem:
-$$\omega^2
-\begin{pmatrix}
-M & 0 \\ 0 & m
-\end{pmatrix}
-\begin{pmatrix}
-u_{1} \\ u_{2}
-\end{pmatrix} = \kappa
-\begin{pmatrix}
-2 & -1 - e^{ika} \\ -1-e^{-ika} & 2
-\end{pmatrix}
-\begin{pmatrix}
-u_{1}\\ u_{2}
-\end{pmatrix},$$
-with eigenfrequencies
-$$\omega^2=\frac{\kappa(M+m)}{Mm}\pm \kappa\left\{\left(\frac{M+m}{Mm}\right)^2-\frac{4}{Mm}{\rm sin}^2\left(\frac{1}{2}ka\right)\right\}^{\frac{1}{2}}$$
-
-Looking at the eigenvectors we see that for every $k$ there are now two values of $\omega$: one corresponding to in-phase motion (–) and anti-phase (+).
-
-```python
-def dispersion_2m(k, kappa=1, M=1.4, m=1, acoustic=True):
-    Mm = M*m
-    m_harm = (M + m) / Mm
-    root = kappa * np.sqrt(m_harm**2 - 4*np.sin(k/2)**2 / Mm)
-    if acoustic:
-        root *= -1
-    return np.sqrt(kappa*m_harm + root)
-
-# TODO: Add panels with eigenvectors
-k = np.linspace(-2*pi, 6*pi, 300)
-fig, ax = pyplot.subplots()
-ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
-ax.plot(k, dispersion_2m(k), label='acoustic')
-ax.set_xlabel('$ka$')
-ax.set_ylabel(r'$\omega$')
-pyplot.xticks([-pi, 0, pi], [r'$-\pi$', '$0$', r'$\pi$'])
-pyplot.yticks([], [])
-pyplot.vlines([-pi, pi], 0, 2.2, linestyles='dashed')
-pyplot.legend()
-pyplot.xlim(-1.75*pi, 3.5*pi)
-pyplot.ylim(bottom=0)
-draw_classic_axes(ax, xlabeloffset=.2)
-```
-
-The lower branch is called _acoustic_ because its linear dispersion near $\omega=0$ matches the behavior of the regular sound waves.
-The upper branch is the _optical branch_ because the high phonon frequency allows whem to efficiently emit and adsorb photons.
-
-Important concept: _density of states_ (DOS): $g(\omega)=\rho_{\rm S}\frac{ {\rm d}k}{ {\rm d}\omega}$
-
-The quantity $\frac{ {\rm d}k}{ {\rm d}\omega}$ can be calculated from the dispersion.
-
-When plotted, the DOS may look something like this:
-
-```python
-matplotlib.rcParams['font.size'] = 24
-k = np.linspace(-.25*pi, 1.5*pi, 300)
-k_dos = np.linspace(0, pi, 20)
-fig, (ax, ax2) = pyplot.subplots(ncols=2, sharey=True, figsize=(10, 5))
-ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
-ax.plot(k, dispersion_2m(k), label='acoustic')
-ax.vlines(k_dos, 0, dispersion_2m(k_dos, acoustic=False),
-          colors=(0.5, 0.5, 0.5, 0.5))
-ax.hlines(
-    np.hstack((dispersion_2m(k_dos, acoustic=False), dispersion_2m(k_dos))),
-    np.hstack((k_dos, k_dos)),
-    1.8*pi,
-    colors=(0.5, 0.5, 0.5, 0.5)
-)
-ax.set_xlabel('$ka$')
-ax.set_ylabel(r'$\omega$')
-ax.set_xticks([0, pi])
-ax.set_xticklabels(['$0$', r'$\pi$'])
-ax.set_yticks([], [])
-ax.set_xlim(-pi/4, 2*pi)
-ax.set_ylim((0, dispersion_2m(0, acoustic=False) + .2))
-draw_classic_axes(ax, xlabeloffset=.2)
-
-k = np.linspace(0, pi, 1000)
-omegas = np.hstack((
-    dispersion_2m(k, acoustic=False), dispersion_2m(k)
-))
-ax2.hist(omegas, orientation='horizontal', bins=75)
-ax2.set_xlabel(r'$g(\omega)$')
-ax2.set_ylabel(r'$\omega$')
-
-# Truncate the singularity in the DOS
-max_x = ax2.get_xlim()[1]
-ax2.set_xlim((0, max_x/2))
-draw_classic_axes(ax2, xlabeloffset=.1)
-matplotlib.rcParams['font.size'] = 16
-```
-
-Note that normally $g(\omega)$ is plotted along the vertical axis and $\omega$ along the horizontal axis – the plot above is just to demonstrate the relation between the dispersion and the DOS. The singularities in $g(\omega)$ at the bottom and top of each branch are called _van Hove singularities_.
-
-### Consistency check with 1 atom per cell
-
-But what if we take $M\rightarrow m$? We should reproduce the previous result with only one band!
-
-The reason why we get two bands is because we have $a\rightarrow 2a$ compared to 1 atom per unit cell!
-
-For $M\rightarrow m$, the dispersion diagram should come back to the diagram obtained for the first example (i.e. with identical masses).
-
-To reconsile the two pictures let's plot two unit cells in the reciprocal space. (Definition: each of these unit cells is called a **Brillouin zone**, a concept that will come up later.)
-
-```python
-k = np.linspace(0, 2*pi, 300)
-k_dos = np.linspace(0, pi, 20)
-fig, ax = pyplot.subplots()
-ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
-ax.plot(k, dispersion_2m(k), label='acoustic')
-omega_max = dispersion_2m(0, acoustic=False)
-ax.plot(k, omega_max * np.sin(k/4), label='equal masses')
-ax.set_xlabel('$ka$')
-ax.set_ylabel(r'$\omega$')
-ax.set_xticks([0, pi, 2*pi])
-ax.set_xticklabels(['$0$', r'$\pi/2a$', r'$\pi/a$'])
-ax.set_yticks([], [])
-ax.set_xlim(-pi/8, 2*pi+.4)
-ax.set_ylim((0, dispersion_2m(0, acoustic=False) + .2))
-ax.legend(loc='lower right')
-pyplot.vlines([pi, 2*pi], 0, 2.2, linestyles='dashed')
-draw_classic_axes(ax, xlabeloffset=.2)
-```
-
-Doubling the lattice constant "folds" the band structure on itself, doubling all the bands.
-
-### Total number of states
-
-Which values of $k$ are allowed for a finite 1D crystal of length $L$?
-
-Assume closed boundary conditions $\rightarrow$ standing waves $\rightarrow k=\frac{\pi}{L},\frac{2\pi}{L},\frac{3\pi}{L},\frac{4\pi}{L}...$
-
-$\rightarrow \rho_{\rm S}(k)=$ (no. of $k$-values / unit of $k$) $=\frac{L}{\pi}$, where only the positive part of $k$-space is filled with allowed $k$-values.
-
-For open boundary conditions $\rightarrow$ running waves $\rightarrow k=$ $...\frac{-4\pi}{L},\frac{-2\pi}{L},0,\frac{2\pi}{L},\frac{4\pi}{L}...$
-
-$\rightarrow \rho_{\rm R}(k)=\frac{L}{2\pi}$, which is lower than for the case of closed boundary conditions, however, in this case the entire $k$-space is filled.
-
-## Summary
-
-* By using plane waves in real space as an Ansatz we found all normal modes and eigenvectors
-* Dispersion relation of a system with period $a$ in real space is periodic with period $2\pi/a$ in $k$-space
-* Computing dispersions explains all the problems we listed before (need for cutoff, lack of scattering with every single atom, existence of insulators).
-* Electrons and phonons have a complicated nonlinear relation between momentum and velocity (**group velocity**), effective mass, and density of states
-* In a system with more than one degree of freedom per unit cell we need to consider independent amplitudes for each degree of freedom, and get multiple bands.
diff --git a/src/8_many_atoms.md b/src/8_many_atoms.md
new file mode 100644
index 00000000..41a08772
--- /dev/null
+++ b/src/8_many_atoms.md
@@ -0,0 +1,224 @@
+```python
+import matplotlib
+from matplotlib import pyplot
+
+import numpy as np
+
+from common import draw_classic_axes, configure_plotting
+
+configure_plotting()
+
+pi = np.pi
+```
+
+## More complex systems
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - formulate equations of motion for electrons or phonons in 1D for systems with more than one degree of freedom per unit cell.
+    - solve these equations to arrive at the dispersion relation.
+    - derive the group velocity, effective mass, and density of states.
+    - explain what happens with the band structure when the periodicity of the lattice is increased or reduced.
+
+### More hoppings in LCAO
+
+Consider electrons in a 1D atomic chain again
+
+![](figures/lattice_potential.svg)
+
+Let's relax the assumption that there is no hopping between atoms that are further than nearest neighbors:
+
+$$ \langle \phi_n | H | \phi_{n+2}\rangle \equiv t' ≠ 0.$$
+
+We now have more terms in the Schrödinger equation:
+
+$$ E c_n = E_0 c_n + t c_{n-1} + t c_{n+1} + t' c_{n-2} + t' c_{n+2},$$
+
+but the same Ansatz as before should work: $ c_n = c_0 \exp(ikna) $.
+
+Indeed, after substituting we get:
+
+$$ E = E_0 + 2t\cos(ka) + 2t'\cos(2ka).$$
+
+*(check that!)*
+
+Therefore including more hopping terms works out of the box using the same Ansatz.
+
+### More degrees of freedom per unit cell:
+
+![](figures/phonons5.svg)
+
+Before we used that all atoms are similar, but this doesn't work anymore. We need to generalize our anzatz. We label all degrees of freedom in each **unit cell** (a repeated element of the system).
+
+For two atoms in a unit cell we have displacements $u_{1,n}$ and $u_{2,n}$. Here the first index is the atom number within the unit cell, second is the unit cell number. The atom with mass $M$ is the atom 1, and the atom with mass $m$ is the atom 2.
+
+Having the degrees of freedom let's write down the equations of motion:
+
+$$
+\begin{aligned}
+M\ddot{u}_{1,n}=\kappa(u_{2,n}-2u_{1,n}+u_{2,n-1})\\
+m\ddot{u}_{2,n}=\kappa(u_{1, n} - 2u_{2,n}+u_{1,n+1}),
+\end{aligned}
+$$
+
+The new Ansatz is conceptually the same as before: all unit cells should be the same up to a phase factor:
+
+$$
+\begin{pmatrix}
+u_{1,n}\\
+u_{2,n}
+\end{pmatrix} =
+e^{i\omega t - ik na}
+\begin{pmatrix}
+u_{1}\\
+u_{2}
+\end{pmatrix}.
+$$
+
+Substituting this ansatz into the equations of motion we end up with an eigenvalue problem:
+$$\omega^2
+\begin{pmatrix}
+M & 0 \\ 0 & m
+\end{pmatrix}
+\begin{pmatrix}
+u_{1} \\ u_{2}
+\end{pmatrix} = \kappa
+\begin{pmatrix}
+2 & -1 - e^{ika} \\ -1-e^{-ika} & 2
+\end{pmatrix}
+\begin{pmatrix}
+u_{1}\\ u_{2}
+\end{pmatrix},$$
+with eigenfrequencies
+$$\omega^2=\frac{\kappa(M+m)}{Mm}\pm \kappa\left\{\left(\frac{M+m}{Mm}\right)^2-\frac{4}{Mm}{\rm sin}^2\left(\frac{1}{2}ka\right)\right\}^{\frac{1}{2}}$$
+
+Looking at the eigenvectors we see that for every $k$ there are now two values of $\omega$: one corresponding to in-phase motion (–) and anti-phase (+).
+
+```python
+def dispersion_2m(k, kappa=1, M=1.4, m=1, acoustic=True):
+    Mm = M*m
+    m_harm = (M + m) / Mm
+    root = kappa * np.sqrt(m_harm**2 - 4*np.sin(k/2)**2 / Mm)
+    if acoustic:
+        root *= -1
+    return np.sqrt(kappa*m_harm + root)
+
+# TODO: Add panels with eigenvectors
+k = np.linspace(-2*pi, 6*pi, 300)
+fig, ax = pyplot.subplots()
+ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
+ax.plot(k, dispersion_2m(k), label='acoustic')
+ax.set_xlabel('$ka$')
+ax.set_ylabel(r'$\omega$')
+pyplot.xticks([-pi, 0, pi], [r'$-\pi$', '$0$', r'$\pi$'])
+pyplot.yticks([], [])
+pyplot.vlines([-pi, pi], 0, 2.2, linestyles='dashed')
+pyplot.legend()
+pyplot.xlim(-1.75*pi, 3.5*pi)
+pyplot.ylim(bottom=0)
+draw_classic_axes(ax, xlabeloffset=.2)
+```
+
+The lower branch is called _acoustic_ because its linear dispersion near $\omega=0$ matches the behavior of the regular sound waves.
+The upper branch is the _optical branch_ because the high phonon frequency allows whem to efficiently emit and adsorb photons.
+
+Important concept: _density of states_ (DOS): $g(\omega)=\rho_{\rm S}\frac{ {\rm d}k}{ {\rm d}\omega}$
+
+The quantity $\frac{ {\rm d}k}{ {\rm d}\omega}$ can be calculated from the dispersion.
+
+When plotted, the DOS may look something like this:
+
+```python
+matplotlib.rcParams['font.size'] = 24
+k = np.linspace(-.25*pi, 1.5*pi, 300)
+k_dos = np.linspace(0, pi, 20)
+fig, (ax, ax2) = pyplot.subplots(ncols=2, sharey=True, figsize=(10, 5))
+ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
+ax.plot(k, dispersion_2m(k), label='acoustic')
+ax.vlines(k_dos, 0, dispersion_2m(k_dos, acoustic=False),
+          colors=(0.5, 0.5, 0.5, 0.5))
+ax.hlines(
+    np.hstack((dispersion_2m(k_dos, acoustic=False), dispersion_2m(k_dos))),
+    np.hstack((k_dos, k_dos)),
+    1.8*pi,
+    colors=(0.5, 0.5, 0.5, 0.5)
+)
+ax.set_xlabel('$ka$')
+ax.set_ylabel(r'$\omega$')
+ax.set_xticks([0, pi])
+ax.set_xticklabels(['$0$', r'$\pi$'])
+ax.set_yticks([], [])
+ax.set_xlim(-pi/4, 2*pi)
+ax.set_ylim((0, dispersion_2m(0, acoustic=False) + .2))
+draw_classic_axes(ax, xlabeloffset=.2)
+
+k = np.linspace(0, pi, 1000)
+omegas = np.hstack((
+    dispersion_2m(k, acoustic=False), dispersion_2m(k)
+))
+ax2.hist(omegas, orientation='horizontal', bins=75)
+ax2.set_xlabel(r'$g(\omega)$')
+ax2.set_ylabel(r'$\omega$')
+
+# Truncate the singularity in the DOS
+max_x = ax2.get_xlim()[1]
+ax2.set_xlim((0, max_x/2))
+draw_classic_axes(ax2, xlabeloffset=.1)
+matplotlib.rcParams['font.size'] = 16
+```
+
+Note that normally $g(\omega)$ is plotted along the vertical axis and $\omega$ along the horizontal axis – the plot above is just to demonstrate the relation between the dispersion and the DOS. The singularities in $g(\omega)$ at the bottom and top of each branch are called _van Hove singularities_.
+
+### Consistency check with 1 atom per cell
+
+But what if we take $M\rightarrow m$? We should reproduce the previous result with only one band!
+
+The reason why we get two bands is because we have $a\rightarrow 2a$ compared to 1 atom per unit cell!
+
+For $M\rightarrow m$, the dispersion diagram should come back to the diagram obtained for the first example (i.e. with identical masses).
+
+To reconsile the two pictures let's plot two unit cells in the reciprocal space. (Definition: each of these unit cells is called a **Brillouin zone**, a concept that will come up later.)
+
+```python
+k = np.linspace(0, 2*pi, 300)
+k_dos = np.linspace(0, pi, 20)
+fig, ax = pyplot.subplots()
+ax.plot(k, dispersion_2m(k, acoustic=False), label='optical')
+ax.plot(k, dispersion_2m(k), label='acoustic')
+omega_max = dispersion_2m(0, acoustic=False)
+ax.plot(k, omega_max * np.sin(k/4), label='equal masses')
+ax.set_xlabel('$ka$')
+ax.set_ylabel(r'$\omega$')
+ax.set_xticks([0, pi, 2*pi])
+ax.set_xticklabels(['$0$', r'$\pi/2a$', r'$\pi/a$'])
+ax.set_yticks([], [])
+ax.set_xlim(-pi/8, 2*pi+.4)
+ax.set_ylim((0, dispersion_2m(0, acoustic=False) + .2))
+ax.legend(loc='lower right')
+pyplot.vlines([pi, 2*pi], 0, 2.2, linestyles='dashed')
+draw_classic_axes(ax, xlabeloffset=.2)
+```
+
+Doubling the lattice constant "folds" the band structure on itself, doubling all the bands.
+
+### Total number of states
+
+Which values of $k$ are allowed for a finite 1D crystal of length $L$?
+
+Assume closed boundary conditions $\rightarrow$ standing waves $\rightarrow k=\frac{\pi}{L},\frac{2\pi}{L},\frac{3\pi}{L},\frac{4\pi}{L}...$
+
+$\rightarrow \rho_{\rm S}(k)=$ (no. of $k$-values / unit of $k$) $=\frac{L}{\pi}$, where only the positive part of $k$-space is filled with allowed $k$-values.
+
+For open boundary conditions $\rightarrow$ running waves $\rightarrow k=$ $...\frac{-4\pi}{L},\frac{-2\pi}{L},0,\frac{2\pi}{L},\frac{4\pi}{L}...$
+
+$\rightarrow \rho_{\rm R}(k)=\frac{L}{2\pi}$, which is lower than for the case of closed boundary conditions, however, in this case the entire $k$-space is filled.
+
+## Summary
+
+* By using plane waves in real space as an Ansatz we found all normal modes and eigenvectors
+* Dispersion relation of a system with period $a$ in real space is periodic with period $2\pi/a$ in $k$-space
+* Computing dispersions explains all the problems we listed before (need for cutoff, lack of scattering with every single atom, existence of insulators).
+* Electrons and phonons have a complicated nonlinear relation between momentum and velocity (**group velocity**), effective mass, and density of states
+* In a system with more than one degree of freedom per unit cell we need to consider independent amplitudes for each degree of freedom, and get multiple bands.
diff --git a/src/9_crystal_structure.md b/src/9_crystal_structure.md
new file mode 100644
index 00000000..df55df04
--- /dev/null
+++ b/src/9_crystal_structure.md
@@ -0,0 +1,103 @@
+# Lecture 5 – Crystal structure and diffraction
+
+_based on chapters 12–14, (up to and including 14.2) of the book_  
+Exercises 12.3, 12.4, 13.3, 13.4, 14.2
+
+!!! summary "Learning goals"
+
+    After this lecture you will be able to:
+
+    - Describe any crystal using crystallographic terminology, and interpret this terminology
+    - Compute the volume filling fraction given a crystal structure
+    - Determine the primitive, conventional, and Wigner-Seitz unit cells of a given lattice
+    - Determine the Miller planes of a given lattice
+
+### Crystal classification
+
+- **_Lattice_**
+    + periodic pattern of *lattice points*, which all have an identical view
+    + lattice points are not necessarily the same as atom positions
+    + there can be multiple atoms per lattice point
+    + freedom of translation
+    + multiple lattices with different point densities possible
+- **_Lattice vectors_**
+    + from lattice point to lattice point
+    + $N$ vectors for $N$ dimensions
+    + multiple combinations possible
+    + not all combinations provide full coverage
+- **_Unit cell_**
+    + spanned by lattice vectors
+    + has 4 corners in 2D, 8 corners in 3D
+    + check if copying unit cell along lattice vectors gives full lattice
+- **_Primitive unit cell_**
+    + smallest possible $\rightarrow$ no identical points skipped
+    + not always most practical choice
+- **_Basis_**
+    + only now we care about the contents (i.e. atoms)
+    + gives element and position of atoms
+    + properly count partial atoms $\rightarrow$ choose which belongs to unit cell
+    + positions in terms of lattice vectors, *not* Cartesian coordinates!
+
+### Example: graphite
+
+![](figures/graphite_mod.svg)
+
+1. Choose origin (can be atom, not necessary)
+2. Find other lattice points that are identical
+3. Choose lattice vectors, either primitive (red) or not primitive (blue)
+    - lengths of lattice vectors and angle(s) between them fully define the crystal lattice
+    - for graphite: $|{\bf a}_1|=|{\bf a}_2|$ = 0.246 nm = 2.46 Ã…, $\gamma$ = 60$^{\circ}$
+4. Specify basis
+    - using ${\bf a}_1$ and ${\bf a}_2$: C$(0,0)$, C$\left(\frac{2}{3},\frac{2}{3}\right)$
+    - using ${\bf a}_1$ and ${\bf a}_{2}'$: C$(0,0)$, C$\left(0,\frac{1}{3}\right)$, C$\left(\frac{1}{2},\frac{1}{2}\right)$, C$\left(\frac{1}{2},\frac{5}{6}\right)$
+
+An alternative type of unit cell is the _Wigner-Seitz cell_: the collection of all points that are closer to one specific lattice point than to any other lattice point. You form this cell by taking all the perpendicular bisectrices or lines connecting a lattice point to its neighboring lattice points.
+
+### Stacking of atoms
+What determines what crystal structure a material adopts? $\rightarrow$ To good approximation, atoms are solid incrompressible spheres that attract each other. How will these organize?
+
+We start with the densest possible packing in 2D:
+![](figures/packing.svg)
+Will the second layer go on sites A, B or C?
+
+ABCABC stacking $\rightarrow$ _cubic close packed_, also known as _face centered cubic_ (fcc):
+
+![](figures/stacking.svg)
+
+- One atom on the center of each side-plane: 'a dice that always throws 1'
+- Conventional unit cell $\neq$ primitive unit cell
+- Cyclic ABC $\rightarrow$ all atoms identical $\rightarrow$ 1 atom per primitive unit cell
+- Conventional cell: $8\times\frac{1}{8}+6\times\frac{1}{2}=1+3=4$ atoms
+
+Examples of fcc crystals: Al, Cu, Ni, Pb, Au and Ag.
+
+### Filling factor
+Filling factor = # of atoms per cell $\times$ volume of 1 atom / volume of cell
+
+$=\frac{4\times\frac{4}{3}\pi R^3}{a^3}=\frac{1}{6}\sqrt{2}\pi\approx 0.74$, where we have used that $\sqrt{2}a=4R$.
+
+Compare this to _body centered cubic_ (bcc), which consists of a cube with atoms on the corners and one atom in the center (Fig. 1.12): filling factor = 0.68. Examples of bcc crystals: Fe, Na, W, Nb.
+
+Question: is 74% the largest possible filling factor? $\rightarrow$ Kepler conjecture (1571 – 1630). Positive proof by Hales _et al._ in 2015!
+
+Crystal structures that are related to fcc:
+
+1. ionic crystals (Fig. 1.13), e.g. NaCl
+2. zincblende (Fig. 1.15), e.g. diamond
+
+ABABAB stacking $\rightarrow$ _hexagonally close-packed_ (hcp), e.g. Co, Zn. In this case there is no cubic symmetry (Fig. 1.11).
+
+### Miller planes
+We start with a simple cubic lattice:
+
+![](figures/cubic_mod.svg)
+
+$|{\bf a}_1|=|{\bf a}_2|=|{\bf a}_3|\equiv a$ (_lattice constant_)
+
+The plane designated by Miller indices $(u,v,w)$ intersects lattice vector ${\bf a}_1$ at $\frac{|{\bf a}_1|}{u}$, ${\bf a}_2$ at $\frac{|{\bf a}_2|}{v}$ and ${\bf a}_3$ at $\frac{|{\bf a}_3|}{w}$.
+
+![](figures/miller.svg)
+
+Miller index 0 means that the plane is parallel to that axis (intersection at "$\frac{|{\bf a}_3|}{0}=\infty$"). A bar above a Miller index means intersection at a negative coordinate.
+
+If a crystal is symmetric under $90^\circ$ rotations, then $(100)$, $(010)$ and $(001)$ are physically indistinguishable. This is indicated with $\{100\}$. $[100]$ is a vector. In a cubic crystal, $[100]$ is perpendicular to $(100)$ $\rightarrow$ proof in problem set.
-- 
GitLab