Implementación eficiente del algoritmo de matriz tridiagonal

12

Estoy resolviendo un problema físico usando un esquema numérico implícito. Esto me lleva a resolver una ecuación lineal con matriz tridiagonal. He codificado este algoritmo de Wikipedia. Me pregunto si hay una biblioteca eficiente que permita resolver este tipo de ecuaciones de manera optimizada. Una nota importante es que la matriz en sí misma solo cambia cuando los parámetros del sistema están cambiando, por lo que tuve la oportunidad de calcular previamente algunos pasos del algoritmo para obtener una buena bonificación de rendimiento. Estoy usando C ++.

linear-algebra libraries c++ gmk
fuente
¿Qué tan grande es un sistema? ¿Necesita ser paralelo?
aterrel
1
El tamaño depende de la precisión requerida (de cientos a decenas de miles de valores). Ahora estoy codificando en una computadora de un solo núcleo, pero es posible obtener acceso a la supercomputadora de la universidad con muchos cpus disponibles, por lo que el soporte de paralelismo sería bueno.
GM

Respuestas:

15

Probablemente debería comenzar con la implementación LAPACK,? Gtsv, por ejemplo, dgtsv . Si desea una versión de memoria distribuida, es posible que desee comenzar con p? Gtsv de ScaLAPACK.

EDITAR: Dado que su matriz no cambia muy a menudo, puede evitar la factorización redundante de la matriz tridiagonal separando la rutina LAPACK? Gtsv en el paso de factorización,? Gttrf, y la etapa de resolución,? Gttrs. Existen rutinas con nombres similares en ScaLAPACK que tienen el mismo propósito.

Jack Poulson
fuente
Gracias, parece lo que necesito. Intentaré ahora ejecutar estas rutinas desde mi código.
GM
1
Como lo está llamando desde C ++, asegúrese de declarar el prototipo dentro de un bloque externo "C" {}. Dependiendo de su sistema, es posible que deba agregar un guión bajo al nombre de la rutina.
Jack Poulson
2

Para sistemas paralelos distribuidos : no he probado ScaLAPACK, que tiene un solucionador tridiagonal paralelo, para el cual hay ejemplos disponibles en línea . He intentado con cierto éxito un método propuesto por David Moulton en una publicación de LANL . Codificar esto puede ser más de lo que quieres hacer, pero al usar LAPACK, es sencillo.

Yann
fuente
1

Encontré un algoritmo recursivo interesante aquí en la página 975. Parece prometedor, me pregunto qué dicen las personas más experimentadas al respecto.

tiam
fuente
Recetas numéricas tiene algunos errores. En términos de una fuente de códigos para usar, no es el mejor, aunque algunos lo consideran un clásico. Me sorprendería si ScaLAPACK no implementara un algoritmo al menos tan eficiente como la reducción cíclica recursiva.
Geoff Oxberry