Two-Level Sketching Alternating Anderson acceleration for Complex Physics Applications

Technische Samenvatting: Twee-niveau schetsing met wisselende Anderson-versnelling voor complexe fysica-toepassingen

Probleemstelling
De wisselende Anderson-Picard (AAP)-methode is een robuuste techniek voor het versnellen van vaste-puntiteraties in single- en multi-fysica-simulaties die worden beheerst door gediscrétiseerde partiële differentiaalvergelijkingen (PDE's). Hoewel AAP de convergentiesnelheden verbetert en convergentie herstelt voor niet-contractieve operatoren, vereist het het oplossen van een kleinste-kwadraten (LS)-probleem bij elke mengstap. In wetenschappelijk rekenen op grote schaal, met name op High-Performance Computing (HPC)-platforms, worden de kosten van het oplossen van deze LS-problemen – die vaak lange en smalle matrices met miljoenen rijen omvatten – een aanzienlijke knelpunt. Bestaande benaderingen om deze kosten te verlagen, vertrouwen ofwel op fysica-specifieke benaderingen (die niet generaliseerbaar zijn) of op adaptieve schetsingstechnieken die instabiliteit of overhead bij geheugentoegang (cache-misses) kunnen introduceren als ze niet zorgvuldig worden beheerd.

Methodologie
De auteurs stellen een nieuw twee-niveau schetsingsextensie van de AAP-methode voor. Deze benadering decomposeert de benadering van het LS-probleem in twee distincte, complementaire fasen:

1. Niveau 1: Statische fysica-gebaseerde projectie ($\Pi_1$):

   • Dit niveau projecteert het residu op een specifiek subset van fysieke velden (bijvoorbeeld druk in Stokes-stroming of extracellulair potentieel in Bidomain-modellen).
   • Motivatie: Multi-fysica-problemen kunnen vaak worden herschreven als systemen die een enkele variabele omvatten (analoog aan de Schur-complement in lineaire systemen). Bovendien opereren verschillende fysieke velden vaak op sterk verschillende schalen; het combineren ervan in een monolithisch LS-probleem kan de conditionering verslechteren.
   • Implementatie: Deze operator is statisch en wordt bepaald voorafgaand aan de uitvoering, wat toestaat voor vooraf toewijzen van geheugen en het verminderen van de dimensie van de LS-matrix van $n$ naar $\ell_1$.

2. Niveau 2: Dynamische algebraïsche schetsing ($\Pi_2$):

   • Dit niveau voert een adaptieve subselectie uit van rijen uit de reeds geprojecteerde matrix.
   • Stabiliteitsanalyse: De auteurs breiden de achterwaartse stabiliteitsanalyse (voorheen toegepast op lineaire systemen) uit naar niet-lineaire vaste-puntiteraties onder de aanname van Lipschitz-continuïteit. Zij leiden voorwaarden af waarbij de verstoring die door schetsing wordt geïntroduceerd, de convergentie van AAP niet in gevaar brengt.
   • Adaptieve Strategie: Het algoritme schat stabiliteitsthrsholds ($\epsilon_{LHS}$ en $\epsilon_{RHS}$) met behulp van goedkope benaderingen van de kleinste singuliere waarde en de Lipschitz-constante. Als deze thresholds worden voldaan, wordt een masker toegepast.
   • Cache-bewustzijn: Om de prestatieboete van willekeurige geheugentoegang (cache-misses) te mitigeren, hebben de auteurs empirische benchmarks uitgevoerd. Zij stelden vast dat voor matrix-vectorproducten het behouden van niet meer dan 10% van de rijen snelheidswinst oplevert, terwijl voor QR-factorisaties het behouden van maximaal 80% optimaal is. De strategie maakt gebruik van ofwel "subselectie" (het behouden van rijen met de grootste residuen) of "gerandomiseerde" selectie, beheerst door convergente reeksparameters ($\eta_j$).

Belangrijkste Bijdragen

• Gestructureerde Decompositie: Het artikel introduceert een gestructureerde twee-niveau decompositie waarbij de statische fysica-gebaseerde restrictie zich richt op conditionering en schaling, terwijl dynamische algebraïsche schetsing gecontroleerde verstoringen introduceert. De auteurs betogen dat deze rollen distinct zijn en niet kunnen worden samengevoegd zonder verlies van robuustheid.
• Niet-lineaire Stabiliteitsanalyse: Het werk breidt de achterwaartse stabiliteitsanalyse uit naar niet-lineaire operatoren, en biedt theoretische rechtvaardiging voor adaptieve rijselectie in niet-lineaire vaste-puntiteraties.
• Geheugen- en Cache-optimalisatie: De implementatie maakt gebruik van in-place-bewerkingen en specifieke richtlijnen voor randomiseerde selectie (10% voor mat-vec, 80% voor QR) om het geheugengebruik en cache-misses te minimaliseren, waarmee praktische HPC-beperkingen worden aangepakt.
• Open-source Implementatie: Het algoritme is geïmplementeerd in Julia, met gebruikmaking van bibliotheken zoals Ferrite en Gridap, en is beschikbaar gesteld via een GitHub-repository.

Resultaten
De methode werd getest op vier benchmarkproblemen: Stokes, $p$-Laplacian, Bidomain-vergelijkingen en stationaire Navier-Stokes, over variërende probleemgroottes (Vrijheidsgraden).

• Prestatiewinst: De twee-niveau schetsing AAP bereikte tot 50% reductie in tijd tot oplossing vergeleken met standaard Anderson-versnelling (AA) zonder de convergentiesnelheden te verslechteren.
• Probleem-specifiek Gedrag:
  • Stokes: Het op druk gebaseerde masker was het meest effectief, met een reductie van 24% in iteraties en 33% in tijd bij grote problemen. De "Randomized Power"-adaptiviteitsstrategie was het meest efficiënt.
  • Bidomain: Het extracellulaire masker werkte robuust, terwijl het intracellulaire masker divergentie veroorzaakte bij grote problemen. Randomiseerde adaptiviteit leverde tot een 43% reductie in CPU-tijd op.
  • Navier-Stokes: Het snelheidsmasker was superieur aan het drukmasker (in tegenstelling tot Stokes), hoewel de winsten bescheiden waren (~9%) omdat de dominante kosten lagen bij de directe lineaire solver en niet bij de LS-oplossing.
  • $p$-Laplacian: Randomiseerde adaptiviteitsstrategieën leverden de beste prestaties voor grotere problemen.
• Robuustheid: De adaptieve strategieën behielden convergentiegedrag dat vergelijkbaar was met niet-adaptieve methoden, wat de stabiliteit van de voorgestelde schetsingsgrenzen bevestigt.

Betekenis en Claims
De auteurs claimen dat hun werk een stabiliteitsgedreven kader biedt voor het schetsen van Anderson-versnelling. Zij benadrukken dat het doel niet is om globale convergentiesnelheden voor schetsing vast te stellen (wat een open probleem blijft, zelfs voor standaard AA), maar om toelaatbare benaderingen te karakteriseren die het convergentiegedrag van de exacte methode behouden.

De betekenis ligt in het vermogen om de rekenkundige overhead van de LS-oplossing te verminderen – een niet-verwaarloosbaar deel van de runtime bij herhaalde AA-toepassingen – zonder de onderliggende niet-lineaire solver of preconditioner te wijzigen. De methode is bijzonder effectief in scenario's waar de residu-evaluatie niet buitensporig dominant is, waardoor de besparingen in de LS-oplossing direct kunnen worden omgezet in reducties van de totale tijd tot oplossing. Het artikel concludeert dat deze aanpak robuust, schaalbaar en klaar is voor integratie in high-performance wetenschappelijke rekenkaders.