用矩阵直接三角分解法求解方程组

//用矩阵直接三角分解法求解方程组
#include<iostream>
#include<stdio.h>
using namespace std;
void main()
{
    int n,i;
    float *a,*x;
    cout<<"请输入方程组元的个数:";
    cin>>n;

    a=new float[n*(n+1)];       //开辟增广矩阵的存储空间
    x=new float[n];             //开辟方程组解的存储空间

    cout<<"请输入方程组的系数矩阵："<<endl;
    for(i=0;i<n*(n+1);i++)
    {
         if ((i+1)%(n+1)==0) continue;
         cin>>*(a+i);
    }

    cout<<"请输入方程组的常数项："<<endl;
    for(i=n;i<n*(n+1);i+=n+1)
      cin>>*(a+i);

    void DirectLU(float*,int,float*);
    DirectLU(a,n,x);
    cout<<"方程组的解是："<<endl;
    for(i=0;i<n;i++)
    cout<<"x["<<i<<"]="<<float(*(x+i))<<endl;
}
void DirectLU(float*u,int n,float*x)
{
    int i,r,k;
    for(r=0;r<=n-1;r++)
    {
        for(i=r;i<=n;i++)
        for(k=0;k<=r-1;k++)
            *(u+r*(n+1)+i)-=*(u+r*(n+1)+k)*(*(u+k*(n+1)+i));
        for(i=r+1;i<=n-1;i++)
        {
            for(k=0;k<=r-1;k++)
                *(u+i*(n+1)+r)-=*(u+i*(n+1)+k)*(*(u+k*(n+1)+r));
            *(u+i*(n+1)+r)/=*(u+r*(n+1)+r);
        }
    }

    for(i=n-1;i>=0;i--)
    {
        for(r=n-1;r>=i+1;r--)
            *(u+i*(n+1)+n)-=*(u+i*(n+1)+r)*x[r];
        x[i]=*(u+i*(n+1)+n)/(*(u+i*(n+1)+i));
    }

}