Invertir matriz 4x4 - Se necesita la solución numérica más estable

Question

Quiero invertir una matriz 4x4. Mis números se almacenan en formato de punto fijo (1.15.16 para ser exactos).

Con aritmética de coma flotante, suelo construir la matriz adjunta y dividir por el determinante (por ejemplo, fuerza bruta la solución). Eso funcionó para mí hasta ahora, pero cuando se trata de números de puntos fijos obtengo una pérdida de precisión inaceptable debido a todas las multiplicaciones utilizadas.

Nota: En la aritmética de punto fijo siempre arrojo algunos de los resultados inmediatos menos significativos.

Entonces, ¿cuál es la forma numérica más estable para invertir una matriz? No me importa mucho el rendimiento, pero simplemente ir a punto flotante sería ralentizar la arquitectura de mi objetivo.

Answer 1

Creo que la respuesta depende de la forma exacta de la matriz. Un método de descomposición estándar (LU, QR, Cholesky, etc.) con pivote (un elemento esencial) es bastante bueno en el punto fijo, especialmente para una matriz pequeña de 4x4. Ver el libro 'Recetas Numéricas' de Press et al. para una descripción de estos métodos.

This paper da algunos algoritmos útiles, pero está detrás de un paywall desafortunadamente. Recomiendan una descomposición de Cholesky (pivotada) con algunas características adicionales demasiado complicadas para enumerarlas aquí.

Answer 2

La simple eliminación de Gauss antigua funcionaría bien.

Depende de qué bibliotecas/clases/estructuras esté utilizando. Puede echar un vistazo al GSL.

Answer 3

Para minimizar los errores de truncamiento y otras irregularidades, utilice "pivoting" - vea el capítulo sobre inversión de matrices en Recetas Numéricas. Tienen la mejor explicación que he encontrado hasta ahora.

Answer 4

Meta-respuesta: ¿Es realmente una matriz 4x4 general? Si su matriz tiene una forma especial, existen fórmulas directas para invertir que serían rápidas y mantener su cuenta atrás.

Por ejemplo, si se trata de un estándar de coordenadas homogéneas transformar de gráficos, como:

[ux vx wx tx] 
[uy vy wy ty] 
[uz vz wz tz] 
[ 0 0 0 1] 

(suponiendo una composición de rotación, escala, matrices de traducción)

entonces hay una easily-derivable direct formula, que es

[ux uy uz -dot(u,t)] 
[vx vy vz -dot(v,t)] 
[wx wy wz -dot(w,t)] 
[ 0 0 0  1 ] 

(matrices ASCII robados de la página enlazada.)

Probablemente no se puede superar eso por la pérdida de precisión en el punto fijo.

Si su matriz proviene de algún dominio donde sabe que tiene más estructura, entonces es probable que haya una respuesta fácil.

Answer 5

Si la matriz representa una transformación afín (muchas veces este es el caso con las matrices 4x4 siempre que no introduzca un componente de escala) la inversa es simplemente la transposición de la parte superior de rotación 3x3 con la última columna negada . Obviamente, si necesita una solución generalizada, es más fácil analizar la eliminación gaussiana.

Answer 6

Puede considerar doblar a 1,31 antes de hacer su algoritmo normal.Duplicará el número de multiplicaciones, pero está haciendo una matriz invertida y todo lo que haga va a estar muy ligado al multiplicador de su procesador.

Para cualquier persona interesada en encontrar las ecuaciones para una inversión de 4x4, puede usar un paquete matemático simbólico para resolverlas por usted. La TI-89 lo hará incluso, aunque tomará varios minutos.

Si nos da una idea de lo que la matriz invertida hace por usted, y cómo encaja con el resto de su procesamiento, podríamos sugerir alternativas.

-Adam

Answer 7

Permítanme hacer una pregunta diferente: hace definitivamente hay que invertir la matriz (lo llaman M), o debe utilizar a la matriz inversa para resolver otras ecuaciones? (por ejemplo, Mx = b para M conocido, b) A menudo hay otras maneras de hacer esto sin necesidad explícita de calcular el inverso. O si la matriz M es una función del tiempo &, cambia lentamente, entonces se puede calcular el inverso completo una vez, &, hay formas iterativas de actualizarlo.

Answer 8

Me gustaría responder a la pregunta planteada por Jason S: ¿está seguro de que necesita invertir su matriz? Esto casi nunca es necesario. No solo eso, a menudo es una mala idea. Si necesita resolver Ax = b, es más numéricamente estable resolver el sistema directamente que multiplicar b por A inverso.

Incluso si usted tiene que resolver Ax = b y otra vez durante muchos valores de b, todavía no es una buena idea para invertir R. Puede factor de A (digamos de factorización LU o factorización de Cholesky) y guardar los factores así que no estás rehaciendo ese trabajo todas las veces, pero igual resolverías el sistema cada vez que uses la factorización.