高斯消元——矩阵分解

高斯消元中的矩阵分解，好用在，可以对系数相同，不同不用重复计算。
矩阵（LU）分解

#include<iostream>
#include<cstdio>
using namespace std;
const int maxn=1000;
int n;
double a[maxn][maxn];
double x[maxn];
//分解系数矩阵A=L*U;
void show()
{
    cout<<endl;
        for(int i=0;i<n;i++)
        {
            for(int j=0;j<n;j++)
            printf("%lf  ",a[i][j]);
            cout<<endl;
        }
        return;
}
void fen(int cur)
{
    if(cur==n)
    {
        show();
        return;
    }
    for(int i=cur;i<n;i++)
    {//U的处理
        double sum=0;
        int d=0;
        while(d<cur)
        {
            sum+=a[cur][d]*a[d][i];
            d++;
        }
        a[cur][i]-=sum;
    }
    //show();
    for(int i=cur+1;i<n;i++)
    {//L的处理
        double sum=0;
        int d=0;
        while(d<cur)
        {
            sum+=a[i][d]*a[d][cur];
            d++;
        }
        a[i][cur]=(a[i][cur]-sum)/a[cur][cur];
    }
    //show();
    fen(cur+1);
}


void f_x()
{
    for(int i=n-1;i>=0;i--)
    {
        double sum=0;
        for(int j=i+1;j<=n-1;j++)
            sum+=a[i][j]*x[j];
        //cout<<sum<<endl;
        x[i]=(x[i]-sum)/a[i][i];
    }
    for(int i=0;i<n;i++)
        cout<<"x"<<i+1<<"= "<<x[i]<<endl;
}
void f_y()
{

    for(int i=1;i<n;i++)
    {
        double sum=0;
        for(int j=0;j<i;j++)
            sum+=a[i][j]*x[j];
        //cout<<sum<<endl;
        x[i]=(x[i]-sum);
    }
    for(int i=0;i<n;i++)
        cout<<"y"<<i+1<<"= "<<x[i]<<endl;
}


int main()
{
    cout<<"Input A"<<endl;
    cin>>n;
    for(int i=0;i<n;i++)
        for(int j=0;j<n;j++)
            cin>>a[i][j];

    fen(0);cout<<endl;
    cout<<"Input the number of X[]"<<endl;
    int m;cin>>m;
    while(m--)
    {
        for(int i=0;i<n;i++)
            cin>>x[i];
        f_y();cout<<endl;
        f_x();
    }

    return 0;
}
/*

4
4 2 1 5
8 7 2 10
4 8 3 6
12 6 11 20
3
-2 -7 -7 -3
*/