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

language-agnostic matrix

Question

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

Con la aritmética de coma flotante, suelo construir la matriz conjunta y dividirla por la 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 bits menos significativos de resultados inmediatos.

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

La simple eliminación Gaussiana funcionaría bien.

Depende de qué bibliotecas / clases / estructuras estés usando. Puedes echarle un vistazo al GSL .

Answer 2

Para minimizar los errores de truncamiento y otras maldades, use "pivoting" - vea el capítulo sobre inversión de matrices en Numerical Recipes. Tienen la mejor explicación que he encontrado hasta ahora.

Answer 3

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 resolverlo 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.

-Adán

Answer 4

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

Answer 5

Déjame hacerte una pregunta diferente: ¿definitivamente necesitas invertir la matriz (llámala M), o necesitas usar la matriz inversa para resolver otras ecuaciones? (por ejemplo, Mx = b para M conocido, b) A menudo hay otras maneras de hacerlo sin necesidad de calcular explícitamente la inversa. O si la matriz M es una función del tiempo y cambia lentamente, entonces podría calcular el inverso completo una vez, y hay formas iterativas de actualizarlo.

Answer 6

Me gustaría secundar la pregunta planteada por Jason S: ¿estás seguro de que necesitas invertir tu 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 tiene que resolver Ax = b una y otra vez para muchos valores de b, aún no es una buena idea invertir A. Puede factorizar A (digamos factorización LU o factorización Cholesky) y guardar los factores para no estar rehaciendo eso funciona todo el tiempo, pero igual resolverías el sistema cada vez que uses la factorización.

Answer 7

Creo que la respuesta a esto 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.

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

Answer 8

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

Por ejemplo, si se trata de una transformación de coordenadas homogénea estándar a partir 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 fórmula directa fácilmente derivable , 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 robadas de la página vinculada).

Probablemente no puedas vencer 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.