2 Ableitungen von Feldern

Stell dir vor, du stehst im dichten Nebel irgendwo an einem Berghang. Du willst zum Gipfel, siehst aber keine zehn Meter weit. Was du hast: einen Höhenmesser und deine Beine. Also machst du einen kleinen Schritt nach Osten und liest ab: 30 cm gewonnen. Zurück, ein Schritt nach Norden: 10 cm verloren. Jetzt weißt du alles: Bergauf geht es nach Osten, leicht südlich — und zwar ohne je die ganze Landschaft gesehen zu haben.

Genau das ist eine Ableitung: Werte an benachbarten Orten vergleichen. Die Mathematik macht die Schritte unendlich klein — der Computer kann das nicht, und das ist die beste Nachricht dieses Buches: Für uns bleibt die Ableitung das, was sie im Nebel war — eine Differenz zweier Messwerte, geteilt durch die Schrittweite. Aus diesem einen Handgriff baut dieses Kapitel alle drei „Ableitungen für Felder”: Gradient, Divergenz und Rotation. Und in Kapitel 3 wirst du sehen: Die Maxwell-Gleichungen sind nichts als vier Sätze über genau diese drei Größen.

Lernziele

Nach diesem Kapitel kannst du …

… Ableitungen numerisch als Nachbar-Differenzen berechnen und am log-log-Plot ablesen, wie schnell der Fehler mit der Schrittweite schrumpft,
… den Gradienten eines Skalarfelds berechnen und deuten: Er zeigt bergauf und steht senkrecht auf den Höhenlinien,
… die Divergenz als Quellstärke deuten („fließt netto etwas aus dem Kästchen heraus?“) und berechnen,
… die Rotation als Wirbelstärke deuten (Paddelrad-Test) und berechnen — auch dort, wo man keinen Wirbel sieht,
… die Sätze von Gauß und Stokes in je einem Satz erklären und numerisch verifizieren.

2.1 Die Ableitung als Differenz

Fangen wir eindimensional an. Die Steigung einer Funktion \(f\) am Punkt \(x\) schätzt man am fairsten symmetrisch: einen halben Schritt zurückblicken, einen halben nach vorn —

\[ f'(x) \;\approx\; \frac{f(x + \tfrac{\Delta x}{2}) - f(x - \tfrac{\Delta x}{2})}{\Delta x}. \]

Das heißt zentrale Differenz. Wie gut ist sie? Wir testen sie an \(f(x) = \sin(x)\), denn da kennen wir die Antwort exakt: \(f'(x) = \cos(x)\).

import numpy as np
import matplotlib.pyplot as plt

x0 = 1.0
# np.logspace(-3, 0, 20): 20 logarithmisch verteilte Werte von 10⁻³
# bis 10⁰ — auf der log-Achse gleichmäßig dicht
schrittweiten = np.logspace(-3, 0, 20)
fehler = []
for dx in schrittweiten:
    geschaetzt = (np.sin(x0 + dx/2) - np.sin(x0 - dx/2)) / dx
    fehler.append(abs(geschaetzt - np.cos(x0)))  # Vergleich mit exakt

fig, ax = plt.subplots(figsize=(5.4, 3.6))
ax.loglog(schrittweiten, fehler, "o-", label="zentrale Differenz")
ax.loglog(schrittweiten, 0.04 * schrittweiten**2, "k--",
          label="Erwartung ~ Δx²")
ax.set_xlabel("Schrittweite Δx")
ax.set_ylabel("Fehler |geschätzt − exakt|")
ax.legend()
plt.tight_layout(); plt.show()

Abbildung 2.1: Der Fehler der zentralen Differenz an f = sin(x), ausgewertet bei x = 1. Gestrichelt die Erwartung ~Δx²: Halbe Schrittweite ergibt ein Viertel des Fehlers.

Die Messpunkte folgen exakt der gestrichelten Linie: Der Fehler fällt quadratisch — Schrittweite halbieren, Fehler viertelt sich. Das ist der Grund, warum wir im ganzen Buch zentrale Differenzen benutzen (die einseitige Variante \(\bigl(f(x+\Delta x) - f(x)\bigr)/\Delta x\) schafft nur einen Faktor 2). Merke dir dieses Bild: Es ist der erste Konvergenztest des Buches, und solche Tests sind unser Wahrheitskriterium — nicht „sieht plausibel aus”, sondern „schrumpft mit der erwarteten Ordnung”.

Für Felder brauchen wir Ableitungen nach einer Richtung: Das Symbol \(\partial f/\partial x\) („partielle Ableitung”) heißt nur: Steigung in x-Richtung, während \(y\) festgehalten wird. Numerisch ist das dieselbe Nachbar-Differenz wie eben — angewendet entlang einer Achse des Arrays. NumPy bringt sie fertig mit: np.gradient(f, dx, axis=...), und im Inneren des Arrays rechnet diese Funktion genau den zentralen Differenzenquotienten von eben:

\[ \frac{\partial f}{\partial x}\bigg|_{\text{Punkt } i} \;\approx\; \frac{f_{i+1} - f_{i-1}}{2\,\Delta x} \qquad \text{(rechter Nachbar minus linker Nachbar,} \\ \text{geteilt durch den doppelten Gitterabstand).} \]

Wann immer in diesem Buch np.gradient auftaucht, steht also nichts Geheimnisvolles dahinter — nur diese eine Zeile Nachbarschafts-Bilanz.

2.2 Der Gradient: bergauf zeigen

Packt man die beiden Steigungen eines Skalarfelds in einen Pfeil, entsteht der Gradient:

\[ \nabla T \;=\; \left( \frac{\partial T}{\partial x},\; \frac{\partial T}{\partial y} \right). \]

Das Dreieck \(\nabla\) („Nabla”) liest du am besten als Pfeil aus Ableitungen. Der Gradient eines Skalarfelds ist also ein Vektorfeld — an jedem Ort der Pfeil, der die Nebel-Frage vom Kapitelanfang beantwortet: In welche Richtung geht es bergauf, und wie steil?

Vorhersage-Punkt

Bevor du weiterliest: Wir zeichnen gleich die Gradientpfeile in die Temperaturlandschaft aus Kapitel 1 (warme Insel links, kalte Senke rechts). In welchem Winkel kreuzen die Pfeile die Konturlinien — und zeigen sie zur warmen Insel hin oder von ihr weg? Lege dich fest.

x, y = np.meshgrid(np.linspace(-2, 2, 200), np.linspace(-2, 2, 200))
dx = x[0, 1] - x[0, 0]

# dieselbe Landschaft wie in Kapitel 1
T = 20 + 6*np.exp(-((x+1)**2 + y**2)) - 4*np.exp(-((x-1)**2 + (y-0.5)**2))

dT_dx = np.gradient(T, dx, axis=1)     # Steigung in x-Richtung
dT_dy = np.gradient(T, dx, axis=0)     # Steigung in y-Richtung

fig, ax = plt.subplots(figsize=(5.6, 4.4))
linien = ax.contour(x, y, T, levels=14, cmap="coolwarm")
ax.clabel(linien, fontsize=7, fmt="%.0f")
# slice(None, None, 14) ist der Indexausdruck ::14 als Objekt —
# x[s] greift also jeden 14. Gitterpunkt in beiden Richtungen heraus
s = (slice(None, None, 14), slice(None, None, 14))
ax.quiver(x[s], y[s], dT_dx[s], dT_dy[s], color="k", scale=60)
ax.set_xlabel("x"); ax.set_ylabel("y")
ax.set_title("∇T: senkrecht auf den Höhenlinien, bergauf")
ax.set_aspect("equal")
plt.tight_layout(); plt.show()

Abbildung 2.2: Der Gradient der Temperaturlandschaft aus Kapitel 1: Jeder Pfeil zeigt bergauf — zur warmen Insel hin, aus der kalten Senke heraus — und kreuzt die Höhenlinien im rechten Winkel.

Beides bestätigt: Die Pfeile stehen senkrecht auf den Konturlinien und zeigen zur warmen Insel hin — bergauf. Das Senkrecht-Stehen ist kein Zufall, sondern Logik: Entlang einer Höhenlinie ändert sich \(T\) per Definition nicht, also kann die Änderung nur quer dazu liegen. Und die Pfeile sind dort am längsten, wo die Linien am engsten stehen — am steilen Hang.

2.3 Die Divergenz: Quellen aufspüren

Der Gradient fragt ein Skalarfeld ab. Die nächsten beiden Werkzeuge fragen ein Vektorfeld ab — und zwar mit den zwei Produkten aus Kapitel 1, angewendet auf den Nabla-Pfeil. Das Skalarprodukt von \(\nabla\) mit dem Feld gibt eine Zahl an jedem Ort:

\[ \operatorname{div} \vec{F} \;=\; \nabla \cdot \vec{F} \;=\; \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y}. \]

Was misst diese Zahl? Denk dir ein winziges Kästchen um den Ort gelegt und mach Bilanz: Fließt netto mehr hinaus als hinein? Wenn ja, muss im Kästchen etwas entspringen — eine Quelle (Divergenz positiv). Fließt netto mehr hinein, verschluckt das Kästchen etwas — eine Senke (negativ). Geht genau so viel hinaus wie hinein, ist das Feld dort quellenfrei (null) — egal wie wild es aussieht.

Halte einen Punkt unbedingt fest, bevor wir weitergehen: Hinein geht ein Vektorfeld — heraus kommt an jedem Ort eine einzige Zahl. Das passt zum Skalarprodukt aus Kapitel 1 (zwei Pfeile rein, Zahl raus): \(\nabla \cdot \vec{F}\) ist „Ableitungspfeil mal Feldpfeil”. Das Ergebnis \(\operatorname{div}\vec{F}\) ist deshalb selbst ein Skalarfeld — eine Zahl an jedem Ort, wie die Temperaturkarte aus Kapitel 1. Man kann es als Farbbild zeichnen, aber nicht mit Pfeilen.

2.3.1 Eine Divergenz komplett von Hand

Damit das nicht abstrakt bleibt, rechnen wir einen einzigen Punkt vollständig durch — mit nichts als den vier Nachbarwerten und dem zentralen Differenzenquotienten. Feld: die Quelle \(\vec{F} = (u, v) = (x, y)\). Punkt: \((1{,}0,\; 0{,}5)\). Schrittweite: \(h = 0{,}1\).

Für \(\partial u/\partial x\) brauchen wir \(u\) beim rechten und linken Nachbarn, für \(\partial v/\partial y\) brauchen wir \(v\) beim oberen und unteren Nachbarn — vier Zahlen insgesamt:

def u(x, y): return x          # x-Komponente des Quellenfelds
def v(x, y): return y          # y-Komponente

px, py, h = 1.0, 0.5, 0.1      # Auswertepunkt und Schrittweite

print(f"u rechts (bei x={px+h}) = {u(px+h, py):.2f}   "
      f"u links (bei x={px-h}) = {u(px-h, py):.2f}")
print(f"v oben   (bei y={py+h}) = {v(px, py+h):.2f}   "
      f"v unten (bei y={py-h}) = {v(px, py-h):.2f}\n")

du_dx = (u(px+h, py) - u(px-h, py)) / (2*h)
dv_dy = (v(px, py+h) - v(px, py-h)) / (2*h)
print(f"∂u/∂x ≈ (1.10 − 0.90) / 0.2 = {du_dx:.2f}")
print(f"∂v/∂y ≈ (0.60 − 0.40) / 0.2 = {dv_dy:.2f}")
print(f"div F am Punkt ({px}, {py}) = {du_dx:.2f} + {dv_dy:.2f} "
      f"= {du_dx + dv_dy:.2f}")

u rechts (bei x=1.1) = 1.10   u links (bei x=0.9) = 0.90
v oben   (bei y=0.6) = 0.60   v unten (bei y=0.4) = 0.40

∂u/∂x ≈ (1.10 − 0.90) / 0.2 = 1.00
∂v/∂y ≈ (0.60 − 0.40) / 0.2 = 1.00
div F am Punkt (1.0, 0.5) = 1.00 + 1.00 = 2.00

Vier Nachbarwerte angeschaut, zwei Differenzenquotienten gebildet, addiert — eine Zahl: 2. Mehr ist eine Divergenz nicht. Wiederhole die Rechnung im Kopf an irgendeinem anderen Punkt: Es kommt wieder 2 heraus (die Differenzen hängen hier nicht vom Ort ab) — das Skalarfeld \(\operatorname{div}\vec{F}\) ist für dieses Feld also überall konstant gleich 2. Und genau diese Rechnung führt np.gradient nachher für alle 40 000 Gitterpunkte gleichzeitig aus, nichts weiter.

2.4 Die Rotation: Wirbel aufspüren

Das Kreuzprodukt mit dem Nabla-Pfeil gibt die zweite Frage. In der Ebene bleibt davon eine Zahl übrig (die z-Komponente — der Pfeil, der aus dem Papier ragt):

\[ \operatorname{rot} \vec{F}\big|_z \;=\; (\nabla \times \vec{F})_z \;=\; \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y}. \]

Die Anschauung: Halte ein winziges Paddelrad in die Strömung. Dreht es sich? Wenn ja, hat das Feld dort Rotation — positiv für Drehung gegen den Uhrzeigersinn. Das Paddelrad reagiert darauf, ob die Strömung auf seiner einen Seite stärker schiebt als auf der anderen.

Auch hier gilt: Vektorfeld rein, an jedem Ort eine Zahl raus — in der Ebene ist \(\operatorname{rot}\vec{F}\) ebenfalls ein Skalarfeld. (Im dreidimensionalen Raum wird daraus ein Pfeil, passend zum Kreuzprodukt; in unserer Ebene bleibt davon nur die Komponente übrig, die aus dem Papier ragt — eine Zahl. Mehr dazu im Kleingedruckten.)

2.4.1 Eine Rotation komplett von Hand

Dieselbe Vier-Nachbarn-Rechnung, jetzt am Wirbelfeld \(\vec{F} = (u, v) = (-y, x)\) — dem Karussell aus Kapitel 1 —, wieder am Punkt \((1{,}0,\; 0{,}5)\) mit \(h = 0{,}1\). Beachte die Überkreuz-Logik der Formel: Für \(\partial v/\partial x\) fragen wir die y-Komponente beim rechten/linken Nachbarn ab, für \(\partial u/\partial y\) die x-Komponente beim oberen/unteren:

def u(x, y): return -y         # x-Komponente des Wirbelfelds
def v(x, y): return x          # y-Komponente

px, py, h = 1.0, 0.5, 0.1      # derselbe Punkt und Schritt wie eben

dv_dx = (v(px+h, py) - v(px-h, py)) / (2*h)
du_dy = (u(px, py+h) - u(px, py-h)) / (2*h)
print(f"v rechts = {v(px+h, py):+.2f}, v links = {v(px-h, py):+.2f}"
      f"   →  ∂v/∂x ≈ {dv_dx:+.2f}")
print(f"u oben   = {u(px, py+h):+.2f}, u unten = {u(px, py-h):+.2f}"
      f"   →  ∂u/∂y ≈ {du_dy:+.2f}")
print(f"rot F am Punkt ({px}, {py}) = {dv_dx:+.2f} − ({du_dy:+.2f}) "
      f"= {dv_dx - du_dy:+.2f}")

v rechts = +1.10, v links = +0.90   →  ∂v/∂x ≈ +1.00
u oben   = -0.60, u unten = -0.40   →  ∂u/∂y ≈ -1.00
rot F am Punkt (1.0, 0.5) = +1.00 − (-1.00) = +2.00

Wieder kommt eine einzige Zahl heraus: \(+2\), Drehung gegen den Uhrzeigersinn — das Paddelrad im Karussell dreht sich, wie erwartet. Sieh dir an, woraus die 2 entsteht: Beide Differenzen tragen \(+1\) bei — die y-Komponente wächst nach rechts (\(\partial v/\partial x = +1\)), und die x-Komponente nimmt nach oben ab (\(\partial u/ \partial y = -1\), mit dem Minus der Formel also auch \(+1\)). Die Formel vergleicht stur Nachbarwerte; ob das Feld dabei „rund aussieht”, ist ihr gleichgültig — das wird gleich noch wichtig.

Vorhersage-Punkt

Bevor du weiterliest: Drei Felder treten an: die Quelle \((x, y)\) (alle Pfeile vom Zentrum weg), der Wirbel \((-y, x)\) (das Karussell aus Kapitel 1) und die Scherung \((y, 0)\) (alle Pfeile waagerecht, oben nach rechts, unten nach links — wie eine Flussströmung, die am Ufer langsamer ist). Sage für jedes der drei Felder voraus: Divergenz — ja oder nein? Rotation — ja oder nein? Sechs Antworten, dann weiterlesen.

# von oben: x, y, dx — das 200er-Gitter aus dem Gradient-Abschnitt
def d_dx(f): return np.gradient(f, dx, axis=1)
def d_dy(f): return np.gradient(f, dx, axis=0)
def divergenz(u, v): return d_dx(u) + d_dy(v)
def rotation(u, v):  return d_dx(v) - d_dy(u)

# np.zeros_like(x): ein Array voller Nullen in der Form von x
felder = [("Quelle (x, y)", x, y),
          ("Wirbel (−y, x)", -y, x),
          ("Scherung (y, 0)", y, np.zeros_like(x))]

fig, achsen = plt.subplots(2, 3, figsize=(9.0, 6.2))
s = (slice(None, None, 20), slice(None, None, 20))
for spalte, (name, u, v) in enumerate(felder):
    for zeile, (wert, art) in enumerate([(divergenz(u, v), "div"),
                                         (rotation(u, v), "rot")]):
        ax = achsen[zeile, spalte]
        bild = ax.imshow(wert, extent=(-2, 2, -2, 2), origin="lower",
                         cmap="RdBu_r", vmin=-2.2, vmax=2.2)
        ax.quiver(x[s], y[s], u[s], v[s], color="k", scale=30, width=0.005)
        ax.set_title(f"{name}: {art} = {wert[130, 130]:+.1f}", fontsize=9)
        ax.set_xticks([]); ax.set_yticks([])
fig.colorbar(bild, ax=list(achsen.flat), shrink=0.75,
             label="div (oben) bzw. rot (unten)")
plt.show()

Abbildung 2.3: Drei Felder unter dem Quellen-Mikroskop (oben, Divergenz als Farbe) und dem Wirbel-Mikroskop (unten, Rotation als Farbe). Die Quelle hat nur Divergenz, der Wirbel nur Rotation — und die Scherung hat Rotation, obwohl keine einzige Feldlinie gekrümmt ist.

Eine Lesehilfe, denn in diesen Bildern stecken zwei Felder übereinander: Die schwarzen Pfeile zeigen die Eingabe — das Vektorfeld \(\vec{F}\), in jeder Spalte oben und unten dasselbe. Die Farbe dahinter zeigt die Ausgabe — das Skalarfeld \(\operatorname{div}\vec{F}\) (obere Zeile) bzw. \(\operatorname{rot}\vec{F}\) (untere Zeile), also an jedem Ort genau die eine Zahl aus den Handrechnungen von eben, dargestellt wie eine Temperaturkarte: rot positiv, blau negativ, weiß null. Kurz: Die Pfeile sind die Frage, die Farbe ist die Antwort.

Und die Antworten: Die Quelle hat Divergenz +2 (genau die Zahl aus unserer Handrechnung — die Farbe bestätigt, dass sie überall gilt) und keine Rotation; der Wirbel hat Rotation +2 und keine Divergenz — so weit, so erwartbar. Die Überraschung ist die Scherung: schnurgerade Feldlinien, und trotzdem \(\operatorname{rot} = -1\). Halte das Paddelrad hinein: Oben schiebt die Strömung nach rechts, unten nach links — es dreht sich (im Uhrzeigersinn, daher negativ), obwohl nirgends etwas „im Kreis fließt”. Rotation sieht man nicht an der Form der Feldlinien, sondern am Quergefälle der Geschwindigkeit. Wenn deine sechste Antwort falsch war, bist du in bester Gesellschaft — und um eine wichtige Erfahrung reicher.

Naheliegende Vermutung

Vermutung: „Das Coulomb-Feld einer Punktladung hat überall positive Divergenz — die Pfeile fächern doch sichtbar auseinander.”

Warum sie naheliegt: Auseinanderlaufende Pfeile sehen exakt so aus wie das Quellenfeld \((x, y)\) oben links, und das hat überall Divergenz.

Was stattdessen stimmt: Die Pfeile fächern auf, ja — aber sie werden mit dem Abstand auch kürzer, und zwar genau im richtigen Maß: Die Kästchen-Bilanz geht überall auf null auf, außer am Ort der Ladung selbst. Das prüfen wir sofort nach:

# von oben: x, y, divergenz()   (Galerie-Abschnitt)
r2 = np.maximum(x**2 + y**2, 1e-12)
ex, ey = x / r2, y / r2          # Coulomb-Feld in der Ebene: r̂/r

div_E = divergenz(ex, ey)
aussen = r2 > 0.5**2             # alles außerhalb von r = 0,5
print(f"|div E| außerhalb der Ladung: höchstens "
      f"{np.max(np.abs(div_E[aussen])):.2e}")
print(f"div E direkt an der Ladung:   {np.max(div_E):.2e}")

|div E| außerhalb der Ladung: höchstens 1.17e-02
div E direkt an der Ladung:   3.96e+03

Fünf Größenordnungen Unterschied: Außerhalb der Ladung ist die Divergenz numerisch null (der winzige Rest ist Diskretisierungsfehler und schrumpft mit feinerem Gitter) — nur an der Ladung selbst explodiert sie. Quellen des elektrischen Feldes sind allein die Ladungen. Diesen Satz wirst du in Kapitel 3 wiedersehen: Er ist die erste Maxwell-Gleichung.

In dieser HTML-Fassung stehen dir beide Mikroskope live im Browser zur Verfügung — die folgende Zelle ist editierbar und läuft per Run-Knopf (oder Strg+Enter). Trage in den Stellschrauben-Zeilen ein eigenes Feld \((u, v)\) ein und sage vorher an, was die beiden Karten zeigen werden. Drei Kandidaten mit wachsendem Aha-Potenzial: das Karussell mit nach außen abfallendem Tempo (u = -y*np.exp(-(x**2+y**2)), v = x*np.exp(-(x**2+y**2)) — wo wechselt die Rotation das Vorzeichen?), der \(1/r\)-Strudel aus Übung 2.4 (u = -y/(x**2+y**2+0.01), v = x/(x**2+y**2+0.01) — dreht das Paddelrad?) und ein selbst erfundenes Feld:

2.5 Zwei Sätze, die Buchhaltung machen

Divergenz und Rotation sind lokale Größen — ein Wert pro winzigem Kästchen. Zwei berühmte Sätze sagen, was passiert, wenn man viele Kästchen zusammenfasst. Der Satz von Gauß:

\[ \oint_{\text{Rand}} \vec{F} \cdot \vec{n} \; \mathrm{d}s \;=\; \iint_{\text{Fläche}} \nabla \cdot \vec{F} \; \mathrm{d}A . \]

Links steht der Fluss: Wir laufen den Rand einer Fläche ab und summieren, wie viel Feld durch ihn hindurch zeigt — das misst das Skalarprodukt mit dem nach außen zeigenden Normalenpfeil \(\vec{n}\) (Kapitel 1 lässt grüßen). Rechts steht die Summe aller Quellstärken im Inneren. Der Satz sagt schlicht: Was netto herauskommt, muss drinnen entspringen. Bei aneinandergrenzenden Kästchen hebt sich alles Innere weg — was Kästchen A nach rechts verlässt, betritt Kästchen B von links —, übrig bleibt der äußere Rand.

Der Satz von Stokes ist derselbe Gedanke für Wirbel: Die Zirkulation (wir laufen den Rand entlang und summieren, wie viel Feld mit uns mitläuft) ist gleich der Summe aller Wirbelstärken auf der Fläche:

\[ \oint_{\text{Rand}} \vec{F} \cdot \mathrm{d}\vec{l} \;=\; \iint_{\text{Fläche}} (\nabla \times \vec{F})_z \; \mathrm{d}A . \]

Glauben müssen wir beides nicht — wir können es nachrechnen. Das Testfeld \(\vec{F} = (x^2 y,\; x - y^3)\) hat ortsabhängige Divergenz und Rotation, ist also ein ehrlicher Prüfstein:

# von oben: x, divergenz(), rotation()   (Galerie-Abschnitt)
U, V = x**2 * y, x - y**3
i = slice(50, 151)               # die Box [-1, 1]² im Gitterausschnitt
rand = x[0, i]                   # Koordinaten entlang einer Boxkante

def flaechenintegral(feld):
    # np.trapezoid(werte, x): numerisches Integral per Trapezregel über
    # Abtastpunkte; zweimal angewendet — erst jede Zeile (axis=1),
    # dann das Ergebnis die Spalte hinunter
    return np.trapezoid(np.trapezoid(feld[i, i], rand, axis=1), rand)

fluss = (np.trapezoid(U[i, 150], rand) - np.trapezoid(U[i, 50], rand)
         + np.trapezoid(V[150, i], rand) - np.trapezoid(V[50, i], rand))
zirkulation = (np.trapezoid(U[50, i], rand) + np.trapezoid(V[i, 150], rand)
               - np.trapezoid(U[150, i], rand) - np.trapezoid(V[i, 50], rand))

print(f"Gauß  :  Fluss        = {fluss:+.5f}   "
      f"Quellsumme  = {flaechenintegral(divergenz(U, V)):+.5f}")
print(f"Stokes:  Zirkulation  = {zirkulation:+.5f}   "
      f"Wirbelsumme = {flaechenintegral(rotation(U, V)):+.5f}")

Gauß  :  Fluss        = -4.08142   Quellsumme  = -4.08387
Stokes:  Zirkulation  = +2.67929   Wirbelsumme = +2.67929

Beide Paare stimmen auf besser als ein Promille überein — der Rest ist Diskretisierungsfehler: verkleinere \(\Delta x\), und er schrumpft mit der Ordnung aus Abbildung 2.1. Das Kapitel-Programm programme/kap02/kap02_quellen_und_wirbel.py enthält Galerie und beide Verifikationen in voller Länge, mit feinerem Gitter und assert-Schranken — führe es aus und vergleiche mit den analytischen Werten \(-4\) und \(8/3\).

Merkkasten

Numerisch ist jede Ableitung eine Nachbar-Differenz; die zentrale Differenz hat Fehler \(\sim \Delta x^2\) — und Konvergenztests sind unser Wahrheitskriterium.
Gradient \(\nabla T\): Pfeil bergauf, senkrecht auf Höhenlinien (Skalarfeld → Vektorfeld).
Divergenz \(\nabla \cdot \vec{F}\): Quellstärke der Kästchen-Bilanz; das Coulomb-Feld ist überall quellenfrei außer an der Ladung.
Rotation \(\nabla \times \vec{F}\): Paddelrad-Test; sie steckt im Quergefälle, nicht in der Krümmung der Feldlinien.
Gauß: Fluss durch den Rand = Quellen im Inneren. Stokes: Zirkulation um den Rand = Wirbel auf der Fläche.

Roter Faden

Die beiden Produkte aus Kapitel 1 sind hier zu Operatoren geworden: Skalarprodukt mit \(\nabla\) → Divergenz, Kreuzprodukt mit \(\nabla\) → Rotation; und die Höhenlinien von dort haben jetzt ihren Gegenspieler, den Gradienten. In Kapitel 3 werden die Maxwell-Gleichungen genau vier Aussagen über \(\nabla \cdot\) und \(\nabla \times\) von \(\vec{E}\) und \(\vec{B}\) sein — die Refutation oben war schon die erste. Die Nachbar-Differenzen auf dem Gitter sind ab Kapitel 5 der Motor des FDTD-Verfahrens, und die Fehlerordnung aus Abbildung 2.1 kehrt dort als Konvergenztest wieder.

Übungen

Ü 2.1 (Verstehen). Das Feld \((x, -y)\) presst von oben und unten zusammen und zieht nach links und rechts auseinander. Divergenz? Rotation? Erst mit der Kästchen-Bilanz und dem Paddelrad argumentieren, dann mit den Formeln nachrechnen.

Musterlösung zu Ü 2.1

Beides null. Bilanz: Was das Kästchen links/rechts verlässt, kommt oben/unten exakt nach (\(\partial u/\partial x = 1\), \(\partial v/\partial y = -1\), Summe 0). Paddelrad: Die Strömung ist spiegelsymmetrisch zu beiden Achsen durchs Rad — kein Drehmoment (\(\partial v/\partial x = \partial u/\partial y = 0\)). Die Nachrechnung auf dem Kapitel-Gitter:

# von oben: x, y, divergenz(), rotation()   (Galerie-Abschnitt)
u_press, v_press = x, -y                  # das Feld (x, −y)
print(f"div, größter Betrag: {np.max(np.abs(divergenz(u_press, v_press))):.1e}")
print(f"rot, größter Betrag: {np.max(np.abs(rotation(u_press, v_press))):.1e}")

div, größter Betrag: 1.1e-14
rot, größter Betrag: 0.0e+00

Ein Feld kann also kräftig verformen, ohne Quelle oder Wirbel zu sein.

Ü 2.2 (Verstehen). Warum steht der Gradient senkrecht auf den Höhenlinien? Argumentiere ohne Formel, nur mit der Bedeutung der beiden Begriffe.

Musterlösung zu Ü 2.2

Entlang einer Höhenlinie ändert sich der Feldwert nicht — in dieser Richtung ist die Steigung null. Der Gradient zeigt in die Richtung der größten Steigung; hätte er eine Komponente entlang der Höhenlinie, würde diese Komponente nichts zur Änderung beitragen und nur Länge verschwenden. Die steilste Richtung ist also die ganz ohne Entlang-Anteil: senkrecht zur Linie.

Ü 2.3 (Verändern). Verschiebe im Kapitel-Programm die Gauß-Box von \([-1, 1]^2\) auf \([0{,}2,\; 1{,}8]^2\) (Indizes anpassen!). Sage vorher, ob Fluss und Quellsumme weiterhin übereinstimmen — und ob sich ihr Wert ändert. Dann ausführen.

Musterlösung zu Ü 2.3

Sie stimmen weiterhin überein — der Satz von Gauß gilt für jede Box. Der Wert ändert sich aber: Die Divergenz \(2xy - 3y^2\) ist ortsabhängig, eine andere Box fängt andere Quellstärke ein. Auf dem 401er-Gitter des Kapitel-Programms entspricht \([0{,}2,\; 1{,}8]\) den Indizes 220 bis 380:

achse4 = np.linspace(-2, 2, 401)
h4 = achse4[1] - achse4[0]
X4, Y4 = np.meshgrid(achse4, achse4)
U4, V4 = X4**2 * Y4, X4 - Y4**3

box = slice(220, 381)                     # x, y ∈ [0,2 ; 1,8]
rand4 = achse4[box]
div4 = np.gradient(U4, h4, axis=1) + np.gradient(V4, h4, axis=0)

quellsumme = np.trapezoid(np.trapezoid(div4[box, box], rand4, axis=1), rand4)
fluss = (np.trapezoid(U4[box, 380], rand4) - np.trapezoid(U4[box, 220], rand4)
         + np.trapezoid(V4[380, box], rand4) - np.trapezoid(V4[220, box], rand4))
print(f"Fluss durch den Rand  : {fluss:+.4f}")
print(f"Divergenz aufsummiert : {quellsumme:+.4f}")

Fluss durch den Rand  : -4.1984
Divergenz aufsummiert : -4.1988

Beide Werte stimmen weiter überein, sind aber andere als bei der zentrierten Box — und negativ, da \(-3y^2\) in der oberen Box dominiert.

Ü 2.4 (Übertragen). Untersuche den Strudel \(\vec{F} = (-y, x)/(x^2 + y^2)\): Berechne seine Rotation numerisch (Vorsicht am Ursprung — schneide ihn aus) und dann die Zirkulation auf einem Kreis oder einer Box um den Ursprung. Du wirst etwas scheinbar Widersprüchliches finden. Wie passt es mit dem Satz von Stokes zusammen?

Musterlösung zu Ü 2.4

Die Rotation ist überall null (außer am Ursprung), die Zirkulation um den Ursprung trotzdem \(2\pi\):

# von oben: x, y, rotation()   (Galerie-Abschnitt)
r2_st = np.maximum(x**2 + y**2, 1e-12)
u_st, v_st = -y/r2_st, x/r2_st             # der Strudel

rot_st = rotation(u_st, v_st)
aussen = r2_st > 0.5**2                    # Ursprung ausgeschnitten
print(f"max |rot| außerhalb r = 0,5 : {np.max(np.abs(rot_st[aussen])):.1e}")

t = np.linspace(0, 2*np.pi, 2000)          # Kreis mit Radius 1
px_, py_ = np.cos(t), np.sin(t)
ub, vb = -py_/(px_**2 + py_**2), px_/(px_**2 + py_**2)
zirk = np.trapezoid(-ub*np.sin(t) + vb*np.cos(t), t)
print(f"Zirkulation um den Ursprung : {zirk:.6f}   (2π = {2*np.pi:.6f})")

max |rot| außerhalb r = 0,5 : 1.2e-02
Zirkulation um den Ursprung : 6.283185   (2π = 6.283185)

Scheinbar ein Widerspruch zu Stokes. Die Auflösung: Die Fläche im Satz von Stokes enthält den Ursprung, und genau dort ist das Feld singulär — die gesamte Wirbelstärke sitzt als unendlich scharfe Spitze in diesem einen Punkt (wie die Divergenz der Punktladung in der Refutation). Jede Schleife, die den Ursprung nicht umschließt, hat Zirkulation null. Merk dir dieses Feld: Es ist das Magnetfeld um einen stromdurchflossenen Draht, und die „Zirkulation = eingeschlossene Stärke”-Regel ist das Ampèresche Gesetz aus Kapitel 3.

Das Kleingedruckte

Ränder: np.gradient benutzt am Arrayrand einseitige Differenzen (Fehler \(\sim \Delta x\) statt \(\Delta x^2\)). Für unsere Bilder egal — für Präzisionsvergleiche legt man die Auswertebox wie oben ins Innere.
Singularitäten: „Divergenz null außer an der Ladung” heißt mathematisch: Die Quelldichte ist eine Delta-Distribution. Numerisch zeigt sie sich als ein einzelner explodierender Gitterwert, dessen Höhe von \(\Delta x\) abhängt — der integrierte Fluss durch eine umschließende Box ist dagegen gitterunabhängig sinnvoll.
2D vs. 3D: In drei Dimensionen ist die Rotation ein voller Vektor mit drei Komponenten; unsere Ebene zeigt nur seine z-Komponente. Divergenz und beide Integralsätze übertragen sich wörtlich (aus der Randkurve wird eine Hüllfläche, aus der Fläche ein Volumen).