//LLUU.cpp
//矩阵的LU分解
#include <iostream>
#include <fstream>
#include <cmath>
using namespace std;
class lluu
{
private:
int n;
double **a,**l,**u;
public:
lluu(int nn)
{
int i;
n = nn;
a = new double*[n];
for(i=0;i<n;i++) a[i] = new double[n+8];
l = new double*[n];
for(i=0;i<n;i++) l[i] = new double[n+8];
u = new double*[n];
for(i=0;i<n;i++) u[i] = new double[n+8];
}
void input();
void lu();
void output();
~lluu()
{
int i;
for(i=0;i<n;i++){delete[] a[i];}
delete[] a;
for(i=0;i<n;i++){delete[] l[i];}
delete[] l;
for(i=0;i<n;i++){delete[] u[i];}
delete[] u;
}
};
void lluu::input()
{
int i,j;
char str1[20];
cout<<"输入文件名:";
cin>>str1;
ifstream fin(str1);
if(!fin)
{cout<<"\n不能打开这个文件"<<str1<<endl;exit(1);}
for(i=0;i<n;i++)
for(j=0;j<n;j++)
fin>>a[i][j];
fin.close();
}
void lluu::lu()
{
int i,j,k;
for(k=0;k<n-1;k++)
{
for(i=k+1;i<n;i++)
{
a[i][k] = a[i][k]/a[k][k];
for(j=k+1;j<n;j++)
a[i][j] = a[i][j] - a[i][k]*a[k][j];
}
}
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
if(j<i) {l[i][j] = a[i][j];u[i][j] = 0;}
if(j==i) {l[i][j] = 1;u[i][j] = a[i][j];}
if(j>i) {l[i][j] = 0;u[i][j] = a[i][j];}
}
}
}
void lluu::output()
{
int i,j;
char str2[20];
cout<<"输出文件名:";
cin>>str2;
ofstream fout (str2);
if(!fout)
{cout<<"\n不能打开这个文件"<<str2<<endl;exit(1);}
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{ fout<<" "<<l[i][j];
cout<<" "<<l[i][j];
}
fout<<endl;cout<<endl;
}
fout<<endl;cout<<endl;
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{ fout<<" "<<u[i][j];
cout<<" "<<u[i][j];
}
fout<<endl;cout<<endl;
}
fout.close();
}
int main()
{
lluu c(4);
c.input();
c.lu();
c.output();
return 0;
}