Algoritmo de Lanczos

Cálculo de valores propios numéricos

El algoritmo de Lanczos es un método iterativo ideado por Cornelius Lanczos que es una adaptación de métodos de potencia para encontrar los valores propios y vectores propios "más útiles" (que tienden hacia los extremos más altos/más bajos) de una matriz hermitiana , donde a menudo, pero no necesariamente, es mucho más pequeño que . ^[1] Aunque en principio es computacionalmente eficiente, el método tal como se formuló inicialmente no fue útil debido a su inestabilidad numérica . ${\displaystyle m}$ ${\displaystyle n\times n}$ ${\displaystyle m}$ ${\displaystyle n}$

En 1970, Ojalvo y Newman demostraron cómo hacer que el método fuera numéricamente estable y lo aplicaron a la solución de estructuras de ingeniería de gran tamaño sometidas a cargas dinámicas. ^[2] Esto se logró usando un método para purificar los vectores de Lanczos (es decir, reortogonalizando repetidamente cada vector recién generado con todos los generados previamente) ^[2] con cualquier grado de precisión, que cuando no se realizaba, producía una serie de vectores que se altamente contaminados por aquellos asociados con las frecuencias naturales más bajas.

En su trabajo original, estos autores también sugirieron cómo seleccionar un vector inicial (es decir, usar un generador de números aleatorios para seleccionar cada elemento del vector inicial) y sugirieron un método determinado empíricamente para determinar el número reducido de vectores (es decir, debería ser seleccionado para que sea aproximadamente 1,5 veces el número de valores propios precisos deseados). Poco después, Paige siguió su trabajo y también proporcionó un análisis de errores. ^{[3] [4]} En 1988, Ojalvo produjo una historia más detallada de este algoritmo y una prueba de error de valor propio eficiente. ^[5] ${\displaystyle m}$

el algoritmo

Ingrese una matriz hermitiana de tamaño y, opcionalmente, una cantidad de iteraciones (por defecto, let ). ${\displaystyle A}$ ${\displaystyle n\times n}$ ${\displaystyle m}$ ${\displaystyle m=n}$

• Estrictamente hablando, el algoritmo no necesita acceso a la matriz explícita, sino sólo una función que calcula el producto de la matriz por un vector arbitrario. Esta función se llama la mayoría de las veces. ${\displaystyle v\mapsto Av}$ ${\displaystyle m}$

Genere una matriz con columnas ortonormales y una matriz de tamaño simétrica real tridiagonal . Si , entonces es unitario y . ${\displaystyle n\times m}$ ${\displaystyle V}$ ${\displaystyle T=V^{*}AV}$ ${\displaystyle m\times m}$ ${\displaystyle m=n}$ ${\displaystyle V}$ ${\displaystyle A=VTV^{*}}$

Advertencia La iteración de Lanczos es propensa a la inestabilidad numérica. Cuando se ejecuta en aritmética no exacta, se deben tomar medidas adicionales (como se describe en secciones posteriores) para garantizar la validez de los resultados.

1. Sea un vector arbitrario con norma euclidiana . ${\displaystyle v_{1}\in \mathbb {C} ^{n}}$ ${\displaystyle 1}$
2. Paso de iteración inicial abreviado:
   1. Dejar . ${\displaystyle w_{1}'=Av_{1}}$
   2. Dejar . ${\displaystyle \alpha _{1}=w_{1}'^{*}v_{1}}$
   3. Dejar . ${\displaystyle w_{1}=w_{1}'-\alpha _{1}v_{1}}$
3. Para hacer: ${\displaystyle j=2,\dots ,m}$
   1. Let (también norma euclidiana ). ${\displaystyle \beta _{j}=\|w_{j-1}\|}$
   2. Si , entonces, deja , ${\displaystyle \beta _{j}\neq 0}$ ${\displaystyle v_{j}=w_{j-1}/\beta _{j}}$

      de lo contrario, elija como un vector arbitrario con norma euclidiana que sea ortogonal a todos . ${\displaystyle v_{j}}$ ${\displaystyle 1}$ ${\displaystyle v_{1},\dots ,v_{j-1}}$

   3. Dejar . ${\displaystyle w_{j}'=Av_{j}}$
   4. Dejar . ${\displaystyle \alpha _{j}=w_{j}'^{*}v_{j}}$
   5. Dejar . ${\displaystyle w_{j}=w_{j}'-\alpha _{j}v_{j}-\beta _{j}v_{j-1}}$
4. Sea la matriz con columnas . Dejar . ${\displaystyle V}$ ${\displaystyle v_{1},\dots ,v_{m}}$ ${\displaystyle T={\begin{pmatrix}\alpha _{1}&\beta _{2}&&&&0\\\beta _{2}&\alpha _{2}&\beta _{3}&&&\\&\beta _{3}&\alpha _{3}&\ddots &&\\&&\ddots &\ddots &\beta _{m-1}&\\&&&\beta _{m-1}&\alpha _{m-1}&\beta _{m}\\0&&&&\beta _{m}&\alpha _{m}\\\end{pmatrix}}}$

Nota para . ${\displaystyle Av_{j}=w_{j}'=\beta _{j+1}v_{j+1}+\alpha _{j}v_{j}+\beta _{j}v_{j-1}}$ ${\displaystyle 2<j<m}$

En principio, existen cuatro formas de escribir el procedimiento de iteración. Paige y otros trabajos muestran que el orden de operaciones anterior es el más estable numéricamente. ^{[6] [7]} En la práctica, el vector inicial puede tomarse como un argumento más del procedimiento, incluyéndose indicadores de imprecisión numérica como condiciones adicionales de terminación del bucle. ${\displaystyle v_{1}}$ ${\displaystyle \beta _{j}=0}$

Sin contar la multiplicación matriz-vector, cada iteración realiza operaciones aritméticas. La multiplicación matriz-vector se puede realizar en operaciones aritméticas donde es el número promedio de elementos distintos de cero en una fila. La complejidad total es así , o si ; El algoritmo de Lanczos puede ser muy rápido para matrices dispersas. Los esquemas para mejorar la estabilidad numérica normalmente se juzgan en función de este alto rendimiento. ${\displaystyle O(n)}$ ${\displaystyle O(dn)}$ ${\displaystyle d}$ ${\displaystyle O(dmn)}$ ${\displaystyle O(dn^{2})}$ ${\displaystyle m=n}$

Los vectores se denominan vectores de Lanczos . El vector no se utiliza después de calcularse y el vector no se utiliza después de calcularse. Por tanto, se puede utilizar el mismo almacenamiento para los tres. Asimismo, si solo se busca la matriz tridiagonal , entonces no se necesita la iteración cruda después de haber calculado , aunque algunos esquemas para mejorar la estabilidad numérica la necesitarían más adelante. A veces, los vectores Lanczos posteriores se vuelven a calcular cuando es necesario. ${\displaystyle v_{j}}$ ${\displaystyle w_{j}'}$ ${\displaystyle w_{j}}$ ${\displaystyle w_{j}}$ ${\displaystyle v_{j+1}}$ ${\displaystyle T}$ ${\displaystyle v_{j-1}}$ ${\displaystyle w_{j}}$ ${\displaystyle v_{1}}$

Aplicación al problema propio

El algoritmo de Lanczos se plantea con mayor frecuencia en el contexto de encontrar los valores propios y los vectores propios de una matriz, pero mientras que una diagonalización ordinaria de una matriz haría que los vectores propios y los valores propios fueran evidentes a partir de la inspección, no ocurre lo mismo con la tridiagonalización realizada por los Lanczos. algoritmo; Se necesitan pasos adicionales no triviales para calcular incluso un solo valor propio o vector propio. No obstante, la aplicación del algoritmo de Lanczos suele ser un importante paso adelante en el cálculo de la descomposición propia.

Si es un valor propio de y su vector propio ( ), entonces es un vector propio correspondiente de con el mismo valor propio: ${\displaystyle \lambda }$ ${\displaystyle T}$ ${\displaystyle x}$ ${\displaystyle Tx=\lambda x}$ ${\displaystyle y=Vx}$ ${\displaystyle A}$

$${\displaystyle {\begin{aligned}Ay&=AVx\\&=VTV^{*}Vx\\&=VTIx\\&=VTx\\&=V(\lambda x)\\&=\lambda Vx\\&=\lambda y.\end{aligned}}}$$

Así, el algoritmo de Lanczos transforma el problema de descomposición propia de en el problema de descomposición propia de . ${\displaystyle A}$ ${\displaystyle T}$

1. Para matrices tridiagonales, existen varios algoritmos especializados, a menudo con mayor complejidad computacional que los algoritmos de propósito general. Por ejemplo, si es una matriz simétrica tridiagonal entonces: ${\displaystyle T}$ ${\displaystyle m\times m}$
   • La recursividad continua permite calcular el polinomio característico en las operaciones y evaluarlo en un punto de las operaciones. ${\displaystyle O(m^{2})}$ ${\displaystyle O(m)}$
   • El algoritmo de valores propios de divide y vencerás se puede utilizar para calcular la descomposición propia completa de las operaciones. ${\displaystyle T}$ ${\displaystyle O(m^{2})}$
   • El método rápido multipolar ^[8] puede calcular todos los valores propios en solo operaciones. ${\displaystyle O(m\log m)}$
2. Se sabe que algunos algoritmos generales de autodescomposición, en particular el algoritmo QR , convergen más rápido para matrices tridiagonales que para matrices generales. La complejidad asintótica del QR tridiagonal es la misma que la del algoritmo de divide y vencerás (aunque el factor constante puede ser diferente); dado que los vectores propios juntos tienen elementos, esto es asintóticamente óptimo. ${\displaystyle O(m^{2})}$ ${\displaystyle m^{2}}$
3. Incluso los algoritmos cuyas tasas de convergencia no se ven afectadas por transformaciones unitarias, como el método de potencia y la iteración inversa , pueden disfrutar de beneficios de rendimiento de bajo nivel al aplicarse a la matriz tridiagonal en lugar de a la matriz original . Dado que es muy escaso y todos los elementos distintos de cero se encuentran en posiciones muy predecibles, permite un almacenamiento compacto con un rendimiento excelente frente al almacenamiento en caché . Asimismo, es una matriz real con todos los vectores propios y valores propios reales, mientras que en general puede tener elementos y vectores propios complejos, por lo que la aritmética real es suficiente para encontrar los vectores propios y valores propios de . ${\displaystyle T}$ ${\displaystyle A}$ ${\displaystyle T}$ ${\displaystyle T}$ ${\displaystyle A}$ ${\displaystyle T}$
4. Si es muy grande, entonces reducirlo para que sea de un tamaño manejable aún permitirá encontrar los valores propios y vectores propios más extremos de ; En la región, el algoritmo de Lanczos puede verse como un esquema de compresión con pérdida para matrices hermitianas, que enfatiza la preservación de los valores propios extremos. ${\displaystyle n}$ ${\displaystyle m}$ ${\displaystyle T}$ ${\displaystyle A}$ ${\displaystyle m\ll n}$

La combinación de buen rendimiento para matrices dispersas y la capacidad de calcular varios (sin calcular todos) valores propios son las principales razones para elegir utilizar el algoritmo de Lanczos.

Aplicación a la tridiagonalización

Aunque el problema propio es a menudo la motivación para aplicar el algoritmo de Lanczos, la operación que realiza principalmente el algoritmo es la tridiagonalización de una matriz, para la cual las transformaciones de Householder numéricamente estables se han favorecido desde la década de 1950. Durante la década de 1960, se hizo caso omiso del algoritmo de Lanczos. El interés en él se vio rejuvenecido por la teoría de la convergencia de Kaniel-Paige y el desarrollo de métodos para prevenir la inestabilidad numérica, pero el algoritmo de Lanczos sigue siendo el algoritmo alternativo que se prueba sólo si Householder no es satisfactorio. ^[9]

Los aspectos en los que difieren los dos algoritmos incluyen:

• Lanczos aprovecha que es una matriz escasa, mientras que Householder no, y generará relleno . ${\displaystyle A}$
• Lanczos trabaja en todo momento con la matriz original (y no tiene ningún problema en que se conozca sólo implícitamente), mientras que Householder en bruto quiere modificar la matriz durante el cálculo (aunque eso se puede evitar). ${\displaystyle A}$
• Cada iteración del algoritmo de Lanczos produce otra columna de la matriz de transformación final , mientras que una iteración de Householder produce otro factor en una factorización unitaria de . Sin embargo, cada factor está determinado por un único vector, por lo que los requisitos de almacenamiento son los mismos para ambos algoritmos y se pueden calcular en el tiempo. ${\displaystyle V}$ ${\displaystyle Q_{1}Q_{2}\dots Q_{n}}$ ${\displaystyle V}$ ${\displaystyle V=Q_{1}Q_{2}\dots Q_{n}}$ ${\displaystyle O(n^{3})}$
• El jefe de familia es numéricamente estable, mientras que el crudo Lanczos no lo es.
• Lanczos es muy paralelo, con sólo puntos de sincronización (los cálculos de y ). El jefe de familia es menos paralelo y tiene una secuencia de cantidades escalares calculadas, cada una de las cuales depende de la cantidad anterior en la secuencia. ${\displaystyle O(n)}$ ${\displaystyle \alpha _{j}}$ ${\displaystyle \beta _{j}}$ ${\displaystyle O(n^{2})}$

Derivación del algoritmo

Hay varias líneas de razonamiento que conducen al algoritmo de Lanczos.

Un método de poder más previsor

El método de potencia para encontrar el valor propio de mayor magnitud y el vector propio correspondiente de una matriz es aproximadamente ${\displaystyle A}$

1. Elija un vector aleatorio . ${\displaystyle u_{1}\neq 0}$
2. Para (hasta que la dirección de haya convergido) haga: ${\displaystyle j\geqslant 1}$ ${\displaystyle u_{j}}$
   1. Dejar ${\displaystyle u_{j+1}'=Au_{j}.}$
   2. Dejar ${\displaystyle u_{j+1}=u_{j+1}'/\|u_{j+1}'\|.}$

• En el límite grande, se aproxima al vector propio normado correspondiente al valor propio de mayor magnitud. ${\displaystyle j}$ ${\displaystyle u_{j}}$

Una crítica que se puede plantear a este método es que es un desperdicio: consume mucho trabajo (los productos matriz-vector en el paso 2.1) extrayendo información de la matriz , pero presta atención sólo al último resultado; Las implementaciones suelen utilizar la misma variable para todos los vectores , y cada nueva iteración sobrescribe los resultados de la anterior. En su lugar, puede ser conveniente conservar todos los resultados intermedios y organizar los datos. ${\displaystyle A}$ ${\displaystyle u_{j}}$

Una información que trivialmente está disponible a partir de los vectores es una cadena de subespacios de Krylov . Una forma de afirmar que sin introducir conjuntos en el algoritmo es afirmar que calcula ${\displaystyle u_{j}}$

un subconjunto de una base tal que para todos y cada uno ${\displaystyle \{v_{j}\}_{j=1}^{m}}$ ${\displaystyle \mathbb {C} ^{n}}$ ${\displaystyle Ax\in \operatorname {span} (v_{1},\dotsc ,v_{j+1})}$ ${\displaystyle x\in \operatorname {span} (v_{1},\dotsc ,v_{j})}$ ${\displaystyle 1\leqslant j<m;}$

esto se satisface trivialmente siempre que sea linealmente independiente de (y en el caso de que exista tal dependencia, entonces se puede continuar la secuencia eligiendo como un vector arbitrario linealmente independiente de ). Sin embargo, es probable que una base que contenga los vectores esté numéricamente mal condicionada , ya que esta secuencia de vectores está diseñada para converger a un vector propio de . Para evitar eso, se puede combinar la iteración de poder con un proceso de Gram-Schmidt , para producir en su lugar una base ortonormal de estos subespacios de Krylov. ${\displaystyle v_{j}=u_{j}}$ ${\displaystyle u_{j}}$ ${\displaystyle u_{1},\dotsc ,u_{j-1}}$ ${\displaystyle v_{j}}$ ${\displaystyle u_{1},\dotsc ,u_{j-1}}$ ${\displaystyle u_{j}}$ ${\displaystyle A}$

1. Elija un vector aleatorio de norma euclidiana . Dejar . ${\displaystyle u_{1}}$ ${\displaystyle 1}$ ${\displaystyle v_{1}=u_{1}}$
2. Para hacer: ${\displaystyle j=1,\dotsc ,m-1}$
   1. Dejar . ${\displaystyle u_{j+1}'=Au_{j}}$
   2. Para todos dejar . (Estas son las coordenadas con respecto a los vectores base ). ${\displaystyle k=1,\dotsc ,j}$ ${\displaystyle g_{k,j}=v_{k}^{*}u_{j+1}'}$ ${\displaystyle Au_{j}=u_{j+1}'}$ ${\displaystyle v_{1},\dotsc ,v_{j}}$
   3. Dejar . (Cancele el componente que está en .) ${\displaystyle w_{j+1}=u_{j+1}'-\sum _{k=1}^{j}g_{k,j}v_{k}}$ ${\displaystyle u_{j+1}'}$ ${\displaystyle \operatorname {span} (v_{1},\dotsc ,v_{j})}$
   4. Si entonces dejamos y , ${\displaystyle w_{j+1}\neq 0}$ ${\displaystyle u_{j+1}=u_{j+1}'/\|u_{j+1}'\|}$ ${\displaystyle v_{j+1}=w_{j+1}/\|w_{j+1}\|}$

      de lo contrario, elija como un vector arbitrario de norma euclidiana que sea ortogonal a todos . ${\displaystyle u_{j+1}=v_{j+1}}$ ${\displaystyle 1}$ ${\displaystyle v_{1},\dotsc ,v_{j}}$

La relación entre los vectores de iteración de potencia y los vectores ortogonales es que ${\displaystyle u_{j}}$ ${\displaystyle v_{j}}$

${\displaystyle Au_{j}=\|u_{j+1}'\|u_{j+1}=u_{j+1}'=w_{j+1}+\sum _{k=1}^{j}g_{k,j}v_{k}=\|w_{j+1}\|v_{j+1}+\sum _{k=1}^{j}g_{k,j}v_{k}}$.

Aquí se puede observar que en realidad no necesitamos los vectores para calcularlos , porque y por lo tanto la diferencia entre y es en , que se cancela mediante el proceso de ortogonalización. Así, la misma base para la cadena de subespacios de Krylov se calcula mediante ${\displaystyle u_{j}}$ ${\displaystyle v_{j}}$ ${\displaystyle u_{j}-v_{j}\in \operatorname {span} (v_{1},\dotsc ,v_{j-1})}$ ${\displaystyle u_{j+1}'=Au_{j}}$ ${\displaystyle w_{j+1}'=Av_{j}}$ ${\displaystyle \operatorname {span} (v_{1},\dotsc ,v_{j})}$

1. Elija un vector aleatorio de norma euclidiana . ${\displaystyle v_{1}}$ ${\displaystyle 1}$
2. Para hacer: ${\displaystyle j=1,\dotsc ,m-1}$
   1. Dejar . ${\displaystyle w_{j+1}'=Av_{j}}$
   2. Para todos dejar . ${\displaystyle k=1,\dotsc ,j}$ ${\displaystyle h_{k,j}=v_{k}^{*}w_{j+1}'}$
   3. Dejar . ${\displaystyle w_{j+1}=w_{j+1}'-\sum _{k=1}^{j}h_{k,j}v_{k}}$
   4. Dejar . ${\displaystyle h_{j+1,j}=\|w_{j+1}\|}$
   5. Si entonces lo dejamos , ${\displaystyle h_{j+1,j}\neq 0}$ ${\displaystyle v_{j+1}=w_{j+1}/h_{j+1,j}}$

      de lo contrario, elija como un vector arbitrario de norma euclidiana que sea ortogonal a todos . ${\displaystyle v_{j+1}}$ ${\displaystyle 1}$ ${\displaystyle v_{1},\dotsc ,v_{j}}$

A priori los coeficientes satisfacen ${\displaystyle h_{k,j}}$

${\displaystyle Av_{j}=\sum _{k=1}^{j+1}h_{k,j}v_{k}}$para todos ; ${\displaystyle j<m}$

La definición puede parecer un poco extraña, pero se ajusta al patrón general ya que ${\displaystyle h_{j+1,j}=\|w_{j+1}\|}$ ${\displaystyle h_{k,j}=v_{k}^{*}w_{j+1}'}$

${\displaystyle v_{j+1}^{*}w_{j+1}'=v_{j+1}^{*}w_{j+1}=\|w_{j+1}\|v_{j+1}^{*}v_{j+1}=\|w_{j+1}\|.}$

Debido a que los vectores de iteración de potencia que se eliminaron de esta recursión satisfacen los vectores y los coeficientes contienen suficiente información para que todo pueda calcularse, no se perdió nada al cambiar los vectores. (De hecho, resulta que los datos recopilados aquí dan aproximaciones significativamente mejores del valor propio más grande que las que se obtienen con un número igual de iteraciones en el método de potencia, aunque eso no es necesariamente obvio en este momento). ${\displaystyle u_{j}}$ ${\displaystyle u_{j}\in \operatorname {span} (v_{1},\ldots ,v_{j}),}$ ${\displaystyle \{v_{j}\}_{j=1}^{m}}$ ${\displaystyle h_{k,j}}$ ${\displaystyle A}$ ${\displaystyle u_{1},\ldots ,u_{m}}$

Este último procedimiento es la iteración de Arnoldi . El algoritmo de Lanczos surge entonces como la simplificación que se obtiene al eliminar pasos de cálculo que resultan triviales cuando es hermitiano; en particular, la mayoría de los coeficientes resultan ser cero. ${\displaystyle A}$ ${\displaystyle h_{k,j}}$

Elementalmente, si es hermitiano entonces ${\displaystyle A}$

${\displaystyle h_{k,j}=v_{k}^{*}w_{j+1}'=v_{k}^{*}Av_{j}=v_{k}^{*}A^{*}v_{j}=(Av_{k})^{*}v_{j}.}$

Porque sabemos que , y dado que por construcción es ortogonal a este subespacio, este producto interno debe ser cero. (Esta es esencialmente también la razón por la cual a las sucesiones de polinomios ortogonales siempre se les puede dar una relación de recurrencia de tres términos ). ${\displaystyle k<j-1}$ ${\displaystyle Av_{k}\in \operatorname {span} (v_{1},\ldots ,v_{j-1})}$ ${\displaystyle v_{j}}$ ${\displaystyle k=j-1}$

${\displaystyle h_{j-1,j}=(Av_{j-1})^{*}v_{j}={\overline {v_{j}^{*}Av_{j-1}}}={\overline {h_{j,j-1}}}=h_{j,j-1}}$

ya que este último es real por ser la norma de un vector. para uno obtiene ${\displaystyle k=j}$

${\displaystyle h_{j,j}=(Av_{j})^{*}v_{j}={\overline {v_{j}^{*}Av_{j}}}={\overline {h_{j,j}}},}$

lo que significa que esto también es real.

De manera más abstracta, si es la matriz con columnas entonces los números pueden identificarse como elementos de la matriz , y para la matriz es Hessenberg superior . Desde ${\displaystyle V}$ ${\displaystyle v_{1},\ldots ,v_{m}}$ ${\displaystyle h_{k,j}}$ ${\displaystyle H=V^{*}AV}$ ${\displaystyle h_{k,j}=0}$ ${\displaystyle k>j+1;}$ ${\displaystyle H}$

${\displaystyle H^{*}=\left(V^{*}AV\right)^{*}=V^{*}A^{*}V=V^{*}AV=H}$

la matriz es hermitiana. Esto implica que también es Hessenberg inferior, por lo que en realidad debe ser tridiagional. Al ser hermitiana, su diagonal principal es real, y dado que su primera subdiagonal es real por construcción, lo mismo ocurre con su primera superdiagonal. Por lo tanto, es una matriz real y simétrica: la matriz de la especificación del algoritmo de Lanczos. ${\displaystyle H}$ ${\displaystyle H}$ ${\displaystyle H}$ ${\displaystyle T}$

Aproximación simultánea de valores propios extremos.

Una forma de caracterizar los vectores propios de una matriz hermitiana es como puntos estacionarios del cociente de Rayleigh. ${\displaystyle A}$

${\displaystyle r(x)={\frac {x^{*}Ax}{x^{*}x}},\qquad x\in \mathbb {C} ^{n}.}$

En particular, el valor propio más grande es el máximo global de y el valor propio más pequeño es el mínimo global de . ${\displaystyle \lambda _{\max }}$ ${\displaystyle r}$ ${\displaystyle \lambda _{\min }}$ ${\displaystyle r}$

Dentro de un subespacio de baja dimensión puede ser factible ubicar el máximo y el mínimo de . Repetir eso para una cadena creciente produce dos secuencias de vectores: y tal que y ${\displaystyle {\mathcal {L}}}$ ${\displaystyle \mathbb {C} ^{n}}$ ${\displaystyle x}$ ${\displaystyle y}$ ${\displaystyle r}$ ${\displaystyle {\mathcal {L}}_{1}\subset {\mathcal {L}}_{2}\subset \cdots }$ ${\displaystyle x_{1},x_{2},\ldots }$ ${\displaystyle y_{1},y_{2},\dotsc }$ ${\displaystyle x_{j},y_{j}\in {\mathcal {L}}_{j}}$

${\displaystyle {\begin{aligned}r(x_{1})&\leqslant r(x_{2})\leqslant \cdots \leqslant \lambda _{\max }\\r(y_{1})&\geqslant r(y_{2})\geqslant \cdots \geqslant \lambda _{\min }\end{aligned}}}$

Entonces surge la pregunta de cómo elegir los subespacios para que estas secuencias converjan a un ritmo óptimo.

Desde , la dirección óptima en la que buscar valores mayores de es la del gradiente , y de la misma manera desde la dirección óptima en la que buscar valores más pequeños es la del gradiente negativo . En general ${\displaystyle x_{j}}$ ${\displaystyle r}$ ${\displaystyle \nabla r(x_{j})}$ ${\displaystyle y_{j}}$ ${\displaystyle r}$ ${\displaystyle -\nabla r(y_{j})}$

${\displaystyle \nabla r(x)={\frac {2}{x^{*}x}}(Ax-r(x)x),}$

por lo tanto, las direcciones de interés son bastante fáciles de calcular en aritmética matricial, pero si uno desea mejorar ambas , entonces hay dos nuevas direcciones a tener en cuenta: y desde y pueden ser vectores linealmente independientes (de hecho, son cercanos a ortogonales) , en general no se puede esperar y ser paralelo. No es necesario aumentar la dimensión de by en cada paso si se consideran subespacios de Krylov, porque entonces para todos , en particular para ambos y . ${\displaystyle x_{j}}$ ${\displaystyle y_{j}}$ ${\displaystyle Ax_{j}}$ ${\displaystyle Ay_{j};}$ ${\displaystyle x_{j}}$ ${\displaystyle y_{j}}$ ${\displaystyle Ax_{j}}$ ${\displaystyle Ay_{j}}$ ${\displaystyle {\mathcal {L}}_{j}}$ ${\displaystyle 2}$ ${\displaystyle \{{\mathcal {L}}_{j}\}_{j=1}^{m}}$ ${\displaystyle Az\in {\mathcal {L}}_{j+1}}$ ${\displaystyle z\in {\mathcal {L}}_{j},}$ ${\displaystyle z=x_{j}}$ ${\displaystyle z=y_{j}}$

En otras palabras, podemos comenzar con algún vector inicial arbitrario para construir los espacios vectoriales. ${\displaystyle x_{1}=y_{1},}$

${\displaystyle {\mathcal {L}}_{j}=\operatorname {span} (x_{1},Ax_{1},\ldots ,A^{j-1}x_{1})}$

y luego buscar tal que ${\displaystyle x_{j},y_{j}\in {\mathcal {L}}_{j}}$

${\displaystyle r(x_{j})=\max _{z\in {\mathcal {L}}_{j}}r(z)\qquad {\text{and}}\qquad r(y_{j})=\min _{z\in {\mathcal {L}}_{j}}r(z).}$

Dado que la iteración del método de potencia pertenece a la enésima, se deduce que una iteración para producir y no puede converger más lentamente que la del método de potencia, y logrará más al aproximar ambos extremos de valores propios. Para el subproblema de optimización en algunos , es conveniente tener una base ortonormal para este espacio vectorial. De esta manera nos encontramos nuevamente con el problema de calcular iterativamente dicha base para la secuencia de subespacios de Krylov. ${\displaystyle j}$ ${\displaystyle u_{j}}$ ${\displaystyle {\mathcal {L}}_{j},}$ ${\displaystyle x_{j}}$ ${\displaystyle y_{j}}$ ${\displaystyle r}$ ${\displaystyle {\mathcal {L}}_{j}}$ ${\displaystyle \{v_{1},\ldots ,v_{j}\}}$

Convergencia y otras dinámicas

Al analizar la dinámica del algoritmo, es conveniente tomar los valores propios y vectores propios de como dados, aunque no sean conocidos explícitamente por el usuario. Para arreglar la notación, sean los valores propios (se sabe que todos son reales y, por lo tanto, es posible ordenarlos) y sea un conjunto ortonormal de vectores propios tales que para todos . ${\displaystyle A}$ ${\displaystyle \lambda _{1}\geqslant \lambda _{2}\geqslant \dotsb \geqslant \lambda _{n}}$ ${\displaystyle z_{1},\dotsc ,z_{n}}$ ${\displaystyle Az_{k}=\lambda _{k}z_{k}}$ ${\displaystyle k=1,\dotsc ,n}$

También es conveniente fijar una notación para los coeficientes del vector de Lanczos inicial con respecto a esta base propia; vamos para todos , para que . Un vector inicial al que se le ha agotado algún componente propio retrasará la convergencia al valor propio correspondiente, y aunque esto simplemente aparece como un factor constante en los límites del error, el agotamiento sigue siendo indeseable. Una técnica común para evitar ser afectado constantemente es seleccionar dibujando primero los elementos aleatoriamente de acuerdo con la misma distribución normal con media y luego reescalar el vector a la norma . Antes del cambio de escala, esto hace que los coeficientes también sean variables estocásticas independientes distribuidas normalmente de la misma distribución normal (ya que el cambio de coordenadas es unitario), y después del cambio de escala el vector tendrá una distribución uniforme en la esfera unitaria en . Esto permite limitar la probabilidad de que, por ejemplo , . ${\displaystyle v_{1}}$ ${\displaystyle d_{k}=z_{k}^{*}v_{1}}$ ${\displaystyle k=1,\dotsc ,n}$ ${\displaystyle \textstyle v_{1}=\sum _{k=1}^{n}d_{k}z_{k}}$ ${\displaystyle v_{1}}$ ${\displaystyle v_{1}}$ ${\displaystyle 0}$ ${\displaystyle 1}$ ${\displaystyle d_{k}}$ ${\displaystyle (d_{1},\dotsc ,d_{n})}$ ${\displaystyle \mathbb {C} ^{n}}$ ${\displaystyle |d_{1}|<\varepsilon }$

El hecho de que el algoritmo de Lanczos sea independiente de las coordenadas (las operaciones solo analizan los productos internos de los vectores, nunca los elementos individuales de los vectores) hace que sea fácil construir ejemplos con una estructura propia conocida para ejecutar el algoritmo: crear una matriz diagonal con los valores propios deseados. en diagonal; Siempre que el vector inicial tenga suficientes elementos distintos de cero, el algoritmo generará una matriz simétrica tridiagonal general como . ${\displaystyle A}$ ${\displaystyle v_{1}}$ ${\displaystyle T}$

Teoría de la convergencia de Kaniel-Paige

Después de los pasos de iteración del algoritmo de Lanczos, es una matriz simétrica real, que de manera similar a la anterior tiene valores propios. Por convergencia se entiende principalmente la convergencia de to (y la convergencia simétrica de to ) a medida que crece, y en segundo lugar la convergencia de algún rango de valores propios. de a sus homólogos de . La convergencia del algoritmo de Lanczos es a menudo órdenes de magnitud más rápida que la del algoritmo de iteración de potencia. ^{[9] : 477} ${\displaystyle m}$ ${\displaystyle T}$ ${\displaystyle m\times m}$ ${\displaystyle m}$ ${\displaystyle \theta _{1}\geqslant \theta _{2}\geqslant \dots \geqslant \theta _{m}.}$ ${\displaystyle \theta _{1}}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle \theta _{m}}$ ${\displaystyle \lambda _{n}}$ ${\displaystyle m}$ ${\displaystyle \theta _{1},\ldots ,\theta _{k}}$ ${\displaystyle T}$ ${\displaystyle \lambda _{1},\ldots ,\lambda _{k}}$ ${\displaystyle A}$

Los límites provienen de la interpretación anterior de los valores propios como valores extremos del cociente de Rayleigh . Dado que es a priori el máximo de en el conjunto de mientras que es simplemente el máximo en un subespacio de Krylov de una dimensión, obtenemos trivialmente . Por el contrario, cualquier punto en ese subespacio de Krylov proporciona un límite inferior para , por lo que si se puede exhibir un punto que es pequeño, entonces esto proporciona un límite estrecho para . ${\displaystyle \theta _{1}}$ ${\displaystyle r(x)}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle r}$ ${\displaystyle \mathbb {C} ^{n},}$ ${\displaystyle \theta _{1}}$ ${\displaystyle m}$ ${\displaystyle \lambda _{1}\geqslant \theta _{1}}$ ${\displaystyle x}$ ${\displaystyle r(x)}$ ${\displaystyle \theta _{1}}$ ${\displaystyle \lambda _{1}-r(x)}$ ${\displaystyle \theta _{1}}$

La dimensión del subespacio de Krylov es ${\displaystyle m}$

${\displaystyle \operatorname {span} \left\{v_{1},Av_{1},A^{2}v_{1},\ldots ,A^{m-1}v_{1}\right\},}$

por lo que cualquier elemento del mismo puede expresarse como para algún polinomio de grado como máximo ; los coeficientes de ese polinomio son simplemente los coeficientes en la combinación lineal de los vectores . El polinomio que queremos tendrá coeficientes reales, pero por el momento también debemos permitir coeficientes complejos, y escribiremos para el polinomio obtenido mediante la conjugación compleja de todos los coeficientes de . En esta parametrización del subespacio de Krylov, tenemos ${\displaystyle p(A)v_{1}}$ ${\displaystyle p}$ ${\displaystyle m-1}$ ${\displaystyle v_{1},Av_{1},A^{2}v_{1},\ldots ,A^{m-1}v_{1}}$ ${\displaystyle p^{*}}$ ${\displaystyle p}$

${\displaystyle r(p(A)v_{1})={\frac {(p(A)v_{1})^{*}Ap(A)v_{1}}{(p(A)v_{1})^{*}p(A)v_{1}}}={\frac {v_{1}^{*}p(A)^{*}Ap(A)v_{1}}{v_{1}^{*}p(A)^{*}p(A)v_{1}}}={\frac {v_{1}^{*}p^{*}(A^{*})Ap(A)v_{1}}{v_{1}^{*}p^{*}(A^{*})p(A)v_{1}}}={\frac {v_{1}^{*}p^{*}(A)Ap(A)v_{1}}{v_{1}^{*}p^{*}(A)p(A)v_{1}}}}$

Usando ahora la expresión para como combinación lineal de vectores propios, obtenemos ${\displaystyle v_{1}}$

${\displaystyle Av_{1}=A\sum _{k=1}^{n}d_{k}z_{k}=\sum _{k=1}^{n}d_{k}\lambda _{k}z_{k}}$

y más en general

${\displaystyle q(A)v_{1}=\sum _{k=1}^{n}d_{k}q(\lambda _{k})z_{k}}$

para cualquier polinomio . ${\displaystyle q}$

De este modo

${\displaystyle \lambda _{1}-r(p(A)v_{1})=\lambda _{1}-{\frac {v_{1}^{*}\sum _{k=1}^{n}d_{k}p^{*}(\lambda _{k})\lambda _{k}p(\lambda _{k})z_{k}}{v_{1}^{*}\sum _{k=1}^{n}d_{k}p^{*}(\lambda _{k})p(\lambda _{k})z_{k}}}=\lambda _{1}-{\frac {\sum _{k=1}^{n}|d_{k}|^{2}\lambda _{k}p(\lambda _{k})^{*}p(\lambda _{k})}{\sum _{k=1}^{n}|d_{k}|^{2}p(\lambda _{k})^{*}p(\lambda _{k})}}={\frac {\sum _{k=1}^{n}|d_{k}|^{2}(\lambda _{1}-\lambda _{k})\left|p(\lambda _{k})\right|^{2}}{\sum _{k=1}^{n}|d_{k}|^{2}\left|p(\lambda _{k})\right|^{2}}}.}$

Una diferencia clave entre numerador y denominador aquí es que el término desaparece en el numerador, pero no en el denominador. Por lo tanto, si se puede elegir que sea grande pero pequeño en todos los demás valores propios, se obtendrá un límite estricto del error . ${\displaystyle k=1}$ ${\displaystyle p}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle \lambda _{1}-\theta _{1}}$

Dado que tiene muchos más valores propios que coeficientes, esto puede parecer una tarea difícil, pero una forma de cumplirlo es utilizar polinomios de Chebyshev . Escribiendo para el polinomio de grado de Chebyshev de primer tipo (el que satisface para todos ), tenemos un polinomio que permanece en el rango del intervalo conocido pero crece rápidamente fuera de él. Con algo de escala del argumento, podemos hacer que mapee todos los valores propios excepto en . Dejar ${\displaystyle A}$ ${\displaystyle p}$ ${\displaystyle c_{k}}$ ${\displaystyle k}$ ${\displaystyle c_{k}(\cos x)=\cos(kx)}$ ${\displaystyle x}$ ${\displaystyle [-1,1]}$ ${\displaystyle [-1,1]}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle [-1,1]}$

${\displaystyle p(x)=c_{m-1}\left({\frac {2x-\lambda _{2}-\lambda _{n}}{\lambda _{2}-\lambda _{n}}}\right)}$

(en caso de que utilice en su lugar el valor propio más grande estrictamente menor que ), entonces el valor máximo de for es y el valor mínimo es , por lo que ${\displaystyle \lambda _{2}=\lambda _{1}}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle |p(\lambda _{k})|^{2}}$ ${\displaystyle k\geqslant 2}$ ${\displaystyle 1}$ ${\displaystyle 0}$

${\displaystyle \lambda _{1}-\theta _{1}\leqslant \lambda _{1}-r(p(A)v_{1})={\frac {\sum _{k=2}^{n}|d_{k}|^{2}(\lambda _{1}-\lambda _{k})|p(\lambda _{k})|^{2}}{\sum _{k=1}^{n}|d_{k}|^{2}|p(\lambda _{k})|^{2}}}\leqslant {\frac {\sum _{k=2}^{n}|d_{k}|^{2}(\lambda _{1}-\lambda _{k})}{|d_{1}|^{2}|p(\lambda _{1})|^{2}}}\leqslant {\frac {(\lambda _{1}-\lambda _{n})\sum _{k=2}^{n}|d_{k}|^{2}}{|p(\lambda _{1})|^{2}|d_{1}|^{2}}}.}$

Además

${\displaystyle p(\lambda _{1})=c_{m-1}\left({\frac {2\lambda _{1}-\lambda _{2}-\lambda _{n}}{\lambda _{2}-\lambda _{n}}}\right)=c_{m-1}\left(2{\frac {\lambda _{1}-\lambda _{2}}{\lambda _{2}-\lambda _{n}}}+1\right);}$

la cantidad

${\displaystyle \rho ={\frac {\lambda _{1}-\lambda _{2}}{\lambda _{2}-\lambda _{n}}}}$

(es decir, la relación entre el primer espacio propio y el diámetro del resto del espectro ) es, por tanto, de importancia clave para la tasa de convergencia. tambien escribiendo

${\displaystyle R=e^{\operatorname {arcosh} (1+2\rho )}=1+2\rho +2{\sqrt {\rho ^{2}+\rho }},}$

podemos concluir que

${\displaystyle {\begin{aligned}\lambda _{1}-\theta _{1}&\leqslant {\frac {(\lambda _{1}-\lambda _{n})\left(1-|d_{1}|^{2}\right)}{c_{m-1}(2\rho +1)^{2}|d_{1}|^{2}}}\\[6pt]&={\frac {1-|d_{1}|^{2}}{|d_{1}|^{2}}}(\lambda _{1}-\lambda _{n}){\frac {1}{\cosh ^{2}((m-1)\operatorname {arcosh} (1+2\rho ))}}\\[6pt]&={\frac {1-|d_{1}|^{2}}{|d_{1}|^{2}}}(\lambda _{1}-\lambda _{n}){\frac {4}{\left(R^{m-1}+R^{-(m-1)}\right)^{2}}}\\[6pt]&\leqslant 4{\frac {1-|d_{1}|^{2}}{|d_{1}|^{2}}}(\lambda _{1}-\lambda _{n})R^{-2(m-1)}\end{aligned}}}$

Por lo tanto, la tasa de convergencia está controlada principalmente por , ya que este límite se reduce en un factor por cada iteración adicional. ${\displaystyle R}$ ${\displaystyle R^{-2}}$

A modo de comparación, se puede considerar cómo depende la tasa de convergencia del método de potencia , pero dado que el método de potencia es principalmente sensible al cociente entre valores absolutos de los valores propios, necesitamos que la brecha propia entre y sea la dominante. Bajo esa restricción, el caso que más favorece el método de potencia es , así que considérelo. Al final del método de potencia, el vector de iteración: ${\displaystyle \rho }$ ${\displaystyle |\lambda _{n}|\leqslant |\lambda _{2}|}$ ${\displaystyle \lambda _{1}}$ ${\displaystyle \lambda _{2}}$ ${\displaystyle \lambda _{n}=-\lambda _{2}}$

${\displaystyle u=(1-t^{2})^{1/2}z_{1}+tz_{2}\approx z_{1}+tz_{2},}$^{[nota 1]}

donde cada nueva iteración multiplica efectivamente la amplitud por ${\displaystyle z_{2}}$ ${\displaystyle t}$

${\displaystyle {\frac {\lambda _{2}}{\lambda _{1}}}={\frac {\lambda _{2}}{\lambda _{2}+(\lambda _{1}-\lambda _{2})}}={\frac {1}{1+{\frac {\lambda _{1}-\lambda _{2}}{\lambda _{2}}}}}={\frac {1}{1+2\rho }}.}$

La estimación del valor propio más grande es entonces

${\displaystyle u^{*}Au=(1-t^{2})\lambda _{1}+t^{2}\lambda _{2},}$

por lo que el límite anterior para la tasa de convergencia del algoritmo de Lanczos debe compararse con

${\displaystyle \lambda _{1}-u^{*}Au=(\lambda _{1}-\lambda _{2})t^{2},}$

que se reduce en un factor de para cada iteración. Por tanto, la diferencia se reduce a la que existe entre y . En la región, este último es más parecido a , y funciona como lo haría el método de potencia con una brecha propia dos veces mayor; una mejora notable. Sin embargo, el caso más desafiante es aquel en el que se produce una mejora aún mayor en la brecha propia; la región es donde el algoritmo de Lanczos en cuanto a convergencia logra la menor mejora con respecto al método de potencia. ${\displaystyle (1+2\rho )^{-2}}$ ${\displaystyle 1+2\rho }$ ${\displaystyle R=1+2\rho +2{\sqrt {\rho ^{2}+\rho }}}$ ${\displaystyle \rho \gg 1}$ ${\displaystyle 1+4\rho }$ ${\displaystyle \rho \ll 1,}$ ${\displaystyle R\approx 1+2{\sqrt {\rho }}}$ ${\displaystyle \rho \gg 1}$

Estabilidad numérica

La estabilidad significa cuánto se verá afectado el algoritmo (es decir, si producirá un resultado aproximado cercano al original) si se introducen y acumulan pequeños errores numéricos. La estabilidad numérica es el criterio central para juzgar la utilidad de implementar un algoritmo en una computadora con redondeo.

Para el algoritmo de Lanczos, se puede demostrar que con aritmética exacta , el conjunto de vectores construye una base ortonormal , y los valores propios/vectores resueltos son buenas aproximaciones a los de la matriz original. Sin embargo, en la práctica (dado que los cálculos se realizan en aritmética de coma flotante donde la inexactitud es inevitable), la ortogonalidad se pierde rápidamente y en algunos casos el nuevo vector podría incluso depender linealmente del conjunto ya construido. Como resultado, algunos de los valores propios de la matriz tridiagonal resultante pueden no ser aproximaciones a la matriz original. Por tanto, el algoritmo de Lanczos no es muy estable. ${\displaystyle v_{1},v_{2},\cdots ,v_{m+1}}$

Los usuarios de este algoritmo deben poder encontrar y eliminar esos valores propios "falsos". Las implementaciones prácticas del algoritmo Lanczos van en tres direcciones para combatir este problema de estabilidad: ^{[6] [7]}

1. Prevenir la pérdida de ortogonalidad,
2. Recuperar la ortogonalidad después de generar la base.
3. Una vez identificados los valores propios buenos y "espurios", elimine los espurios.

Variaciones

Existen variaciones del algoritmo de Lanczos donde los vectores involucrados son matrices altas y estrechas en lugar de vectores y las constantes de normalización son matrices cuadradas pequeñas. Estos se denominan algoritmos Lanczos de "bloque" y pueden ser mucho más rápidos en computadoras con una gran cantidad de registros y tiempos de recuperación de memoria prolongados.

Muchas implementaciones del algoritmo Lanczos se reinician después de un cierto número de iteraciones. Una de las variaciones reiniciadas más influyentes es el método Lanczos reiniciado implícitamente, ^[10] que se implementa en ARPACK . ^[11] Esto ha llevado a una serie de otras variaciones reiniciadas, como la bidiagonalización reiniciada de Lanczos. ^[12] Otra variación reiniciada exitosa es el método Thick-Restart Lanczos, ^[13] que se ha implementado en un paquete de software llamado TRLan. ^[14]

Espacio nulo sobre un campo finito

En 1995, Peter Montgomery publicó un algoritmo, basado en el algoritmo de Lanczos, para encontrar elementos del espacio nulo de una matriz dispersa grande sobre GF(2) ; Dado que el conjunto de personas interesadas en grandes matrices dispersas sobre campos finitos y el conjunto de personas interesadas en grandes problemas de valores propios apenas se superponen, esto a menudo también se denomina algoritmo de Lanczos en bloque sin causar una confusión irrazonable. ^{[ cita necesaria ]}

Aplicaciones

Los algoritmos de Lanczos son muy atractivos porque la multiplicación por es la única operación lineal a gran escala. Dado que los motores de recuperación de texto de términos ponderados implementan precisamente esta operación, el algoritmo Lanczos se puede aplicar de manera eficiente a documentos de texto (consulte indexación semántica latente ). Los vectores propios también son importantes para métodos de clasificación a gran escala, como el algoritmo HITS desarrollado por Jon Kleinberg o el algoritmo PageRank utilizado por Google. ${\displaystyle A\,}$

Los algoritmos de Lanczos también se utilizan en física de la materia condensada como método para resolver hamiltonianos de sistemas de electrones fuertemente correlacionados , ^[15] así como en códigos de modelos de capas en física nuclear . ^[16]

Implementaciones

La biblioteca NAG contiene varias rutinas ^[17] para la solución de sistemas lineales a gran escala y problemas propios que utilizan el algoritmo de Lanczos.

MATLAB y GNU Octave vienen con ARPACK integrado. Tanto las matrices almacenadas como las implícitas se pueden analizar mediante la función eigs() (Matlab/Octave).

De manera similar, en Python , el paquete SciPy tiene scipy.sparse.linalg.eigsh, que también es un contenedor para las funciones SSEUPD y DSEUPD de ARPACK que utilizan el método Lanczos reiniciado implícitamente.

Una implementación en Matlab del algoritmo Lanczos (tenga en cuenta los problemas de precisión) está disponible como parte del paquete Matlab de propagación de creencias gaussianas. La biblioteca de filtrado colaborativo GraphLab ^[18] incorpora una implementación paralela a gran escala del algoritmo Lanczos (en C++) para multinúcleo.

La biblioteca PRIMME también implementa un algoritmo similar a Lanczos.

Notas

1. No es necesario que ambos coeficientes sean reales, pero la fase tiene poca importancia. Tampoco es necesario que los componentes de otros vectores propios hayan desaparecido por completo, pero se reducen al menos tan rápido como los de , así describe el peor de los casos. ${\displaystyle z_{2}}$ ${\displaystyle u\approx z_{1}+tz_{2}}$

Referencias

1. Lanczos, C. (1950). "Un método de iteración para la solución del problema de valores propios de operadores integrales y diferenciales lineales" (PDF) . Revista de Investigación de la Oficina Nacional de Normas . 45 (4): 255–282. doi : 10.6028/jres.045.026 .
2. Ojalvo, IU; Newman, M. (1970). "Modos de vibración de grandes estructuras mediante un método automático de reducción matricial". Revista AIAA . 8 (7): 1234-1239. Código bibliográfico : 1970AIAAJ...8.1234N. doi : 10.2514/3.5878.
3. Paige, CC (1971). El cálculo de valores propios y vectores propios de matrices dispersas muy grandes (tesis doctoral). U. de Londres. OCLC  654214109.
4. Paige, CC (1972). "Variantes computacionales del método Lanczos para el problema propio". J.Inst. Aplicaciones de matemáticas . 10 (3): 373–381. doi :10.1093/imamat/10.3.373.
5. Ojalvo, IU (1988). "Orígenes y ventajas de los vectores Lanczos para grandes sistemas dinámicos". Proc. Sexta Conferencia de Análisis Modal (IMAC), Kissimmee, FL . págs. 489–494.
6. Cullum; Willoughby (1985). Algoritmos de Lanczos para cálculos de valores propios simétricos grandes . vol. 1.ISBN​ 0-8176-3058-9.
7. Yousef Saad (22 de junio de 1992). Métodos numéricos para problemas de valores propios grandes. ISBN 0-470-21820-7.
8. Coakley, Ed S.; Rokhlin, Vladimir (2013). "Un algoritmo rápido de divide y vencerás para calcular los espectros de matrices tridiagonales simétricas reales". Análisis Armónico Aplicado y Computacional . 34 (3): 379–414. doi :10.1016/j.acha.2012.06.003.
9. Golub, Gene H.; Préstamo de Van, Charles F. (1996). Cálculos matriciales (3. ed.). Baltimore: Universidad Johns Hopkins. Prensa. ISBN 0-8018-5413-X.
10. D.Calvetti ; L. Reichel; DC Sorensen (1994). "Un método Lanczos reiniciado implícitamente para grandes problemas de valores propios simétricos". Transacciones Electrónicas sobre Análisis Numérico . 2 : 1–21.
11. RB Lehoucq; DC Sorensen; C. Yang (1998). Guía del usuario de ARPACK: Solución de problemas de valores propios a gran escala con métodos Arnoldi reiniciados implícitamente . SIAM. doi :10.1137/1.9780898719628. ISBN 978-0-89871-407-4.
12. E. Kokiopoulou; C. Bekas; E. Gallopoulos (2004). "Cálculo de los tripletes singulares más pequeños con bidiagonalización de Lanczos reiniciada implícitamente" (PDF) . Aplica. Número. Matemáticas . 49 : 39–61. doi :10.1016/j.apnum.2003.11.011.
13. Kesheng Wu; Horst Simon (2000). "Método Lanczos de reinicio completo para grandes problemas de valores propios simétricos". Revista SIAM sobre Análisis y Aplicaciones de Matrices . 22 (2). SIAM: 602–616. doi :10.1137/S0895479898334605.
14. Kesheng Wu; Horst Simon (2001). "Paquete de software TRLan". Archivado desde el original el 1 de julio de 2007 . Consultado el 30 de junio de 2007 .
15. Chen, HY; Atkinson, Washington; Wortis, R. (julio de 2011). "Anomalía de sesgo cero inducida por trastornos en el modelo de Anderson-Hubbard: cálculos numéricos y analíticos". Revisión física B. 84 (4): 045113. arXiv : 1012.1031 . Código bibliográfico : 2011PhRvB..84d5113C. doi : 10.1103/PhysRevB.84.045113. S2CID  118722138.
16. Shimizu, Noritaka (21 de octubre de 2013). "Código de modelo de capa nuclear para computación paralela masiva, "KSHELL"". arXiv : 1310.5431 [nucl-th].
17. El grupo de algoritmos numéricos. "Índice de palabras clave: Lanczos". Manual de la biblioteca NAG, Marcos 23 . Consultado el 9 de febrero de 2012 .
18. GraphLab Archivado el 14 de marzo de 2011 en la Wayback Machine.

Lectura adicional

• Golub, Gene H .; Préstamo de Van, Charles F. (1996). "Métodos Lanczos". Cálculos matriciales . Baltimore: Prensa de la Universidad Johns Hopkins. págs. 470–507. ISBN 0-8018-5414-8.
• Ng, Andrew Y .; Zheng, Alice X.; Jordania, Michael I. (2001). "Análisis de enlaces, vectores propios y estabilidad" (PDF) . IJCAI'01 Actas de la 17ª Conferencia Internacional Conjunta sobre Inteligencia Artificial . 2 : 903–910.
• Erik Koch (2019). «Diagonalización Exacta y Método Lanczos» (PDF) . En E. Pavarini; E. Koch; S. Zhang (eds.). "Métodos de muchos cuerpos para materiales reales" . Jülich. ISBN 978-3-95806-400-3.