¿Cuándo superan las transformaciones ortogonales la eliminación gaussiana?

22

Como sabemos, los métodos de transformaciones ortogonales (rotaciones de Givens y reflexiones de Housholder) para sistemas de ecuaciones lineales son más caros que la eliminación gaussiana, pero en teoría tienen mejores propiedades de estabilidad en el sentido de que no cambian el número de condición del sistema. Aunque conozco solo un ejemplo académico de una matriz que se ve afectada por la eliminación gaussiana con pivote parcial. Y existe la opinión común de que es muy poco probable que se encuentre con este tipo de comportamiento en la práctica (consulte las notas de esta conferencia [pdf] ).

Entonces, ¿dónde buscaremos la respuesta sobre el tema? Implementaciones paralelas? Actualizando? ..

faleichik
fuente

Respuestas:

24

Exactitud

Trefethen y Schreiber escribieron un excelente artículo, Estabilidad en el caso promedio de la eliminación gaussiana , que analiza el lado de la precisión de su pregunta. Estas son algunas de sus conclusiones:

  1. "Para factorización QR con o sin pivotante columna, el elemento media máxima de la matriz residual es , mientras que para la eliminación de Gauss es O ( n ) . Esta comparación revela que la eliminación gaussiana es ligeramente inestable, pero la la inestabilidad solo sería detectable para problemas de matriz muy grandes resueltos con baja precisión. Para la mayoría de los problemas prácticos, la eliminación gaussiana es altamente estable en promedio " . (El énfasis es mío)O(n1/2)O(n)

  2. "Después de los primeros pasos de la eliminación gaussiana, los elementos restantes de la matriz se distribuyen aproximadamente de manera normal, independientemente de si comenzaron de esa manera".

Hay mucho más en el documento que no puedo capturar aquí, incluida la discusión de la matriz del peor de los casos que mencionó, por lo que le recomiendo encarecidamente que la lea.

Actuación

Para matrices reales cuadrados, LU con pivoteo parcial requiere aproximadamente fracasos, mientras que QR basado en Householder requiere aproximadamente 4 / 3 n 3 flops. Por lo tanto, para matrices cuadradas razonablemente grandes, la factorización QR solo será aproximadamente el doble de cara que la factorización LU.2/3n34/3n3

Para matrices, donde m n , LU con pivoteo parcial requiere m n 2 - n 3 / 3 flop, frente de QR 2 m n 2 - 2 n 3 / 3 (que todavía es el doble de factorización LU). Sin embargo , es sorprendentemente común que las aplicaciones produzcan matrices flacas muy altas ( m n ), y Demmel et al. tenga un buen papel, Factorización QR secuencial y paralela que evite la comunicaciónm×nmnmn2n3/32mn22n3/3mn, que (en la sección 4) analiza un algoritmo inteligente que solo requiere el envío de mensajes cuando se utilizan procesadores p , frente a los mensajes n log p de los enfoques tradicionales. El gasto es que se realizan flops adicionales O ( n 3 log p ) , pero para n muy pequeños, esto a menudo se prefiere al costo de latencia de enviar más mensajes (al menos cuando solo se necesita realizar una factorización QR).logppnlogpO(n3logp)n

Jack Poulson
fuente
10

Me sorprende que nadie haya mencionado problemas de mínimos cuadrados lineales , que ocurren con frecuencia en la informática científica. Si desea utilizar la eliminación gaussiana, debe formar y resolver las ecuaciones normales, que se ven así:

ATAx=ATb,

Axb

ATAA

Geoff Oxberry
fuente
2
n3AHA2/3n36n3
1
Además de la estabilidad garantizada por el uso de transformaciones ortogonales, la gran ventaja de SVD es que la descomposición proporciona su propia verificación de condición, ya que la relación del valor singular más grande al más pequeño es precisamente el número de condición (norma 2). Para las otras descomposiciones, el uso de un estimador de condición (por ejemplo, Hager-Higham) es, aunque no tan costoso como la descomposición propiamente dicha, algo "agregado".
JM
1
@JackPoulson Solo por curiosidad, ¿tienes una referencia para tu conteo de flops para SVD? Por lo que puedo decir de un vistazo rápido en Golub & Van Loan (p. 254 3a edición), la constante parecería más alta para usar el SVD en la resolución de problemas de mínimos cuadrados, pero podría estar equivocado. Gracias por adelantado.
OscarB
1
8/3n3A=FBGHCB=UΣVHx:=(G(V(inv(Σ)(UH(FHb)))))O(n2)CO(n2)
1
σ1σn
3

¿Cómo se mide el rendimiento? ¿Velocidad? ¿Exactitud? ¿Estabilidad? Una prueba rápida en Matlab da lo siguiente:

>> N = 100;
>> A = randn(N); b = randn(N,1);
>> tic, for k=1:10000, [L,U,p] = lu(A,'vector'); x = U\(L\b(p)); end; norm(A*x-b), toc
ans =
   1.4303e-13
Elapsed time is 2.232487 seconds.
>> tic, for k=1:10000, [Q,R] = qr(A); x = R\(Q'*b); end; norm(A*x-b), toc             
ans =
   5.0311e-14
Elapsed time is 7.563242 seconds.

Por lo tanto, resolver un solo sistema con una descomposición LU es aproximadamente tres veces más rápido que resolverlo con una descomposición QR, a costa de medio dígito decimal de precisión (¡este ejemplo!).

Pedro
fuente
Cualquiera de los méritos que sugirió son bienvenidos.
faleichik
3

El artículo que cita defiende la Eliminación Gaussiana al decir que, aunque es numéricamente inestable, tiende a funcionar bien en matrices aleatorias y, dado que la mayoría de las matrices que uno puede pensar son como matrices aleatorias, deberíamos estar bien. Esta misma afirmación se puede decir de muchos métodos numéricamente inestables.

Considere el espacio de todas las matrices. Estos métodos funcionan bien en casi todas partes. Eso es 99.999 ...% de todas las matrices que uno podría crear no tendrá problemas con los métodos inestables. Solo hay una fracción muy pequeña de matrices para las que GE y otros tendrán dificultades.

Los problemas que preocupan a los investigadores tienden a estar en esa pequeña fracción.

No construimos matrices al azar. Construimos matrices con propiedades muy especiales que corresponden a sistemas muy especiales, no aleatorios. Estas matrices a menudo están mal acondicionadas.

Geométricamente puedes considerar el espacio lineal de todas las matrices. Hay un subespacio de volumen / medida cero de matrices singulares que atraviesa este espacio. Muchos problemas que construimos se agrupan alrededor de este subespacio. No se distribuyen al azar.

Como ejemplo, considere la ecuación de calor o la dispersión. Estos sistemas tienden a eliminar información del sistema (todos los estados iniciales gravitan a un solo estado final) y, como resultado, las matrices que describen estas ecuaciones son enormemente singulares. Este proceso es muy poco probable en una situación aleatoria pero ubicua en los sistemas físicos.

MRocklin
fuente
2
Si el sistema lineal está inicialmente mal acondicionado, no importa qué método utilice: tanto la descomposición de LU como la de QR darán resultados inexactos. QR puede ganar solo en los casos en que el proceso de eliminación gaussiana "arruine" una buena matriz. El problema principal es que no se conocen casos prácticos de tal comportamiento.
faleichik
Para la mayoría de las aplicaciones científicas, generalmente obtenemos matrices que son escasas, simétricas, positivas definidas y / o diagonalmente dominantes. Con muy pocas excepciones, hay una estructura en la matriz que nos permite explotar ciertas técnicas sobre la eliminación gaussiana tradicional.
Paul
@Paul: Por otro lado, la eliminación de Gauss densa es donde se pasa la mayor parte del tiempo en el método multifrontal para matrices no asimétricas dispersas.
Jack Poulson
66
@Paul Simplemente no es cierto que "la mayoría de las aplicaciones producen SPD / matrices diagonalmente dominantes". Sí, generalmente hay una estructura explotable de algún tipo, pero los problemas no simétricos e indefinidos son extremadamente comunes.
Jed Brown
44
"En cincuenta años de computación, no se sabe que surjan problemas matriciales que provoquen una inestabilidad explosiva en circunstancias naturales". - LN Trefethen y D. Bau Dan un interesante análisis probabilístico en su libro.
JM