LU分解法求解线性方程组

LU分解法是求解线性方程组的一种算法。

先将系数矩阵分别转换成上三角矩阵U和下三角矩阵L，其中U(k,j)=a(k,j) - ∑L(k,m)*U(m,j)    m<k，L(i,k)=(a(i,k) - ∑L(i,m)*U(m,j))/(u(k,k))   m<k;

求解Ly=b，b为因变量矩阵；

求解Ux=y；

c++代码如下：

#include<stdio.h>

using namespace std;

#define  n  3    //矩阵阶数
int main()
{
	int i, j, k, m;
	float t, a[n][n], u[n][n], l[n][n], b[n], x[n], y[n];
	for (i = 0; i<n; i++)
	for (j = 0; j<n; j++)
	{
		scanf("%f", &t);
		a[i][j] = t;
	}
	for (i = 0; i<n; i++)
		scanf("%f", &b[i]);
	for (k = 0; k<n; k++)
	{
		for (j = k; j<n; j++)
		{
			t = 0;
			for (m = 0; m <= k - 1; m++)
				t = t + l[k][m] * u[m][j];
			u[k][j] = a[k][j] - t;
		}
		for (i = k + 1; i<n; i++)
		{
			t = 0;
			for (m = 0; m <= k - 1; m++)
				t = t + l[i][m] * u[m][k];
			l[i][k] = (a[i][k] - t) / u[k][k];
		}
	}
	for (i = 0; i<n; i++)
	{
		t = 0;
		for (j = 0; j <= i - 1; j++)
			t = t + l[i][j] * y[j];
		y[i] = b[i] - t;
	}
	for (i = n - 1; i >= 0; i--)
	{
		t = 0;
		for (k = i + 1; k<n; k++)
			t = t + u[i][k] * x[k];
		x[i] = (y[i] - t) / u[i][i];
	}
	for (i = 0; i<n; i++)
		printf("%f\n", x[i]);
	return 0;
}