15 Das große Finale
Erinnerst du dich an Kapitel 1? Ein Stahlträger unter einer Betondecke, links im kühlen Auflager der Außenwand, und unten rechts ein Brandherd. Drei Regler — die Temperatur der Brandgase, die Netzfeinheit, die Zeit —, ein Knopfdruck, und ein Farbfilm zeigte, wie die Wärme vom Feuer her durch den Träger kroch und wie er sich dabei krümmte. Damals war das eine Blackbox: Du durftest daran drehen, aber der Motorraum war verschweißt. „Jede der rund vierhundert Zeilen dahinter”, hieß es, „öffnen wir bis zum Schluss.”
Jetzt ist der Schluss da. In diesem Kapitel steht dasselbe Programm noch einmal vor dir — mit offener Motorhaube. Wir bauen es aus Teilen zusammen, die du alle schon geschrieben hast: den Löser aus 1, die Dreiecke aus Kapitel 10, die Rechtecke aus Kapitel 11, die Zeitschleife aus Kapitel 9, die Thermoelastik aus Kapitel 14. Kein einziger neuer Baustein kommt hinzu. Das Finale erklärt nichts Neues — es ist eine geführte Werksbesichtigung und danach eine Experiment-Werkstatt, in der du mit allem, was du gelernt hast, gezielte Versuche fährst.
Lernziele
Nach diesem Kapitel kannst du …
- … die Architektur eines vollständigen FEM-Programms aus dem Gedächtnis zeichnen (Netz → Assemblierung → Randbedingungen → Zeitschleife → Löser → Kopplung → Ausgabe) und jedes Modul einem Buchkapitel zuordnen,
- … ein gezieltes Experiment entwerfen: eine Frage stellen, die Regler wählen, das Ergebnis gegen eine Vorhersage prüfen,
- … den Einfluss von Netzdichte, Zeitschritt und Randbedingungen auf ein gekoppeltes Ergebnis trennen — welcher Regler wirkt worauf,
- … die Grenzen dieses Programms benennen (2D, linear, isotrop, Einweg-Kopplung, kleine Netze) und sagen, welches Werkzeug jenseits davon weiterhilft,
- … jedes Ergebnis mit der Simulations-Checkliste aus Kapitel 12 absichern, bevor du ihm glaubst.
Vermutung: „Ein ‚richtiges’ FEM-Programm — OpenGeoSys, Abaqus, ANSYS — ist etwas grundlegend anderes als unsere Kapitel-Programme. Das hier war didaktisches Spielzeug; die Profis machen es völlig anders.”
Warum sie naheliegt: Die Profi-Programme haben Millionen Codezeilen, Handbücher wie Telefonbücher und lösen Netze mit Milliarden Elementen. Dagegen wirken unsere paar hundert Zeilen reines Python wie ein Bastelmodell neben einem Flugzeug.
Was stattdessen stimmt: Sie haben exakt dieselben Module — Geometrie und Netz, Elementmatrizen, Assemblierung, Randbedingungen, Zeitschritte, Löser, Kopplung. Was sie unterscheidet, ist nicht die Architektur, sondern die Ausstattung: mehr Elementtypen, schlauere Löser, krumme Ränder, Jahrzehnte Feinschliff. Wer dieses Finale versteht, liest ihre Handbücher als Landkarte statt als Fremdsprache — jeder Fachbegriff darin hat in diesem Buch seinen Platz. Genau das prüfen wir am Ende an OpenGeoSys nach.
15.1 Die Werksbesichtigung: das Programm als Fließband
Das vollständige Programm liegt als programme/kap15/fem_balken.py neben diesem Kapitel — rund vierhundert Zeilen reines Python, jede Funktion mit einem Kapitel-Rückverweis in ihrem Docstring. Bevor wir einzelne Module aufklappen, hier das ganze Werk aus der Vogelperspektive.
Abbildung 15.1 zeigt es als zwei ineinandergreifende Fließbänder. Oben läuft die Zeitschleife der Wärme (aus Kapitel 9): Sie assembliert einmal die Leitungsmatrix und die Wärmeeimer und lässt dann die Wärmefront Schritt für Schritt weiterkriechen. An den gewünschten Ausgabeschritten reicht sie das aktuelle Temperaturfeld nach unten an das Kopplungs-Fließband der Mechanik (aus Kapitel 14), das daraus Verformung und Spannung rechnet. Daneben läuft, ohne je einzugreifen, der Energiebilanz-Wachhund aus Kapitel 12 mit.
Code
import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(11.0, 4.6))
def kasten(x, y, w, text, farbe, kapitel):
ax.add_patch(plt.Rectangle((x, y), w, 0.9, facecolor=farbe,
edgecolor="black", lw=1.2))
ax.text(x + w / 2, y + 0.58, text, ha="center", va="center", fontsize=8.5)
ax.text(x + w / 2, y + 0.20, kapitel, ha="center", va="center",
fontsize=7.5, style="italic", color="0.25")
def pfeil(x1, y1, x2, y2):
ax.annotate("", xy=(x2, y2), xytext=(x1, y1),
arrowprops=dict(arrowstyle="-|>", color="black", lw=1.6))
blau = "#dfe8f0"; orange = "#fde8d0"; gruen = "#d9f0d3"; gelb = "#fff2c8"
# Obere Reihe: die Zeitschleife der Waerme.
kasten(0.0, 3.2, 1.7, "Netz bauen", blau, "Kap. 10/11")
kasten(2.1, 3.2, 1.9, "K und Wärme-\neimer C", blau, "Kap. 9/10")
kasten(4.4, 3.2, 2.0, "Ränder +\nimpliziter Schritt", orange, "Kap. 8/9")
kasten(6.8, 3.2, 1.7, "lösen", blau, "Kap. 4")
pfeil(1.7, 3.65, 2.1, 3.65); pfeil(4.0, 3.65, 4.4, 3.65)
pfeil(6.4, 3.65, 6.8, 3.65)
# Zeitschleifen-Ruecklauf.
ax.annotate("", xy=(5.4, 4.1), xytext=(7.65, 4.1),
arrowprops=dict(arrowstyle="-|>", color="0.4", lw=1.3,
connectionstyle="arc3,rad=-0.35"))
ax.text(6.5, 4.6, "nächster Zeitschritt", ha="center", fontsize=8,
style="italic", color="0.4")
# Uebergabe nach unten (Ausgabeschritt).
pfeil(5.4, 3.2, 5.4, 2.5)
ax.text(5.7, 2.85, "Ausgabeschritt:\nTemperaturfeld", ha="left", va="center",
fontsize=7.5, style="italic", color="0.3")
# Untere Reihe: die Kopplungskette der Mechanik.
kasten(0.0, 1.3, 1.7, "Element-\nΔT", orange, "Kap. 14")
kasten(2.1, 1.3, 1.9, "thermische\nKräfte F", blau, "Kap. 13")
kasten(4.4, 1.3, 1.7, "lösen\nK u = F", blau, "Kap. 4")
kasten(6.5, 1.3, 1.7, "u → ε → σ", gruen, "Kap. 14")
kasten(8.6, 1.3, 1.7, "Farbfilm +\nverformtes Netz", gruen, "Kap. 15")
pfeil(1.7, 1.75, 2.1, 1.75); pfeil(4.0, 1.75, 4.4, 1.75)
pfeil(6.1, 1.75, 6.5, 1.75); pfeil(8.2, 1.75, 8.6, 1.75)
# Der Wachhund.
kasten(8.6, 3.2, 1.7, "Energie-\nbilanz", gelb, "Kap. 12")
ax.text(9.45, 2.9, "läuft immer mit", ha="center", fontsize=7.5,
style="italic", color="0.4")
ax.set_xlim(-0.2, 10.5)
ax.set_ylim(1.0, 5.0)
ax.axis("off")
plt.tight_layout()
plt.show()
Was man hier sieht: die ganze Reise dieses Buches auf einem Bild. Von links nach rechts, von oben nach unten fließt eine gekoppelte Simulation — und jeder Kasten trägt die Nummer des Kapitels, das ihn gebaut hat. Es gibt keinen Kasten ohne Kapitel. Genau das ist die Botschaft des Finales: Du kennst jedes Glied.
15.1.1 Fünf Module, fünf Ausschnitte
Jetzt klappen wir den Motorraum auf. Für fünf Kernmodule zeigen wir den entscheidenden Ausschnitt aus fem_balken.py und sagen in einem Absatz, warum er genau so aussieht. Dies ist die einzige Stelle im ganzen Buch, an der Code im Mittelpunkt steht — nicht, weil wir hier etwas erklären (das ist alles längst geschehen), sondern weil der Code selbst hier der Lerngegenstand ist: der Beweis, dass alles zusammenpasst.
Modul 1 — die Wärmeeimer (gelumpte_kapazitaeten, ab Zeile 191). Jedes Element verteilt seine Wärmekapazität \(\rho\,c\cdot\text{Fläche}
\cdot\text{Tiefe}\) zu gleichen Teilen auf seine Knoten — beim Dreieck auf drei, beim Rechteck auf vier. Das ist das Lumping aus Kapitel 9, nur eben in der Fläche:
speicher = material["dichte"] * material["kapazitaet"]
for element in elemente:
gesamt = speicher * elementflaeche(knoten, element) * TIEFE_M
anteil = gesamt / len(element) # Dreieck: /3, Rechteck: /4
for ecke in element:
kapazitaeten[ecke] = kapazitaeten[ecke] + anteilWarum genau so: len(element) ist drei oder vier — dieselbe eine Zeile bedient beide Elementfamilien, weil die verteilte Wärmeträgheit in beiden Fällen einfach zu gleichen Eimern zusammengeklumpt wird. Kein Sonderfall, kein if.
Modul 2 — der implizite Zeitschritt mit allen Rändern (ein_waermeschritt, ab Zeile 217). Das Herzstück. Für jeden Knoten steht die Gleichung des impliziten Euler aus Kapitel 9,
\[ \frac{C_i}{\Delta t}\,T^{\text{neu}}_i \;+\;\sum_j K_{ij}\,T^{\text{neu}}_j \;=\;\frac{C_i}{\Delta t}\,T^{\text{alt}}_i \;+\; \text{Last}_i, \]
nur ist die linke Seite jetzt kein tridiagonales Band mehr, sondern der volle 9-Punkt-Stern der Fläche — also lösen wir sie mit Gauß-Seidel aus 1 statt mit Thomas. Die vier Randtypen aus Kapitel 8 greifen genau dort an, wo sie dort hergeleitet wurden:
if typ == "leistung": # Neumann: Watt auf die rechte Seite
rechte_seite[knoten_nr] += bedingung[1] * rand_flaeche / gesamt_flaeche
elif typ == "robin": # Brandgas/Luft: Zwitter (Kapitel 8)
matrix[knoten_nr][knoten_nr] += alpha * rand_flaeche # -> Diagonale
rechte_seite[knoten_nr] += alpha * rand_flaeche * t_umgebung # -> rechte Seite
elif typ == "temperatur": # Dirichlet: festhalten und streichen
festwerte[knoten_nr] = bedingung[1]
# "isoliert": nichts tun -- der Fluss ist null.Warum genau so: Der Robin-Rand ist im Finale kein Nebendarsteller mehr, sondern die tragende Randbedingung — mit ihm treffen die Brandgase auf den Träger. Er ist der Zwitter, den Kapitel 8 zwar beschrieben, aber nicht ausgerechnet hat: Der Übergangsstrom \(Q = \alpha_c\,A\,(T_{\text{Gas}} - T)\) hängt von der gesuchten Temperatur \(T\) selbst ab. Deshalb wandert sein einer Teil (\(\alpha_c A\)) auf die Diagonale der Matrix, sein anderer (\(\alpha_c A\,T_{\text{Gas}}\)) auf die rechte Seite. Genau dieselbe Aufteilung trägt auch der brandherd-Rand der Unterseite — er ist ein Robin-Rand mit ortsabhängiger Umgebung: unter dem rechten Bereich das heiße Brandgas, links die kühle Luft. Damit ist das Versprechen aus Kapitel 8 eingelöst. Die feste Temperatur des Auflagers (Dirichlet) dagegen greift nicht in die Matrix, sondern wird per Streich-Rezept aus Kapitel 8 herausgelöst — dieselbe Funktion reduziere_system, die schon die Dreiecke aus Kapitel 10 benutzten.
Modul 3 — die Einweg-Kopplung (koppel_mechanik, ab Zeile 367). An den Ausgabeschritten übergibt die Wärme ihr Temperaturfeld an die Mechanik — und nur in diese eine Richtung (Kapitel 14):
delta_t = element_uebertemperaturen(temperaturen, dreiecke) # Kap. 14
matrix = assembliere_elastik(knoten, dreiecke, e_modul, nu, TIEFE_M)
kraefte = assembliere_thermokraefte(knoten, dreiecke, e_modul, nu, TIEFE_M,
alpha, delta_t)
verschiebungen = loese_elastik(matrix, kraefte, lager) # Kap. 4/14
spannungen = elementspannungen(knoten, dreiecke, verschiebungen, ...)Warum genau so: Es ist Wort für Wort die Rechenkette aus Kapitel 14 — assembliere_elastik, assembliere_thermokraefte, loese_elastik, elementspannungen sind dieselben Funktionen aus programme/gemeinsam/elastik2d.py, die dort das Kapitel getragen haben. Das Finale importiert sie, es kopiert sie nicht. Die Einbahnstraße (Wärme treibt Mechanik, nicht umgekehrt) steckt schon darin, dass wir erst die ganze Wärmerechnung fertig machen und dann ihr Ergebnis als feste Last einsetzen.
Modul 4 — der Wachhund (nodale_leistungen, ab Zeile 400). In jedem Schritt sammelt er die Leistung ein, die über die Ränder herein- oder hinausfließt — als Rest-Term der impliziten Gleichung:
leitungsfluss = 0.0
for j in range(knotenzahl):
leitungsfluss = leitungsfluss + systemmatrix[i][j] * temperaturen_neu[j]
traegheit = kapazitaeten[i] / dt * (temperaturen_neu[i] - temperaturen_alt[i])
leistungen.append(traegheit + leitungsfluss) # >0 herein, <0 hinausWarum genau so: Was Trägheit (Aufheizen) und Leitung zusammen verlangen, muss von außen über den Rand geliefert worden sein. Im Inneren und an isolierten Rändern ist dieser Rest null — nur wo ein Rand Wärme ein- oder ausspeist, steht eine Zahl. Am Ende muss hereingeflossen = gespeichert + hinausgeflossen aufgehen (Kapitel 12). Schließt die Bilanz nicht, steckt ein Fehler in der Assemblierung — der Wachhund bellt, ganz ohne Vergleichslösung.
Modul 5 — die Zeitschleife selbst (simuliere, ab Zeile 449). Sie hält alles zusammen und ist kurz genug, um sie ganz zu zeigen — das ist das ganze Fließband aus Abbildung 15.1 in Python:
for schritt in range(1, zeitschritte + 1):
alt = temperaturen
temperaturen = ein_waermeschritt(alt, systemmatrix, kapazitaeten, dt, ...)
leistungen = nodale_leistungen(temperaturen, alt, systemmatrix, kapazitaeten, dt)
# ... Wachhund-Summen aktualisieren ...
if schritt in kopplungsschritte:
verschiebungen, spannungen = koppel_mechanik(knoten, dreiecke,
temperaturen, material, lager)
momentaufnahmen.append(...)Warum genau so: Ein Zeitschritt der Wärme, dann der Wachhund, dann — nur an den Ausgabeschritten — ein Sprung in die Mechanik. Mehr ist es nicht. Vierhundert Zeilen, aber ihr Kern passt in zehn.
15.2 Die Blackbox aus Kapitel 1 — jetzt mit offenem Motorraum
In Kapitel 1 lief die Blackbox mit ihren drei Reglern. Was damals nicht zu sehen war: Ihr Maschinenraum war schon genau dieses Programm — die volle 2D-Kette aus Dreiecksnetz, Brand-Rändern und Thermoelastik, nur in einer kompakten Fassung, die das Rechenpaket NumPy als Rechenknecht benutzt, damit der Browser flott bleibt. fem_balken.py rechnet dieselbe Kette in reinem Python nach — Zeile für Zeile aus den Kapiteln dieses Buches. Stellen wir beide für dieselben Regler nebeneinander — Brandgase \(500\,°\mathrm{C}\), mittleres Netz, 30 Minuten.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
from matplotlib.collections import PolyCollection
import fem_balken as fb
erg = fb.simuliere(nx=16, ny=4, elementtyp="dreiecke", material_name="stahl",
zeit_minuten=30.0, zeitschritte=10, kopplungsschritte=[10])
knoten = erg["knoten"]; dreiecke = erg["dreiecke"]
temp = erg["verlauf"][-1]; u = erg["momentaufnahmen"][-1]["verschiebungen"]
sig = erg["momentaufnahmen"][-1]["spannungen"]
fig, (a1, a2) = plt.subplots(2, 1, figsize=(9.0, 4.4))
tp = []; tw = []
for d in dreiecke:
tp.append([(knoten[d[k]][0], knoten[d[k]][1]) for k in range(3)])
tw.append((temp[d[0]] + temp[d[1]] + temp[d[2]]) / 3.0)
c1 = PolyCollection(tp, array=tw, cmap="inferno", edgecolors="0.7", linewidths=0.2)
c1.set_clim(20, max(tw)); a1.add_collection(c1); a1.autoscale_view()
a1.set_aspect("equal"); a1.set_title("Temperaturfeld nach 30 min (°C)")
fig.colorbar(c1, ax=a1, fraction=0.03, pad=0.02)
ueber = 60.0
sp = [s[0] / 1e6 for s in sig]; g = max(1.0, max(abs(v) for v in sp))
vp = []
for d in dreiecke:
vp.append([(knoten[d[k]][0] + ueber * u[2 * d[k]],
knoten[d[k]][1] + ueber * u[2 * d[k] + 1]) for k in range(3)])
c2 = PolyCollection(vp, array=sp, cmap="RdBu_r", edgecolors="0.6", linewidths=0.2)
c2.set_clim(-g, g); a2.add_collection(c2); a2.autoscale_view()
a2.set_aspect("equal"); a2.set_title("Verformung (× 60), Farbe σ_xx (N/mm²)")
a2.set_xlabel("Länge x (m)")
fig.colorbar(c2, ax=a2, fraction=0.03, pad=0.02)
plt.tight_layout(); plt.show()
Was man sieht: Die Brandgase heizen den Träger von unten und rechts, das Auflager und die Decke halten dagegen — nach 30 Minuten ist die Mitte erst rund \(52\,°\mathrm{C}\) warm, das heißeste Eck knapp über \(110\,°\mathrm{C}\). Weil der Untergurt wärmer ist als der Obergurt, wölbt sich der Träger leicht nach oben. Aber die Wölbung ist winzig — bei realistischer Konvektion braucht es Bruchteile eines Millimeters, und das Bild muss sechzigfach überhöhen, um sie zu zeigen.
| Größe | Blackbox (Kapitel 1, NumPy) | fem_balken.py (reines Python) |
|---|---|---|
| Mittentemperatur | \(52\,°\mathrm{C}\) | \(52\,°\mathrm{C}\) |
| Heißester Punkt | \(112\,°\mathrm{C}\) | \(112\,°\mathrm{C}\) |
| Hebung des freien Endes | \(0{,}15\ \mathrm{mm}\) | \(0{,}15\ \mathrm{mm}\) |
Identisch — und genau das ist die Pointe. Zwei unabhängig geschriebene Fassungen desselben Verfahrens, eine kompakt mit NumPy, eine ausführlich in reinem Python, liefern bis auf Rundung dieselben Zahlen. Das ist nebenbei ein waschechter Verifikationstest im Geist von Kapitel 12: Wer dasselbe Modell zweimal verschieden implementiert und zweimal dasselbe herausbekommt, hat eine der besten Proben bestanden, die es gibt. Der Unterschied zwischen Kapitel 1 und heute liegt also nicht im Ergebnis — er liegt darin, dass du jetzt jede Zeile verstehst.
Und die Physik dahinter bleibt bemerkenswert: Ein realistisches Konvektionsfeuer wölbt den Träger in einer halben Stunde erst um gut einen Zehntelmillimeter. Die sichtbare, kräftige Krümmung — bis zu den kanonischen \(12\ \mathrm{mm}\) — entsteht erst, wenn der Brand voll entwickelt ist und die Unterseite selbst richtig heiß wird (Kapitel 14). Genau diesen Fall bauen wir gleich in Reihe A und in der Finale-Animation.
15.3 Die Experiment-Werkstatt
Von hier an ist das Kapitel eine Werkstatt. Ein FEM-Programm ist kein Orakel, dem man eine Frage stellt und eine Wahrheit erhält — es ist ein Instrument, mit dem man Experimente fährt. Und Experimente haben eine feste Form, damit am Ende ein Erkenntnisgewinn steht und nicht bloß ein weiteres buntes Bild. Abbildung 15.3 zeigt den Rahmen, den wir für jede Versuchsreihe benutzen.
Code
import matplotlib.pyplot as plt
schritte = [("1. Frage", "Was will ich wissen?", "#e8eef5"),
("2. Vorhersage", "Was erwarte ich — vor dem Rechnen?", "#ffe9cc"),
("3. Regler", "An welchem einen Regler drehe ich?", "#e8eef5"),
("4. Ergebnis", "Was liefert das Programm?", "#e0efe0"),
("5. Schluss", "Traf die Vorhersage? Checkliste!", "#e8dff0")]
fig, ax = plt.subplots(figsize=(10.5, 2.0))
x = 0.0
for titel, text, farbe in schritte:
ax.add_patch(plt.Rectangle((x, 0), 2.0, 1.0, facecolor=farbe,
edgecolor="black", lw=1.1))
ax.text(x + 1.0, 0.66, titel, ha="center", va="center", fontsize=9.5,
weight="bold")
ax.text(x + 1.0, 0.28, text, ha="center", va="center", fontsize=7.3)
if x < 8.0:
ax.annotate("", xy=(x + 2.3, 0.5), xytext=(x + 2.0, 0.5),
arrowprops=dict(arrowstyle="-|>", color="black", lw=1.4))
x += 2.3
ax.set_xlim(-0.1, x); ax.set_ylim(-0.15, 1.15); ax.axis("off")
plt.tight_layout(); plt.show()
Was man sieht: fünf Schritte, und der zweite ist der wichtigste. Die Vorhersage vor dem Rechnen ist es, die aus Klicken ein Experiment macht — nur wer sich festlegt, kann überrascht werden, und nur wer überrascht wird, lernt. Wir fahren jetzt fünf Reihen: Netz, Zeit, Ränder, Material, Elemente. (Damit lösen wir zugleich Übung 1.2 ein — dort solltest du Regler-Experimente durchführen und protokollieren; hier ist die Vorlage.)
15.3.1 Reihe A — Das Netz: wann darf ich der Zahl glauben?
Frage. Wenn ich das Netz feiner mache, ändert sich das Ergebnis. Ab wann ändert es sich nicht mehr — wann ist das Netz fein genug? Wir nehmen den voll entwickelten Brand: Die Flammen halten die Unterseite auf \(220\,°\mathrm{C}\), die Betondecke die Oberseite auf \(20\,°\mathrm{C}\) (Gefälle \(\Delta T = 200\ \mathrm{K}\) über die Höhe, links eingespannt), und messen die Auslenkung des freien Endes auf immer feineren Netzen.
Vorhersage. Aus Kapitel 14 kennen wir die Handformel für dieses Gefälle: Krümmung \(\kappa = \alpha\,\Delta T/H = 12\cdot 10^{-6}\cdot 200/0{,}1 = 0{,}024\ \mathrm{1/m}\), Kragarm-Auslenkung \(w = L^2\kappa/2 = 12{,}0\ \mathrm{mm}\) nach oben. Und wir wissen: Grobe Dreiecksnetze biegen zu steif — die Auslenkung nähert sich der \(12\ \mathrm{mm}\) von unten. Wir erwarten also eine Zahlenreihe, die wächst und sich langsam einpendelt.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
import fem_balken as fb
gefaelle = {"links": ("isoliert",), "rechts": ("isoliert",),
"oben": ("temperatur", 20.0), "unten": ("temperatur", 220.0)}
netze = [(4, 1), (8, 2), (16, 4), (24, 6)]
knotenzahlen = []; biegungen = []; spannungen = []
for nx, ny in netze:
e = fb.simuliere(nx=nx, ny=ny, raender=gefaelle, zeit_minuten=60.0,
zeitschritte=4, kopplungsschritte=[4])
knotenzahlen.append(len(e["knoten"]))
biegungen.append(fb.enddurchbiegung(e))
spannungen.append(fb.groesste_vergleichsspannung(e))
fig, ax = plt.subplots(figsize=(8.0, 3.4))
ax.plot(knotenzahlen, biegungen, "o-", color="tab:blue", ms=8,
label="Auslenkung (mm)")
ax.axhline(12.0, ls="--", color="tab:blue", alpha=0.6,
label="Handformel 12 mm")
ax.set_xlabel("Knotenzahl (feineres Netz →)")
ax.set_ylabel("Auslenkung (mm)", color="tab:blue")
ax.tick_params(axis="y", labelcolor="tab:blue")
ax2 = ax.twinx()
ax2.plot(knotenzahlen, spannungen, "s--", color="tab:orange", ms=7,
label="max |σ_xx| (N/mm²)")
ax2.set_ylabel("max |σ_xx| (N/mm²)", color="tab:orange")
ax2.tick_params(axis="y", labelcolor="tab:orange")
ax.set_title("Reihe A: Netzkonvergenz")
ax.legend(loc="lower right", fontsize=8)
plt.tight_layout(); plt.show()
Ergebnis. Die Auslenkung wächst brav — von rund \(1\ \mathrm{mm}\) auf dem gröbsten Netz (4×1) über \(3{,}7\) und \(7{,}7\) bis auf \(9{,}6\ \mathrm{mm}\) auf dem feinsten (24×6), immer näher an die \(12\ \mathrm{mm}\). Die Maximalspannung dagegen pendelt (\(82 \to 150
\to 180 \to 159\ \mathrm{N/mm^2}\)) und will sich nicht einpendeln.
Schluss. Für die Auslenkung ist das Netz 16×4 brauchbar (rund \(64\ \%\) der Handformel), 24×6 besser (\(80\ \%\)); ganz konvergiert ist sie noch nicht, für den qualitativen Eindruck aber gut. Die Spannung ist die härtere Nuss: Ihr Höchstwert sitzt an der Einspann-Ecke, wo die Lagerung einen scharfen Knick erzwingt — eine Spannungsspitze, die ein Dreiecksnetz mit konstanter Spannung je Element nie ganz einfängt (das ist dieselbe Kerbwirkung, die in Übung 15.2 das Loch erzeugt). Merke: Verschiebungen konvergieren schneller als Spannungen, und Spannungsspitzen an Ecken konvergieren gar nicht — genau davor warnt die Konvergenzprobe aus Kapitel 12.
15.3.2 Reihe B — Die Zeit: warum implizit, und wie groß darf Δt sein?
Frage. In Kapitel 9 explodierte der explizite Zeitschritt, sobald \(\Delta t\) zu groß wurde. Unser Finale rechnet implizit. Was passiert, wenn wir hier \(\Delta t\) vergrößern — explodiert es auch?
Vorhersage. Der implizite Euler war immer stabil (Kapitel 9). Wir erwarten also keine Explosion, egal wie groß \(\Delta t\) ist — aber grobe Schritte sollten die Wärmefront verschmieren, die Mittentemperatur also unterschätzen.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
import fem_balken as fb
schrittzahlen = [1, 2, 5, 10, 20, 40]
dts = []; mitten = []
for zs in schrittzahlen:
e = fb.simuliere(nx=16, ny=4, zeit_minuten=30.0, zeitschritte=zs,
kopplungsschritte=[zs])
dts.append(30.0 * 60.0 / zs)
mitten.append(fb.mittentemperatur(e))
fig, ax = plt.subplots(figsize=(8.0, 3.2))
ax.plot(dts, mitten, "o-", color="tab:green", ms=8)
for x, y in zip(dts, mitten):
ax.annotate("%.0f °C" % y, (x, y), textcoords="offset points",
xytext=(0, 8), fontsize=7.5, ha="center")
ax.set_xlabel("Zeitschritt Δt (s) — größer →")
ax.set_ylabel("Mittentemperatur nach 30 min (°C)")
ax.set_title("Reihe B: implizit ist immer stabil — grobe Δt verschmieren nur")
ax.invert_xaxis()
plt.tight_layout(); plt.show()
Ergebnis. Von rund \(52\,°\mathrm{C}\) bei feinen Schritten fällt die Mittentemperatur auf \(49\,°\mathrm{C}\) beim einen Riesenschritt (\(\Delta t = 1800\ \mathrm{s}\)) — aber nichts explodiert. Der grobe Schritt ist ungenau, nicht instabil.
Schluss. Genau das ist der Grund, warum das Finale (und jedes ernsthafte transiente Programm) implizit rechnet: Man darf große Zeitschritte nehmen, ohne dass die Rechnung entgleist — man zahlt nur mit Genauigkeit, nicht mit einer Katastrophe. Die Faustregel bleibt die Konvergenzprobe aus Kapitel 12: \(\Delta t\) halbieren; ändert sich das Ergebnis kaum noch, war es fein genug. (Beim impliziten Euler halbiert sich der Zeitfehler beim Halbieren von \(\Delta t\) — „Zeit ist gröber als Raum”.)
15.3.3 Reihe C — Die Ränder: drei Flammen, drei Geschichten
Frage. Die Brandgase an der rechten Stirn lassen sich auf drei Arten modellieren (Kapitel 8): als feste Temperatur (Dirichlet), als eingespeiste Leistung (Neumann) oder als Übergang an ein Gas fester Temperatur (Robin). Wie unterscheiden sich die drei?
Vorhersage. Die feste Temperatur deckelt das Ende bei \(500\,°\mathrm{C}\). Die feste Leistung kennt keine solche Grenze — das Ende wird so heiß, wie die eingespeisten Watt es treiben, und kann \(500\,°\mathrm{C}\) überschreiten. Der Robin-Rand dagegen kann den Träger höchstens auf die Gastemperatur bringen, nie darüber: Sobald die Oberfläche das Gas eingeholt hat, versiegt der antreibende Temperaturunterschied.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
import fem_balken as fb
def profil_halbe_hoehe(e):
kn = e["knoten"]; t = e["verlauf"][-1]
paare = []
for i in range(len(kn)):
if abs(kn[i][1] - 0.05) < 1e-9:
paare.append((kn[i][0], t[i]))
paare.sort()
return [p[0] for p in paare], [p[1] for p in paare]
szenen = [
("Temperatur 500 °C", {"links": ("temperatur", 20.0),
"rechts": ("temperatur", 500.0), "oben": ("isoliert",),
"unten": ("isoliert",)}, "tab:red"),
("Leistung 1500 W", {"links": ("temperatur", 20.0),
"rechts": ("leistung", 1500.0), "oben": ("isoliert",),
"unten": ("isoliert",)}, "tab:orange"),
("Brandgas 500 °C (Robin)", {"links": ("temperatur", 20.0),
"rechts": ("robin", 25.0, 500.0), "oben": ("isoliert",),
"unten": ("isoliert",)}, "tab:blue"),
]
fig, ax = plt.subplots(figsize=(8.0, 3.4))
for name, raender, farbe in szenen:
e = fb.simuliere(nx=16, ny=4, raender=raender, zeit_minuten=30.0,
zeitschritte=10, kopplungsschritte=[10])
xs, ts = profil_halbe_hoehe(e)
ax.plot(xs, ts, "o-", color=farbe, ms=5, label=name)
ax.set_xlabel("Ort x (m)"); ax.set_ylabel("Temperatur (°C)")
ax.set_title("Reihe C: drei Ränder für die Brandgase")
ax.legend(fontsize=8); plt.tight_layout(); plt.show()
Ergebnis. Genau wie vorhergesagt: Die Temperatur-Flamme klebt bei \(500\,°\mathrm{C}\); die Leistungs-Flamme schießt darüber hinaus (bei \(1500\ \mathrm{W}\) auf \(542\,°\mathrm{C}\), bei \(3000\ \mathrm{W}\) auf über \(1000\,°\mathrm{C}\) — sie kennt keinen Deckel); der Robin-Rand erwärmt die kleine Stirnfläche in 30 Minuten erst auf \(59\,°\mathrm{C}\).
In Kapitel 8 blieb der Übergangsrand ein Versprechen: Der Strom \(Q = \alpha_c\,A\,(T_{\text{Gas}} - T)\) hängt von der gesuchten Temperatur \(T\) selbst ab, deshalb — so hieß es — „rechnen wir das in Kapitel 15 aus”. Hier ist es. Der Wärmeübergangskoeffizient \(\alpha_c\) beschreibt, wie gut die Brandgase ihre Wärme an den Stahl abgeben; die Norm für den Brandfall setzt \(\alpha_c = 25\ \mathrm{W/(m^2 K)}\) (CEN (Europäisches Komitee für Normung) 2005). Im Code wandert sein temperaturabhängiger Teil auf die Diagonale der Systemmatrix, sein bekannter Teil auf die rechte Seite (Modul 2 der Werksbesichtigung) — genau die „Zwitter”-Aufteilung aus Kapitel 8.
Der entscheidende Unterschied zur Leistung: Die Brandgase können den Träger höchstens auf ihre eigene Temperatur bringen, nie darüber. Lässt man dieselbe rechte Stirn stundenlang weiterheizen, kriecht das Ende beim Robin-Rand von \(59\) über \(94\) und \(147\) langsam gegen die \(500\,°\mathrm{C}\) der Gase — und bleibt darunter, weil der antreibende Unterschied \(T_{\text{Gas}} - T\) mit jedem Grad kleiner wird. Die feste Leistung dagegen kennt diese Bremse nicht: \(3000\ \mathrm{W}\) treiben das Ende in denselben Stunden auf \(1000\), \(2000\), \(4000\,°\mathrm{C}\) — längst geschmolzen. Ein echtes Feuer ist ein Gas fester Temperatur, kein Brenner fester Leistung; darum ist der Robin-Rand das ehrlichere Brandmodell. Für sehr großes \(\alpha_c\) presst er das Ende sofort auf Gastemperatur (das ist Dirichlet), für \(\alpha_c = 0\) fließt gar nichts (das ist der isolierte Rand) — alle wirklichen Übergänge liegen dazwischen.
Schluss. Welchen Rand man wählt, ist keine Rechenfrage, sondern eine Modellierungsfrage (Kapitel 12): Sind die Brandgase ein Thermostat (Temperatur), ein Brenner fester Leistung (Neumann) oder ein Gas, das über einen Übergang heizt (Robin)? Dieselbe Geometrie liefert je nach Randwahl verschiedene Wahrheiten — und nur eine passt zur Wirklichkeit.
15.3.4 Reihe D — Das Material: Stahl gegen Aluminium
Frage. Unser Träger ist aus Stahl. Wir stellen ihm einen Aluminium-Herausforderer aus derselben Geometrie gegenüber und fahren ein Renn-Duell mit drei Fragen: Wer wird schneller heiß? Wer krümmt sich stärker? Wer trägt mehr Spannung?
Vorhersage. Aluminium leitet Wärme viermal besser und speichert weniger — es sollte schneller heiß werden. Es dehnt sich fast doppelt so stark pro Grad (\(\alpha_{\text{Alu}} = 23\) gegen \(\alpha_{\text{Stahl}} = 12\cdot 10^{-6}\)) — es sollte sich stärker krümmen. Aber Stahl ist dreimal steifer (\(E = 210\) gegen \(70\ \mathrm{GPa}\)), und für die Einspann-Spannung zählt das Produkt \(E\,\alpha\) (Kapitel 13) — da liegt Stahl vorn.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
import fem_balken as fb
# Links: Aufheizkurve der Mitte, wenn die rechte Stirn heiß gehalten wird.
heiss = {"links": ("temperatur", 20.0), "rechts": ("temperatur", 500.0),
"oben": ("isoliert",), "unten": ("isoliert",)}
fig, (a1, a2) = plt.subplots(1, 2, figsize=(10.0, 3.4))
for mat, farbe in [("stahl", "tab:gray"), ("alu", "tab:blue")]:
zeiten = []; temps = []
for minute in [2, 5, 10, 20, 30]:
e = fb.simuliere(nx=12, ny=3, material_name=mat, raender=heiss,
zeit_minuten=float(minute), zeitschritte=8,
kopplungsschritte=[8])
zeiten.append(minute); temps.append(fb.mittentemperatur(e))
a1.plot(zeiten, temps, "o-", color=farbe, ms=6,
label=fb.MATERIALIEN[mat]["name"])
a1.set_xlabel("Zeit (min)"); a1.set_ylabel("Mittentemperatur (°C)")
a1.set_title("Wer wird schneller heiß?"); a1.legend(fontsize=8)
# Rechts: Krümmung und Spannung (voll entwickelter Brand, Gefälle 200 K).
gefaelle = {"links": ("isoliert",), "rechts": ("isoliert",),
"oben": ("temperatur", 20.0), "unten": ("temperatur", 220.0)}
namen = []; bieg = []; spann = []
for mat in ["stahl", "alu"]:
e = fb.simuliere(nx=16, ny=4, material_name=mat, raender=gefaelle,
zeit_minuten=60.0, zeitschritte=6, kopplungsschritte=[6])
namen.append(fb.MATERIALIEN[mat]["name"])
bieg.append(fb.enddurchbiegung(e)); spann.append(fb.groesste_vergleichsspannung(e))
x = [0, 1]
a2.bar([p - 0.18 for p in x], bieg, 0.35, color="tab:blue",
label="Auslenkung (mm)")
a2.bar([p + 0.18 for p in x], spann, 0.35, color="tab:red",
label="max |σ_xx| (N/mm²)")
a2.set_xticks(x); a2.set_xticklabels(namen)
a2.set_title("Wer krümmt sich, wer trägt Spannung?"); a2.legend(fontsize=8)
plt.tight_layout(); plt.show()
Ergebnis. Alle drei Vorhersagen treffen zu. Aluminium wird dramatisch schneller heiß (nach 10 Minuten \(74\) gegen kaum \(21\,°\mathrm{C}\) in der Mitte). Aluminium krümmt sich stärker (\(14{,}8\) gegen \(7{,}7\ \mathrm{mm}\)). Aber Stahl trägt mehr Spannung (\(180\) gegen \(122\ \mathrm{N/mm^2}\)).
Schluss. Es gibt keinen „besseren” Werkstoff — es kommt auf die Frage an. Wer schnelle, gleichmäßige Wärme will (Kühlkörper, Bratpfanne), nimmt Aluminium. Wer Kräfte tragen und Form halten will (Träger, Welle), nimmt Stahl. Und wer die Wärmespannung fürchtet, erinnert sich: Nicht die Dehnung allein spannt, sondern \(E\,\alpha\) — der steife Stahl bezahlt seine Steifigkeit mit höherer Spannung (das ist Übung 13.2, jetzt am 2D-Träger bestätigt).
15.3.5 Reihe E — Die Elemente: Dreiecke gegen Rechtecke
Frage. Wir haben zwei Elementfamilien gebaut: P1-Dreiecke (Kapitel 10) und Q1-Rechtecke (Kapitel 11). Auf demselben Gitter, demselben Problem — liefern sie dasselbe?
Vorhersage. Beide sind lineare Elemente zweiter Konvergenzordnung (Kapitel 12), beide erfüllen die Energiebilanz exakt. Auf einem glatten Problem sollten sie fast dasselbe liefern; kleine Unterschiede sind möglich, weil das Rechteck auch die diagonalen Nachbarn koppelt (der 9-Punkt-Stern aus Kapitel 11).
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import matplotlib.pyplot as plt
import fem_balken as fb
fig, ax = plt.subplots(figsize=(8.0, 3.2))
stile = {"dreiecke": ("tab:blue", "o", "Dreiecke (P1)"),
"rechtecke": ("tab:orange", "s", "Rechtecke (Q1)")}
for et in ["dreiecke", "rechtecke"]:
e = fb.simuliere(nx=16, ny=4, elementtyp=et, zeit_minuten=30.0,
zeitschritte=10, kopplungsschritte=[10])
kn = e["knoten"]; t = e["verlauf"][-1]
paare = []
for i in range(len(kn)):
if abs(kn[i][1] - 0.05) < 1e-9:
paare.append((kn[i][0], t[i]))
paare.sort()
farbe, marker, name = stile[et]
ax.plot([p[0] for p in paare], [p[1] for p in paare], marker=marker,
color=farbe, ms=6, alpha=0.75,
label="%s: Mitte %.1f °C" % (name, fb.mittentemperatur(e)))
ax.set_xlabel("Ort x (m)"); ax.set_ylabel("Temperatur (°C)")
ax.set_title("Reihe E: Dreiecke und Rechtecke — praktisch deckungsgleich")
ax.legend(fontsize=8); plt.tight_layout(); plt.show()
Ergebnis. Die beiden Kurven liegen übereinander — beide treffen rund \(51{,}5\,°\mathrm{C}\) in der Mitte (Dreiecke \(51{,}5\), Rechtecke \(51{,}3\)), beide schließen die Energiebilanz perfekt.
Schluss. Für dieses glatte, achsparallele Problem ist der Elementtyp zweitrangig — ein beruhigendes Ergebnis, denn es heißt, dass die Physik und nicht die Elementwahl das Bild bestimmt. Unterschiede treten erst auf, wo es darauf ankommt: bei krummen Rändern (dort können Rechtecke gar nicht, und man braucht die isoparametrischen Elemente aus der Grenzen-Tafel) oder bei sehr groben Netzen. Die Wahl der Elementfamilie ist eine Frage der Geometrie, nicht der Genauigkeit.
15.4 Das Finale in Bewegung
Die fünf Reihen haben je einen Regler isoliert. Jetzt lassen wir alles zusammenlaufen — die Auflösung des Kapitel-1-Trailers, Glied für Glied verstanden: der voll entwickelte Brand hält die Unterseite heiß, die Wärme kriecht von unten herein, das Temperaturgefälle über die Höhe wächst, und der Träger wölbt sich synchron nach oben.
Code
import os, sys
sys.path.insert(0, os.path.join("..", "programme", "kap15"))
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.collections import PolyCollection
from matplotlib.animation import FuncAnimation
from IPython.display import HTML
import fem_balken as fb
gefaelle = {"links": ("isoliert",), "rechts": ("isoliert",),
"oben": ("temperatur", 20.0), "unten": ("temperatur", 220.0)}
schritte = list(range(1, 13))
erg = fb.simuliere(nx=16, ny=4, raender=gefaelle, zeit_minuten=30.0,
zeitschritte=12, kopplungsschritte=schritte)
knoten = erg["knoten"]; dreiecke = erg["dreiecke"]
momente = erg["momentaufnahmen"]; ueber = 8.0
def dreieck_polys(verformt):
polys = []
for d in dreiecke:
polys.append([verformt[d[k]] for k in range(3)])
return polys
fig, (axl, axr) = plt.subplots(1, 2, figsize=(11.0, 3.0),
gridspec_kw={"width_ratios": [1.4, 1.6]})
def zeichne(k):
axl.clear(); axr.clear()
m = momente[k]; temp = m["temperaturen"]; u = m["verschiebungen"]; sig = m["spannungen"]
tp = [[(knoten[d[j]][0], knoten[d[j]][1]) for j in range(3)] for d in dreiecke]
tw = [(temp[d[0]] + temp[d[1]] + temp[d[2]]) / 3.0 for d in dreiecke]
cl = PolyCollection(tp, array=np.array(tw), cmap="inferno", edgecolors="none")
cl.set_clim(20, 220); axl.add_collection(cl); axl.autoscale_view()
axl.set_aspect("equal"); axl.set_xlabel("x (m)")
axl.set_title("Temperatur (%.0f min)" % (m["zeit"] / 60.0), fontsize=9)
verf = [(knoten[i][0] + ueber * u[2 * i], knoten[i][1] + ueber * u[2 * i + 1])
for i in range(len(knoten))]
sp = [s[0] / 1e6 for s in sig]
cr = PolyCollection(dreieck_polys(verf), array=np.array(sp), cmap="RdBu_r",
edgecolors="0.6", linewidths=0.2)
cr.set_clim(-150, 150); axr.add_collection(cr)
axr.plot([0, 1], [0, 0], "--", color="0.7", lw=0.7)
axr.set_xlim(-0.05, 1.15); axr.set_ylim(-0.05, 0.20); axr.set_aspect("equal")
axr.set_title("Träger wölbt sich (× 8), σ_xx", fontsize=9); axr.axis("off")
return []
ani = FuncAnimation(fig, zeichne, frames=len(momente), interval=260, blit=False)
plt.close(fig)
HTML(ani.to_jshtml())Was man sieht: Genau der Trailer aus Kapitel 1 — nur ist jetzt jede Farbe, jede Krümmung, jeder Millimeter das Ergebnis einer Rechnung, die du Zeile für Zeile kennst. Das ist das Versprechen des ersten Kapitels, eingelöst.
15.5 Die interaktive Einheit: dein Finale
Jetzt gehört das Programm dir. Die folgende Einheit ist das ganze Finale im Browser — alle Regler offen. Der Code darin ist eine kompakte Inline-Fassung; das vollständige, ausführlich kommentierte Programm ist programme/kap15/fem_balken.py, an dem sich dieses Kapitel entlanghangelt.
Bevor du startest, eine Meister-Vorhersage über zwei Regler auf einmal:
Auflösung der Meister-Vorhersage. Überraschung: Aluminium mit halber Flamme wölbt sich fast genau gleich wie Stahl mit voller — \(5{,}8\) gegen \(6{,}0\ \mathrm{mm}\) (auf diesem groben 12×3-Netz). Der Grund steht in der Krümmungsformel \(\kappa = \alpha\,\Delta T / H\): Aluminiums doppeltes \(\alpha\) mal dem halben \(\Delta T\) landet wieder bei Stahls Krümmung. Aber schau auf die Spannung — der Stahl trägt \(184\ \mathrm{N/mm^2}\) gegen Aluminiums \(62\): rund dreimal so viel. Zwei Regler, die sich in der Krümmung aufheben, aber in der Spannung weit auseinanderlaufen. Genau solche Wechselwirkungen macht ein Experiment sichtbar, das eine bloße Vorhersage nie fände.
15.6 Die Landkarte, alle Häkchen gesetzt
Wir sind am Ziel. Abbildung 15.11 ist die Landkarte aus Kapitel 1, zum letzten Mal — jetzt mit allen Häkchen. Jeder der fünf Bausteine ist gefüllt, und die drei Erweiterungen (Verformung, Prüfen, Finale) stehen daneben, ebenfalls abgehakt.
Code
import matplotlib.pyplot as plt
bausteine = [
("Geometrie\n& Netz", "Kap. 10–11"),
("Physik", "Kap. 5, 6, 9"),
("Rand-\nbedingungen", "Kap. 8"),
("Gleichungs-\nsystem", "Kap. 3, 7"),
("Löser", "Kap. 4"),
]
farben = ["#cfe8f3", "#ffe0b3", "#ffd6d6", "#d9f0d3", "#e6dcf0"]
fig, ax = plt.subplots(figsize=(9.5, 3.4))
breite, luecke = 1.5, 0.4
for k, ((name, kapitel), farbe) in enumerate(zip(bausteine, farben)):
x0 = k * (breite + luecke)
ax.add_patch(plt.Rectangle((x0, 0), breite, 1.0, facecolor=farbe,
edgecolor="black", lw=1.3))
ax.text(x0 + breite / 2, 0.66, name, ha="center", va="center",
fontsize=10, weight="bold")
ax.text(x0 + breite / 2, 0.30, kapitel, ha="center", va="center", fontsize=9)
ax.text(x0 + breite - 0.14, 0.86, "✓", ha="center", va="center",
fontsize=15, color="tab:green", weight="bold")
if k < len(bausteine) - 1:
ax.annotate("", xy=(x0 + breite + luecke, 0.5),
xytext=(x0 + breite, 0.5),
arrowprops=dict(arrowstyle="-|>", color="black", lw=1.3))
mitte = (len(bausteine) * (breite + luecke) - luecke) / 2
ax.text(mitte, 1.62,
"✓ Verformung: Kap. 13–14 ✓ Prüfen & Vertrauen: Kap. 12 "
"✓ das Finale: Kap. 15",
ha="center", va="center", fontsize=9.5, weight="bold", color="tab:green")
ax.text(mitte, -0.25, "Python-Grundlagen: Kap. 2", ha="center", va="center",
fontsize=9, style="italic")
ax.set_xlim(-0.2, len(bausteine) * (breite + luecke) - luecke + 0.2)
ax.set_ylim(-0.45, 1.85)
ax.axis("off")
plt.tight_layout()
plt.show()
Was man sieht: dasselbe Bild wie am Anfang, aber jeder Kasten trägt jetzt sein grünes Häkchen. Das ist kein Schmuck — es ist die Buchhaltung deines Wissens. Du hast jeden dieser Kästen selbst gebaut.
15.7 Grenzen und Ausblick: die Landkarte des Weiterlernens
Ehrlichkeit gehört zum Handwerk (Kapitel 12). Dieses Programm hat bewusste Grenzen, und jede davon ist die Tür zu einem größeren Werkzeug. Eine dieser Grenzen wiegt im Brandfall schwerer als alle anderen und gehört deshalb an die Spitze der Tafel.
Abbildung 15.12 sammelt die übrigen Grenzen als Landkarte des Weiterlernens.
Code
import matplotlib.pyplot as plt
tueren = [
("2D → 3D", "Tetraeder statt Dreiecke,\ndrei Verschiebungen je Knoten"),
("gerade → krumme Ränder", "isoparametrische Elemente,\nJacobi & Gauß-Quadratur (Kap. 11)"),
("linear → nichtlinear", "E, λ von Temperatur abhängig;\nPlastizität, großes Verformen"),
("Einweg → Rückkopplung", "Verformung ändert Wärme;\niterativ gekoppelt (Kap. 14)"),
("fest → adaptiv", "Netz dort verfeinern,\nwo der Fehler groß ist (Kap. 12)"),
("Gauß-Seidel → schnell", "Mehrgitter, direkte Bandlöser,\nvorkonditionierte Verfahren (Kap. 4)"),
]
fig, ax = plt.subplots(figsize=(10.0, 4.2))
for k, (titel, text) in enumerate(tueren):
zeile = k // 3
spalte = k % 3
x0 = spalte * 3.3
y0 = (1 - zeile) * 2.0
ax.add_patch(plt.Rectangle((x0, y0), 3.0, 1.6, facecolor="#eef1f5",
edgecolor="black", lw=1.2))
ax.text(x0 + 1.5, y0 + 1.2, titel, ha="center", va="center",
fontsize=10, weight="bold")
ax.text(x0 + 1.5, y0 + 0.5, text, ha="center", va="center", fontsize=8)
ax.set_xlim(-0.2, 9.8)
ax.set_ylim(-0.3, 3.9)
ax.axis("off")
plt.tight_layout()
plt.show()
Was man sieht: sechs Türen, und hinter jeder wartet ein Stück Profi-FEM. Der Reihe nach:
- Von 2D zu 3D. Unser Träger lebt in der Ebene. Ein echtes Bauteil ist räumlich — dann treten an die Stelle der Dreiecke Tetraeder, jeder Knoten trägt drei Verschiebungen, und die Matrizen wachsen. Die Idee bleibt Wort für Wort dieselbe; es wird nur größer.
- Von geraden zu krummen Rändern. Unsere Rechtecke müssen achsparallel sein (Kapitel 11). Für runde Löcher und geschwungene Formen braucht man isoparametrische Elemente mit Jacobi-Matrix und Gauß-Quadratur — der Ausblick, den Kapitel 11 bewusst offenließ.
- Von linearem zu nichtlinearem Material. Wir rechnen mit festen Kennwerten (Kapitel 13). In Wahrheit ändern sich \(E\) und \(\lambda\) mit der Temperatur — die Erweichung aus dem Kasten oben —, und jenseits der Streckgrenze fließt das Metall (EN 1993-1-2, (CEN (Europäisches Komitee für Normung) 2005)). Dann wird aus einem Gleichungssystem eine Folge von Gleichungssystemen.
- Von der Einweg- zur echten Rückkopplung. Bei uns treibt die Wärme die Verformung, nicht umgekehrt (Kapitel 14). Echte Kopplung lässt auch die Verformung auf die Wärme zurückwirken — man rechnet dann iterativ hin und her, bis sich beide einig sind.
- Von festen zu adaptiven Netzen. Wir setzen das Netz von Hand (Kapitel 12). Profi-Programme verfeinern selbst dort, wo der Fehler groß ist — etwa an der Einspann-Ecke oder am Loch, wo unsere Spannung in Reihe A nicht konvergieren wollte.
- Von Gauß-Seidel zu schnellen Lösern. Unser Entspannungsverfahren
- ist ehrlich, aber langsam. Für Millionen Unbekannte nimmt man Mehrgitter- und vorkonditionierte Verfahren, die das Einpendeln um Größenordnungen beschleunigen.
Schlusswort
Am Anfang stand eine Blackbox mit drei Reglern und ein Versprechen. Jetzt steht dasselbe Programm offen vor dir, und du weißt, was jede Zeile tut.
Du kannst ein reales Ding — einen Träger, eine Pfanne, eine Brücke — in ein Netz aus Elementen zerlegen. Du kannst für jedes Element die Physik in eine kleine Matrix gießen und die vielen kleinen Matrizen zu einem großen Gleichungssystem assemblieren. Du kannst die Randbedingungen richtig anbringen — fest, beflammt, isoliert, luftgekühlt. Du kannst das System lösen, in der Zeit voranschreiten und die Wärme an die Verformung koppeln. Und, am wichtigsten: Du kannst dem Ergebnis misstrauen, bis es die Checkliste bestanden hat. Das ist die ganze Finite-Elemente-Methode. Es gibt nichts Geheimes mehr.
Die schönste Art, dieses Wissen zu behalten, ist, es umzubauen. Nimm fem_balken.py und mach es zu deinem Programm: ein anderes Bauteil, ein anderes Material, eine andere Frage aus deinem Alltag. Die drei Abschlussprojekte unten sind Startpunkte — aber der beste Umbau ist der, den du selbst erfindest.
Roter Faden
Zurück: überall — das war der Zweck dieses Kapitels. Jedes Modul trug seine Kapitelnummer: das Netz aus Kapitel 10 und Kapitel 11, die Zeitschleife aus Kapitel 9, die Ränder aus Kapitel 8, der Löser aus 1, die Kopplung aus Kapitel 13 und Kapitel 14, der Wachhund und die Checkliste aus Kapitel 12. Der Bimetall-Cliffhanger aus Übung 13.3 wird in Übung 15.3 endgültig eingelöst, und der Robin-Ausblick aus Kapitel 8 ist in Reihe C ausgerechnet.
Vor: hinaus aus dem Buch. Die Grenzen-Tafel zeigt sechs Türen, und OpenGeoSys wartet dahinter. Von hier an ist die Landkarte deine.
Übungen
Ü 15.1 (Verstehen). Führe zwei der fünf Versuchsreihen selbst durch und schreibe für jede ein Protokoll nach dem Rahmen aus Abbildung 15.3 (Frage → Vorhersage → Regler → Ergebnis → Schluss). Wähle etwa Reihe A (Netz) und Reihe C (Ränder). Wichtig ist der zweite Schritt: Leg deine Vorhersage schriftlich fest, bevor du rechnest. (Dies ist die ausformulierte Fassung von Übung 1.2 — jetzt mit dem offenen Programm statt der Blackbox.)
Ein Beispiel-Protokoll für Reihe A: Frage — reicht 16×4 für die Auslenkung? Vorhersage — grobe Netze biegen zu steif (Kapitel 14), also wächst die Zahl mit feinerem Netz und pendelt sich unter der Handformel \(12\ \mathrm{mm}\) ein. Regler — nur NX, NY (Elementtyp, Material, Zeit fest). Ergebnis — \(1{,}2 \to 3{,}7 \to
7{,}7 \to 9{,}6\ \mathrm{mm}\) für 4×1, 8×2, 16×4, 24×6. Schluss — Vorhersage bestätigt; 16×4 trifft rund 64 %, 24×6 rund 80 % der Handformel. Die Auslenkung ist noch nicht ganz konvergiert, für einen qualitativen Eindruck aber gut. Der Schluss hält der Checkliste stand (Konvergenz geprüft, Nulltest wäre der gleichmäßig warme freie Träger mit Spannung null).
Ü 15.2 (Verändern). Schneide ein Loch in den Träger. Das Programm bringt dafür die Funktion netz_mit_loch mit (ab Zeile 621 in fem_balken.py), die ein rechteckiges Loch aus dem Dreiecksnetz stanzt. Belaste den Träger auf gleichmäßige Erwärmung (\(\Delta T = 100\
\mathrm{K}\)) und halte ihn dabei längs fest (beide_enden_laengs_fest), sodass überall Druckspannung herrscht. Wo wird die Spannung am größten? Vergleiche mit dem vollen Träger ohne Loch.
Im vollen Träger steht überall dieselbe Längsspannung \(E\,\alpha\,\Delta T = 210\,000 \cdot 12\cdot 10^{-6}\cdot 100 = 252\ \mathrm{N/mm^2}\) Druck (Kapitel 13) — schon das liegt über der S235-Streckgrenze, der voll eingespannte Träger fließt also bereits. Mit Loch drängen sich die Kraftlinien um das Loch herum und rücken an seinen Ecken eng zusammen: An den Loch-Ecken steigt die Spannung auf rund \(468\ \mathrm{N/mm^2}\) — ein Konzentrationsfaktor von etwa \(1{,}9\) auf diesem groben Netz (ein feineres Netz und ein rundes Loch zeigen den bekannten Faktor \(\approx 3\)). Das ist die Kerbwirkung.
Die Kerbwirkung ist keine Rechenspielerei. Im Zweiten Weltkrieg baute die USA rund \(2700\) Liberty-Frachter in Rekordzeit — durchgehend geschweißt statt genietet, mit rechteckigen Ladeluken im Deck. Etliche brachen im kalten Nordatlantik plötzlich auseinander, manche mitten entzwei. Die Rissanalyse zeigte: Fast alle Brüche begannen an den scharfen Ecken der Luken, wo sich die Spannung konzentrierte — verschlimmert durch Schweißnähte genau in der Ecke und durch den bei Kälte spröden Stahl (Sridhar 2016). Eine runde Ecke hätte die Spitze entschärft. Seither werden Luken, Fenster und Ausschnitte gerundet — dieselbe Physik, die dein Loch-Experiment zeigt, nur mit tödlichem Ausgang. Genau deshalb warnt Kapitel 12: Ein grobes Netz unterschätzt solche Spitzen — man muss gezielt verfeinern, um sie zu sehen.
Ü 15.3 (Übertragen — Abschlussprojekt, eines wählen). Bau das Programm zu deinem um. Jedes Projekt braucht ein Verifikations-Protokoll: Womit prüfst du, dass das Ergebnis stimmt (Nulltest, Handformel, Größenordnung)?
(a) Der echte Bimetall-Streifen. Schichte zwei Materialzonen übereinander — untere Hälfte Aluminium, obere Hälfte Stahl — und erwärme gleichmäßig. Jetzt krümmt sich der Streifen wirklich, denn die beiden Metalle liegen nebeneinander und halten einander fest (das löst den Cliffhanger aus Übung 13.3 ein). Prüfe die Durchbiegung gegen die Bimetall-Handformel von Timoshenko. Lösungsskript: loesungen/kap15_ue3_bimetall.py.
(b) Die Kühlrippe. Ersetze den Träger durch eine schmale, hohe Rippe, die unten an einem heißen Körper sitzt und oben in die Luft ragt: unten ein Temperatur-Rand, ringsum Robin-Luftkühlung. Wie weit oben ist die Rippe noch warm? (Verifikation: die Wärmestrombilanz — was unten hereinkommt, muss seitlich als Kühlung wieder hinaus.)
(c) Die zweite Flamme von unten. Setze eine zweite Flamme an die Unterkante des Trägers, versetzt zur ersten. Wo treffen sich die beiden Wärmefronten, und wie verbiegt sich der Träger, wenn oben und unten verschieden geheizt wird? (Verifikation: Schalte beide Flammen gleich heiß — dann muss das Höhengefälle und damit die Biegung verschwinden.)
Das Skript loesungen/kap15_ue3_bimetall.py weist jedem Dreieck sein eigenes Material zu (unten Alu, oben Stahl) und assembliert damit — sonst ist es die Rechenkette aus Kapitel 14. Ergebnis: Der Streifen krümmt sich, die Alu-Seite außen (sie will sich stärker dehnen). Die Timoshenko-Handformel liefert \(\kappa = 24\,(\alpha_{\text{oben}} -
\alpha_{\text{unten}})\,\Delta T / [H\,(14 + n + 1/n)]\) mit \(n =
E_{\text{unten}}/E_{\text{oben}}\), also eine Durchbiegung von \(7{,}6\
\mathrm{mm}\). Die FEM nähert sich ihr von unten (\(1{,}5 \to 4{,}2 \to
5{,}7 \to 6{,}5\ \mathrm{mm}\) für 8×2 bis 32×8) — wieder die CST-Biegesteifheit aus Kapitel 14. Der Bimetall-Cliffhanger aus Übung 13.3 ist damit endgültig eingelöst: Nebeneinander liegende, verschieden dehnende Fasern krümmen den Verbund und erzeugen Spannung.
Das Kleingedruckte
Auch dieses Finale bleibt bei den bewussten Vereinfachungen des Buches (Strategie §7), und es benennt sie ehrlich: Die Mechanik rechnet nur auf Dreiecken (P1) — löst die Wärme auf Rechtecken, triangulieren wir dasselbe Gitter für die Kopplung (die Knoten sind ja identisch). Die Kopplung ist Einweg (Wärme → Verformung, nicht zurück). Das Material ist linear und isotrop, seine Kennwerte temperaturunabhängig — die folgenreichste Grenze im Brandfall, wie der Kasten an der Spitze der Grenzen-Tafel zeigt. Die Netze sind klein genug fürs reine Python im Browser. Und die Brandgase koppeln über einen Robin-Übergang mit dem Normwert \(\alpha_c = 25\ \mathrm{W/(m^2 K)}\): realistisch, aber schwach — ein Konvektionsfeuer wölbt den Träger in einer halben Stunde nur um Bruchteile eines Millimeters; die kräftige, sichtbare Krümmung entsteht erst im voll entwickelten Brand, wenn die Unterseite selbst richtig heiß wird. Jede dieser Grenzen steht in der Tafel des Weiterlernens, und jede ist eine Einladung, weiterzugehen.