-π/d

π/d

k

E(k)

Drag red dot to trace E(k).
   Ctrl-F to clear trace.
Be patient - can be sluggish.

V(x)

x

d

a

b

Kronig-Penney potential

Energy E

This sheet shows the mathematical construction used
to find the energies of the eigenstates of the Hamiltonian 
as a function of wavenumber k, for the Kronig-Penney
periodic potential. 

The energy eigenvalues are determined by the points of intersection
of the purple curve with the red line given by cos(kd), where d is the 
periodicity of the potential. For values of E where the purple curve lies
above or below the dotted lines at +1 and -1, no such intersection is 
possible, so no eigenstates with these energies can exist -- these are
'forbidden bands'.  For any given value of k, there are multiple possible 
energies, representing different bands. The first three are shown.

axoff = 7
axoffx = 30
α = 0.61 (slider from 0.1 to 0.95, step 0.01)
k1a(x) = sqrt(x) * α
V = 19.8 (slider from 1 to 50, step 0.2)
k2b(x) = sqrt(x - V) * (1 - α)
Kb(x) = sqrt(V - x) * (1 - α)
Coeff2(x) = abs((x - V / 2) / sqrt(x * (x - V)))
LHS2(x) = cos(k1a(x)) * cos(k2b(x)) - Coeff2(x) * sin(k1a(x)) * sin(k2b(x))
Coeff(x) = (x - V / 2) / sqrt(x * (V - x))
LHS1(x) = cos(k1a(x)) * cosh(Kb(x)) - Coeff(x) * sin(k1a(x)) * sinh(Kb(x))
LHS = If(0 < x < V, LHS1, LHS2)
bottom = If(x ≥ 0, -1)
K = Intersect(bottom, LHS, (96.21333853128843, -1))
Orig: (x(K) + axoffx, -axoff)
KMax: (x(Orig) + 20, y(Orig))
KMin: (x(Orig) - 20, y(Orig))
kax = Segment(KMin, KMax)
KP = Point(kax)
kd1: 2 * π * (x(KP) - x(Orig)) / Length[kax]
top = If(x ≥ 0, 1)
J = Intersect(LHS, top, (50.58448417704378, 1))
I = Intersect(top, LHS, (44.74889994655909, 1))
H = Intersect(LHS, bottom, (22.631411922605803, -1.0000000002080458))
LHSb2 = If(x(H) ≤ x ≤ x(I), LHS)
LHSb3 = If(x(J) ≤ x ≤ x(K), LHS)
Scalefactor = 0.2
Ymax: (x(Orig), y(Orig) + 0.8 * axoff)
yax = Segment(Orig, Ymax)
EScale: Scalefactor * yax / V
F = Intersect(LHS, top, (6.0536779794700655, 1.0000000000225735))
kScale: kax / (2 * π)
RHS = If(x ≥ 0, cos(kd1))
A = Intersect(RHS, LHS, F)
B = Intersect(RHS, LHSb2, ((x(H) + x(I)) / 2, LHSb2((x(H) + x(I)) / 2)))
G = Intersect(LHS, bottom, (11.44848929600421, -0.9999999999999953))
band1l = Segment((x(Orig) - 20, y(Orig) + EScale * x(F)), (x(Orig) + 20, y(Orig) + EScale * x(F)))
band1t = Segment((x(Orig) - 20, y(Orig) + EScale * x(G)), (x(Orig) + 20, y(Orig) + EScale * x(G)))
Band1, f, a, g, b = Polygon(F, (x(F), y(G)), G, (x(G), y(F)))
Band2, h, c, i, d = Polygon(H, (x(H), y(I)), I, (x(I), y(H)))
band2l = Segment((x(Orig) - 20, y(Orig) + EScale * x(H)), (x(Orig) + 20, y(Orig) + EScale * x(H)))
band2t = Segment((x(Orig) - 20, y(Orig) + EScale * x(I)), (x(Orig) + 20, y(Orig) + EScale * x(I)))
C = Intersect(RHS, LHS, I)
L = Intersect(bottom, LHS, (97.4048316369187, -1))
Band3, j, e, k, l = Polygon(J, (x(J), y(K)), K, (x(K), y(J)))
E = Intersect(RHS, LHSb3, K)
D = Intersect(RHS, LHSb3, ((x(J) + x(K)) / 2, LHSb3((x(J) + x(K)) / 2)))
OriginV: (10, -axoff)
VXmax: (x(KMin) - 20, -axoff)
vxaxis = Segment(OriginV, VXmax)
VYmax: (x(OriginV), y(OriginV) + 0.8 * axoff)
vyaxis = Segment(OriginV, VYmax)
M: (x(OriginV) + Length[vxaxis] / 3.2, y(OriginV))
vd = Segment(OriginV, M)
Alphad: (x(OriginV) + α * Length[vd], y(OriginV))
M2: (x(OriginV) + 2 * Length[vd], y(OriginV))
M3: (x(OriginV) + 3 * Length[vd], y(OriginV))
Alphad2: (x(OriginV) + (1 + α) * Length[vd], y(OriginV))
Alphad3: (x(OriginV) + (2 + α) * Length[vd], y(OriginV))
v1 = Segment(OriginV, Alphad)
v2 = Segment(Alphad, (x(Alphad), y(Alphad) + 0.6 * axoff))
v3 = Segment((x(Alphad), y(Alphad) + 0.6 * axoff), (x(M), y(M) + 0.6 * axoff))
v4 = Segment((x(M), y(M) + 0.6 * axoff), M)
v1b = Segment(M, Alphad2)
v2b = Segment(Alphad2, (x(Alphad2), y(Alphad2) + 0.6 * axoff))
v3b = Segment((x(Alphad2), y(Alphad2) + 0.6 * axoff), (x(M2), y(M2) + 0.6 * axoff))
v4b = Segment((x(M2), y(M2) + 0.6 * axoff), M2)
v1c = Segment(M2, Alphad3)
v2c = Segment(Alphad3, (x(Alphad3), y(Alphad3) + 0.6 * axoff))
v3c = Segment((x(Alphad3), y(Alphad3) + 0.6 * axoff), (x(M3), y(M3) + 0.6 * axoff))
v4c = Segment((x(M3), y(M3) + 0.6 * axoff), M3)
B2: (x(Orig) + kd1 * kScale, y(Orig) + x(B) * EScale)
A2: (x(Orig) + kd1 * kScale, y(Orig) + x(A) * EScale)
C2: (x(Orig) + kd1 * kScale, y(Orig) + x(C) * EScale)
D2: (x(Orig) + kd1 * kScale, y(Orig) + x(D) * EScale)
band3l = Segment((x(Orig) - 20, y(Orig) + EScale * x(J)), (x(Orig) + 20, y(Orig) + EScale * x(J)))
band3t = Segment((x(Orig) - 20, y(Orig) + EScale * x(K)), (x(Orig) + 20, y(Orig) + EScale * x(K)))