Articulo de referencia

Algoritmo QR

En álgebra lineal numérica , el algoritmo QR o iteración QR es un algoritmo de valores propios : es decir, un procedimiento para calcular los valores propios y los vectores prop...

En álgebra lineal numérica , el algoritmo QR o iteración QR es un algoritmo de valores propios : es decir, un procedimiento para calcular los valores propios y los vectores propios de una matriz . El algoritmo QR fue desarrollado a finales de la década de 1950 por John GF Francis y por Vera N. Kublanovskaya , trabajando de forma independiente. [ 1 ] [ 2 ] [ 3 ] La idea básica es realizar una descomposición QR , escribiendo la matriz como un producto de una matriz ortogonal y una matriz triangular superior , multiplicar los factores en orden inverso e iterar.

El algoritmo QR práctico

Formalmente, sea A una matriz real de la cual queremos calcular los valores propios, y sea A 0  := A . En el k -ésimo paso (comenzando con k = 0 ), calculamos la descomposición QR A k = Q k R k donde Q k es una matriz ortogonal (es decir, Q T = Q −1 ) y R k es una matriz triangular superior. Luego formamos A k +1 = R k Q k . Nótese que Ak+1=RkQk=Qk1QkRkQk=Qk1AkQk=QkTAkQk,{\displaystyle A_{k+1}=R_{k}Q_{k}=Q_{k}^{-1}Q_{k}R_{k}Q_{k}=Q_{k}^{-1}A_{k}Q_{k}=Q_{k}^{\mathsf {T}}A_{k}Q_{k},} Por lo tanto, todos los A k son similares y, en consecuencia, tienen los mismos valores propios. El algoritmo es numéricamente estable porque procede mediante transformaciones de similitud ortogonales .

Bajo ciertas condiciones, [ 4 ] las matrices A k convergen a una matriz triangular, la forma de Schur de A. Los autovalores de una matriz triangular se enumeran en la diagonal y se resuelve el problema de autovalores. Al comprobar la convergencia, no es práctico exigir ceros exactos, pero el teorema del círculo de Gershgorin proporciona una cota para el error.

Si las matrices convergen, entonces los valores propios a lo largo de la diagonal aparecerán según su multiplicidad geométrica. Para garantizar la convergencia, A debe ser una matriz simétrica, y para todos los valores propios distintos de ceroλ{\displaystyle \lambda }No debe existir un valor propio correspondiente.λ{\displaystyle -\lambda }. [ 5 ] Debido a que una sola iteración QR tiene un costo deO(norte3){\displaystyle {\mathcal {O}}(n^{3})}y la convergencia es lineal, el algoritmo QR estándar es extremadamente costoso de calcular, especialmente considerando que no se garantiza su convergencia. [ 6 ]

Utilizando la forma de Hessenberg

En la forma cruda anterior, las iteraciones son relativamente costosas. Esto se puede mitigar llevando primero la matriz A a la forma de Hessenberg superior (lo que cuesta103norte3+O(norte2){\textstyle {\tfrac {10}{3}}n^{3}+{\mathcal {O}}(n^{2})}operaciones aritméticas utilizando una técnica basada en la reducción de Householder ), con una secuencia finita de transformaciones de similitud ortogonales, algo parecido a una descomposición QR bilateral. [ 7 ] [ 8 ] (Para la descomposición QR, los reflectores de Householder se multiplican solo por la izquierda, pero para el caso de Hessenberg se multiplican tanto por la izquierda como por la derecha). Determinar la descomposición QR de una matriz de Hessenberg superior cuesta6norte2+O(norte){\textstyle 6n^{2}+{\mathcal {O}}(n)}operaciones aritméticas. Además, como la forma de Hessenberg ya es casi triangular superior (tiene solo una entrada distinta de cero debajo de cada diagonal), usarla como punto de partida reduce el número de pasos necesarios para la convergencia del algoritmo QR.

Si la matriz original es simétrica , entonces la matriz de Hessenberg superior también es simétrica y, por lo tanto, tridiagonal , al igual que todos los A k . En este caso, alcanzar la forma de Hessenberg cuesta43norte3+O(norte2){\textstyle {\tfrac {4}{3}}n^{3}+{\mathcal {O}}(n^{2})}operaciones aritméticas utilizando una técnica basada en la reducción de Householder. [ 7 ] [ 8 ] Determinar la descomposición QR de una matriz tridiagonal simétrica cuestaO(norte){\displaystyle {\mathcal {O}}(n)}operaciones. [ 9 ]

Fase de iteración

Si una matriz de HessenbergA{\displaystyle A}tiene elementoak,k1=0{\displaystyle a_{k,k-1}=0}para algunosk{\displaystyle k}, es decir, si uno de los elementos justo debajo de la diagonal es de hecho cero, entonces se descompone en bloques cuyos problemas de autovalores pueden resolverse por separado; un autovalor es o bien un autovalor de la submatriz de la primerak1{\displaystyle k-1}filas y columnas, o un valor propio de la submatriz de filas y columnas restantes. El propósito del paso de iteración QR es reducir una de estasak,k1{\displaystyle a_{k,k-1}}elementos de modo que efectivamente un pequeño bloque a lo largo de la diagonal se separa del resto de la matriz. En el caso de un valor propio real, eso suele ser el1×1{\displaystyle 1\times 1}bloque en la esquina inferior derecha (en cuyo caso elementoanortenorte{\displaystyle a_{nn}}sostiene ese autovalor), mientras que en el caso de un par de autovalores complejos conjugados es el2×2{\displaystyle 2\times 2}bloque en la esquina inferior derecha.

La tasa de convergencia depende de la separación entre los autovalores, por lo que un algoritmo práctico utilizará desplazamientos, ya sean explícitos o implícitos, para aumentar la separación y acelerar la convergencia. Un algoritmo QR simétrico típico aísla cada autovalor (y luego reduce el tamaño de la matriz) con solo una o dos iteraciones, lo que lo hace eficiente y robusto.

Una sola iteración con cambio explícito

Los pasos de una iteración QR con desplazamiento explícito en una matriz de Hessenberg realA{\displaystyle A}son:

  1. Elige un turnoμ{\displaystyle \mu }y restarlo de todos los elementos diagonales, produciendo la matrizAμI{\displaystyle A-\mu I}Una estrategia básica es utilizarμ=anorte,norte{\displaystyle \mu =a_{n,n}}, pero existen estrategias más refinadas que acelerarían aún más la convergencia. La idea es queμ{\displaystyle \mu }debería estar cerca de un valor propio, ya que realizar este cambio acelerará la convergencia hacia ese valor propio.
  2. Realizar una secuencia de rotaciones de GivensGRAMO1,GRAMO2,,GRAMOnorte1{\displaystyle G_{1},G_{2},\dots ,G_{n-1}}enAμI{\displaystyle A-\mu I}, dóndeGRAMOi{\displaystyle G_{i}}actos en filasi{\displaystyle i}yi+1{\displaystyle i+1}, yGRAMOi{\displaystyle G_{i}}se elige para poner la posición a cero(i+1,i){\displaystyle (i+1,i)}deGRAMOi1GRAMO1(AμI){\displaystyle G_{i-1}\dotsb G_{1}(A-\mu I)}Esto produce la matriz triangular superior.R=GRAMOnorte1GRAMO1(AμI){\displaystyle R=G_{n-1}\dotsb G_{1}(A-\mu I)}. El factor ortogonalQ{\displaystyle Q}seríaGRAMO1TGRAMO2TGRAMOnorte1T{\displaystyle G_{1}^{\mathrm {T} }G_{2}^{\mathrm {T} }\dotsb G_{n-1}^{\mathrm {T} }}, pero no es necesario ni eficiente producirlo explícitamente.
  3. Ahora multiplicaR{\displaystyle R}mediante las matrices de GivensGRAMO1T{\displaystyle G_{1}^{\mathrm {T} }},GRAMO2T{\displaystyle G_{2}^{\mathrm {T} }}, ...,GRAMOnorte1T{\displaystyle G_{n-1}^{\mathrm {T} }}a la derecha, dondeGRAMOiT{\displaystyle G_{i}^{\mathrm {T} }}en cambio actúa sobre columnasi{\displaystyle i}yi+1{\displaystyle i+1}Esto produce la matrizRQ=RGRAMO1TGRAMO2TGRAMOnorte1T{\displaystyle RQ=RG_{1}^{\mathrm {T} }G_{2}^{\mathrm {T} }\dotsb G_{n-1}^{\mathrm {T} }}, que de nuevo está en plena forma, al estilo Hessenberg.
  4. Finalmente, deshaga el cambio agregandoμ{\displaystyle \mu }a todas las entradas diagonales. El resultado esA=RQ+μI{\displaystyle A'=RQ+\mu I}. DesdeQ{\displaystyle Q}se desplaza conI{\displaystyle I}, tenemos esoA=QT(AμI)Q+μI=QTAQ{\displaystyle A'=Q^{\mathrm {T} }(A-\mu I)Q+\mu I=Q^{\mathrm {T} }AQ}.

El objetivo del cambio es modificar qué rotaciones de Givens se seleccionan.

En más detalle, la estructura de uno de estosGRAMOi{\displaystyle G_{i}}Las matrices son GRAMOi=[I0000dos00sdo0000I]{\displaystyle G_{i}={\begin{bmatrix}I&0&0&0\\0&c&-s&0\\0&s&c&0\\0&0&0&I\end{bmatrix}}} donde elI{\displaystyle I}en la esquina superior izquierda hay un(norte1)×(norte1){\displaystyle (n-1)\times (n-1)}matriz identidad y los dos escalaresdo=porqueθ{\displaystyle c=\cos \theta }ys=pecadoθ{\displaystyle s=\sin \theta }están determinados por qué ángulo de rotaciónθ{\displaystyle \theta }es apropiado para la posición cero(i+1,i){\displaystyle (i+1,i)}No es necesario exhibirθ{\displaystyle \theta }; los factoresdo{\displaystyle c}ys{\displaystyle s}se puede determinar directamente a partir de los elementos de la matrizGRAMOi{\displaystyle G_{i}}debe actuar sobre. Tampoco es necesario producir la matriz completa; multiplicación (desde la izquierda) porGRAMOi{\displaystyle G_{i}}solo afecta a las filasi{\displaystyle i}yi+1{\displaystyle i+1}, por lo que es más fácil simplemente actualizar esas dos filas en su lugar. Del mismo modo, para la multiplicación del paso 3 porGRAMOiT{\displaystyle G_{i}^{\mathrm {T} }}Desde la derecha, basta con recordari{\displaystyle i},do{\displaystyle c}, ys{\displaystyle s}.

Si se utiliza el método simpleμ=anorte,norte{\displaystyle \mu =a_{n,n}}estrategia, luego al comienzo del Paso 2 tenemos una matriz Aanorte,norteI=(××××××××××0××××00×××000×0){\displaystyle A-a_{n,n}I={\begin{pmatrix}\times &\times &\times &\times &\times \\\times &\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&\times &0\end{pmatrix}}} donde el×{\displaystyle \times }denota “podría ser cualquier cosa”. La primera rotación de GivensGRAMO1{\displaystyle G_{1}}pone a cero el(i+1,i){\displaystyle (i+1,i)}posición de esto, produciendo GRAMO1(Aanorte,norteI)=(×××××0××××0××××00×××000×0).{\displaystyle G_{1}(A-a_{n,n}I)={\begin{pmatrix}\times &\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&\times &0\end{pmatrix}}{\text{.}}} Cada nueva rotación anula otro elemento subdiagonal, aumentando así el número de ceros conocidos hasta que llegamos a H=GRAMOnorte2GRAMO1(Aanorte,norteI)=(×××××0××××00×××000hnorte1,norte1hnorte1,norte000hnorte,norte10).{\displaystyle {\begin{aligned}H&=G_{n-2}\dotsb G_{1}(A-a_{n,n}I)\\[1ex]&={\begin{pmatrix}\times &\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&h_{n-1,n-1}&h_{n-1,n}\\0&0&0&h_{n,n-1}&0\end{pmatrix}}.\end{aligned}}} La rotación finalGRAMOnorte1{\displaystyle G_{n-1}}tiene(do,s){\displaystyle (c,s)}elegido de modo queshnorte1,norte1+dohnorte,norte1=0{\displaystyle sh_{n-1,n-1}+ch_{n,n-1}=0}. Si|hnorte1,norte1||hnorte,norte1|{\displaystyle |h_{n-1,n-1}|\gg |h_{n,n-1}|}, como suele ocurrir cuando nos acercamos a la convergencia, entoncesdo1{\displaystyle c\approx 1}y|s|1{\displaystyle |s|\ll 1}Al realizar esta rotación se produce R=GRAMOnorte1GRAMOnorte2GRAMO1(Aanorte,norteI)=(×××××0××××00×××000×dohnorte1,norte0000shnorte1,norte),{\displaystyle {\begin{aligned}R&=G_{n-1}G_{n-2}\dotsb G_{1}(A-a_{n,n}I)\\[1ex]&={\begin{pmatrix}\times &\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&\times &ch_{n-1,n}\\0&0&0&0&sh_{n-1,n}\end{pmatrix}},\end{aligned}}} que es nuestra matriz triangular superior. Pero ahora llegamos al Paso 3, y necesitamos comenzar a rotar los datos entre columnas. La primera rotación actúa sobre las columnas.1{\displaystyle 1}y2{\displaystyle 2}, produciendo RGRAMO1T=(××××××××××00×××000×dohnorte1,norte0000shnorte1,norte).{\displaystyle RG_{1}^{\mathrm {T} }={\begin{pmatrix}\times &\times &\times &\times &\times \\\times &\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&\times &ch_{n-1,n}\\0&0&0&0&sh_{n-1,n}\end{pmatrix}}{\text{.}}} El patrón esperado es que cada rotación mueve algún valor distinto de cero de la diagonal hacia la subdiagonal, devolviendo la matriz a la forma de Hessenberg. Esto termina en RGRAMO1TGRAMOnorte1T=(××××××××××0××××00×××000s2hnorte1,nortedoshnorte1,norte).{\displaystyle RG_{1}^{\mathrm {T} }\dotsb G_{n-1}^{\mathrm {T} }={\begin{pmatrix}\times &\times &\times &\times &\times \\\times &\times &\times &\times &\times \\0&\times &\times &\times &\times \\0&0&\times &\times &\times \\0&0&0&-s^{2}h_{n-1,n}&csh_{n-1,n}\end{pmatrix}}{\text{.}}} Algebraicamente la forma no cambia, pero numéricamente el elemento en posición(norte,norte1){\displaystyle (n,n{-}1)}se ha acercado mucho a cero: solía haber un factors{\displaystyle s}espacio entre él y el elemento diagonal de arriba, pero ahora el espacio es más como un factors2{\displaystyle s^{2}}y otra iteración lo convertiría en factor.s4{\displaystyle s^{4}}; tenemos convergencia cuadrática. En la práctica eso significaO(1){\displaystyle O(1)}Las iteraciones por valor propio son suficientes para la convergencia y, por lo tanto, en general podemos completar enO(norte){\displaystyle O(n)}pasos QR, cada uno de los cuales hace un simpleO(norte2){\displaystyle O(n^{2})}operaciones aritméticas (o tan poco comoO(norte){\displaystyle O(n)}operaciones, en el caso de queA{\displaystyle A}es simétrico).

Visualización

Figura 1: Cómo varía la salida de una sola iteración del algoritmo QR o LR en función de su entrada.

El algoritmo QR básico se puede visualizar cuando A es una matriz simétrica definida positiva. En ese caso, A se puede representar como una elipse en dos dimensiones o un elipsoide en dimensiones superiores. La relación entre la entrada del algoritmo y una sola iteración se puede representar como en la Figura 1 (haga clic para ver una animación). Cabe destacar que el algoritmo LR se muestra junto al algoritmo QR.

Una sola iteración hace que la elipse se incline o "caiga" hacia el eje x. En el caso de que el semieje mayor de la elipse sea paralelo al eje x, una iteración de QR no hace nada. Otra situación en la que el algoritmo "no hace nada" es cuando el semieje mayor es paralelo al eje y en lugar del eje x. En ese caso, se puede pensar que la elipse se equilibra precariamente sin poder caer en ninguna dirección. En ambas situaciones, la matriz es diagonal. Una situación en la que una iteración del algoritmo "no hace nada" se denomina punto fijo . La estrategia empleada por el algoritmo es la iteración hacia un punto fijo . Observe que un punto fijo es estable, mientras que el otro es inestable. Si la elipse se inclinara ligeramente alejándose del punto fijo inestable, una iteración de QR haría que la elipse se inclinara alejándose del punto fijo en lugar de acercándose. Sin embargo, eventualmente, el algoritmo convergería a un punto fijo diferente, pero tardaría mucho tiempo.

Hallar valores propios frente a hallar vectores propios

Figura 2: Cómo se ve afectado el resultado de una sola iteración de QR o LR cuando dos valores propios se aproximan entre sí.

Cabe destacar que hallar incluso un solo vector propio de una matriz simétrica no es computable (en aritmética real exacta según las definiciones del análisis computable ). [ 10 ] Esta dificultad surge cuando se desconocen las multiplicidades de los valores propios de una matriz. Por otro lado, este problema no se presenta al hallar los valores propios. Los valores propios de una matriz siempre son computables.

Ahora analizaremos cómo se manifiestan estas dificultades en el algoritmo QR básico. Esto se ilustra en la Figura 2. Recordemos que las elipses representan matrices simétricas definidas positivas. A medida que los dos autovalores de la matriz de entrada se aproximan, la elipse de entrada se transforma en un círculo. Un círculo corresponde a un múltiplo de la matriz identidad. Un casi círculo corresponde a un múltiplo cercano de la matriz identidad cuyos autovalores son casi iguales a los elementos de la diagonal de la matriz. Por lo tanto, el problema de encontrar aproximadamente los autovalores resulta sencillo en este caso. Sin embargo, observemos qué sucede con los semiejes de las elipses. Una iteración de QR (o LR) inclina los semiejes cada vez menos a medida que la elipse de entrada se aproxima a un círculo. Los autovectores solo se pueden conocer cuando los semiejes son paralelos a los ejes x e y. El número de iteraciones necesarias para lograr un paralelismo casi perfecto aumenta indefinidamente a medida que la elipse de entrada se vuelve más circular.

Si bien puede resultar imposible calcular la descomposición en valores propios de una matriz simétrica arbitraria, siempre es posible perturbarla con una pequeña cantidad y calcular la descomposición en valores propios de la matriz resultante. En el caso de que la matriz se represente como un círculo casi perfecto, se puede reemplazar por una matriz cuya representación sea un círculo perfecto. En ese caso, la matriz resultante es un múltiplo de la matriz identidad y su descomposición en valores propios es inmediata. Sin embargo, cabe destacar que la base de valores propios resultante puede estar bastante alejada de la base de valores propios original.

Aceleración: Cambios y desinflación

La ralentización que se produce cuando la elipse se vuelve más circular tiene una recíproca: resulta que cuando la elipse se estira más (y se vuelve menos circular), la rotación de la elipse se acelera. Este estiramiento puede inducirse cuando la matrizMETRO{\displaystyle M}que la elipse representa se reemplaza porMETROλI{\displaystyle M-\lambda I}dóndeλ{\displaystyle \lambda }es aproximadamente el valor propio más pequeño deMETRO{\displaystyle M}En este caso, la relación de los dos semiejes de la elipse se aproxima{\displaystyle \infty }En dimensiones superiores, este tipo de desplazamiento reduce la longitud del semieje menor de un elipsoide en relación con los demás semiejes, lo que acelera la convergencia al autovalor menor, pero no la convergencia a los demás autovalores. Esto resulta inútil cuando el autovalor menor está completamente determinado, por lo que la matriz debe ser desinflada , lo que simplemente significa eliminar su última fila y columna.

También es necesario abordar el problema del punto fijo inestable. La heurística de desplazamiento suele diseñarse para solucionar este problema: los desplazamientos prácticos suelen ser discontinuos y aleatorios. El desplazamiento de Wilkinson, muy adecuado para matrices simétricas como las que estamos visualizando, es particularmente discontinuo.

El algoritmo QR implícito

En la práctica computacional moderna, el algoritmo QR se realiza en una versión implícita que facilita la introducción del uso de múltiples desplazamientos. [ 4 ] La matriz se transforma primero a la forma de Hessenberg superior.A0=QAQT{\displaystyle A_{0}=QAQ^{\mathsf {T}}}como en la versión explícita; luego, en cada paso, la primera columna deAk{\displaystyle A_{k}}se transforma mediante una transformación de similitud de Householder de tamaño pequeño a la primera columna depag(Ak){\displaystyle p(A_{k})}(opag(Ak)mi1{\displaystyle p(A_{k})e_{1}}), dóndepag(Ak){\displaystyle p(A_{k})}, de grador{\displaystyle r}, es el polinomio que define la estrategia de cambio (a menudopag(incógnita)=(incógnitaλ)(incógnitaλ¯){\displaystyle p(x)=(x-\lambda )(x-{\bar {\lambda }})}, dóndeλ{\displaystyle \lambda }yλ¯{\displaystyle {\bar {\lambda }}}son los dos autovalores de la cola2×2{\displaystyle 2\times 2}submatriz principal deAk{\displaystyle A_{k}}, el llamado doble desplazamiento implícito ). Luego, transformaciones sucesivas de Householder de tamañor+1{\displaystyle r+1}se realizan para devolver la matriz de trabajoAk{\displaystyle A_{k}}a la forma de Hessenberg superior. Esta operación se conoce como persecución de abultamiento , debido a la peculiar forma de las entradas no nulas de la matriz a lo largo de los pasos del algoritmo. Como en la primera versión, la deflación se realiza tan pronto como una de las entradas subdiagonales deAk{\displaystyle A_{k}}es suficientemente pequeño.

Propuesta de cambio de nombre

Dado que en la versión implícita moderna del procedimiento no se realizan explícitamente descomposiciones QR , algunos autores, por ejemplo Watkins, [ 11 ] sugirieron cambiar su nombre a algoritmo de Francis . Golub y Van Loan utilizan el término paso QR de Francis .

Interpretación y convergencia

El algoritmo QR puede considerarse una variación más sofisticada del algoritmo básico de autovalores de "potencia" . Recordemos que el algoritmo de potencia multiplica repetidamente A por un único vector, normalizando después de cada iteración. El vector converge a un autovector del autovalor más grande. En cambio, el algoritmo QR trabaja con una base completa de vectores, utilizando la descomposición QR para renormalizar (y ortogonalizar). Para una matriz simétrica A , al converger, AQ = Q Λ , donde Λ es la matriz diagonal de autovalores a la que convergió A , y donde Q es una composición de todas las transformaciones de similitud ortogonales necesarias para llegar allí. Por lo tanto, las columnas de Q son los autovectores.

Historia

El algoritmo QR fue precedido por el algoritmo LR , que utiliza la descomposición LU en lugar de la descomposición QR. El algoritmo QR es más estable, por lo que el algoritmo LR se usa poco hoy en día. Sin embargo, representa un paso importante en el desarrollo del algoritmo QR.

El algoritmo LR fue desarrollado a principios de la década de 1950 por Heinz Rutishauser , quien trabajaba en ese entonces como asistente de investigación de Eduard Stiefel en la ETH Zúrich . Stiefel sugirió que Rutishauser utilizara la secuencia de momentos y 0 T A k x 0 , k = 0, 1, ... (donde x 0 e y 0 son vectores arbitrarios) para encontrar los valores propios de A. Rutishauser tomó un algoritmo de Alexander Aitken para esta tarea y lo desarrolló en el algoritmo de cociente - diferencia o algoritmo qd . Después de organizar el cálculo en una forma adecuada, descubrió que el algoritmo qd es en realidad la iteración A k = L k U k (descomposición LU), A k +1 = U k L k , aplicada a una matriz tridiagonal, de la cual se deriva el algoritmo LR. [ 12 ]

Otras variantes

Una variante del algoritmo QR , el algoritmo de Golub-Kahan-Reinsch, comienza reduciendo una matriz general a una bidiagonal. [ 13 ] Esta variante del algoritmo QR para el cálculo de valores singulares fue descrita por primera vez por Golub y Kahan (1965) . La subrutina DBDSQR de LAPACK implementa este método iterativo , con algunas modificaciones para cubrir el caso en que los valores singulares son muy pequeños ( Demmel y Kahan 1990 ) . Junto con un primer paso que utiliza reflexiones de Householder y, si es apropiado, la descomposición QR , esto forma la rutina DGESVD para el cálculo de la descomposición de valores singulares . El algoritmo QR también puede implementarse en dimensiones infinitas con los correspondientes resultados de convergencia. [ 14 ] [ 15 ]

Referencias

  1. JGF Francis, "La transformación QR, I", The Computer Journal , 4 (3), páginas 265 271 (1961, recibido en octubre de 1959). doi:10.1093/comjnl/4.3.265
  2. Francis, JGF (1962). "La transformación QR, II" . The Computer Journal . 4 (4): 332– 345. doi : 10.1093/comjnl/4.4.332 .
  3. Vera N. Kublanovskaya, «Sobre algunos algoritmos para la solución del problema completo de valores propios», Matemáticas Computacionales y Física Matemática de la URSS , vol. 1, n.º 3, páginas 637-657 (1963, recibido en febrero de 1961). También publicado en: Zhurnal Vychislitel'noi Matematiki i Matematicheskoi Fiziki , vol. 1, n.º 4, páginas 555-570 (1961). doi:10.1016/0041-5553(63)90168-X
  4. 1 2 Golub, GH; Van Loan, CF (1996). Cálculos matriciales (3.ª ed.). Baltimore: Johns Hopkins University Press. ISBN  0-8018-5414-8.
  5. Holmes, Mark H. (2023). Introducción a la computación científica y al análisis de datos . Textos en ciencia e ingeniería computacional (Segunda edición). Cham: Springer. ISBN  978-3-031-22429-4.
  6. Golub, Gene H.; Van Loan, Charles F. (2013). Cálculos matriciales . Estudios de Johns Hopkins en ciencias matemáticas (Cuarta ed.). Baltimore: The Johns Hopkins University Press. ISBN  978-1-4214-0794-4.
  7. 1 2 Demmel, James W. (1997). Álgebra lineal numérica aplicada . SIAM.
  8. ^ Trefethen , Lloyd N .; Bau, David (1997). Álgebra lineal numérica . SIAM.
  9. Ortega, James M.; Kaiser, Henry F. (1963). "Los métodos LL T y QR para matrices tridiagonales simétricas" . The Computer Journal . 6 (1): 99– 101. doi : 10.1093/comjnl/6.1.99 .
  10. "Álgebra lineal: ¿Por qué la imposibilidad de computar la descomposición espectral no es un problema?" . MathOverflow . Consultado el 9 de agosto de 2021 .
  11. Watkins, David S. (2007). El problema de los valores propios de la matriz: métodos de GR y subespacio de Krylov . Filadelfia, PA: SIAM. ISBN 978-0-89871-641-2.
  12. Parlett, Beresford N.; Gutknecht, Martin H. (2011), "De qd a LR, o ¿cómo se descubrieron los algoritmos qd y LR?" (PDF) , IMA Journal of Numerical Analysis , 31 (3): 741–754 , doi : 10.1093/imanum/drq003 , hdl : 20.500.11850/159536 , ISSN 0272-4979 
  13. Bochkanov Sergey Anatolyevich. Guía del usuario de ALGLIB - Operaciones generales con matrices - Descomposición en valores singulares. Proyecto ALGLIB. 11/12/2010. URL:Consultado: 11/12/2010. (Archivado por WebCite en
  14. Deift, Percy; Li, Luenchau C.; Tomei, Carlos (1985). "Toda fluye con infinitas variables" . Journal of Functional Analysis . 64 (3): 358– 402. doi : 10.1016/0022-1236(85)90065-5 .
  15. ^ Colbrook, Mateo J.; Hansen, Anders C. (2019). "Sobre el algoritmo QR de dimensión infinita" . Matemática numérica . 143 (1): 17– 83. arXiv : 2011.08172 . doi : 10.1007/s00211-019-01047-5 .

Fuentes

  • Demmel, James ; Kahan, William (1990). "Valores singulares precisos de matrices bidiagonales". SIAM Journal on Scientific and Statistical Computing . 11 (5): 873– 912. CiteSeerX 10.1.1.48.3740 . doi : 10.1137/0911052 . 
  • Golub, Gene H.; Kahan , William (1965). "Cálculo de los valores singulares y la pseudoinversa de una matriz". Journal of the Society for Industrial and Applied Mathematics, Series B: Numerical Analysis . 2 (2): 205– 224. Bibcode : 1965SJNA....2..205G . doi : 10.1137/0702016 . JSTOR 2949777 . 
  • Problema de valores propios en PlanetMath .
  • Notas sobre bases ortogonales y el funcionamiento del algoritmo QR por Peter J. Olver
  • Módulo para el método QR
  • Biblioteca C++