Randomized Strong Recursive Skeletonization: Simultaneous Compression and LU Factorization of Hierarchical Matrices using Matrix-Vector Products

Resumen Técnico: Esqueletización Recursiva Fuerte Aleatorizada

Planteamiento del Problema
Las matrices jerárquicas ($H^2$-matrices) bajo la condición de admisibilidad fuerte ofrecen un marco para comprimir matrices densas derivadas de ecuaciones diferenciales parciales (PDE) y ecuaciones integrales. Estas matrices particionan las interacciones en componentes de campo cercano (almacenados densamente) y de campo lejano (aproximados por bloques de bajo rango), utilizando una estructura de base anidada para lograr una complejidad de almacenamiento lineal y una multiplicación matriz-vector eficiente. Si bien los productos rápidos de matriz-vector están bien establecidos, la inversión directa (factorización LU) de matrices $H^2$ sigue siendo un desafío, particularmente para problemas tridimensionales. Los métodos de inversión existentes suelen requerir el acceso explícito a las entradas de la matriz y dependen de la recompresión repetida de los bloques fuera de la diagonal durante el proceso de factorización recursiva. Además, obtener una representación $H^2$ para matrices disponibles solo a través de su acción sobre vectores (entornos libres de matriz o matrix-free) ha sido un problema abierto, ya que los enfoques aleatorizados previos sufrían de prefactores grandes o requerían matrices de prueba estructuradas que son difíciles de diseñar dinámicamente.

Metodología
El artículo presenta la Esqueletización Recursiva Fuerte Aleatorizada (RSRS), un algoritmo que comprime y factoriza simultáneamente matrices $H^2$ generales utilizando únicamente la acción de la matriz $A$ y su adjunta $A^*$ sobre vectores. El método opera bajo la condición de admisibilidad fuerte y evita el acceso directo a las entradas de la matriz.

La metodología central se basa en dos componentes principales:

1. Esquematización Aleatorizada con Matrices Gaussianas Densas: En lugar de diseñar matrices de prueba estructuradas con patrones de ceros específicos (lo que aumenta la complejidad de muestreo y requiere conocimiento previo de los índices de los esqueletos), RSRS emplea matrices de prueba gaussianas densas $\Omega$ y $\Psi$. El algoritmo genera esquemas iniciales $Y = A\Omega$ y $Z = A^*\Psi$.
2. Imposición Estructural Dinámica mediante Transformaciones Lineales: A medida que la factorización progresa recursivamente a través de un árbol jerárquico (desde los niveles más finos a los más gruesos), el algoritmo introduce la estructura necesaria en las matrices de prueba mediante post-procesamiento. Se aplican dos operaciones clave a los esquemas densos:
   • Anulación de Bloques: Calcula una base de nulidad para anular las contribuciones de los índices de campo cercano, aislando efectivamente las interacciones de campo lejano para la compresión de bajo rango (Descomposición Interpolativa).
   • Extracción de Bloques: Utiliza transformaciones de pseudoinversa para extraer los subbloques densos de campo cercano específicos requeridos para los pasos de eliminación.

El algoritmo aplica recursivamente un procedimiento de "esparcificar-luego-diagonalizar". En cada caja de la jerarquía, calcula una descomposición interpolativa (ID) de las interacciones de campo lejano para identificar índices redundantes (esqueletos). Luego, aplica matrices de eliminación de bloques para desacoplar estos índices del resto del sistema. Crucialmente, las matrices de esquema ($Y, Z$) y las matrices de prueba ($\Omega, \Psi$) se actualizan en cada paso utilizando las matrices de eliminación y esparcimiento derivadas ($L, U, E, F$). Esto permite al algoritmo mantener esquemas del sistema diagonalizado en evolución sin acceder nunca a las entradas individuales de la matriz ni volver a generar muestras aleatorias.

Contribuciones Clave

• Inversión Libre de Matriz (Matrix-Free): RSRS extiende la inversión de $H^2$ a entornos donde las entradas de la matriz son inaccesibles, dependiendo únicamente de los productos matriz-vector. Esto amplía la aplicabilidad a operadores integrales resolubles mediante Métodos de Multipolo Rápido (FMM) o FFTs, y a bloques densos derivados de resolvedores de PDE dispersos.
• Compresión y Factorización Simultáneas: A diferencia de los métodos tradicionales que separan la compresión y la inversión, RSRS integra estos pasos. Esto evita las estructuras de datos de almacenamiento intensivo y las recompresiones repetidas requeridas por enfoques previos, manteniendo los bloques fuera de la diagonal actualizados en forma comprimida.
• Eficiencia de Muestreo: El número de productos matriz-vector requeridos ($s$) es independiente del tamaño de la matriz $N$ y depende linealmente solo del rango máximo de los bloques fuera de la diagonal y la geometría del problema. El método logra esto utilizando matrices de prueba gaussianas densas genéricas, evitando los altos costos de muestreo asociados con las matrices de prueba estructuradas.
• Robustez: El algoritmo está diseñado para manejar matrices indefinidas e mal condicionadas, como las que surgen de las ecuaciones de Helmholtz 3D, sin requerir que la matriz sea definida positiva.

Resultados Numéricos
Los autores validan RSRS a través de tres experimentos numéricos:

1. Ecuación Integral de Volumen 2D: Una ecuación integral de Laplace de primera clase en una rejilla uniforme. A pesar de que la matriz se vuelve cada vez más mal condicionada a medida que $N$ crece (número de condición $\sim 10^6$), RSRS mantuvo errores de compresión relativa bajos ($10^{-7}$ a $10^{-6}$) y redujo las iteraciones de GMRES de miles a menos de 5 cuando se utilizó como precondicionador.
2. Ecuación Integral de Contorno 3D: Una ecuación de Fredholm de segunda clase en una esfera. El método demostró una precisión estable y escalamiento lineal en el tiempo de reconstrucción, memoria y tiempo de resolución con respecto a $N$.
3. Complemento de Schur de Helmholtz 3D: Un complemento de Schur denso derivado de un resolvedor directo disperso para la ecuación de Helmholtz 3D. La matriz es indefinida. RSRS comprimió con éxito el complemento de Schur, reduciendo las iteraciones de GMRES de la no convergencia (más de 10,000) a dígitos únicos (3–9) para varios tamaños de problema.

En todos los casos, el tiempo de reconstrucción ($T_{rec}$), el tiempo de resolución ($T_{sol}$) y el uso de memoria escalaron linealmente con $N$, consistentes con el análisis de complejidad teórica de $O(k^2 N)$ y $O(k N)$.

Significancia y Reivindicaciones
El artículo afirma que RSRS representa un cambio significativo respecto a las técnicas tradicionales de inversión de $H^2$ al permitir la factorización directa en entornos libres de matriz. Al combinar la compresión y la inversión en un único proceso recursivo impulsado por la esquematización aleatorizada, el método ofrece una eficiencia computacional competitiva y requisitos de almacenamiento reducidos. Los autores enfatizan que el uso de matrices de prueba gaussianas densas, junto con las transformaciones lineales dinámicas, elimina la necesidad de diseños de matrices de prueba complejos y específicos del problema. El método se presenta como un resolvedor directo aproximado robusto adecuado para una amplia gama de aplicaciones, incluyendo ecuaciones integrales de contorno, discretizaciones de PDE y cuantificación de la incertidumbre, particularmente donde las entradas de la matriz son costosas o no están disponibles. El trabajo sugiere que el efecto de retroalimentación de usar esquemas para definir las transformaciones que actualizan dichos esquemas es benigno en la práctica, permitiendo una factorización precisa sin volver a generar muestras.