Este artículo proviene de un motor de traducción automática.
C#
Descomposición de una matriz
Descargar el ejemplo de código
Descomposición de la matriz, una técnica que descompone una matriz cuadrada numérica en dos diferentes matrices cuadradas, es la base para resolver eficientemente un sistema de ecuaciones, que a su vez es la base para invertir una matriz.E invertir una matriz es una parte de muchos algoritmos importantes.Este artículo presenta y explica el código C# que realiza la descomposición de la matriz, inversión de la matriz, un sistema de solución de ecuaciones y las operaciones conexas.
En efecto, descomposición de la matriz no es un tema llamativo, pero una colección de métodos de matriz puede ser un complemento importante a la biblioteca de código personal.Los métodos se explican por lo que puede modificar el código fuente para satisfacer sus propias necesidades.Además, algunas de las técnicas utilizadas en los métodos de la matriz pueden ser reutilizado en otros escenarios de codificación.
La mejor manera para que usted consiga una sensación para el tipo de información presentada en este artículo es echar un vistazo a la captura de pantalla en figura 1.El programa de demostración comienza creando una matriz 4 x 4 y mostrar sus valores.A continuación, la matriz se descompone en lo que denomina una matriz LUP.L significa inferior y U para la parte superior.La parte de P (P significa permutación) es una matriz de valores {3.120} e indica que filas 0 y 3 fueron intercambiados durante el proceso de descomposición.La descomposición también genera un valor de alternar de -1, que indica que se produjo un número impar de intercambios de fila.El programa de demostración muestra la descomposición de dos maneras: primero como una matriz combinada de LU y matrices L y U luego por separado.A continuación, el programa calcula y muestra la inversa de la matriz original, utilizando la matriz LUP detrás de las escenas.El programa de demostración calcula el determinante de la matriz original, mediante la descomposición.A continuación, la inversa de la matriz se utiliza para resolver un sistema de ecuaciones lineales y concluye mediante la combinación de las matrices L y U en la matriz original.
.png)
Figura 1 matriz descomposición Demo
Pero ¿por qué ir a todo el problema de la creación de un método de descomposición de la matriz personalizada y una biblioteca de métodos relacionados?Si bien se dispone de muchas herramientas de matriz de independiente, a veces puede ser difíciles de integrar en un sistema o aplicación.Y a pesar de la importancia fundamental de la descomposición de la matriz, hay pocas implementaciones de código .NET gratis, sin derechos de autor disponibles; las que existen a menudo no se explican con suficiente detalle para modificar el código fuente para adaptarse a sus escenarios de codificación.
Este artículo asume que tienes conocimientos de programación de C# intermedios y al menos un conocimiento básico de la terminología y las operaciones de la matriz.En este artículo se presenta todo el código clave de C#.El código también está disponible en el sitio de descarga de código MSDN en archve.msdn.microsoft.com/mag201212matrix.
Definición de matriz
Hay varias maneras de implementar una matriz en C#.El enfoque tradicional y el que se usa en este artículo, es utilizar una matriz de matrices, a veces llamado una matriz escalonada.Por ejemplo, este código define una matriz con tres filas y dos columnas:
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
A diferencia de la mayoría de los lenguajes programación, C# tiene un tipo de matriz multidimensional integrado, que proporciona un enfoque alternativo. Por ejemplo:
double[,] m = new double[3,2];
m[2,1] = 5.0;
Un tercer enfoque para implementar matrices en C# es utilizar una sola matriz combinada con la manipulación del í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;
Sin importar el esquema de almacenamiento utilizado, matrices pueden implementarse usando un enfoque de método estático o una programación orientada a objetos. Por ejemplo, podría parecerse un enfoque de programación orientada a objetos:
public class MyMatrix
{
int m; // number rows
int n; // number columns
data[][]; // the values
...
}
No hay sola mejor opción para el diseño de la matriz; el mejor diseño depende de la situación particular de codificación en que está operando y sobre sus preferencias de codificación. Este artículo utiliza un método estático porque es el más fácil de entender y refactorizar.
Cuando se utiliza un diseño de matriz de matrices de matrices, porque cada fila debe asignarse por separado, a menudo es conveniente definir un método auxiliar para realizar la asignación de memoria. Por ejemplo:
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;
}
El método puede ser llamado como tal:
double[][] m = MatrixCreate(3,2);
m[2][1] = 5.0;
Este método muestra una de las ventajas de crear su propia biblioteca de métodos de matriz: Si desea mejorar el rendimiento, se puede omitir la comprobación de los parámetros de entrada — a expensas de aumentar el riesgo de que el código de llamada causando una excepción. (Para mantener este artículo corto, la mayoría de la comprobación de errores se ha eliminado.) Otra ventaja es que puede personalizar su biblioteca para optimizar su situación exacta. La principal desventaja de crear su propia biblioteca es que puede tomar más tiempo que usando una biblioteca existente.
Otro método conveniente para agregar a la biblioteca de la matriz es uno que puede utilizarse para mostrar una matriz como una cadena. Aquí es una posibilidad:
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;
}
Puede que desee parametrizar el número de decimales a mostrar, o el relleno de ancho de columna o ambos.
Multiplicación de matrices
Microsoft .NET Framework 4 y posterior proporcionan una buena manera para mejorar significativamente el rendimiento de un método de multiplicación de la matriz. Multiplicación de matrices se ilustra en la figura 2.
.png)
Figura 2 multiplicación de matrices
Tenga en cuenta que el cálculo de cada valor de la celda del resultado no depende de ningún otro valor de la celda en el resultado, así que cada cálculo es independiente y potencialmente se pueden realizar en paralelo en un equipo con varios procesadores. Aquí es un enfoque estándar para la multiplicación de matrices:
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 la multiplicación de matrices es una operación de O(n^3), el rendimiento puede ser un problema. Por ejemplo, si la matriz A tiene tamaño 300 x 200 y la matriz B tiene tamaño 200 x 400, calcular el producto de la A y B requiere 300 * 200 * 400 = 24.000.000 multiplicaciones. La biblioteca paralelo de tarea (TPL) en el espacio de nombres de System.Threading.Tasks en el .NET Framework 4 y más tarde facilita el código una simple versión paralelizada de multiplicación de matrices. Una posibilidad es:
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 versión chuletas hasta los cálculos por filas. Detrás de las escenas, la TPL genera todo el código de plomería de sincronización complejo para realizar los cálculos de varios procesadores.
Pruebas de consistencia
Un aspecto interesante de una biblioteca de métodos que están relacionados entre sí es que a menudo es posible probarlos comprobando si producen resultados consistentes. Por ejemplo, supongamos que tiene un método que crea una matriz al azar:
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;
}
Además, supongamos que tiene una matriz que crea la matriz de identidad — es decir, una matriz cuadrada con 1.0s de la diagonal principal, 0.0s en otros lugares:
static double[][] MatrixIdentity(int n)
{
double[][] result = MatrixCreate(n, n);
for (int i = 0; i < n; ++i)
result[i][i] = 1.0;
return result;
}
Y supongamos que tiene un método que compara dos matrices para la igualdad:
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 el método MatrixAreEqual no comparar los valores de celda de igualdad exacta porque los valores son de tipo double. En cambio, el método comprueba si todos los valores de celda son muy cerca (en epsilon) mutuamente.
Porque el producto de cualquier matriz cuadrada m y la matriz de identidad de la misma dimensión son iguales a la original m de matriz, puede probar el método de producto de la matriz a lo largo de las líneas de este 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.
Descomposición de una matriz
Descomposición de la matriz tiene una matriz cuadrada M y computa dos nuevas matrices cuadradas que al multiplicar juntos dar la matriz original. La idea es similar a factoring número ordinario: el número 6 puede ser factorizado en 2 y 3 porque 2 * 3 = 6. Al principio puede parecer que hay poco sentido en descomposición de una matriz, pero resulta que la descomposición de la matriz hace la muy difícil tarea de inversión de la matriz bastante un poco más simple. Hay muchos tipos diferentes de descomposición de la matriz, y cada tipo puede calcularse mediante varios algoritmos diferentes. La técnica presentada en este artículo se llama descomposición LUP y utiliza el método de Doolittle pivotante parcial.
Para entender la descomposición LUP, resulta útil comprender primero la descomposición LU más simple, que fue presentada por el famoso matemático Alan Turing. Supongamos que tiene esta matriz 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
Una posible descomposición de LU de M es 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
Y 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
Esto funciona porque L * U = M. Observe que la matriz L inferior tiene 1.0s en diagonal y 0.0s en la parte superior derecha. En otras palabras, los valores de celda significativa de la matriz inferior están en la parte inferior izquierda. Del mismo modo, los valores de celda significativa de la matriz superior son de la diagonal principal y en la parte superior derecha.
Observe también que no hay ninguna superposición de las ubicaciones de los valores de celda significativa en L y U. Así, en lugar de generar dos matrices de resultados, L y U, descomposición de la matriz normalmente almacena los resultados en una sola matriz que sostiene L y U para ahorrar espacio de memoria superiores e inferiores.
LUP matriz descomposición es una variación leve pero importante en la descomposición LU. Descomposición de LUP toma una matriz M y produce matrices L y U, pero también una matriz P. El producto de L y U en LUP no es exactamente la matriz original M pero en cambio es una versión de M donde algunas de las filas han sido reordenados. La forma en que se han reorganizado las filas se almacena en la matriz de P; Esta información puede utilizarse para reconstruir la matriz original.
Un primo cercano a la descomposición de Doolittle presentado en este artículo se denomina descomposición de Crout. La principal diferencia entre Doolittle y Crout es que Doolittle coloca 1.0s de la diagonal principal de la matriz de L y Crout 1.0s de la diagonal principal de la matriz de la U.
La razón por la descomposición de LUP se utiliza más a menudo que la descomposición LU es sutil. Como verás poco, descomposición de la matriz se utiliza para calcular la inversa de una matriz. Sin embargo, cuando la descomposición de la matriz se utiliza como un ayudante para inversión de la matriz, resulta que la inversión se producirá un error si hay un valor 0.0 de la diagonal principal de la matriz de LU. Tan en descomposición LUP, cuando un 0.0 termina el valor de la diagonal principal, el algoritmo de intercambio dos filas para mover el valor 0.0 off la diagonal y realiza un seguimiento de qué filas se intercambiaron en la matriz P.
Figura 3 muestra un método de descomposición de la matriz.
Figura 3 método de descomposición de la matriz
static double[][] MatrixDecompose(double[][] matrix,
out int[] perm, out int toggle)
{
// Doolittle LUP decomposition.</span>
<span class="sentence">// assumes matrix is square.</span>
<span class="sentence">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;
}
El método podría llamarse así:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0);
int[] perm;
int toggle;
double luMatrix = MatrixDecompose(m, out perm, out toggle);
El método MatrixDecompose acepta como entrada una matriz cuadrada. El método tiene tres valores de retorno. El retorno explícito es una matriz de LU permutada. El método devuelve dos valores como parámetros out. Uno es una matriz de permutación que contiene información sobre cómo fueron permutadas las filas. El segundo parámetro es un valor de palanca que es + 1 o -1 según el número de intercambios de fila incluso (+ 1) o impar (-1). El valor de la palanca no se utiliza para la inversión de la matriz, pero es necesario si la descomposición de la matriz se utiliza para calcular el determinante de una matriz.
El método de MatrixDecompose es bastante complicado, pero siendo realistas, hay sólo algunos detalles, que es necesario comprender para modificar el código. La versión presentada aquí asigna memoria nueva para la matriz de retorno de LU utilizando un método auxiliar MatrixDuplicate:
static double[][] MatrixDuplicate(double[][] matrix)
{
// assumes matrix is not null.</span>
<span class="sentence">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;
}
Un enfoque alternativo es calcular el resultado en la matriz de entrada con el fin de ahorrar memoria. Con C# semántica, esto haría que el parámetro de matriz en un parámetro ref porque se utiliza para entrada y salida. Con este método, la firma del método sería:
static void MatrixDecompose(ref double[][] matrix, out int[] perm,
out int toggle)
O bien, porque la explícita devuelve valor ha sido eliminada, se puede utilizar para la matriz permutación o alternar los intercambios. Por ejemplo:
static int[] MatrixDecompose(ref double[][] matrix, out int toggle)
Quizá desee eliminar el parámetro toggle para simplificar la firma del método si no desea utilizar la descomposición de matriz para calcular un determinante.
Otra área de MatrixDecompose puede modificar es esta declaración:
if (Math.Abs(result[j][j]) < 1.0E-20)
return null;
En palabras, este código significa: "Si, incluso después de intercambiar dos filas porque había un valor 0.0 de la diagonal principal, todavía hay un 0.0 allí, entonces devolver null." Puede modificar el valor de epsilon arbitraria para la igualdad de cero cheque de 1.0E-20 a algún otro valor se basa en las características de la precisión del sistema. Y en lugar de devolver null, puede que desee iniciar una excepción; Si el método se llama como parte de la inversión de la matriz, no la inversión. Por último, si está utilizando la descomposición de la matriz para algún propósito que no sea la inversión de la matriz, puede eliminar por completo esta declaración.
Inversión de la matriz
La clave para utilizar la descomposición de la matriz para invertir una matriz es escribir un método auxiliar que resuelve un sistema de ecuaciones. Este método de clave auxiliar se presenta en figura 4.
Figura 4 el 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;
}
El método de HelperSolve encuentra una matriz x que multiplicado por una matriz de LU da matriz b. El método es muy inteligente, y sólo usted puede entender por seguimiento a través de algunos ejemplos. Hay dos bucles. El primer bucle utiliza sustitución hacia delante en la parte inferior de la matriz de LU. El segundo bucle utiliza sustitución hacia atrás en la parte superior de la matriz de LU. Algunas implementaciones de descomposición de la matriz diferentes llaman a su método análogo algo como luBackSub.
Aunque el código es corto, es complicado, pero no debería haber ningún escenario donde tendría que modificar el código. Aviso que ayudanteSolve acepta una matriz de LU de MatrixDecompose pero no utiliza el parámetro out de perm. Esto implica que helpersolve es en realidad un código de contenedor adicional ayudante necesidades y método para resolver un sistema de ecuaciones. Si refactorizar el código de este artículo a un diseño de programación orientada a objetos, puede hacer HelperSolve un método privado.
Con el método HelperSolve en su lugar, se puede implementar el método de inversión de la matriz, como se muestra en la figura 5.
Figura 5 el 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;
}
Una vez más, el código es complicado. El algoritmo de inversión se basa en la idea de que el producto de una matriz M y su inverso es la matriz identidad. Método MatrixInverse esencialmente resuelve un sistema de ecuaciones Ax = b donde A es una descomposición de la matriz de LU y las constantes b 1.0 sean 0,0 y corresponden a la matriz identidad. Observe que MatrixInverse utiliza la matriz de perm en la llamada a MatrixDecompose.
Llamar al MatrixInverse método podría tener este aspecto:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0);
double[][] inverse = MatrixInverse(m);
Console.WriteLine(MatrixAsString(inverse));
Para resumir, una operación importante de la matriz es la inversión de la matriz, que es bastante difícil. Un enfoque consiste en descomponer la matriz original, escribir un ayudante resolver el método que realiza la sustitución hacia adelante y hacia atrás y luego uso la descomposición junto con la matriz de permutación de LU y el ayudante resolver el método para encontrar la inversa. Este enfoque puede parecer complicado, pero suele ser más eficiente y más fácil que una matriz inversa de computación directamente.
Determinante de la matriz
Con un método de descomposición de la matriz en la mano, es fácil de código de un método para calcular el determinante de una 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;
}
Resulta sorprendente que, un poco, el determinante de una matriz es sólo más o menos (dependiendo del valor de la palanca) el producto de los valores de la diagonal principal de la matriz de descomposición LUP. Observe que hay una conversión de tipo implícito del valor de la palanca de int a doble. Además de agregar a MatrixDeterminant de comprobación de errores, puede que desee agregar un cortocircuito regreso en situaciones donde la matriz de entrada tiene tamaño 1 (se vuelve el único valor) o tamaño 2 x 2 (volver un * d - b * c). Llamar al método determinante podría tener este aspecto:
double[][] m = MatrixRandom(4, 4, -9.0, 9.0, 0);
double det = MatrixDeterminant(m);
Console.WriteLine("Determinant = " + det.ToString("F2"));
Resolver sistemas de ecuaciones
El método HelperSolve puede ser fácilmente adaptado para resolver un sistema de ecuaciones lineales:
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;
}
Aquí está el código que produjo la captura de pantalla en figura 1 para resolver el siguiente 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));
Aviso que SystemSolve reorganiza su parámetro de entrada b utilizando la matriz de la ondulación permanente de MatrixDecompose antes de llamar al HelperSolve.
Entender la matriz permutación
Las últimas pocas líneas de salida en la captura de pantalla en figura 1 indican que es posible multiplicar las matrices L y U de manera que se obtenga la matriz original. Saber cómo hacer esto no le permitirá resolver problemas prácticos de la matriz, pero le ayudará a entender la parte P de descomposición LUP. Regeneración de la matriz original de sus componentes L y U también puede ser útil para probar los métodos de biblioteca de matriz de consistencia.
Una forma de generar una matriz original después de la descomposición de LUP es multiplicar L y U, y luego permutar las filas del producto utilizando la matriz 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");
El método UnPermute puede codificarse 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;
}
Un segundo método para regenerar la matriz original de su descomposición LUP es convertir la matriz de perm a una matriz de perm y luego multiplicar la matriz de perm y la LU combinado:
// 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");
Una matriz de perm es una matriz cuadrada con un 1.0 valor en cada fila y cada columna. El método que crea una matriz de ondulación permanente de una matriz de ondulación permanente puede ser codificado como tal:
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;
}
En resumen
Hay muchos algoritmos que requieren para resolver un sistema de ecuaciones lineales, encontrar la inversa de una matriz o encontrar el determinante de una matriz. Usando la descomposición de la matriz es una técnica eficaz para realizar todas estas operaciones. El código presentado aquí puede utilizarse en situaciones donde no desee dependencias externas en su base de código o necesita la capacidad de personalizar las operaciones con el fin de mejorar el rendimiento o modificar la funcionalidad. Tomar la píldora roja!
Dr.James McCaffrey trabaja para voltios Information Sciences Inc., donde dirige la formación técnica para ingenieros de software trabajando en el Microsoft Redmond, Wash., campus. Ha colaborado en el desarrollo de varios productos de Microsoft como, por ejemplo, Internet Explorer y MSN Search. McCaffrey es el autor de .NET Test Automation Recipes (Apress, 2006). Él puede ser contactado en jammc@microsoft.com.
Gracias a los siguientes expertos técnicos por su ayuda en la revisión de este artículo: Paul Koch y Dan Liebling