implement numerical doped semiconductor figures

docs/14_doping_and_devices.md

@@ -2,8 +2,7 @@
 from matplotlib import pyplot
 
 import numpy as np
-from scipy.optimize import curve_fit
-from scipy.integrate import quad
+from scipy.optimize import root_scalar
 
 import plotly.graph_objs as go
 
@@ -11,9 +10,11 @@ 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
@@ -92,7 +93,7 @@ That allows us to use the previous results and to conclude that an acceptor crea
 
 ```python
 E = np.linspace(-3, 3, 1000)
-fig, ax = pyplot.subplots()
+fig_dos, ax = pyplot.subplots()
 
 n_F = 1/(np.exp(2*(E - E_F)) + 1)
 g_e = m_e * sqrt_plus(E - E_C)
@@ -201,14 +202,58 @@ $$E_F = kT\ln[N_V/(N_A-N_D)], \textrm{ for } N_A > N_D.$$
 ## Temperature dependence of the carrier density and Fermi level
 
 It is instructive to consider how $E_F$, $n_e$ and $n_h$ depend on carrier concentrations.
-In this case, we consider an n-doped semiconductor, however, the same logic applies to p-doped semiconductors.
+Below we show what happens in an n-doped semiconductor, however, the same logic applies to p-doped semiconductors.
+
+```python
+def n(E_F, T, N_C, N_V, N_D, N_A, E_C, E_V, E_D, E_A, /):
+    """Compute total charge density of a doped semiconductor."""
+    n_e = N_C * np.exp(-(E_C - E_F) / T)
+    n_h = N_V * np.exp((E_V - E_F) / T)
+    n_ionized_D = N_D / (1 + np.exp((E_F - E_D) / T))
+    n_ionized_A = N_A / (1 + np.exp((E_A - E_F) / T))
+    return n_e - n_h - n_ionized_D + n_ionized_A
+
+
+def equilibrium_E_F(*args):
+    """Compute the Fermi energy in equilibrium"""
+    E_C, E_V = args[5], args[6]
+    return root_scalar(n, args, bracket=(E_C, E_V)).root
+
 
-![](figures/E_F_and_carrier_density.svg)
+N_C, N_V, N_D, N_A, E_C, E_V, E_D, E_A = 1, 4, 0.01, 0.003, 1, -1, 0.8, -0.9
+Ts = np.linspace(0.001, .4, 100)
+E_Fs = np.array(
+    [equilibrium_E_F(T, N_C, N_V, N_D, N_A, E_C, E_V, E_D, E_A) for T in Ts]
+)
+
+fig = pyplot.figure(figsize=(10, 10))
+ax1, ax2 = fig.subplots(nrows=2)
+ax1.plot(Ts, E_Fs)
+ax1.hlines([E_C, (E_C + E_V) / 2, E_D], Ts[0], Ts[-1], 'k', linestyles="--")
+ax1.set_ylim(-.1, 1.1)
+ax1.set_ylabel("$E_F$")
+ax1.set_xlabel("$T$")
+ax1.set_yticks([E_C, (E_C + E_V) / 2, E_D])
+ax1.set_yticklabels(["$E_C$", "$(E_C + E_V) / 2$", "$E_D$"])
+draw_classic_axes(ax1, y=-.05, ylabeloffset=0.01)
+
+ax2.plot(Ts, N_V * np.exp((E_V - E_Fs) / Ts), label="$n_h$")
+ax2.plot(Ts, N_C * np.exp(-(E_C - E_Fs) / Ts), label="$n_e$")
+ax2.hlines((N_D - N_A), Ts[0], Ts[-1], 'k', linestyles="--")
+ax2.set_ylim(0, (N_D - N_A) * 3)
+ax2.legend()
+ax2.set_ylabel("$n$")
+ax2.set_xlabel("$T$")
+ax2.set_yticks([(N_D - N_A)])
+ax2.set_yticklabels(["$N_D - N_A$"])
+draw_classic_axes(ax2, ylabeloffset=0.01)
+```
 
-There are several relevant temperature limits:
+<!-- #region tags=[] -->
+As we go from highest to lowest temperature, we observe several regimes:
 
-* **Intrinsic limit** . If the temperature is sufficiently large, then $n_i \gg |N_D-N_A|$ and therefore $n_e = n_h = n_i$. Additionally, if holes are heavier than electrons, then $E_F$ has an upturn in this limit.
-* **Extrinsic limit**. If we decrease the temperature, we decrease the number of intrinsic carriers to the point where most of the charge carriers come form the fully ionized donors. As a result, the number of carriers stays approximately constant in this temperature range.
+* **Intrinsic limit** . If the temperature is sufficiently large, then $n_i \gg |N_D-N_A|$ and therefore $n_e \approx n_h \approx n_i$. Additionally, if holes are heavier than electrons, then $E_F$ grows with temperature in this limit.
+* **Extrinsic limit**. As we decrease the temperature, we decrease the number of intrinsic carriers to the point where most of the charge carriers come form the fully ionized donors. As a result, the number of carriers stays approximately constant in this temperature range.
 * **Freeze-out limit**. Once the temperature is sufficiently low $kT \ll E_G - E_D$, we expect the electrons to "freeze away" from the conduction band to the donor band. The charge carriers still come from the donors, however, not all donors are ionized now.
 * **Zero temperature**. There are no charge carriers in neither conduction nor valance bands. The highest energy electrons are in the donor band and therefore $E_F$ should match the donor band. 
 
@@ -226,8 +271,9 @@ We represent the properties of inhomogeneous materials using the **band diagram*
 The main idea is to plot the dependence of various energies ($E_F$, bottom of conduction band $E_C$, top of the valence band $E_V$) as a function of position.
 
 Let us build up the band diagram step by step:
+<!-- #endregion -->
 
-```python
+```python tags=[]
 
 
 def trans(x, a, b):
@@ -453,9 +499,6 @@ annot = [
     
     *base_annot
         ]
-         
-    
-         
 
 updatemenus=list([ 
     dict(
@@ -616,7 +659,7 @@ See the book for details.
 Density of states in a doped semiconductor:
 
 ```python
-fig
+fig_dos
 ```
 
 Charge balance determines the number of electrons and holes as well as the position of the Fermi level.