C語言 例題6-5 LU分解法 LU Decomposition
Let A be a square matrix. An LU factorization refers to the factorization of A, with proper row and/or column orderings or permutations, into two factors, a lower triangular matrix L and an upper triangular matrix U, A=LU.
For a general n×n matrix A, we assume that an LU decomposition exists, and write the form of L and U explicitly. We then systematically solve for the entries in L and U from the equations that result from the multiplications necessary for A=LU.
// C Program to decompose a matrix into
// lower and upper traingular matrix
#include <stdio.h>
#include <stdlib.h>
const int MAX = 3;
void luDecomposition(int mat[][MAX], int n)
{
int lower[n][n], upper[n][n];
for (int i=0;i<=n;i++)
{
for (int j=0;j<=n;j++)
{
lower[i][j]=0;
upper[i][j]=0;
}
}
//memset(lower, 0, sizeof(lower));
//memset(upper, 0, sizeof(upper));
// Decomposing matrix into Upper and Lower
// triangular matrix
for (int i = 0; i < n; i++) {
// Upper Triangular
for (int k = i; k < n; k++) {
// Summation of L(i, j) * U(j, k)
int sum = 0;
for (int j = 0; j < i; j++)
sum += (lower[i][j] * upper[j][k]);
// Evaluating U(i, k)
upper[i][k] = mat[i][k] - sum;
}
// Lower Triangular
for (int k = i; k < n; k++) {
if (i == k)
lower[i][i] = 1; // Diagonal as 1
else {
// Summation of L(k, j) * U(j, i)
int sum = 0;
for (int j = 0; j < i; j++)
sum += (lower[k][j] * upper[j][i]);
// Evaluating L(k, i)
lower[k][i] = (mat[k][i] - sum) / upper[i][i];
}
}
}
// setw is for displaying nicely
printf(" Lower Triangular Upper Triangular\n");
// Displaying the result :
for (int i = 0; i < n; i++) {
// Lower
for (int j = 0; j < n; j++)
printf("\t%2d\t", lower[i][j]);
printf("\t");
// Upper
for (int j = 0; j < n; j++)
printf("\t%2d\t", upper[i][j]);
printf("\n");
}
}
// Driver code
int main()
{
int mat[3][3] = { { 2, -1, -2 },
{ -4, 6, 3 },
{ -4, -2, 8 } };
luDecomposition(mat, 3);
return 0;
}
輸出畫面
Lower Triangular Upper Triangular
1 0 0 2 -1 -2
-2 1 0 0 4 -1
-2 -1 1 0 0 3
沒有留言:
張貼留言