Inverter matriz 4x4 - numérico mais solução estável necessário

StackOverflow https://stackoverflow.com/questions/155670

Pergunta

Eu quero inverter uma matriz 4x4. Meus números são armazenados em formato de ponto fixo (1.15.16 para ser exato).

Com de ponto flutuante I aritmética geralmente apenas construir a matriz adjunta e dividir pelo determinante (por exemplo bruta forçar a solução). Isso funcionou para mim até agora, mas quando se lida com números de ponto fixo recebo uma perda de precisão inaceitável devido a todas as multiplicações usados.

Nota:. No ponto fixo aritmética eu sempre jogar fora alguns dos bits menos significativos de resultados imediatos

Então - Qual é a forma mais estável numérica para inverter uma matriz? Eu não me importo muito sobre o desempenho, mas simplesmente vai de ponto flutuante seria lento no meu arquitetura alvo.

Foi útil?

Solução

Eu acho que a resposta a esta depende da forma exata da matriz. Um método de decomposição padrão (LU, QR, Cholesky etc.) com pivotante (essencial) é bastante boa no ponto fixo, especialmente para uma matriz 4x4 pequeno. Veja o livro 'Numerical Recipes' por Press et al. para uma descrição destes métodos.

Este artigo dá alguns algoritmos úteis, mas está atrás de um paywall infelizmente. Eles recomendam um (articulado) decomposição de Cholesky com alguns recursos adicionais muito complicado para a lista aqui.

Outras dicas

Meta-resposta: É realmente uma matriz geral 4x4? Se a sua matriz tem uma forma especial, então há fórmulas diretas para invertendo que seria rápido e manter a sua operação de contagem regressiva.

Por exemplo, se é uma homogênea coordenadas padrão transformar a partir de gráficos, como: 

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

(assumindo uma composição de rotação, escala, matrizes de tradução)

, em seguida, há um facilmente derivável fórmula direta , que é 

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

(matrizes ASCII roubado da página vinculada).

Você provavelmente não pode bater isso para perda de precisão no ponto fixo.

Se a sua matriz vem de algum domínio em que você sabe que tem mais estrutura, então não é provável que seja uma resposta fácil.

Eu gostaria de segunda a pergunta Jason S levantada: você está certo de que você precisa para inverter a sua matriz? Isso quase nunca é necessário. Não só isso, muitas vezes é uma má idéia. Se você precisa resolver Ax = b, é mais numericamente estável para resolver o sistema diretamente do que a b multiplicar por A inversa.

Mesmo que você tem que resolver Ax = b repetidamente para muitos valores de b, ainda não é uma boa idéia para inverter A. Você pode fator A (dizem LU fatoração ou Cholesky fatoração) e salvar os fatores que você não está refazendo que o trabalho de cada vez, mas você ainda resolver o sistema cada vez usando a fatoração.

Deixe-me fazer uma pergunta diferente: você definitivamente precisa de inverter a matriz (chamemos-lhe M), ou você precisa usar a matriz inversa para resolver outras equações? (Por exemplo Mx = b para conhecida H, b) Muitas vezes existem outras maneiras de fazer isso w / o necessitando explicitamente para calcular o inverso. Ou se a matriz M é uma função do tempo e muda lentamente, então você pode calcular o inverso completo uma vez, e existem maneiras iterativos para atualizá-lo.

Para minimizar erros de truncamento e outros maldade, use "pivotantes" - veja o capítulo sobre invertendo matrizes no Numerical Recipes. Eles têm a melhor explicação que eu encontrei até agora.

Você pode considerar dobrando para 1,31 antes de fazer o seu algoritmo normal. Ele vai dobrar o número de multiplicações, mas você está fazendo uma inversão da matriz e qualquer coisa que você vai ser muito ligada ao multiplicador em seu processador.

Para qualquer um interessado em encontrar as equações para um invertido 4x4, você pode usar um pacote de matemática simbólica para resolvê-los para você. A TI-89 vai fazer isso mesmo, embora isso vai demorar vários minutos.

Se você nos dá uma idéia do que a inversão da matriz faz para você, e como ele se encaixa com o resto do seu processamento que pode ser capaz de sugerir alternativas.

-Adam

Se a matriz representa uma transformação afim (muitas vezes este é o caso com 4x4 matrizes, desde que você não introduzir um componente de escalonamento) o inverso é simplesmente a transposição da parte de rotação 3x3 superior com a última coluna negada. Obviamente, se você precisar de uma solução generalizada, em seguida, olhando para eliminação de Gauss é provavelmente o mais fácil.

Licenciado em: CC-BY-SA com atribuição
Não afiliado a StackOverflow
scroll top