Fast solvers for Tokamak fluid models with PETSC

Dieser Beitrag stellt einen neuartigen semi-vergröbernden geometrischen Mehrgitterlöser vor, der in PETSc für den Tokamak-Code M3D-C1 implementiert ist und die Konvergenzbeschränkungen des bestehenden Block-Jacobi-Vorkonditionierers durch Ausnutzung der toroidalen Gitterstruktur überwindet, um eine überlegene Robustheit und Leistung bei komplexen magnetohydrodynamischen Modellen zu erzielen.

Das Wesentliche

Stellen Sie sich vor, Sie versuchen vorherzusagen, wie sich eine wirbelnde, superheiße Suppe aus geladenen Teilchen (Plasma) in einer riesigen, donutförmigen Maschine namens Tokamak verhält. Diese Maschine ist darauf ausgelegt, Fusionsenergie zu erzeugen, ähnlich der Kraft der Sonne. Diese „Suppe" ist jedoch unglaublich chaotisch. Wenn Sie versuchen, ihre Bewegung schrittweise mit einem Computer zu berechnen, wird die Mathematik so kompliziert, dass der Computer stecken bleibt oder so lange braucht, dass die Antwort bis zu ihrem Eintreffen nutzlos ist.

Dieser Artikel handelt vom Bau eines intelligenteren, schnelleren Rechners für genau diese Art von Problem.

Hier ist die Aufschlüsselung dessen, was die Autoren getan haben, unter Verwendung einfacher Analogien:

1. Das Problem: Der „Stau" der Mathematik

Der von ihnen verwendete Computercode (genannt M3D-C1) versucht, Gleichungen zu lösen, die beschreiben, wie sich das Plasma bewegt. Um dies zu tun, muss er ein riesiges Puzzle Millionen Mal lösen.

• Der alte Weg (Block-Jacobi): Stellen Sie sich vor, Sie haben eine riesige Stadtkarte mit Staus. Die alte Methode war so, als würde man eine andere Person bitten, den Verkehr nur auf einer Straße zu beheben, während der Rest der Stadt ignoriert wird. Wenn die Stadt klein ist, funktioniert dies. Aber je größer die Stadt wird (mehr „Ebenen" oder Scheiben der Donutform), desto weniger können die Personen, die die Straßen reparieren, schnell miteinander kommunizieren. Die Staus werden schlimmer, und die Lösung verlangsamt sich oder funktioniert gar nicht mehr.
• Die spezifische Herausforderung: Das Plasma in diesen Maschinen ist „anisotrop". Stellen Sie sich einen Stapel Papier vor. Es ist sehr einfach, ein Blatt Papier über die Oberfläche zu schieben (leichte Richtung), aber sehr schwer, es durch den Stapel zu drücken (schwere Richtung). Der alte mathematische Löser verstand diese „Papierstapel"-Struktur nicht und versuchte daher, die schwere und die leichte Richtung mit derselben ungeschickten Methode zu lösen.

2. Die Lösung: Der „Multigitter"-Aufzug

Die Autoren bauten einen neuen Löser mit einer Methode namens Multigitter (MG).

• Die Analogie: Stellen Sie sich vor, Sie versuchen, ein verlorenes Spielzeug in einem riesigen, mehrstöckigen Herrenhaus zu finden.
  • Der alte Weg: Sie prüfen jeden einzelnen Raum, jede einzelne Schublade und jede einzelne Ecke im Erdgeschoss, bevor Sie nach oben gehen. Das dauert ewig.
  • Der Multigitter-Weg: Sie schauen sich zuerst ein Miniaturmodell des gesamten Herrenhauses aus der Vogelperspektive an. Sie erkennen schnell den allgemeinen Bereich, in dem das Spielzeug fehlt (das „grobe" Gitter). Dann zoomen Sie auf eine mittlere Karte, um den Bereich einzugrenzen. Schließlich gehen Sie in den tatsächlichen Raum (das „feine" Gitter), um das Spielzeug aufzuheben.
  • Indem das Problem zunächst auf den „Gesamt"-Ebenen gelöst wird, weiß der Löser genau, wo er suchen muss, wenn er zu den winzigen Details hinabsteigt. Dies macht ihn unglaublich schnell.

3. Das „Geheimrezept": Semi-Vergröberung

Die Autoren erkannten, dass der Tokamak wie ein Stapel 2D-Scheiben (poloidale Ebenen) ist, die zu einem 3D-Donut verdreht sind.

• Sie wandten ihren „Multigitter-Aufzug" spezifisch auf die Stapelrichtung (die toroidale Richtung) an.
• Anstatt zu versuchen, das gesamte 3D-Chaos auf einmal zu vereinfachen, behielten sie die 2D-Scheiben detailliert (da die Form der Tokamak-Wand komplex ist), machten aber den Stapel der Scheiben einfacher, je höher sie stiegen.
• Dies ist wie bei einem dicken Buch, bei dem man nur die Anzahl der Seiten reduziert, während der Text auf jeder Seite klar bleibt. Es ist eine perfekte Passform für die Form der Maschine.

4. Die Ergebnisse: Geschwindigkeit und Zuverlässigkeit

Das Team testete diesen neuen Löser an zwei sehr schwierigen Szenarien:

• Szenario A: Der „Runaway-Elektron" (SPARC): Dies simuliert ein gefährliches Ereignis, bei dem sich Teilchen unkontrolliert beschleunigen.
  • Ergebnis: Der neue Löser war bei kleineren Aufbauten wettbewerbsfähig mit dem alten und bei den größten, komplexesten Aufbauten viel schneller. Er löste das Problem in weniger Schritten und sparte Zeit.
• Szenario B: Der „Stellator" (eine andere, stärker verdrehte Maschine): Diese Geometrie ist noch stärker verdreht und unregelmäßiger als ein Standard-Donut.
  • Ergebnis: Der alte Löser versagte vollständig und konnte keine Antwort finden. Der neue Multigitter-Löser gelangte zum Erfolg. Er war robust genug, um die verdrehte Geometrie zu bewältigen, die das alte Werkzeug zum Scheitern brachte.

5. Die Hardware: Nutzung von Supercomputern

Sie führten diese Tests auf Perlmutter durch, einem der schnellsten Supercomputer der Welt, der sowohl leistungsstarke CPUs als auch GPUs (Grafikkarten) verwendet.

• Sie stellten fest, dass zwar das „Setup" (der Bau der Miniaturmodelle) teuer war, die eigentliche Lösung auf den GPUs jedoch unglaublich schnell erfolgte.
• Sie entdeckten, dass sie für die schwierigsten Probleme einen „schweren" Glätter (ein spezieller mathematischer Trick) verwenden mussten, um zu verhindern, dass der Löser stecken bleibt, was etwas mehr Rechenleistung erforderte, sich aber in der Geschwindigkeit auszahlte.

Zusammenfassung

Der Artikel behauptet, dass sie durch das Verständnis der spezifischen „Papierstapel"-Form von Fusionsplasma-Maschinen ein neues mathematisches Werkzeug (Multigitter) geschaffen haben, das:

1. Probleme schneller löst als die aktuelle Standardmethode bei großen, komplexen Simulationen.
2. Nicht abstürzt bei verdrehten, komplexen Formen, bei denen die alte Methode versagt.
3. Ein entscheidender erster Schritt ist, um Fusionsenergie-Simulationen praktisch und schnell genug zu machen, um beim Entwurf echter Kraftwerke zu helfen.

Sie behaupteten nicht, dass dies die Fusionsenergie selbst löst, sondern dass es den schnellen, zuverlässigen Rechner bereitstellt, der benötigt wird, um die Physik zu simulieren, die schließlich zu Fusionsenergie führen wird.

Technische Zusammenfassung

Technischer Zusammenfassung: Schnelle Löser für Tokamak-Fluidmodelle mit PETSc

Problemstellung
Der M3D-C1-Code ist ein Hauptwerkzeug zur Simulation von erweiterten magnetohydrodynamischen (MHD) Plasmen in der Fusionsenergieforschung, das sich speziell mit großskaligen Instabilitäten und Störungen in Tokamaks und Stellaratoren befasst. Der Code verwendet einen semi-impliziten Zeitintegrator, der die Lösung zahlreicher linearer Systeme pro Zeitschritt erfordert. Unter diesen wird der implizite Fortschritt der Impulsgleichung als die herausforderndste und am wenigsten robuste Komponente identifiziert.

Der aktuelle Produktionslöser verwendet eine Krylov-Methode (FGMRES/GMRES), die durch ein Block-Jacobi-Verfahren (BJ) vorkonditioniert wird, wobei Blöcke Freiheitsgrade auf Ebenen mit konstantem toroidalem Koordinatenwert gruppieren. Während der BJ-Löser für kleinere Konfigurationen effektiv ist, verschlechtert sich seine Konvergenz mit zunehmender Anzahl toroidaler Ebenen aufgrund der spektralen Eigenschaften der vorkonditionierten Matrix erheblich. Darüber hinaus versagt der BJ-Löser bei nicht-achsensymmetrischen Stellarator-Konfigurationen vollständig, was die Anwendbarkeit des Codes auf komplexe Fusionsgeräte einschränkt.

Methodik
Diese Arbeit entwickelt einen geometrischen Mehrgitterlöser (MG), der auf die spezifische Diskretisierung und Gitterstruktur von M3D-C1 zugeschnitten ist. Die Methodik nutzt folgende strukturelle Merkmale:

• Gitterstruktur: M3D-C1 verwendet unstrukturierte 2D-Gitter in der poloidalen Ebene (um den Torus extrudiert), kombiniert mit einem regulären 1D-Gitter in der toroidalen Richtung.
• Semi-Vergröberung: Der Löser verwendet eine geometrische Semi-Vergröberung speziell in der toroidalen Richtung. Dieser Ansatz wird gewählt, weil das toroidale Gitter strukturiert ist, was Standard-geometrische MG-Techniken (2:1-Vergröberung) ermöglicht, während die komplexe poloidale Ebene für zukünftige Arbeiten unvergröbert bleibt.
• Algebraische Schnittstelle: Obwohl die Anwendung Jacobimatrizen auf feinen Gittern berechnet, werden die Operatoren des groben Gitters algebraisch unter Verwendung eines Galerkin-Prozesses konstruiert ($A_{i+1} = P^T A_i P$), wobei $P$ ein 3-Punkt-Stencil-Verlängerungsoperator ist. Dies erzeugt einen geometrisch/algebraischen Hybridlöser mit einer rein algebraischen Schnittstelle, die mit PETSc kompatibel ist.
• Glätter: Zwei Glätter werden untersucht:
  1. Punkt-Block-Jacobi (PBJ): Berechnet explizite dichte Inverse der Diagonalblöcke (Größe 36 pro Vertex). Dies ist schnell und speichereffizient und eignet sich für GPU-Beschleunigung.
  2. Block-Jacobi (BJ): Der bestehende Produktionslöser wird als Glätter wiederverwendet. Dies beinhaltet vollständige LU-Zerlegungen auf poloidalen Ebenen. Obwohl dies teuer ist, wird es als Nachglätter in einem $V(0,1)$-Zyklus für die schwierigsten Probleme verwendet, bei denen PBJ stagniert.
• Implementierung: Der Löser ist innerhalb der PETSc-Bibliothek implementiert und nutzt deren Krylov-Löser sowie abstrakte Schnittstellen für lineare Algebra-Back-Ends (einschließlich MUMPS und SuperLU für direkte Lösungen).

Hauptergebnisse
Der Löser wurde auf zwei Testfällen unter Verwendung des Perlmutter-Supercomputers (CPU- und GPU-Knoten) evaluiert: einem Runaway-Elektronen-Modell einer SPARC-Tokamak-Störung und einem quasi-achsensymmetrischen Stellarator-Modell.

1. SPARC-Tokamak (Skalierbarkeit und Geschwindigkeit):

   • Bei Konfigurationen mit 4 bis 16 toroidalen Ebenen konvergierte der MG-Löser (unter Verwendung von PBJ-Glättung) bei erneuter Aktualisierung des Vorkonditionierers in einer einzigen Iteration und zeigte eine hervorragende algorithmische Skalierbarkeit.
   • Bei der größten Konfiguration (64 Ebenen) nahm die Problemkomplexität zu, was eine strengere Lösertoleranz ($rtol = 10^{-7}$) und einen aggressiveren $V(6,6)$-Zyklus mit häufigen Aktualisierungen erforderte.
   • Leistung: Der MG-Löser erzielte eine wettbewerbsfähige Leistung im Vergleich zum BJ-Löser. Im 64-Ebenen-Fall reduzierte der MG-Löser die Gesamtzahl der Krylov-Iterationen um einen Faktor von etwa 16 im Vergleich zum BJ-Löser. Während die Gesamt-Lösezeit auf GPUs ungefähr vergleichbar war, zeigte der MG-Löser auf CPUs signifikant schnellere Lösezeiten.
   • Komplexität: Theoretische Analysen deuten darauf hin, dass der MG-Löser eine asymptotische Arbeitskomplexität von etwa dem 2-fachen des BJ-Lösers aufweist, bedingt durch die geometrische Summe der Arbeit über die Gitterebenen hinweg. Die beobachtete Leistung zeigt, dass die Lösungen auf groben Gittern derzeit schlecht performen, was Optimierungspotenzial aufzeigt.

2. Stellarator (Robustheit):

   • Bei einem nicht-achsensymmetrischen Stellarator-Testfall versagte der Standard-Block-Jacobi-Löser vollständig bei der Konvergenz.
   • Im Gegensatz dazu konvergierte der MG-Löser (unter Verwendung von BJ als Glätter) erfolgreich. Dies demonstriert einen kritischen Robustheitsvorteil des MG-Ansatzes für nicht-achsensymmetrische Geometrien, bei denen die spektralen Eigenschaften der Matrix den BJ-Vorkonditionierer ausschalten.

Bedeutung und Behauptungen
Die Arbeit behauptet, den ersten Schritt zur Anwendung von Mehrgitterverfahren auf komplexe, wissenschaftlich relevante MHD-Tokamak-Modelle innerhalb des M3D-C1-Codes darzustellen. Die Hauptbeiträge sind:

• Nachweis der Machbarkeit: Der Beweis, dass Mehrgitterverfahren einen Weg zu schnelleren, skalierbaren und robusteren Lösern für hochgradig anisotrope, schlecht konditionierte magnetisierte Plasmamodelle bieten können.
• Toroidale Semi-Vergröberung: Die Feststellung, dass eine Semi-Vergröberung in der toroidalen Richtung ein nützlicher und effektiver erster Schritt für gängige Tokamak-Codes mit extrudierten Gitterstrukturen ist.
• Robustheit: Der Nachweis, dass MG-Löser Probleme lösen können (insbesondere Stellarator-Störungen), bei denen der aktuelle Produktions-Block-Jacobi-Löser versagt.

Zukünftige Arbeiten und Einschränkungen
Die Autoren bleiben bezüglich des aktuellen Zustands des Lösers bescheiden und weisen auf mehrere Verbesserungsbereiche hin:

• Leistung auf groben Gittern: Die aktuelle Leistung wird durch ineffiziente Lösungen auf groben Gittern begrenzt, teilweise aufgrund redundanter Parallelisierungsanforderungen und MPI-Fehlschläge bei direkten Lösern (MUMPS/SuperLU) auf groben Gittern.
• Vollständiges 3D-Mehrgitter: Die wichtigste zukünftige Optimierung ist die Hinzufügung von 2D-Mehrgitter in der poloidalen Ebene, um einen vollständigen 3D-MG-Löser zu erstellen. Dies würde die Handhabung unstrukturierter geometrischer grober Gitter erfordern.
• Glätter: Zukünftige Arbeiten umfassen die Optimierung der Blockgrößen für PBJ-Glätter, um besser zu GPU-Thread-Gruppen zu passen, sowie die Untersuchung von FieldSplit-Vorkonditionierern, die die drei entkoppelten MHD-Wellen separat behandeln.
• Aktualisierungsstrategien: Das aktuelle Aktualisierungsmodell (periodische LU-Zerlegung) ist einfach; intelligentere Algorithmen sind erforderlich, um quasi-stationäre Zustände effizienter zu handhaben.

Zusammenfassend etabliert diese Arbeit einen viable, robusten Mehrgitter-Rahmen für M3D-C1, der den bestehenden Löser in Bezug auf Iterationszahlen und Robustheit bei schwierigen Geometrien übertrifft, und identifiziert spezifische algorithmische und Implementierungshürden, die angegangen werden müssen, um eine volle textbook-Mehrgittereffizienz zu erreichen.