Local cubic spline interpolation for Vlasov-type equations on a multi-patch geometry

Sintesi Tecnica: Interpolazione Spline Cubica Locale per Equazioni di tipo Vlasov su una Geometria Multi-patch

Enunciato del Problema
Il documento affronta la risoluzione numerica di equazioni di tipo Vlasov, specificamente modelli di centro di guida (guiding-center), su geometrie multi-patch complesse tipiche delle simulazioni tokamak. In queste simulazioni, le linee di mesh sono allineate con le superfici di flusso magnetico per minimizzare le variazioni lungo le linee di campo. Questa geometria crea sfide specifiche: il nucleo presenta un punto O dove le linee di mesh convergono, e il bordo presenta un punto X (separatrice) dove le linee divergono. Sebbene la suddivisione del dominio globale in sottodomini (patch) permetta il raffinamento locale e l'ottimizzazione di mesh uniformi all'interno di ogni patch, ciò introduce difficoltà nel mantenere la regolarità della soluzione alle interfacce tra i patch. Studi precedenti che utilizzavano la continuità $C^0$ (corrispondenza solo dei valori della funzione) alle interfacce sono stati giudicati potenzialmente instabili per schemi di advezione semi-Lagrangiani. Il problema centrale è come imporre una continuità di ordine superiore ($C^1$) e calcolare accuratamente le derivate d'interfaccia per interpolazioni spline cubiche B locali su mesh multi-patch arbitrarie (uniformi e non uniformi) senza sacrificare l'efficienza dei calcoli locali.

Metodologia
Gli autori propongono un metodo semi-Lagrangiano all'indietro (BSL) utilizzando spline B cubiche locali definite su ogni patch. Per garantire stabilità e accuratezza, il metodo impone la regolarità $C^1$ alle interfacce dei patch facendo corrispondere sia i valori della funzione che le loro prime derivate.

1. Ricostruzione della Derivata: L'innovazione centrale è una formula generalizzata per il calcolo delle derivate d'interfaccia. Basandosi su una formula di Crouseilles et al. [11] progettata per mesh uniformi, gli autori estendono l'approccio a mesh non uniformi arbitrarie utilizzando polinomi di Hermite.

   • Formulazione Esatta: Imponendo la continuità $C^2$ tra i polinomi di Hermite locali definiti su celle adiacenti, gli autori derivano una relazione ricorsiva che lega la derivata all'interfaccia alle derivate alle interfacce vicine e ai valori della funzione. Ciò risulta in un piccolo sistema lineare globale (di dimensione $(N_p-1) \times (N_p-1)$ per $N_p$ patch) che deve essere risolto per ottenere derivate equivalenti a quelle di una spline globale.
   • Formulazione Approssimata: Gli autori analizzano il comportamento asintotico dei coefficienti del sistema globale. Dimostrano che l'influenza delle interfacce distanti decade esponenzialmente all'aumentare del numero di celle per patch. Di conseguenza, per patch sufficientemente grandi (ad esempio, $\ge 30$ celle), la dipendenza dalle derivate distanti diventa trascurabile (raggiungendo la precisione di macchina). Ciò consente un'approssimazione locale troncata, dove le derivate vengono calcolate utilizzando solo un numero limitato di celle vicine, evitando la risoluzione globale ma mantenendo un'alta precisione.
   • Formule Esplicite: Per il caso specifico di patch uniformi, gli autori forniscono espressioni chiuse esplicite per i coefficienti, ottimizzando ulteriormente il calcolo.

2. Estensione 2D: Il metodo è esteso a geometrie multi-patch 2D, gestendo tre casi:

   • Mesh Conformes: Un'applicazione diretta del prodotto tensoriale delle formule 1D per calcolare le derivate $r$, le derivate $\theta$ e le derivate incrociate.
   • Mesh Non Conformes: Casi in cui i punti di rottura (break points) non sono allineati tra i patch. Gli autori propongono un algoritmo che calcola progressivamente le derivate su interfacce più grossolane e le interpola su patch più raffinati per imporre la continuità $C^1$.
   • T-Joints: Connessioni che coinvolgono tre o più patch dove un bordo è connesso solo parzialmente. Il metodo adatta l'ordine di calcolo delle derivate in base alla geometria per garantire che tutte le derivate necessarie siano disponibili.

3. Implementazione: L'algoritmo è implementato nella libreria Gyselalib++. L'equazione di Poisson è risolta con un metodo a Elementi Finiti (FEM) con spline B cubiche su una spline globale equivalente, mentre l'equazione di advezione è risolta con il proposto metodo BSL su spline locali.

Contributi Chiave

• Generalizzazione delle Formule delle Derivate d'Interfaccia: Il lavoro generalizza la formula del calcolo delle derivate da mesh uniformi a mesh non uniformi arbitrarie utilizzando i polinomi di Hermite, consentendo la gestione di situazioni non conformi e il raffinamento locale.
• Analisi della Stabilità: Gli autori dimostrano che la regolarità $C^0$ è insufficiente per la stabilità negli schemi BSL su domini multi-patch e che è richiesta la regolarità $C^1$.
• Compromesso Locale vs. Globale: Il lavoro stabilisce che, sebbene un sistema globale esatto garantisca la continuità $C^2$, un'approssimazione locale utilizzando un stencil troncato (ad esempio, 30 celle) raggiunge una precisione di macchina equivalente alla soluzione globale, offrendo un'alternativa computazionalmente efficiente.
• Gestione di Geometrie Complesse: Il metodo gestisce con successo i T-joint e le mesh non conformi, che sono necessari per modellare le sezioni trasversali complesse dei tokamak che coinvolgono sia punti O che punti X.

Risultati Numerici
Il metodo è stato validato utilizzando la libreria Gyselalib++ con i seguenti test:

1. Interpolazione Spline: Test su patch concentrici (uniformi e non uniformi) hanno confermato che le rappresentazioni spline locali convergono alla spline globale equivalente con ordine 4. L'errore nelle derivate d'interfaccia utilizzando l'approssimazione troncata ha eguagliato la formula esatta una volta che il numero di celle ha superato le $\approx 30$, raggiungendo la precisione di macchina ($10^{-13}$).
2. Test di Advezione:
   • Patch non conformi: Le spline locali su una mesh non conforme hanno raggiunto un'accuratezza paragonabile a una mesh globale raffinata, nonostante l'uso di un numero totale di celle inferiore.
   • T-joint: In una configurazione di T-joint 2D, i risultati della spline locale hanno concordato con i risultati della spline globale fino alla precisione di macchina ($\approx 6 \times 10^{-13}$).
3. Modello di Centro di Guida (Instabilità Diocotrona):
   • Le simulazioni su griglie multi-patch conformi hanno riprodotto i risultati delle simulazioni spline globali (tassi di crescita e conservazione della massa) con errori coerenti con la precisione dell'interpolazione.
   • Una simulazione con una mesh non ottimale (patch centrale depopolata e interfacce spostate) ha mostrato che, sebbene la scelta della mesh influenzi l'ordine di grandezza dell'errore assoluto, il metodo della spline locale ha comunque tracciato correttamente il comportamento della soluzione globale, validando la robustezza dell'approccio.

Significato e Rivendicazioni
Il documento sostiene di fornire un framework robusto per risolvere equazioni di tipo Vlasov su geometrie multi-patch, essenziali per simulazioni realistiche dei tokamak. La significatività primaria risiede nel consentire l'uso di interpolazione spline cubica locale con condizioni al contorno di Hermite per mantenere la stabilità e l'accuratezza di alto ordine attraverso le interfacce dei patch.

Gli autori sottolineano che il loro metodo rimuove il vincolo riscontrato negli studi precedenti riguardante il numero di punti per patch nei casi conformi, sostituendolo con la necessità di risolvere un piccolo sistema globale (o un'approssimazione locale). Notano che il metodo è attualmente validato per geometrie con punti O. Sebbene l'introduzione dei T-joint sia un passo verso la connessione di geometrie con punti O e punti X, gli autori dichiarano esplicitamente che la gestione delle singolarità dei punti X richiede ulteriori studi e non è ancora stata completamente risolta in questo lavoro. Il metodo è presentato come un passo significativo verso solver semi-Lagrangiani efficienti e ad alto ordine per geometrie di confinamento magnetico complesse.