Fast solvers for Tokamak fluid models with PETSC

Questo articolo presenta un nuovo risolutore geometrico multigriglia a semi-avvolgimento implementato in PETSc per il codice tokamak M3D-C1, che affronta i limiti di convergenza del precondizionatore Jacobi a blocchi esistente sfruttando la struttura della griglia toroidale per ottenere una robustezza e prestazioni superiori su modelli complessi di magnetoidrodinamica.

L'essenziale

Immagina di cercare di prevedere il comportamento di una zuppa supercalda e vorticosa di particelle cariche (plasma) all'interno di una macchina gigantesca a forma di ciambella chiamata tokamak. Questa macchina è progettata per creare energia da fusione, simile alla potenza del sole. Tuttavia, questa "zuppa" è incredibilmente caotica. Se provi a calcolare il suo movimento passo dopo passo utilizzando un computer, la matematica diventa così complicata che il computer si blocca, o impiega così tanto tempo che la risposta è inutile al momento in cui arriva.

Questo articolo riguarda la costruzione di una calcolatrice più intelligente e veloce per questo specifico tipo di problema.

Ecco la spiegazione di ciò che gli autori hanno fatto, utilizzando semplici analogie:

1. Il Problema: Il "Traffico" della Matematica

Il codice informatico che utilizzano (chiamato M3D-C1) cerca di risolvere equazioni che descrivono come si muove il plasma. Per fare questo, deve risolvere un enorme puzzle milioni di volte.

• Il Vecchio Metodo (Block Jacobi): Immagina di avere una vasta mappa di una città con ingorghi stradali. Il vecchio metodo era come chiedere a una persona diversa di risolvere il traffico su una sola strada alla volta, ignorando il resto della città. Se la città è piccola, questo funziona. Ma man mano che la città diventa più grande (più "piani" o fette della forma a ciambella), le persone che risolvono i problemi stradali non riescono a comunicare abbastanza velocemente tra loro. Gli ingorghi peggiorano e la soluzione rallenta o smette di funzionare completamente.
• La Sfida Specifica: Il plasma in queste macchine è "anisotropo". Pensalo come un mazzo di fogli di carta. È molto facile far scivolare un foglio lungo la superficie (direzione facile), ma molto difficile spingerlo attraverso il mazzo (direzione difficile). Il vecchio risolutore matematico non capiva questa struttura a "mazzo di fogli", quindi cercava di risolvere la direzione difficile e quella facile con lo stesso metodo goffo.

2. La Soluzione: L'Ascensore "Multigrid"

Gli autori hanno costruito un nuovo risolutore utilizzando un metodo chiamato Multigrid (MG).

• L'Analogia: Immagina di cercare un giocattolo smarrito in un'enorme villa a più piani.
  • Il Vecchio Metodo: Controlli ogni singola stanza, ogni singolo cassetto e ogni singolo angolo al piano terra prima di salire. Ci vuole un'eternità.
  • Il Metodo Multigrid: Prima guardi un modello in miniatura dell'intera villa da una vista a volo d'uccello. Individui rapidamente l'area generale dove manca il giocattolo (la "griglia" grossolana). Poi, ingrandisci su una mappa di dimensioni medie per restringere il campo. Infine, vai nella stanza reale (la "griglia" fine) per prendere il giocattolo.
  • Risolvendo il problema prima sui livelli della "grande immagine", il risolutore sa esattamente dove guardare quando scende ai dettagli minuscoli. Questo lo rende incredibilmente veloce.

3. La "Salsa Segreta": Semi-Coarsening

Gli autori hanno realizzato che il tokamak è come un mazzo di fette 2D (piani poloidali) attorcigliate in una ciambella 3D.

• Hanno applicato il loro "Ascensore Multigrid" specificamente alla direzione di impilamento (la direzione toroidale).
• Invece di cercare di semplificare tutto il caos 3D in una volta sola, hanno mantenuto le fette 2D dettagliate (perché la forma della parete del tokamak è complessa) ma hanno reso il mazzo di fette più semplice man mano che salivano di livello.
• È come prendere un libro spesso e ridurre solo il numero di pagine, mantenendo il testo su ogni pagina chiaro. È una soluzione perfetta per la forma della macchina.

4. I Risultati: Velocità e Affidabilità

Il team ha testato questo nuovo risolutore su due scenari molto difficili:

• Scenario A: L'"Elettrone Fuggitivo" (SPARC): Questo simula un evento pericoloso in cui le particelle accelerano in modo incontrollabile.
  • Risultato: Il nuovo risolutore è stato competitivo con il vecchio su configurazioni più piccole e molto più veloce sulle configurazioni più grandi e complesse. Ha risolto il problema in meno passaggi, risparmiando tempo.
• Scenario B: Lo "Stellarator" (Una macchina diversa, più attorcigliata): Questa geometria è ancora più attorcigliata e irregolare di una ciambella standard.
  • Risultato: Il vecchio risolutore ha fallito completamente e non è riuscito a trovare una risposta. Il nuovo risolutore Multigrid ha avuto successo. È stato abbastanza robusto da gestire la geometria attorcigliata che aveva fatto fallire il vecchio strumento.

5. L'Hardware: Utilizzo dei Supercomputer

Hanno eseguito questi test su Perlmutter, uno dei supercomputer più veloci al mondo, che utilizza sia potenti CPU che GPU (schede grafiche).

• Hanno scoperto che mentre la "configurazione" (costruire i modelli in miniatura) era costosa, la risoluzione effettiva era incredibilmente veloce sulle GPU.
• Hanno scoperto che per i problemi più difficili, dovevano utilizzare un "smussatore" (un trucco matematico specifico) "pesante" per evitare che il risolutore si bloccasse, il quale richiedeva un po' più di potenza di calcolo ma ripagava in velocità.

Riassunto

L'articolo afferma che, comprendendo la specifica forma a "mazzo di fogli" delle macchine a fusione al plasma, hanno creato un nuovo strumento matematico (Multigrid) che:

1. Risolve i problemi più velocemente del metodo standard attuale su simulazioni grandi e complesse.
2. Non si blocca su forme attorcigliate e complesse dove il vecchio metodo fallisce.
3. È un primo passo cruciale verso la realizzazione di simulazioni di energia da fusione pratiche e abbastanza veloci da aiutare a progettare vere centrali elettriche.

Non hanno affermato che questo risolva l'energia da fusione in sé, ma piuttosto che fornisce la calcolatrice veloce e affidabile necessaria per simulare la fisica che porterà eventualmente all'energia da fusione.

Sommario tecnico

Riepilogo Tecnico: Solutori Veloci per Modelli Fluidi dei Tokamak con PETSc

Enunciato del Problema
Il codice M3D-C1 è uno strumento primario per la simulazione di plasmi magnetoidrodinamici (MHD) estesi nella ricerca sull'energia da fusione, affrontando specificamente instabilità su larga scala e interruzioni (disruptions) nei tokamak e negli stellarator. Il codice impiega un integratore temporale semi-implicito che richiede la risoluzione di numerosi sistemi lineari per ogni passo temporale. Tra questi, l'avanzamento implicito dell'equazione della quantità di moto è identificato come il componente più impegnativo e meno robusto.

Il solver di produzione attuale utilizza un metodo di Krylov (FGMRES/GMRES) precondizionato da un metodo di Jacobi a Blocchi (BJ), dove i blocchi raggruppano i gradi di libertà su piani di coordinata toroidale costante. Sebbene efficace per configurazioni più piccole, la convergenza del solver BJ degrada significativamente all'aumentare del numero di piani toroidali a causa delle proprietà spettrali della matrice precondizionata. Inoltre, il solver BJ non converge affatto su configurazioni di stellarator non assialsimmetrici, limitando l'applicabilità del codice a dispositivi di fusione complessi.

Metodologia
Questo lavoro sviluppa un solver multigrid (MG) geometrico adattato alla specifica discretizzazione e struttura della griglia di M3D-C1. La metodologia sfrutta le seguenti caratteristiche strutturali:

• Struttura della Griglia: M3D-C1 utilizza griglie 2D non strutturate nel piano poloidale (estruse attorno al toro) combinate con una griglia 1D regolare nella direzione toroidale.
• Semi-Raffinamento: Il solver impiega un semi-raffinamento geometrico specificamente nella direzione toroidale. Questo approccio è scelto perché la griglia toroidale è strutturata, permettendo l'uso di tecniche MG geometriche standard (raffinamento 2:1) lasciando il complesso piano poloidale non raffinato per lavori futuri.
• Interfaccia Algebrica: Sebbene l'applicazione calcoli le Jacobiane su griglie fini, gli operatori sulla griglia grossa sono costruiti algebricamente utilizzando un processo di Galerkin ($A_{i+1} = P^T A_i P$), dove $P$ è un operatore di prolungamento a 3 punti. Ciò crea un solver ibrido geometrico/algebrico con un'interfaccia puramente algebrica compatibile con PETSc.
• Smussatori: Vengono investigati due smussatori:
  1. Jacobi a Blocchi Puntuali (PBJ): Calcola inverse dense esplicite dei blocchi diagonali (dimensione 36 per vertice). È veloce ed efficiente in termini di memoria, adatto all'accelerazione GPU.
  2. Jacobi a Blocchi (BJ): Il solver di produzione esistente è riutilizzato come smussatore. Coinvolge fattorizzazioni LU complete sui piani poloidali. Sebbene costoso, è utilizzato come post-smussatore in un ciclo $V(0,1)$ per i problemi più difficili dove il PBJ stagna.
• Implementazione: Il solver è implementato all'interno della libreria PETSc, sfruttando i suoi solver di Krylov e le interfacce astratte per i back-end di algebra lineare (inclusi MUMPS e SuperLU per risoluzioni dirette).

Risultati Chiave
Il solver è stato valutato su due casi di test utilizzando il supercomputer Perlmutter (nodi CPU e GPU): un modello di elettroni runaway di un'interruzione del tokamak SPARC e un modello di stellarator quasi-assialsimmetrico.

1. Tokamak SPARC (Scalabilità e Velocità):

   • Su configurazioni con 4-16 piani toroidali, il solver MG (utilizzando lo smussamento PBJ) ha converguto in una singola iterazione quando il precondizionatore è stato aggiornato, dimostrando un'eccellente scalabilità algoritmica.
   • Sulla configurazione più grande (64 piani), la difficoltà del problema è aumentata, richiedendo una tolleranza del solver più stretta ($rtol = 10^{-7}$) e un ciclo $V(6,6)$ più aggressivo con aggiornamenti frequenti.
   • Prestazioni: Il solver MG ha raggiunto prestazioni competitive rispetto al solver BJ. Sul caso a 64 piani, il solver MG ha ridotto il numero totale di iterazioni di Krylov di un fattore di circa 16x rispetto al solver BJ. Mentre il tempo totale di risoluzione era grossomodo comparabile sulle GPU, il solver MG ha mostrato tempi di risoluzione significativamente più veloci sulle CPU.
   • Complessità: L'analisi teorica suggerisce che il solver MG ha una complessità di lavoro asintotica di circa il doppio rispetto al solver BJ a causa della somma geometrica del lavoro attraverso i livelli della griglia. Le prestazioni osservate indicano che le risoluzioni sulla griglia grossa stanno attualmente performando male, suggerendo spazio per ottimizzazioni.

2. Stellarator (Robustezza):

   • Su un caso di test di stellarator non assialsimmetrico, il solver standard Jacobi a Blocchi non è riuscito a convergere completamente.
   • Al contrario, il solver MG (utilizzando BJ come smussatore) ha converguto con successo. Ciò dimostra un vantaggio critico di robustezza dell'approccio MG per geometrie non assialsimmetriche dove le proprietà spettrali della matrice sconfiggono il precondizionatore BJ.

Significato e Affermazioni
Il documento afferma di presentare il primo passo nell'applicazione di metodi multigrid ai complessi modelli MHD dei tokamak scientificamente rilevanti all'interno del codice M3D-C1. I contributi principali sono:

• Dimostrazione di Fattibilità: Dimostrare che i metodi multigrid possono fornire una via verso solver più veloci, scalabili e robusti per modelli di plasma magnetizzato altamente anisotropi e mal condizionati.
• Semi-Raffinamento Toroidale: Stabilire che il semi-raffinamento nella direzione toroidale è un primo passo utile ed efficace per i comuni codici tokamak con strutture di griglia estruse.
• Robustezza: Mostrare che i solver MG possono risolvere problemi (specificamente interruzioni di stellarator) dove l'attuale solver di produzione Jacobi a Blocchi fallisce.

Lavori Futuri e Limitazioni
Gli autori rimangono modesti riguardo allo stato attuale del solver, notando diverse aree di miglioramento:

• Prestazioni della Griglia Grossa: Le prestazioni attuali sono limitate da risoluzioni inefficienti sulla griglia grossa, in parte a causa di requisiti di parallelismo ridondanti e fallimenti MPI nei solver diretti (MUMPS/SuperLU) sulle griglie grosse.
• Multigrid 3D Completo: L'ottimizzazione futura più impattante è l'aggiunta di un multigrid 2D nel piano poloidale per creare un solver MG 3D completo. Ciò richiederebbe la gestione di griglie geometriche grosse non strutturate.
• Smussatori: I lavori futuri includono l'ottimizzazione delle dimensioni dei blocchi per gli smussatori PBJ per adattarsi meglio ai gruppi di thread GPU e l'investigazione di precondizionatori FieldSplit che trattano separatamente le tre onde MHD disaccoppiate.
• Strategie di Aggiornamento: Il modello di aggiornamento corrente (fattorizzazione LU periodica) è semplice; sono necessari algoritmi più intelligenti per gestire gli stati quasi-stazionari in modo più efficiente.

In sintesi, questo lavoro stabilisce un framework multigrid vitale e robusto per M3D-C1 che supera il solver esistente in termini di conteggi di iterazioni e robustezza su geometrie difficili, identificando al contempo specifiche ostacoli algoritmici e di implementazione da affrontare per raggiungere l'efficienza multigrid completa da manuale.