Este artigo foi traduzido por máquina.
C#
Decomposição de matriz
Decomposição da matriz, uma técnica que quebra em uma matriz numérica quadrada em duas diferentes matrizes de Praça, é a base para resolver com eficiência um sistema de equações, que por sua vez é a base para a inversão de uma matriz.E inversão de uma matriz é uma parte de muitos algoritmos importantes.Este artigo apresenta e explica o código c# que realiza a decomposição da matriz, inversão de matriz, um sistema de solução de equações e operações conexas.
Reconhecidamente, decomposição da matriz não é um tema chamativo, mas uma coleção de métodos de matriz pode ser um complemento importante à sua biblioteca de código pessoal.Os métodos são explicados para que você pode modificar o código-fonte para atender às suas necessidades.Além disso, algumas das técnicas usadas em métodos de matriz podem ser reutilizadas em outros cenários de codificação.
A melhor maneira para você começar uma sensação para o tipo de informação apresentada neste artigo é dar uma olhada no screenshot no Figura 1.O programa de demonstração começa criando uma matriz quadrada de 4 x 4 e exibindo seus valores.Em seguida, a matriz é decomposta no que é chamado de uma matriz LUP.L significa inferior e U representa o superior.A parte de P (P significa permutação) é uma matriz com valores {3.120} e indica que as linhas 0 e 3 foram trocadas durante o processo de decomposição.A decomposição também gerou um valor de alternância de -1, indicando que ocorreu um número ímpar de trocas de linha.O programa de demonstração mostra a decomposição de duas maneiras: primeiro como uma matriz de LU combinada e matrizes L e U, em seguida, separadamente.Em seguida, o programa calcula e exibe a inversa da matriz original, usando a matriz LUP nos bastidores.O programa de demonstração calcula o determinante da matriz original, novamente usando a decomposição.Ele usa o inverso da matriz para resolver um sistema de equações lineares e conclui combinando as matrizes L e U volta para a matriz original.
.png)
Figura 1 matriz decomposição Demo
Mas por que ir para todos os problemas da criação de um método de decomposição de matriz personalizada e uma biblioteca de métodos relacionados?Embora existam muitas ferramentas de matriz independente disponível, às vezes pode ser difícil de integrar em um aplicativo ou sistema.E apesar da importância fundamental da decomposição da matriz, há poucas implementações de código de .NET de livre, sem direitos autorais disponíveis; aqueles que existem muitas vezes não são explicados em detalhes suficientes para modificar o código-fonte para atender seus cenários de codificação.
Este artigo pressupõe que você tenha habilidades intermediárias c# programação e pelo menos um conhecimento básico de operações de matriz e terminologia.Todo o código c# a chave é apresentada neste artigo.O código também está disponível a partir do site de download de código do MSDN em archve.msdn.microsoft.com/mag201212matrix.
Definição de matriz
Existem várias maneiras de implementar uma matriz em C#.A abordagem tradicional e usada neste artigo, é usar uma matriz de matrizes, às vezes chamado de uma matriz denteada.Por exemplo, este código define uma matriz com três linhas e duas colunas:
double[][] m = new double[3][]; m[0] = new double[2]; m[1] = new double[2]; m[2] = new double[2]; m[2][1] = 5.0; // set row 2, col 1
Ao contrário da maioria das linguagens de programação, C# tem um tipo de matriz multidimensional incorporado, que fornece uma abordagem alternativa. Por exemplo:
double[,] m = new double[3,2]; m[2,1] = 5.0;
Uma terceira abordagem para implementar matrizes em c# é usar uma única matriz combinada com manipulação de índice de matriz, como este:
int rows = 3; int cols = 2; double[] m = new double[rows * cols]; int i = 2; int j = 1; m[i * cols + j] = 5.0;
Independentemente do esquema de armazenamento usado, matrizes podem ser implementados usando uma abordagem do método estático ou um OOP. Por exemplo, uma abordagem OOP poderia ser semelhante:
public class MyMatrix { int m; // number rows int n; // number columns data[][]; // the values ...
}
Não há nenhum única melhor escolha para a concepção de matriz; o melhor projeto depende o cenário de codificação específico em que você está operando suas preferências pessoais de codificação. Este artigo usa uma abordagem estática-método porque é mais fácil de entender e refatorar.
Quando usando um design de matriz de matrizes de matrizes, porque cada linha deve ser alocada separadamente, muitas vezes é conveniente definir um método auxiliar para realizar a alocação de memória. Por exemplo:
static double[][] MatrixCreate(int rows, int cols) { // creates a matrix initialized to all 0.0s // do error checking here?
double[][] result = new double[rows][]; for (int i = 0; i < rows; ++i) result[i] = new double[cols]; // auto init to 0.0 return result; }
O método pode ser chamado assim:
double[][] m = MatrixCreate(3,2); m[2][1] = 5.0;
Este método demonstra uma das vantagens de criar sua própria biblioteca de métodos de matriz: Se você deseja melhorar o desempenho, você pode omitir os parâmetros de entrada a verificação de erros — à custa de aumentar o risco do código de chamada, fazendo com que uma exceção. (Para manter este artigo curto, verificação de erro mais foi removido.) Outra vantagem é que você pode personalizar a sua biblioteca para otimizar seu cenário exato. A principal desvantagem de criar sua própria biblioteca é que ele pode demorar mais do que usar uma biblioteca existente.
Um outro método conveniente para adicionar à sua biblioteca de matriz é aquele que pode ser usado para exibir uma matriz como uma seqüência de caracteres. Aqui está uma possibilidade:
static string MatrixAsString(double[][] matrix) { string s = ""; for (int i = 0; i < matrix.Length; ++i) { for (int j = 0; j < matrix[i].Length; ++j) s += matrix[i][j].ToString("F3").PadLeft(8) + " "; s += Environment.NewLine; } return s; }
Você pode querer parametrizar o número de casas decimais para exibir, ou o preenchimento de largura da coluna ou ambos.
Multiplicação de matrizes
O Microsoft .NET Framework 4 e posterior fornecem uma forma elegante para melhorar significativamente o desempenho de um método de multiplicação de matriz. Multiplicação de matrizes é ilustrada em Figura 2.
.png)
Figura 2 multiplicação de matrizes
Observe que o cálculo de cada valor de célula do resultado não depende em todos os outros valores de célula no resultado, então cada cálculo é independente e potencialmente podem ser executadas em paralelo em uma máquina com vários processadores. Aqui está uma abordagem padrão para a multiplicação de matriz:
static double[][] MatrixProduct(double[][] matrixA, double[][] matrixB) { int aRows = matrixA.Length; int aCols = matrixA[0].Length; int bRows = matrixB.Length; int bCols = matrixB[0].Length; if (aCols != bRows) throw new Exception("Non-conformable matrices in MatrixProduct"); double[][] result = MatrixCreate(aRows, bCols); for (int i = 0; i < aRows; ++i) // each row of A for (int j = 0; j < bCols; ++j) // each col of B for (int k = 0; k < aCols; ++k) result[i][j] += matrixA[i][k] * matrixB[k][j]; return result; }
Porque a multiplicação de matrizes é uma operação de O(n^3), o desempenho pode ser um problema. Por exemplo, se a matriz A tem tamanho 300 x 200 e a matriz B tem tamanho 200 x 400, computação produto da e B requer 300 * 200 * 400 = 24,000,000 multiplicações. A tarefa paralela TPL (biblioteca) no namespace System.Threading.Tasks no .NET Framework 4 e mais tarde torna-se fácil codificar uma versão paralelizada simples de multiplicação de matrizes. Uma possibilidade é:
static double[][] MatrixProduct(double[][] matrixA, double[][] matrixB) { // error check, compute aRows, aCols, bCols double[][] result = MatrixCreate(aRows, bCols); Parallel.For(0, aRows, i => { for (int j = 0; j < bCols; ++j) for (int k = 0; k < aCols; ++k) result[i][j] += matrixA[i][k] * matrixB[k][j]; } ); return result; }
Esta versão costeletas até os cálculos por linhas. Nos bastidores, a TPL gera todo o código de encanamento de complexos de sincronização para executar os cálculos em vários processadores.
Testes de consistência
Um aspecto interessante de uma biblioteca de métodos que relacionam-se uns aos outros é que muitas vezes é possível testá-los, verificando para ver se eles produzem resultados consistentes. Por exemplo, suponha que você tem um método que cria uma matriz aleatória:
static double[][] MatrixRandom(int rows, int cols, double minVal, double maxVal, int seed) { // return matrix with values between minVal and maxVal Random ran = new Random(seed); double[][] result = MatrixCreate(rows, cols); for (int i = 0; i < rows; ++i) for (int j = 0; j < cols; ++j) result[i][j] = (maxVal - minVal) * ran.NextDouble() + minVal; return result; }
Além disso, suponha que você tem uma matriz que cria a matriz de identidade — ou seja, uma matriz quadrada com 1.0s na diagonal principal, 0.0s em outro lugar:
static double[][] MatrixIdentity(int n) { double[][] result = MatrixCreate(n, n); for (int i = 0; i < n; ++i) result[i][i] = 1.0; return result; }
E suponha que você tenha um método que compara duas matrizes de igualdade:
static bool MatrixAreEqual(double[][] matrixA, double[][] matrixB, double epsilon) { // true if all values in A == corresponding values in B int aRows = matrixA.Length; int bCols = matrixB[0].Length; for (int i = 0; i < aRows; ++i) // each row of A and B for (int j = 0; j < aCols; ++j) // each col of A and B if (Math.Abs(matrixA[i][j] - matrixB[i][j]) > epsilon) return false; return true; }
Observe que o método de MatrixAreEqual não se compara valores de célula de igualdade exata porque os valores são de tipo double. Em vez disso, o método verifica se todos os valores de célula são muito perto (dentro de epsilon) uns aos outros.
Porque o produto de qualquer matriz quadrada m e a matriz identidade de mesma dimensão são iguais ao original matriz m, você pode testar o método do produto de matriz ao longo das linhas deste código:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0); double[][] i = MatrixIdentity(4); double[][] mi = MatrixProduct(m, i); if (MatrixAreEqual(m, mi, 0.00000001) == true) Console.WriteLine("Consistent result"); else Console.WriteLine("Something is wrong"); Consistency checking lends itself well to random input testing.
Decomposição de matriz
Decomposição da matriz é uma matriz quadrada M e computa duas novas matrizes quadradas que quando multiplicados juntos dar a matriz original M. A idéia é semelhante ao factoring número ordinário: o número 6 pode ser fatorado em 2 e 3, porque 2 * 3 = 6. A princípio pode parecer que há pouco ponto em decomposição de uma matriz, mas acontece que decomposição da matriz torna a tarefa muito difícil de inversão de matriz bastante um pouco mais simples. Há muitos tipos diferentes de decomposição da matriz, e cada tipo pode ser calculado usando vários algoritmos diferentes. A técnica apresentada neste artigo é chamada de decomposição LUP e usa o método de em Doolittle com giro parcial.
Para compreender a decomposição de LUP, é útil compreender primeiro a mais simples decomposição LU, que foi introduzida pelo famoso matemático Alan Turing. Suponha que você tenha essa matriz de 4 x 4 m:
9.000 5.000 3.000 4.000 4.000 8.000 2.000 5.000 3.000 5.000 7.000 1.000 2.000 6.000 0.000 8.000
Uma possível decomposição de LU m é L =
1.000 0.000 0.000 0.000 0.444 1.000 0.000 0.000 0.333 0.577 1.000 0.000 0.222 0.846 -0.219 1.000
E U =
9.000 5.000 3.000 4.000 0.000 5.778 0.667 3.222 0.000 0.000 5.615 -2.192 0.000 0.000 0.000 3.904
Isso funciona porque L * U = M. Observe que a matriz de L inferior tem 1.0s na diagonal e 0.0s no canto superior direito. Em outras palavras, os valores de célula significativos da matriz inferior estão no canto inferior esquerdo. Da mesma forma, os valores de célula significativos da matriz superior estão na diagonal principal e no canto superior direito.
Observe também que não há nenhuma sobreposição das localizações dos valores de célula significativa em L e U. Assim, em vez de gerar duas matrizes de resultados, L e U, decomposição da matriz normalmente armazena os resultados superiores e inferiores em uma única matriz que contém L e U para economizar espaço de memória.
Decomposição de matriz LUP é uma variação ligeira mas importante na decomposição LU. Decomposição de LUP leva uma matriz M e produz matrizes L e U, mas também uma matriz de P. O produto de L e U em LUP não é exatamente a matriz original M, mas em vez disso, é uma versão do M onde algumas linhas têm sido reorganizadas. A maneira em que as linhas tenham sido reorganizadas é armazenada na matriz P; Esta informação pode ser usada para reconstruir a matriz original M.
Um primo próximo a decomposição de Doolittle apresentada neste artigo é chamado de decomposição de na Crout. A principal diferença entre Doolittle e Crout é que Doolittle coloca 1.0s na diagonal principal da matriz L e Crout coloca 1.0s na diagonal principal da matriz U.
A razão de decomposição LUP é usada mais frequentemente do que a decomposição LU é sutil. Como você verá em breve, decomposição da matriz é usada para calcular a inversa de uma matriz. No entanto, quando a decomposição da matriz é usada como um auxiliar para inversão de matriz, verifica-se que a inversão falhará se não houver um valor 0.0 na diagonal principal da matriz LU. Assim, na decomposição de LUP, quando um 0,0 valor acaba na diagonal principal, o algoritmo troca duas linhas para mover o valor 0.0 fora da diagonal e mantém o controle de quais linhas foram trocadas na matriz P.
Figura 3 listas de um método de decomposição da matriz.
Figura 3 método de decomposição de Matrix
static double[][] MatrixDecompose(double[][] matrix, out int[] perm, out int toggle) { // Doolittle LUP decomposition.
// assumes matrix is square.
int n = matrix.Length; // convenience double[][] result = MatrixDuplicate(matrix); perm = new int[n]; for (int i = 0; i < n; ++i) { perm[i] = i; } toggle = 1; for (int j = 0; j < n - 1; ++j) // each column { double colMax = Math.Abs(result[j][j]); // largest val in col j int pRow = j; for (int i = j + 1; i < n; ++i) { if (result[i][j] > colMax) { colMax = result[i][j]; pRow = i; } } if (pRow != j) // swap rows { double[] rowPtr = result[pRow]; result[pRow] = result[j]; result[j] = rowPtr; int tmp = perm[pRow]; // and swap perm info perm[pRow] = perm[j]; perm[j] = tmp; toggle = -toggle; // row-swap toggle } if (Math.Abs(result[j][j]) < 1.0E-20) return null; // consider a throw for (int i = j + 1; i < n; ++i) { result[i][j] /= result[j][j]; for (int k = j + 1; k < n; ++k) result[i][k] -= result[i][j] * result[j][k]; } } // main j column loop return result; }
O método pode ser chamado como este:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0); int[] perm; int toggle; double luMatrix = MatrixDecompose(m, out perm, out toggle);
O método MatrixDecompose aceita como entrada uma matriz quadrada. O método tem três valores de retorno. O retorno explícito é uma matriz de LU permutada. O método retorna dois valores como parâmetros de saída. Um é uma matriz de permutação que contém informações sobre como as linhas foram permutadas. O segundo parâmetro out é um valor de alternância que é + 1 ou -1, dependendo se o número de trocas de linha foi mesmo (+ 1) ou ímpar (-1). O valor de alternância não é usado para inversão de matriz, mas é necessário se a decomposição da matriz é usada para calcular o determinante de uma matriz.
O método MatrixDecompose é bastante complicado, mas na verdade, existem apenas alguns detalhes que você precisa entender para modificar o código. A versão aqui apresentada aloca memória nova para a matriz de retorno de LU usando um método auxiliar MatrixDuplicate:
static double[][] MatrixDuplicate(double[][] matrix) { // assumes matrix is not null.
double[][] result = MatrixCreate(matrix.Length, matrix[0].Length); for (int i = 0; i < matrix.Length; ++i) // copy the values for (int j = 0; j < matrix[i].Length; ++j) result[i][j] = matrix[i][j]; return result; }
Uma abordagem alternativa é para calcular o resultado na matriz de entrada a fim de economizar memória. Com c# semântica, isso tornaria o parâmetro matriz um parâmetro ref porque ele é usado para entrada e saída. Usando essa abordagem, a assinatura de método seria:
static void MatrixDecompose(ref double[][] matrix, out int[] perm, out int toggle)
Ou, como explícita retorna valor foi eliminada, você poderia usá-lo para a matriz de permutação ou alternar as trocas. Por exemplo:
static int[] MatrixDecompose(ref double[][] matrix, out int toggle)
Você pode querer eliminar o parâmetro toggle para simplificar a assinatura do método, se você não pretende usar a decomposição da matriz para calcular um determinante.
Outra área de MatrixDecompose, você pode querer modificar é esta declaração:
if (Math.Abs(result[j][j]) < 1.0E-20) return null;
Em palavras, significa que este código: "Se, mesmo depois de trocar duas linhas porque não havia um valor 0.0 na diagonal principal, há ainda um 0,0 lá, em seguida, retornar nulo." Você pode desejar modificar o valor de epsilon arbitrário para a igualdade de zero seleção de 1.0 e-20 para algum outro valor baseado nas características de precisão do seu sistema. E em vez de retornar null, você pode querer lançar uma exceção; se o método a ser chamado como parte da inversão de matriz, a inversão falharia. Finalmente, se você estiver usando decomposição da matriz para alguma finalidade diferente de inversão de matriz, você pode querer eliminar completamente esta declaração.
Inversão de matriz
A chave para usar a decomposição da matriz para inverter uma matriz é escrever um método auxiliar que resolve um sistema de equações. Esse método auxiliar chave é apresentado em Figura 4.
Figura 4 o método de HelperSolve
static double[] HelperSolve(double[][] luMatrix, double[] b) { // solve luMatrix * x = b int n = luMatrix.Length; double[] x = new double[n]; b.CopyTo(x, 0); for (int i = 1; i < n; ++i) { double sum = x[i]; for (int j = 0; j < i; ++j) sum -= luMatrix[i][j] * x[j]; x[i] = sum; } x[n - 1] /= luMatrix[n - 1][n - 1]; for (int i = n - 2; i >= 0; --i) { double sum = x[i]; for (int j = i + 1; j < n; ++j) sum -= luMatrix[i][j] * x[j]; x[i] = sum / luMatrix[i][i]; } return x; }
O método de HelperSolve encontra uma matriz x, que, quando multiplicado por uma matriz de LU dá a matriz b. O método é bastante inteligente, e você pode apenas compreender pelo rastreamento através de alguns exemplos. Existem dois loops. O primeiro circuito usa a substituição para frente na parte inferior da matriz LU. O segundo loop usa a substituição para trás na parte superior da matriz LU. Algumas implementações de decomposição de diferentes matriz chamam seu método análogo algo como luBackSub.
Embora o código seja curto, é complicado, mas não deve haver qualquer cenários onde você precisa modificar o código. Aviso que auxiliarSolve aceita uma matriz de LU do MatrixDecompose, mas não usa o parâmetro out do perm. Isso implica que helpersolve é na verdade um código auxiliar método e necessidades adicionais de wrapper para resolver um sistema de equações. Se você refatorar o código neste artigo para um design OOP, talvez queira fazer HelperSolve um método particular.
Com o método de HelperSolve no lugar, o método de inversão de matriz pode ser implementado, conforme mostrado no Figura 5.
Figura 5 o método de MatrixInverse
static double[][] MatrixInverse(double[][] matrix) { int n = matrix.Length; double[][] result = MatrixDuplicate(matrix); int[] perm; int toggle; double[][] lum = MatrixDecompose(matrix, out perm, out toggle); if (lum == null) throw new Exception("Unable to compute inverse"); double[] b = new double[n]; for (int i = 0; i < n; ++i) { for (int j = 0; j < n; ++j) { if (i == perm[j]) b[j] = 1.0; else b[j] = 0.0; } double[] x = HelperSolve(lum, b); for (int j = 0; j < n; ++j) result[j][i] = x[j]; } return result; }
Novamente, o código é complicado. O algoritmo de inversão baseia-se na idéia de que o produto de uma matriz M e sua inversa é a matriz identidade. Método MatrixInverse essencialmente resolve um sistema de equações Ax = b onde A é uma decomposição da matriz de LU e as constantes de b são 0.0 ou 1.0 e correspondem à matriz de identidade. Observe que o MatrixInverse usa a matriz de perm da chamada para MatrixDecompose.
Chamando o MatrixInverse método poderia ser assim:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0); double[][] inverse = MatrixInverse(m); Console.WriteLine(MatrixAsString(inverse));
Para resumir, uma operação de matriz importante é a inversão de matriz, que é bastante difícil. Uma abordagem é decompor a matriz original, escrever um ajudante resolver método que realiza a substituição para frente e para trás, e então usar a decomposição em conjunto com a matriz de permutação de LU e o auxiliar resolver método para encontrar o inverso. Essa abordagem pode parecer complicada, mas é geralmente mais eficiente e mais fácil do que uma matriz inversa de computação diretamente.
Determinante da matriz
Com um método de decomposição de matriz na mão, é fácil de codificar um método para calcular o determinante de uma matriz:
static double MatrixDeterminant(double[][] matrix) { int[] perm; int toggle; double[][] lum = MatrixDecompose(matrix, out perm, out toggle); if (lum == null) throw new Exception("Unable to compute MatrixDeterminant"); double result = toggle; for (int i = 0; i < lum.Length; ++i) result *= lum[i][i]; return result; }
Como se constata, um pouco surpreendentemente, o determinante de uma matriz é apenas mais ou menos (dependendo do valor de alternar) o produto dos valores na diagonal principal da matriz decomposição LUP. Observe que há uma conversão de tipo implícito do valor da alternância de int para duplo. Além de adicionar verificação para MatrixDeterminant de erros, você pode querer adicionar um curto-circuito de retorno em situações onde a matriz de entrada tem tamanho 1 (depois voltar o valor único) ou tamanho 2 x 2 (em seguida, retornar um * d - b * c). Chamando o método determinante poderia ficar assim:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0); double det = MatrixDeterminant(m); Console.WriteLine("Determinant = " + det.ToString("F2"));
Resolução de sistemas de equações
O método HelperSolve pode ser facilmente adaptado para resolver um sistema de equações lineares:
static double[] SystemSolve(double[][] A, double[] b) { // Solve Ax = b int n = A.Length; int[] perm; int toggle; double[][] luMatrix = MatrixDecompose(A, out perm, out toggle); if (luMatrix == null) return null; // or throw double[] bp = new double[b.Length]; for (int i = 0; i < n; ++i) bp[i] = b[perm[i]]; double[] x = HelperSolve(luMatrix, bp); return x; }
Aqui está o código que produziu a captura de tela em Figura 1 para resolver o seguinte sistema:
3x0 + 7x1 + 2x2 + 5x3 = 49 x0 + 8x1 + 4x2 + 2x3 = 30 2x0 + x1 + 9x2 + 3x3 = 43 5x0 + 4x1 + 7x2 + x3 = 52 double[][] m = MatrixCreate(4, 4); m[0][0] = 3.0; m[0][1] = 7.0; m[0][2] = 2.0; m[0][3] = 5.0; m[1][0] = 1.0; m[1][1] = 8.0; m[1][2] = 4.0; m[1][3] = 2.0; m[2][0] = 2.0; m[2][1] = 1.0; m[2][2] = 9.0; m[2][3] = 3.0; m[3][0] = 5.0; m[3][1] = 4.0; m[3][2] = 7.0; m[3][3] = 1.0; double[] b = new double[] { 49.0, 30.0, 43.0, 52.0 }; double[] x = SystemSolve(m, b); Console.WriteLine("\nSolution is x = \n" + VectorAsString(x));
Observe que SystemSolve reorganiza seu parâmetro de entrada b usando a matriz de perm de MatrixDecompose antes de chamar HelperSolve.
Compreendendo a matriz permutação
As últimas algumas linhas de saída da imagem em Figura 1 indicam que é possível multiplicar as matrizes L e U, de modo a obter a matriz original. Saber como fazer isso não permitirá que você resolver problemas práticos de matriz, mas ele vai te ajudar a entender a parte de P da decomposição de LUP. Regenerar um matrix original de seus componentes L e U também pode ser útil para testar seus métodos de biblioteca de matriz para a consistência.
Uma maneira para regenerar um matrix original após a decomposição de LUP é multiplicar L e U e, em seguida, permutar as linhas de produto usando a matriz de P:
// create matrix m // call MatrixDecompose, returning lu and perm // extract the lower part of lu as matrix 'lower' // extract the upper part of lu as matrix 'upper' double[][] lu = MatrixProduct(lower, upper); double[][] orig = UnPermute(lu, perm); if (MatrixAreEqual(orig, m, 0.000000001) == true) Console.WriteLine("L and U unpermuted using perm array");
O método UnPermute pode ser codificado como este:
static double[][] UnPermute(double[][] luProduct, int[] perm) { double[][] result = MatrixDuplicate(luProduct); int[] unperm = new int[perm.Length]; for (int i = 0; i < perm.Length; ++i) unperm[perm[i]] = i; // create un-perm array for (int r = 0; r < luProduct.Length; ++r) // each row result[r] = luProduct[unperm[r]]; return result; }
Uma segunda abordagem para regenerar um matrix original da sua decomposição LUP é converter a matriz de perm para uma matriz de perm e depois multiplicar a matriz de perm e a matriz de LU combinada:
// create matrix m // call MatrixDecompose, returning lu and perm // extract the lower part of lu as matrix 'lower' // extract the upper part of lu as matrix 'upper' double[][] permMatrix = PermArrayToMatrix(perm); double[][] orig = MatrixProduct(permMatrix, lu); if (MatrixAreEqual(orig, m, 0.000000001) == true) Console.WriteLine("L and U unpermuted using perm matrix");
Uma matriz de perm é uma matriz quadrada com valor de um 1,0 em cada linha e cada coluna. O método que cria uma matriz de perm de uma matriz de perm pode ser codificado da seguinte forma:
static double[][] PermArrayToMatrix(int[] perm) { // Doolittle perm array to corresponding matrix int n = perm.Length; double[][] result = MatrixCreate(n, n); for (int i = 0; i < n; ++i) result[i][perm[i]] = 1.0; return result; }
Conclusão
Existem muitos algoritmos que exigem a resolução de um sistema de equações lineares, encontrar o inverso de uma matriz ou encontrar o determinante de uma matriz.Utilizando a decomposição da matriz é uma técnica eficaz para a realização de todas essas operações.O código apresentado aqui pode ser usado em situações em que você quer sem dependências externas em sua base de código, ou você precisa a capacidade de personalizar as operações a fim de melhorar o desempenho ou modificar a funcionalidade.Tome a pílula vermelha!
Dr. James McCaffrey* trabalha para a Volt Information Sciences Inc., onde gerencia o treinamento técnico para engenheiros de software trabalham no Microsoft Redmond, Wash., campus. Ele trabalhou em vários produtos da Microsoft, como o Internet Explorer e o MSN Busca. McCaffrey é o autor de .NET Test Automation Recipes (Apress, 2006). Ele pode ser contatado em jammc@microsoft.com.*
Agradecemos aos seguintes especialistas técnicos pela revisão deste artigo: Paul Koch e Dan Liebling