import numpy as np
import matplotlib.pyplot as plt
from scipy.constants import c, epsilon_0, mu_0
Z0 = np.sqrt(mu_0 / epsilon_0) # Wellenwiderstand des Vakuums, ≈ 377 Ω
N = 1201 # E-Punkte: 12 m bei Δx = 1 cm
DX = 0.01
DT = DX / c # S = 1, wie in Kapitel 5
X = np.arange(N) * DX
def simuliere(schritte, signal, quellort=100, hart=False, eps_r=None,
sigma=None, raender=("mur", "mur"), merke=(), sonden=()):
"""Das 1D-Wellenlabor: Kern aus Kap. 5 + Quelle, Material, Ränder."""
e, b = np.zeros(N), np.zeros(N - 1)
# np.ones(n): Array voller Einsen — hier „überall Vakuum"
eps = np.ones(N) if eps_r is None else eps_r
sig = np.zeros(N) if sigma is None else sigma
# Verlust-Koeffizienten: ca staucht das alte Feld bei jedem Schritt,
# cb wiegt die Nachbar-Differenz (σ = 0 → ca = 1: Kern aus Kap. 5)
v = sig * DT / (2 * epsilon_0 * eps)
ca = (1 - v) / (1 + v)
cb = (c**2 * DT / DX / eps) / (1 + v)
bilder, pegel = {}, {s: 0.0 for s in sonden}
huellkurve = np.zeros(N)
for n in range(1, schritte + 1):
b -= DT / DX * np.diff(e) # B springt (halber Takt)
alt_links, alt_rechts = e[1], e[-2] # für die Mur-Ränder
# E springt — nur innen; Anatomie der Zeile direkt unter der Zelle
e[1:-1] = ca[1:-1] * e[1:-1] - cb[1:-1] * np.diff(b)
if hart:
e[quellort] = signal(n * DT) # überschreibt!
else:
e[quellort] += signal(n * DT) # addiert
if raender[0] == "mur":
e[0] = alt_links # Welle tritt links aus
if raender[1] == "mur":
e[-1] = alt_rechts # … und rechts
if n in merke:
bilder[n] = e.copy()
for s in sonden:
pegel[s] = max(pegel[s], abs(e[s]))
huellkurve = np.maximum(huellkurve, np.abs(e))
return e, bilder, pegel, huellkurve
def gauss(t, t0=4.2e-9, tau=0.7e-9):
"""Glockenpuls; t0 = 6·tau, damit er bei t = 0 praktisch null ist."""
return np.exp(-((t - t0) / tau)**2)6 Das 1D-Wellenlabor
Wenn es draußen dunkel ist und drinnen hell, zeigt dir die Fensterscheibe dein Spiegelbild — blass, aber deutlich. Tagsüber ist es auch da, nur überstrahlt. Offenbar wirft Glas einen Teil des Lichts zurück, obwohl es doch durchsichtig ist. Wie viel Prozent sind es? Und woher „weiß” das Licht an der Glasoberfläche überhaupt, dass es sich teilen soll?
Am Ende dieses Kapitels hast du die Antwort selbst gemessen — in einem Wellenlabor aus etwa vierzig Zeilen Python. Der Zwei-Zeilen-Kern aus Kapitel 5 konnte bisher nur zusehen, wie ein Anfangspuls davonläuft. Ein Labor braucht mehr: eine Quelle, die auf Befehl sendet; Wände, die Wellen schlucken statt sie zurückzuwerfen; Materialien, die dämpfen oder brechen. Genau diese drei Anbauten bauen wir jetzt ein — und mit jedem führen wir ein Experiment durch.
Lernziele
Nach diesem Kapitel kannst du …
- … eine Quelle einbauen und erklären, warum man ihren Wert addiert statt setzt (weiche vs. harte Quelle),
- … den Mur-Rand erklären (bei \(S = 1\): der Randpunkt übernimmt den alten Wert des Nachbarn) und seine Restreflexion messen,
- … mit Spiegel und Dauerquelle eine stehende Welle erzeugen und den Knotenabstand \(\lambda/2\) nachmessen,
- … Verluste über die Leitfähigkeit \(\sigma\) einbauen (ca/cb-Koeffizienten), das exponentielle Abklingen gegen die Theorie prüfen — und erklären, warum eine \(\sigma\)-Rampe kaum reflektiert (die PML-Idee),
- … Reflexion und Transmission an einer Glasgrenze messen und mit den Fresnel-Formeln vergleichen.
6.1 Das Labor
Hier ist das komplette Labor — der Kern aus Kapitel 5, erweitert um Quelle, Material und Ränder. Bevor du den Code liest, der Bauplan in Worten. Die Bühne ist für alle Experimente dieses Kapitels dieselbe: eine 12 m lange Strecke aus \(N = 1201\) E-Punkten im Abstand \(\Delta x = 1\) cm, gerechnet mit dem magischen Zeitschritt \(\Delta t = \Delta x/c\) aus Kapitel 5. Die Stellschrauben der Funktion simuliere() sind die Zutaten der kommenden Abschnitte: eine Quelle (Signalform, Ort, weich oder hart), wahlweise Material (\(\varepsilon_r\)- und \(\sigma\)-Profile entlang der Strecke) und pro Rand die Wahl zwischen schluckend ("mur") und spiegelnd ("pec", der festgenagelte Rand aus Kapitel 5). Die Messgeräte: Über merke lassen sich Momentaufnahmen des Feldes zu gewählten Zeitschritten fotografieren, über sonden der größte je gesehene Feldbetrag an einzelnen Orten protokollieren, und die huellkurve hält für jeden Ort den größten Betrag über die gesamte Laufzeit fest. Lies die Funktion einmal ganz; jede Zutat wird anschließend einzeln vorgeführt:
Die E-Update-Zeile ist wieder die dichteste der Funktion — sie ist die Kapitel-5-Zeile, nur mit Material-Koeffizienten davor. Auseinandergefaltet passiert dies:
rest = ca[1:-1] * e[1:-1] # was vom alten Feld übrig bleibt
antrieb = cb[1:-1] * np.diff(b) # Schub durch die B-Nachbarn
e[1:-1] = rest - antrieb # nur die inneren N−2 PunkteDie Längen-Buchhaltung ist die aus Kapitel 5: e hat \(N\) Einträge, das dazwischen wohnende b hat \(N-1\), also liefert np.diff(b) (Nachbar-Differenzen, eine pro Lücke) genau \(N-2\) Werte — einen für jeden inneren E-Punkt, elementweise verheiratet mit den passend zugeschnittenen ca[1:-1] und e[1:-1]. Im Vakuum ist ca = 1 und cb = c²Δt/Δx, und die Zeile schrumpft exakt auf den Kapitel-5-Kern zusammen. Und wie dort gilt: Die Randpunkte e[0] und e[-1] lässt die Zeile absichtlich in Ruhe — ihnen fehlt der symmetrische B-Partner, sie werden von den Randzeilen darunter regiert (Mur übernimmt den alten Nachbarwert, PEC lässt die Null stehen).
Drei Dinge sind neu, und jedes bekommt jetzt sein Experiment: die Quellzeile (+= oder =), die Mur-Randzeilen — und die Koeffizienten ca/cb, die Material ins Spiel bringen.
6.2 Wände, die schlucken: der Mur-Rand
Die festen Ränder aus Kapitel 5 (\(E = 0\)) sind Metallwände — jede Welle kommt zurück. Für ein Labor, das den freien Raum nachbilden soll, brauchen wir das Gegenteil: einen Rand, durch den Wellen hinauslaufen, als ginge das Gitter ewig weiter. Bei \(S = 1\) ist die Lösung verblüffend kurz. Eine linkslaufende Welle rückt pro Zeitschritt exakt eine Zelle nach links — der Randpunkt muss also einfach den alten Wert seines Nachbarn übernehmen:
\[ E_0^{\,n+1} \;=\; E_1^{\,n}. \]
Das ist die Randbedingung von Mur (in ihrer 1D-Bestform). Ob sie hält, was sie verspricht, sagt der einfachste Test, den das Labor hergibt: Die weiche Standardquelle bei \(x = 1\) m gibt den Gauß-Puls ab (er teilt sich in zwei Hälften der Höhe 0,5, eine läuft nach links, eine nach rechts), beide Ränder stehen auf "mur", und wir lassen 2200 Schritte verstreichen — mehr als genug Zeit, dass beide Hälften ihre Ränder erreicht haben. Erfolgskriterium: Was danach noch im Gitter steht, wurde reflektiert; bei einem perfekten Rand wäre es exakt null.
# von oben: simuliere(), gauss() (Labor-Zelle)
e, _, _, _ = simuliere(2200, gauss)
print(f"Rest im Gitter, nachdem der Puls beide Ränder passiert hat: "
f"{np.max(np.abs(e)):.1e}")Rest im Gitter, nachdem der Puls beide Ränder passiert hat: 7.1e-16Vom Puls der Höhe 0,5 bleibt nur Maschinenrauschen — der Rand verschluckt die Welle praktisch perfekt. (Dass dieser Trick genau so nur bei \(S = 1\) und nur in 1D funktioniert, gehört ins Kleingedruckte; die allgemeine Fassung kommt in Kapitel 9.)
6.3 Senden, ohne zu stören: weiche und harte Quellen
Bisher starteten unsere Wellen als Anfangsbedingung. Eine Quelle sendet stattdessen mitten im Lauf: An einem festen Punkt wird das Feld im Takt eines Signals bewegt — wie ein Finger, der auf der Wasseroberfläche tippt. Im Code ist das eine Zeile, aber es gibt zwei Versionen, und der Unterschied ist eine Lehrstunde:
- weich:
e[quellort] += signal(t)— das Signal wird zum Feld addiert, - hart:
e[quellort] = signal(t)— das Signal überschreibt das Feld.
Solange nur gesendet wird, sehen beide gleich aus. Der Unterschied zeigt sich erst, wenn eine Welle zur Quelle zurückkommt:
# von oben: simuliere(), gauss(), X (Labor-Zelle)
ERG = {}
for art, hart in (("weich", False), ("hart", True)):
_, bilder, _, _ = simuliere(2600, gauss, quellort=300, hart=hart,
raender=("mur", "pec"), merke=(2500,))
ERG[art] = bilder[2500]
fig, ax = plt.subplots(figsize=(7.0, 3.0))
ax.plot(X, ERG["weich"], label="weiche Quelle")
ax.plot(X, ERG["hart"], label="harte Quelle")
ax.axvline(3.0, color="gray", lw=0.8, ls=":")
ax.set_xlabel("x (m)"); ax.set_ylabel("E (V/m)")
ax.legend(); plt.tight_layout(); plt.show()
print(f"max|E| nach Quell-Durchgang — weich: "
f"{np.max(np.abs(ERG['weich'])):.1e}, hart: "
f"{np.max(np.abs(ERG['hart'])):.2f}")
max|E| nach Quell-Durchgang — weich: 1.2e-33, hart: 1.00Das Standbild zeigt nur das Ende der Geschichte — der eigentliche Unterschied passiert mittendrin, und in dieser HTML-Fassung kannst du ihm zusehen. Beide Läufe laufen unten als Film übereinander, und bis etwa Schritt 1900 erzählen sie dieselbe Geschichte: Der Puls verlässt die Quelle nach beiden Seiten (die weiche Welle ist halb so hoch — das Addieren verteilt das Signal auf beide Laufrichtungen, das Überschreiben erzwingt die volle Höhe am Quellpunkt), links schluckt der Mur-Rand seine Hälfte, rechts wirft die Metallwand ihre zurück. Dann erreicht der Rückläufer die gepunktete Quellort-Linie — jetzt lohnen die Einzelschritt-Knöpfe: Oben gleitet er ungerührt hindurch und verschwindet im linken Mur-Rand, unten prallt er ab, als stünde dort eine zweite Metallwand, und pendelt fortan als Geisterpuls zwischen Quelle und Wand.
Code der Animation (nur in der HTML-Fassung)
# von oben: np, plt, simuliere(), gauss(), X
from matplotlib import animation
from IPython.display import HTML
TAKTE = range(40, 2601, 40)
filme = {}
for art, hart in (("weich", False), ("hart", True)):
_, bilder, _, _ = simuliere(2600, gauss, quellort=300, hart=hart,
raender=("mur", "pec"), merke=TAKTE)
filme[art] = bilder
fig_a, (oben, unten) = plt.subplots(2, 1, figsize=(7.0, 4.6),
sharex=True, sharey=True)
linien = {}
for ax, art, farbe in ((oben, "weich", "C0"), (unten, "hart", "C1")):
linien[art], = ax.plot(X, filme[art][40], color=farbe, lw=1.4)
ax.axvline(3.0, color="gray", lw=0.8, ls=":") # Quellort
ax.axvline(12.0, color="black", lw=2.0) # Metallwand
ax.set_ylim(-1.15, 1.15)
ax.set_ylabel("E (V/m)")
ax.text(0.15, 0.88, f"{art}e Quelle", color=farbe, fontsize=10)
unten.set_xlabel("x (m)")
def zeichne(j):
n = TAKTE[j]
for art in ("weich", "hart"):
linien[art].set_ydata(filme[art][n])
oben.set_title(f"Schritt {n} von 2600", fontsize=10)
return list(linien.values())
anim = animation.FuncAnimation(fig_a, zeichne, frames=len(TAKTE),
interval=80)
plt.close(fig_a)
HTML(anim.to_jshtml(default_mode="loop"))
WarnungNaheliegende Vermutung
Vermutung: „Eine Quelle gibt ihren Wert vor — natürlich überschreibt man das Feld am Quellpunkt.”
Warum sie naheliegt: So steuert man in jedem Programm eine Variable: Man weist ihr den gewünschten Wert zu. „Hart” wirkt wie die ehrlichere, direktere Quelle.
Was stattdessen stimmt: Wer das Feld überschreibt, erklärt den Quellpunkt zur Randbedingung: Dort gilt ab sofort „E ist, was ich sage” — egal was ankommt. Ein zurückkehrender Puls trifft auf einen Punkt, der sich keinen Millimeter bewegen lässt — also auf eine Metallwand, und wird vollständig reflektiert (im Bild: der Geisterpuls der harten Quelle). Die weiche Quelle dagegen addiert nur ihren Beitrag und lässt alles Übrige ungestört passieren — Superposition statt Diktat. Deshalb gilt: senden mit +=.
6.4 Interferenz nach Plan: die stehende Welle
Quelle plus Spiegel ergibt das erste klassische Experiment: Eine Dauerquelle (Sinus statt Puls) sendet gegen die Metallwand; hinlaufende und reflektierte Welle überlagern sich. Aus Kapitel 4 kennst du die Zutaten — was dabei herauskommen muss, sagt eine Zeile Trigonometrie: \(\sin(kx - \omega t) + \text{Rücklauf}\) ergibt ein ortsfestes Muster mal Takt, mit Nullstellen (Knoten) im Abstand \(\lambda/2\) und einem Knoten direkt auf dem Spiegel (dort erzwingt das Metall \(E = 0\)).
Der Aufbau: Der rechte Rand wird auf "pec" gestellt — die Metallwand bei 12 m ist unser Spiegel —, der linke bleibt schluckend. Die Dauerquelle sendet mit \(f = 300\) MHz (\(\lambda = 1\) m) und sanft anschwellender Amplitude; sie sitzt bei \(x = 1{,}13\) m, und die krumme Zahl ist Absicht: Stünde sie exakt auf dem \(\lambda/2\)-Raster der Wand (etwa bei glatt 1 m, elf Wellenlängen vor dem Spiegel), käme die reflektierte Welle am Quellort genau gegenphasig an und würde links der Quelle alles auslöschen — ein hübscher Spezialfall, der aber das Bild verfälschen würde. Gemessen wird im eingeschwungenen Zustand: Nach 2900 Schritten (die Welle ist dann mehrfach hin- und hergelaufen) fotografieren wir das Feld über eine volle Schwingungsperiode — 100 Schnappschüsse — und bilden daraus je Ort den Maximalausschlag, die Einhüllende. Sie zeigt auf einen Blick, wo das Feld schwingt und wo dauerhaft Ruhe herrscht.
# von oben: simuliere(), X; c (Import) (Labor-Zelle)
LAMBDA = 1.0
F = c / LAMBDA # ≈ 300 MHz
def sinus(t):
return np.sin(2*np.pi*F*t) * min(t / 20e-9, 1.0) # sanfte Rampe
# eingeschwungen nach 2 Hin- und Rückwegen; eine Periode = 100 Schritte
_, bilder, _, _ = simuliere(5000, sinus, quellort=113,
raender=("mur", "pec"), merke=range(2901, 3001))
einhuellende = np.max(np.abs(np.array(list(bilder.values()))), axis=0)
fig, ax = plt.subplots(figsize=(7.4, 3.2))
ax.plot(X, einhuellende, lw=1.2)
for k in range(1, 12): # Theorie: Knoten
ax.axvline(12.0 - k*LAMBDA/2, color="gray", lw=0.6, ls="--", alpha=0.7)
ax.axvline(1.13, color="tab:red", lw=0.8, ls=":")
ax.annotate("Quelle", (1.18, 1.6), fontsize=8, color="tab:red")
ax.set_xlabel("x (m)"); ax.set_ylabel("max |E| über eine Periode")
plt.tight_layout(); plt.show()
Die Knoten sitzen exakt auf dem \(\lambda/2\)-Raster — der erste eine halbe Wellenlänge vor der Wand, der „nullte” auf ihr. An den Knoten ist dauerhaft Ruhe, obwohl zwei Wellen voller Stärke hindurchlaufen; an den Bäuchen schaukelt sich das Feld auf die doppelte Amplitude. Links der Quelle dagegen: glatte Einhüllende, keine Knoten — dort laufen alle Wellenanteile in dieselbe Richtung hinaus, und ohne Gegenläufer steht nichts. Dieses Experiment kehrt wieder: als Resonator in Kapitel 15 und als Stehwellenverhältnis auf Leitungen in Kapitel 22.
Die Einhüllende ist nur das Endergebnis — wie das stehende Muster entsteht, zeigt in dieser HTML-Fassung der Film. Die graue Hülle im Hintergrund ist die gemessene Einhüllende von oben, das blaue Feld schwingt darin. Sieh dir drei Phasen an: Zuerst läuft die anschwellende Welle einfach nach rechts (überall gleiche Höhe). Dann kommt die Reflexion von der Wand zurück, und von rechts nach links rastet das Feld Stück für Stück in die Hülle ein. Im eingeschwungenen Zustand schließlich lohnen die Einzelschritt-Knöpfe: Die Knoten rühren sich nie — dauerhafte Ruhe, während beiderseits alles pumpt —, die Bäuche schlagen bis zur doppelten Sendeamplitude aus, und links der Quelle wogt die Welle ohne Knoten einfach hinaus.
Code der Animation (nur in der HTML-Fassung)
# von oben: np, plt, simuliere(), sinus(), X, einhuellende
from matplotlib import animation
from IPython.display import HTML
TAKTE = range(40, 3001, 40)
_, bilder_film, _, _ = simuliere(5000, sinus, quellort=113,
raender=("mur", "pec"), merke=TAKTE)
fig_a, ax_a = plt.subplots(figsize=(7.2, 3.4))
ax_a.fill_between(X, -einhuellende, einhuellende,
color="gray", alpha=0.25, lw=0)
linie, = ax_a.plot(X, bilder_film[40], color="C0", lw=1.2)
ax_a.axvline(1.13, color="tab:red", lw=0.8, ls=":") # Quellort
ax_a.axvline(12.0, color="black", lw=2.0) # Spiegelwand
ax_a.set_ylim(-2.3, 2.3)
ax_a.set_xlabel("x (m)"); ax_a.set_ylabel("E (V/m)")
def zeichne(j):
n = TAKTE[j]
linie.set_ydata(bilder_film[n])
ax_a.set_title(f"Schritt {n} von 5000 — grau: die Einhüllende",
fontsize=10)
return [linie]
anim = animation.FuncAnimation(fig_a, zeichne, frames=len(TAKTE),
interval=80)
plt.close(fig_a)
HTML(anim.to_jshtml(default_mode="loop"))6.5 Verluste: Wellen ausbremsen — und die PML-Idee
Bisher ist unser Material nur Vakuum. Der erste Fremdstoff ist ein verlustbehaftetes Medium: eines mit Leitfähigkeit \(\sigma\), in dem das E-Feld Ströme treibt, deren Reibung Feldenergie in Wärme verwandelt (das Ohmsche \(\vec{J} = \sigma\vec{E}\) wird in Kapitel 18 genauer beleuchtet). In der Update-Regel landet das als zusätzlicher Term, und nach Auflösen stehen die zwei Koeffizienten aus dem Labor:
\[ E_{\text{neu}} = c_a\,E_{\text{alt}} - c_b\,\frac{\Delta B}{\Delta x}, \qquad c_a = \frac{1 - v}{1 + v},\quad v = \frac{\sigma\,\Delta t}{2\,\varepsilon_0\varepsilon_r}. \]
\(c_a < 1\) staucht das Feld bei jedem Schritt ein wenig — Reibung als Faktor. Eine Handrechnung für \(\sigma = 0{,}001\) S/m im Vakuum, mit einer Feinheit, die leicht übersehen wird:
# von oben: DT, Z0; epsilon_0 (Import) (Labor-Zelle)
sigma0 = 0.001
v = sigma0 * DT / (2 * epsilon_0)
ca = (1 - v) / (1 + v)
print(f"v = {v:.2e} → ca = {ca:.5f}")
print(f"ca^100 (würde nur E gedämpft) : {ca**100:.3f}")
print(f"Theorie fürs Abklingen der WELLE: e^(−αL) = "
f"{np.exp(-sigma0*Z0/2*1.0):.3f} (α = σZ₀/2, L = 1 m)")v = 1.88e-03 → ca = 0.99624
ca^100 (würde nur E gedämpft) : 0.686
Theorie fürs Abklingen der WELLE: e^(−αL) = 0.828 (α = σZ₀/2, L = 1 m)Pro Schritt fast nichts — aber Vorsicht: Die Welle klingt nicht mit \(c_a^{100} = 0{,}686\) pro Meter ab, sondern nur mit \(e^{-\alpha L} = 0{,}829\). Der Faktor 2 im Exponenten hat eine schöne physikalische Erklärung: \(c_a\) staucht allein das E-Feld — doch eine Welle trägt die Hälfte ihrer Energie im B-Feld (Kapitel 4), und das lässt die Reibung ungeschoren. Jeder Energieanteil verbringt nur die halbe Zeit „angreifbar” im E-Feld, also wirkt effektiv der halbe Exponent: \(\alpha = \sigma Z_0 / 2\). Ob das stimmt, entscheidet das Experiment — ein 600-MHz-Wellenzug durchquert vier Meter Verlustzone:
# von oben: simuliere(), X, sigma0, Z0 (Labor- und Handrechnungs-Zelle)
def traeger(t, f0=600e6, t0=15e-9, tau=2.5e-9):
return np.sin(2*np.pi*f0*t) * np.exp(-((t - t0)/tau)**2)
ZONE = (X >= 6.0) & (X <= 10.0)
sigma_sprung = np.where(ZONE, sigma0, 0.0)
_, _, pegel, huelle = simuliere(4200, traeger, sigma=sigma_sprung,
sonden=(580, 1020))
fig, ax = plt.subplots(figsize=(7.0, 3.0))
ax.plot(X, huelle, label="Einhüllende (gemessen)")
maske = (X >= 5.8) & (X <= 10.0)
ax.plot(X[maske], pegel[580]*np.exp(-sigma0*Z0/2*(X[maske] - X[580])),
"k--", label="Theorie e^(−αx)")
ax.axvspan(6, 10, color="orange", alpha=0.15)
ax.set_xlabel("x (m)"); ax.set_ylabel("max |E|")
ax.legend(); plt.tight_layout(); plt.show()
print(f"Dämpfung über die Zone: gemessen {pegel[1020]/pegel[580]:.3f},"
f" Theorie {np.exp(-sigma0*Z0/2*4.0):.3f}")
Dämpfung über die Zone: gemessen 0.472, Theorie 0.471Aber sieh genau hin: Am Eintritt der Zone (bei 6 m) hat die gemessene Einhüllende eine kleine Beule — dort wird ein Teil des Wellenzugs reflektiert. Ein Sprung in den Materialeigenschaften wirft immer etwas zurück, auch ein Sprung in \(\sigma\). Was hilft? Dasselbe wie beim Betreten kalten Wassers: langsam hineingehen. Wir ersetzen den Sprung durch eine Rampe gleicher Gesamtdämpfung und messen, was aus der Zone zurückkommt:
# von oben: simuliere(), traeger(), ZONE, sigma0, X
sigma_rampe = np.where(ZONE, sigma0 * 2*(X - 6.0)/4.0, 0.0)
echo = {}
for name, sig in (("σ-Sprung", sigma_sprung), ("σ-Rampe ", sigma_rampe)):
_, bilder, _, _ = simuliere(1400, traeger, sigma=sig, merke=(1280,))
echo[name] = np.max(np.abs(bilder[1280][:450]))
print(f"Echo aus der Zone, {name}: {echo[name]:.1e}")Echo aus der Zone, σ-Sprung: 3.9e-03
Echo aus der Zone, σ-Rampe : 7.8e-05Die Rampe reflektiert fünfzigmal weniger als der Sprung. Genau das ist die Idee der PML (perfectly matched layer), des Standard-Randes aller professionellen FDTD-Programme: ein Verlustmaterial, dessen Dämpfung sanft von null anwächst — die Welle findet keine Kante, an der sie sich spiegeln könnte, läuft hinein und wird aufgezehrt. Unsere \(\sigma\)-Rampe ist die halbe PML; was ihr zur ganzen fehlt (und warum es erst in 2D nötig wird), steht im Kleingedruckten.
In dieser HTML-Fassung kannst du dem Echo bei der Entstehung zusehen — durch eine Lupe: Die y-Achse ist fest auf \(\pm 0{,}008\) gekappt, der Hauptpuls (Amplitude \(\approx 0{,}5\)) ragt also sechzigfach über den Bildrand hinaus. Das ist Absicht, denn uns interessiert nur, was zurückkommt. Oben läuft der Wellenzug gegen den \(\sigma\)-Sprung, unten gegen die Rampe gleicher Gesamtdämpfung. Beobachte den Moment des Zoneneintritts (Einzelschritt!): Oben löst sich an der Kante ein deutliches Echo und läuft nach links davon — unten findet die Welle keine Kante, und es bleibt still. Der Titel misst das bisher größte Echo links der Zone in beiden Läufen mit — er beginnt damit erst ab Schritt 1100, wenn der Hinläufer das Messfenster geräumt hat (vorher stünde dort schlicht der Puls selbst).
Code der Animation (nur in der HTML-Fassung)
# von oben: np, plt, simuliere(), traeger(), sigma_sprung, sigma_rampe, X
from matplotlib import animation
from IPython.display import HTML
TAKTE = range(20, 1401, 20)
filme = {}
for name, sig in (("σ-Sprung", sigma_sprung), ("σ-Rampe", sigma_rampe)):
_, bilder, _, _ = simuliere(1400, traeger, sigma=sig, merke=TAKTE)
filme[name] = bilder
fig_a, achsen = plt.subplots(2, 1, figsize=(7.2, 4.6),
sharex=True, sharey=True)
linien, echo_max = {}, {}
for ax, name, farbe in ((achsen[0], "σ-Sprung", "C0"),
(achsen[1], "σ-Rampe", "C2")):
linien[name], = ax.plot(X, filme[name][20], color=farbe, lw=1.0)
echo_max[name] = 0.0
ax.axvspan(6, 10, color="orange", alpha=0.15) # Verlustzone
ax.set_ylim(-0.008, 0.008) # die Lupe!
ax.set_ylabel("E (V/m)")
ax.text(0.15, 0.0062, name, color=farbe, fontsize=10)
achsen[1].set_xlabel("x (m)")
def zeichne(j):
n = TAKTE[j]
for name in filme:
linien[name].set_ydata(filme[name][n])
if n >= 1100: # Hinläufer hat das Fenster geräumt
echo_max[name] = max(echo_max[name],
np.abs(filme[name][n][:450]).max())
if n >= 1100:
stand = (f"größtes Echo bisher: Sprung "
f"{echo_max['σ-Sprung']:.1e}, "
f"Rampe {echo_max['σ-Rampe']:.1e}")
else:
stand = "Echo-Messung startet bei Schritt 1100"
achsen[0].set_title(f"Schritt {n} — {stand}", fontsize=10)
return list(linien.values())
anim = animation.FuncAnimation(fig_a, zeichne, frames=len(TAKTE),
interval=80)
plt.close(fig_a)
HTML(anim.to_jshtml(default_mode="loop"))6.6 Das Spiegelbild vermessen: Fresnel an der Glasgrenze
Jetzt zur Fensterscheibe. Glas ist verlustfrei, aber dichter als Luft: relative Permittivität \(\varepsilon_r = n^2\) mit Brechungsindex \(n = 1{,}5\) — im Labor genügt dafür das eps_r-Array. Was an der Grenze passieren muss, verraten zwei Stetigkeitsforderungen: Das E-Feld und das B-Feld können an einer (unmagnetischen, ladungsfreien) Grenzfläche nicht springen. Mit den Amplituden \(1\) (einfallend), \(r\) (reflektiert) und \(t\) (durchgelassen) heißt das
\[ \underbrace{1 + r = t}_{\text{E stetig}}, \qquad \underbrace{1 - r = n\,t}_{\text{B stetig, da } B = nE/c \text{ im Glas}} \;\;\Rightarrow\;\; r = \frac{1 - n}{1 + n}, \quad t = \frac{2}{1 + n}. \]
Die Handrechnung für Glas: \(r = (1 - 1{,}5)/(1 + 1{,}5) = -0{,}2\) und \(t = 2/2{,}5 = 0{,}8\). Das Minuszeichen heißt: Der reflektierte Puls kommt invertiert zurück. Und die reflektierte Leistung ist \(r^2 = 0{,}04\).
# von oben: simuliere(), gauss(), X (Labor-Zelle)
N_GLAS = 1.5
eps_glas = np.where(X >= 7.0, N_GLAS**2, 1.0)
_, bilder, _, _ = simuliere(1250, gauss, eps_r=eps_glas, merke=(600, 1124))
einfall = np.max(bilder[600])
r_mess = np.min(bilder[1124][:600]) / einfall
t_mess = np.max(bilder[1124][700:]) / einfall
fig, ax = plt.subplots(figsize=(7.2, 3.2))
ax.plot(X, bilder[600], "gray", label="vorher (einfallend)")
ax.plot(X, bilder[1124], "tab:red", label="nachher")
ax.axvspan(7, 12, color="tab:blue", alpha=0.12)
ax.set_xlabel("x (m)"); ax.set_ylabel("E (V/m)")
ax.legend(); plt.tight_layout(); plt.show()
print(f"gemessen: r = {r_mess:+.4f} (Fresnel: {(1-N_GLAS)/(1+N_GLAS):+.1f})")
print(f" t = {t_mess:+.4f} (Fresnel: {2/(1+N_GLAS):+.1f})")
print(f"reflektierte Leistung r² = {r_mess**2*100:.1f} %")
gemessen: r = -0.2003 (Fresnel: -0.2)
t = +0.7999 (Fresnel: +0.8)
reflektierte Leistung r² = 4.0 %Das Labor liefert \(r = -0{,}2003\) und \(t = 0{,}7999\) — Fresnel auf drei Nachkommastellen, aus nichts als Nachbar-Differenzen. Die 4 % sind dein Spiegelbild in der Fensterscheibe. (Genau genommen hat die Scheibe zwei Oberflächen, die je 4 % zurückwerfen — was das für Echos bedeutet, untersucht Übung 6.4.) Und der durchgelassene Puls zeigt beide vorhergesagten Effekte: volle Form, aber räumlich um den Faktor \(n\) gestaucht, weil er im Glas nur noch \(c/n\) schafft.
Im Standbild sieht man die Stauchung — das Tempo sieht man nur im Film. Die beiden gepunkteten Senkrechten unten sind reine Theorie: Sie wandern mit exakt \(c\) (grau) beziehungsweise \(c/n\) (blau), ohne je auf die Simulation zu schauen — und die Pulsspitzen kleben an ihnen. Drei Momente lohnen den Einzelschritt: Ganz am Anfang teilt sich der Quellpuls, und seine linke Hälfte verlässt die Bühne durch den Mur-Rand (Abschnitt 6.2 in Aktion). Beim Aufprall auf die Grenze existieren für einige Schritte alle drei Pulse gleichzeitig. Und danach sieht man dem Glas-Puls förmlich an, wie er gegen den grauen Takt zurückfällt: ein Drittel langsamer, exakt der Faktor \(n = 1{,}5\).
Code der Animation (nur in der HTML-Fassung)
# von oben: np, plt, simuliere(), gauss(), X, DT, N_GLAS; c (Import)
from matplotlib import animation
from IPython.display import HTML
TAKTE = range(20, 1251, 20)
_, bilder_film, _, _ = simuliere(1250, gauss, eps_r=eps_glas,
merke=TAKTE)
T0G = 4.2e-9 # Signal-Maximum der Quelle (gauss)
T_AUFPRALL = T0G + 6.0 / c # Spitze legt 6 m bis zur Grenze zurück
fig_a, ax_a = plt.subplots(figsize=(7.2, 3.4))
linie, = ax_a.plot(X, bilder_film[20], color="tab:red", lw=1.3)
ax_a.axvspan(7, 12, color="tab:blue", alpha=0.12)
m_vakuum, = ax_a.plot([np.nan] * 2, [-0.2, 0.6], ":", color="gray",
lw=1.2)
m_glas, = ax_a.plot([np.nan] * 2, [-0.2, 0.6], ":", color="tab:blue",
lw=1.2)
ax_a.set_ylim(-0.2, 0.6)
ax_a.set_xlabel("x (m)"); ax_a.set_ylabel("E (V/m)")
def zeichne(j):
n = TAKTE[j]
t = n * DT
linie.set_ydata(bilder_film[n])
if t < T_AUFPRALL: # Anlauf: Spitze bei 1 m + c·t
x_hin = 1.0 + c * (t - T0G)
m_vakuum.set_xdata([x_hin] * 2 if 1.0 <= x_hin else [np.nan] * 2)
m_glas.set_xdata([np.nan] * 2)
else: # geteilte Pulse
x_r = 7.0 - c * (t - T_AUFPRALL) # reflektiert: Tempo c
x_t = 7.0 + c / N_GLAS * (t - T_AUFPRALL) # im Glas: Tempo c/n
m_vakuum.set_xdata([x_r] * 2 if x_r > 0 else [np.nan] * 2)
m_glas.set_xdata([x_t] * 2 if x_t < 12 else [np.nan] * 2)
ax_a.set_title(f"Schritt {n} — gepunktet: Theorie-Tempo c (grau) "
f"und c/n (blau)", fontsize=10)
return [linie, m_vakuum, m_glas]
anim = animation.FuncAnimation(fig_a, zeichne, frames=len(TAKTE),
interval=80)
plt.close(fig_a)
HTML(anim.to_jshtml(default_mode="loop"))Und das komplette Labor steht in dieser HTML-Fassung live in deinem Browser bereit — die folgende Zelle ist editierbar und läuft per Run-Knopf (oder Strg+Enter), ganz ohne Installation. Vorschläge: Stelle Wasser ein (N_MEDIUM = 1.33 — wie viel Prozent Spiegelbild erwartest du vorher?), mach das Medium dichter (N_MEDIUM = 3), oder ersetze den harten Sprung durch eine Rampe (der auskommentierte Block) und sieh nach, was vom Spiegelbild übrig bleibt — die PML-Idee aus diesem Kapitel, zum Selbermessen:
Das Kapitel-Programm programme/kap06/kap06_wellenlabor_1d.py enthält das Labor samt aller vier Experimente mit assert-Schranken (Mur-Rest \(< 10^{-12}\), Dämpfung und Fresnel auf Prozentbruchteile).
Roter Faden
Das Labor verbindet fast alles Bisherige: den Kern aus Kapitel 5, die stehende Welle, die in Kapitel 4 als Formel auftauchte (und hier als Messung), die Energiebuchhaltung aus Kapitel 3 (das Ohmsche \(\sigma\) wird in Kapitel 18 zur Wärmeleistung \(\vec{J}\cdot\vec{E}\)). Nach vorn: Kapitel 7 zeigt, dass unsere Update-Schleife eine Matrix-Multiplikation ist; Kapitel 8 klärt, welcher Zeitschritt erlaubt ist; in Kapitel 12 wird aus der 1D-Glasgrenze die 2D-Brechung; und in Kapitel 22 kehrt Fresnel als Reflexion an Leitungs-Stoßstellen wieder — dort heißt \(n\) dann Impedanzverhältnis, und „Anpassung” wird zur Ingenieurskunst.
Übungen
Ü 6.1 (Verstehen). Die harte Quelle reflektiert den Rückläufer vollständig — wie eine Metallwand. Eine Metallwand erzwingt \(E = 0\); die harte Quelle erzwingt \(E = s(t) \neq 0\). Warum reflektiert sie trotzdem genauso total?
HinweisMusterlösung zu Ü 6.1
Entscheidend ist nicht, welcher Wert erzwungen wird, sondern dass einer erzwungen wird: Der Punkt reagiert nicht auf das ankommende Feld. Dank Superposition kann man das Gesamtfeld als „vorgeschriebener Anteil + Abweichung” lesen — und für die Abweichung gilt am Quellpunkt exakt die Bedingung \(= 0\), also die Metallwand-Bedingung. Der Rückläufer ist eine solche Abweichung und wird daher wie an einer Wand (invertiert und vollständig) reflektiert.
Ü 6.2 (Verstehen). Berechne \(r\) und \(t\) für den umgekehrten Weg — einen Puls, der aus dem Glas in die Luft läuft (\(n_1 = 1{,}5 \to n_2 = 1\); allgemein gilt \(r = (n_1 - n_2)/(n_1 + n_2)\), \(t = 2n_1/(n_1+n_2)\)). Was fällt am Vorzeichen von \(r\) auf — und wieso ist \(t > 1\) kein Verstoß gegen die Energieerhaltung?
HinweisMusterlösung zu Ü 6.2
\(r = (1{,}5 - 1)/(2{,}5) = +0{,}2\): gleicher Betrag, aber nicht invertiert — die Inversion gibt es nur beim Übergang ins dichtere Medium. Und \(t = 3/2{,}5 = 1{,}2 > 1\): Die Feldamplitude hinter der Grenze ist größer als davor! Kein Widerspruch: Energie fließt mit \(S \propto E^2/Z\), und in der dünneren Luft gehört zur selben Leistung ein größeres \(E\). Die Bilanz stimmt exakt: \(r^2 + (1/n)\,t^2 = 0{,}04 + 0{,}96 = 1\). Leistung, nicht Amplitude, ist die Erhaltungsgröße.
Ü 6.3 (Verändern). Miss im Labor die Reflexion an Wasser (\(n = 1{,}33\)). Sage den Wert vorher mit der Fresnel-Formel voraus.
HinweisMusterlösung zu Ü 6.3
Vorhersage: \(r = (1 - 1{,}33)/(2{,}33) = -0{,}1416\), also \(r^2 \approx 2\,\%\). Die Messung — nur eine Zeile anders:
# von oben: simuliere(), gauss(), X (Labor-Zelle)
eps_wasser = np.where(X >= 7.0, 1.33**2, 1.0)
_, bilder, _, _ = simuliere(1250, gauss, eps_r=eps_wasser,
merke=(600, 1124))
r_w = np.min(bilder[1124][:600]) / np.max(bilder[600])
print(f"r gemessen: {r_w:+.4f} (Fresnel: {(1-1.33)/(1+1.33):+.4f})")r gemessen: -0.1418 (Fresnel: -0.1416)Wasser spiegelt also nur halb so stark wie Glas — deshalb ist das Spiegelbild in einem stillen, dunklen Teich blasser als im Fenster.
Ü 6.4 (Übertragen). Eine echte Fensterscheibe hat zwei Oberflächen. Baue eine 50 cm dicke Glasschicht (von 7 m bis 7,5 m) und sieh dir an, was zurückkommt: Es ist nicht ein Echo, sondern eine Echofolge. Miss die Amplituden der ersten beiden Echos und erkläre sie mit den Ergebnissen von Ü 6.2.
HinweisMusterlösung zu Ü 6.4
# von oben: simuliere(), gauss(), X, N_GLAS
eps_scheibe = np.where((X >= 7.0) & (X < 7.5), N_GLAS**2, 1.0)
_, bilder, _, _ = simuliere(1350, gauss, eps_r=eps_scheibe,
merke=(600, 1300))
spur = bilder[1300][:650] / np.max(bilder[600])
fig, ax = plt.subplots(figsize=(7.0, 2.6))
ax.plot(X[:650], spur)
ax.set_xlabel("x (m)"); ax.set_ylabel("E / E_einfall")
ax.set_title("Was die Scheibe zurückschickt")
plt.tight_layout(); plt.show()
print(f"Echo 1: {np.min(spur):+.4f} Echo 2: {np.max(spur):+.4f}")
Echo 1: -0.2003 Echo 2: +0.1923Echo 1 (\(-0{,}200\)) ist die bekannte Reflexion an der Vorderseite — invertiert. Echo 2 (\(+0{,}192\)) erzählt seine Reiseroute: hinein (\(t = 0{,}8\)), an der Rückseite reflektiert (\(r' = +0{,}2\), nicht invertiert — Ü 6.2!), wieder hinaus (\(t' = 1{,}2\), ebenfalls Ü 6.2): \(0{,}8 \cdot 0{,}2 \cdot 1{,}2 = +0{,}192\) — exakt der Messwert, inklusive des „unmöglichen” Faktors \(t' > 1\). Jeder weitere Umlauf in der Scheibe multipliziert mit \(r'^2 = 0{,}04\): Echo 3 wäre \(+0{,}0077\), kaum noch sichtbar. Die Echos werden schnell leise, aber sie hören nie ganz auf — dieses unendliche Echospiel heißt Fabry-Pérot und bekommt in Kapitel 27 sein eigenes Kapitel, denn richtig eingesetzt macht es aus zwei Glasflächen einen Spiegel, ein Filter oder einen Laser-Resonator.
Das Kleingedruckte
- Mur allgemein: Für \(S < 1\) lautet der Rand \(E_0^{n+1} = E_1^n + \frac{c\Delta t - \Delta x}{c\Delta t + \Delta x} (E_1^{n+1} - E_0^n)\) — bei \(S = 1\) kollabiert er auf unsere Kurzform. In 2D/3D bleibt vom Mur-Rand nur eine Näherung, die schräg einfallende Wellen teilweise reflektiert; deshalb dominiert dort die PML.
- Die ganze PML dämpft nicht nur über \(\sigma\), sondern zugleich über einen magnetischen Verlustterm, abgestimmt auf \(\sigma_m/\mu_0 = \sigma/\varepsilon_0\) — dann ist die Schicht für jede Frequenz und (in 2D/3D) jeden Einfallswinkel reflexionsfrei „angepasst”. Unsere Rampe zeigt das Prinzip; die volle Konstruktion kommt in Kapitel 9.
- Amplituden in verschiedenen Medien vergleicht man nur über die Leistung \(\propto E^2/Z\) (siehe Ü 6.2: \(t = 1{,}2\) und trotzdem Energieerhaltung). Systematisch wird diese Buchführung mit dem Impedanzbegriff in Kapitel 22.
- Harte Quellen haben ihren Platz: Wo man ein exakt vorgeschriebenes Feldprofil einspeisen will (z. B. als Rand eines Wellenleiters), ist Überschreiben genau richtig — man muss nur wissen, dass man damit eine reflektierende Randbedingung gesetzt hat.