Optimización de la factorización Cholesky para múltiples matrices dispersas con el mismo patrón no nulo

Question

Estoy utilizando una factorización Cholesky para resolver el paso lineal de un sistema de ecuaciones no lineal (análisis de elementos finitos no lineales). En el PETSc se puede especificar un parámetro para SAME_NONZERO_PATTERN durante las sucesivas resoluciones.

El aumento de velocidad es bueno, casi demasiado bueno (100 veces más rápido para la segunda solución). Esto me hace preguntarme qué tipo de optimización se puede hacer para las factorizaciones Cholesky cuando se resuelven muchas matrices dispersas con el mismo patrón de no ceros (pero con valores diferentes, posiblemente similares).

Answer 1

Casi todos los paquetes de solución directa dispersa dividen la factorización en dos etapas:

• La factorización simbólica calcula una ordenación, a menudo utilizando la disección anidada o un algoritmo de grado mínimo aproximado, de forma que los factores sean lo más dispersos posible y asigna espacio para mantener el resultado. El valor de las entradas de la matriz no se utiliza o sólo se utiliza para hacer estimaciones sobre el pivote.

• La factorización numérica calcula la factorización dada la ordenación y la dispersión calculada por la factorización simbólica.

En la mayoría de las circunstancias, la factorización simbólica es mucho menos costosa que la factorización numérica, por lo que SAME_NONZERO_PATTERN ofrece un beneficio limitado. Esto cambia cuando se ejecuta en paralelo con muchos procesos, ya que la factorización simbólica no escala tan bien (y algunos paquetes populares, incluyendo MUMPS, la calculan en serie).

Es muy poco probable que su factor de 100 provenga de la etapa de factorización simbólica. ¿Quizás estés reutilizando los factores de la iteración anterior? Para cuestiones de rendimiento como ésta, suele ser útil ejecutar con -log_summary y comparar los tiempos y el equilibrio de carga para el Mat*FactorSym y Mat*FactorNum eventos. Envíe la salida a petsc-users@mcs.anl.gov o petsc-maint@mcs.anl.gov si desea ayuda para interpretar la salida.

En cuanto a la "actualización" de los factores tras pequeños cambios en las entradas de la matriz, los intentos de hacerlo han sido normalmente infructuosos. Una alternativa que recomiendo es utilizar la antigua factorización como precondicionador de la nueva matriz. Para los problemas no lineales usando SNES, se puede experimentar usando -snes_lag_preconditioner o -snes_lag_jacobian . Si utiliza KSP directamente, entonces pase SAME_PRECONDITIONER a KSPSetOperators() .

Answer 2

Un buen ejemplo es el caso de las señales de suavizado mediante splines cúbicos. En ese caso, surgen sistemas de ecuaciones definidos positivamente que tienen una estructura rica, lo que permite calcular la factorización de Cholesky con una complejidad lineal (en lugar de una complejidad cúbica). En el ejemplo de los splines cúbicos, la matriz de coeficientes $A$ es pentadiagonal, simétrica, persimétrica y de Toeplitz. Entonces, considerando la ecuación de la factorización de Cholesky $A = L L^T$ , donde digamos $L$ es triangular inferior, se puede demostrar que $L$ tiene también una estructura muy rica y que sus elementos pueden calcularse resolviendo un sistema acoplado de tres ecuaciones en diferencia de segundo orden. Sin embargo, la posibilidad de optimizar la factorización Cholesky de una matriz dispersa definida positivamente depende en gran medida de la estructura particular de la dispersidad.

Answer 3

No tengo mucha experiencia con factorizaciones dispersas y no estoy seguro de haber entendido bien tu pregunta. Pero hasta donde yo sé, los métodos de factorización dispersa primero calculan un reordenamiento de los elementos de la matriz para que el resultado de la factorización (la matriz triangular $\mathbf{L}$ para lo cual $\mathbf{A}=\mathbf{LL}^T$ o $\mathbf{L}^T\mathbf{L}$ ) tiene una estructura muy dispersa. Así que asumo que la optimización simple es, que el reordenamiento no tiene que ser recalculado cada vez y simplemente puede ser reutilizado para el mismo patrón de dispersión de $\mathbf{A}$ ya que debería ser independiente de los valores reales.