C++ Exercises（十四）--解线性方程组

帮一个朋友做的作业题，现在真是什么专业都要学编程了呀。。

//向量类
/**//*
**author:phinecos
**date:7/17/2008
*/
classCVector
{
public:
CVector(unsignedintd=0);//由向量的维数创建向量，向量元素值初始化为
CVector(unsignedintd,double*pe);//由向量的维数和向量元素数组创建数组
CVector(doublex,doubley);//由两个元素生成二维向量
CVector(doublex,doubley,doublez);//由三个元素生成三维向量
CVector(CVector&v);//复制构造函数

~CVector();//析构函数

CVector&operator=(constCVector&v);//重载赋值运算符
CVector&operator+(void);//重载一元运算符+
CVectoroperator+(constCVector&v)const;//重载二元运算符+
CVector&operator-(void);//重载一元运算符-
CVectoroperator-(constCVector&v)const;//重载二元运算符-
CVectoroperator*(constCVector&v)const;//重载二元运算符*,表示向量的叉乘（向量积）
doubleoperator%(constCVector&v)const;//重载二元运算符%,表示向量的点乘（数量积）
CVectoroperator*(constdouble&d)const;//重载二元运算符*，表示向量的数乘

double&operator[](constunsignedinti);//重载操作符[]，对指定向量元素进行操作

unsignedintGetDegree()const;//获取向量维数
voidprintElements();
friendclassCLinearEquation;

private:
double*pElement;//向量元素存储地址
unsignedintnDegree;//向量的维数
voidinitByZero();//初始化为
voidinitByArray(double*pe);//用数组初始化
voidalloc(unsignedintn=0);//分配空间
voidreverse();//变号
voiddestory();
};

向量类实现
<!--<br /><br />Code highlighting produced by Actipro CodeHighlighter (freeware)<br />http://www.CodeHighlighter.com/<br /><br />-->#include"stdafx.h"
#include"Vector.h"
#include<cassert>
#include<memory.h>
#include<iostream>
usingnamespacestd;
/**//*
**author:phinecos
**date:7/17/2008
*/
CVector::CVector(unsignedintd):nDegree(d)
{//由向量的维数创建向量，向量元素值初始化为
if(d>0)
{
this->alloc(this->nDegree);
this->initByZero();//向量元素值初始化为

}
}

CVector::CVector(doublex,doubley):nDegree(2)
{
this->alloc(nDegree);
this->pElement[0]=x;
this->pElement[1]=y;
}

CVector::CVector(doublex,doubley,doublez):nDegree(3)
{//由三个元素生成三维向量
this->alloc(nDegree);
this->pElement[0]=x;
this->pElement[1]=y;
this->pElement[2]=z;
}
CVector::CVector(unsignedintd,double*pe):nDegree(d)
{//由向量的维数和向量元素数组创建数组
if(d>0)
{
this->alloc(this->nDegree);
this->initByArray(pe);
}
}
CVector::CVector(CVector&v)
{//复制构造函数
this->nDegree=v.GetDegree();
if(this->nDegree>0)
{
this->alloc(this->nDegree);
this->initByArray(v.pElement);
}
}

CVector&CVector::operator=(constCVector&v)
{//重载赋值运算符
if(v.GetDegree()>0)
{
this->destory();//销毁原数据
this->nDegree=v.GetDegree();//新向量大小
this->alloc(this->nDegree);//分配空间
this->initByArray(v.pElement);//复制数据
}
return*this;
}
CVector&CVector::operator+(void)
{//重载一元运算符+
return*this;
}

CVectorCVector::operator+(constCVector&v)const
{//重载二元运算符+
assert(this->nDegree==v.GetDegree());//两个向量的维数应该相等
CVectorresult(this->nDegree,this->pElement);
for(unsignedinti=0;i<this->nDegree;++i)
{
result.pElement[i]+=v.pElement[i];
}
returnresult;
}
CVector&CVector::operator-(void)
{//重载一元运算符-
this->reverse();
return*this;
}
CVectorCVector::operator-(constCVector&v)const
{//重载二元运算符-
assert(this->nDegree==v.GetDegree());//两个向量的维数应该相等
CVectorresult(this->nDegree,this->pElement);
for(unsignedinti=0;i<this->nDegree;++i)
{
result.pElement[i]-=v.pElement[i];
}
returnresult;
}
CVectorCVector::operator*(constCVector&v)const
{//重载二元运算符*,表示向量的叉乘（向量积）
assert(this->nDegree==v.GetDegree());
assert(this->nDegree==3);
doublea1=this->pElement[0];
doubleb1=this->pElement[1];
doublec1=this->pElement[2];

doublea2=v.pElement[0];
doubleb2=v.pElement[1];
doublec2=v.pElement[2];

returnCVector(b1*c2-b2*c1,c1*a2-a1*c2,a1*b2-a2*b1);
}

doubleCVector::operator%(constCVector&v)const
{//重载二元运算符%,表示向量的点乘（数量积）
assert(this->nDegree==v.GetDegree());
doubleresult=0.0f;
for(unsignedinti=0;i<this->nDegree;++i)
{
result+=this->pElement[i]*v.pElement[i];
}
returnresult;
}
double&CVector::operator[](constunsignedinti)
{//重载操作符[]，对指定向量元素进行操作
assert(i>=0&&i<this->nDegree);
returnthis->pElement[i];
}
voidCVector::reverse()
{
for(unsignedinti=0;i<this->nDegree;++i)
{
this->pElement[i]=-this->pElement[i];
}
}
voidCVector::alloc(unsignedintn)
{//分配大小为n的存储区
if(n>0)
{
this->pElement=newdouble[n];
}
}
voidCVector::initByZero()
{
for(unsignedinti=0;i<nDegree;++i)
{
this->pElement[i]=0;
}
}
voidCVector::initByArray(double*pe)
{
for(unsignedinti=0;i<nDegree;++i)
{
this->pElement[i]=pe[i];
}
}

voidCVector::destory()
{
if(this->pElement!=NULL)
{
delete[]pElement;
this->pElement=NULL;
this->nDegree=0;
}
}

CVector::~CVector(void)
{
this->destory();
}

unsignedintCVector::GetDegree()const
{//获取向量维数
returnthis->nDegree;
}

voidCVector::printElements()
{
for(unsignedinti=0;i<this->nDegree;++i)
{
cout<<"x["<<i<<"]="<<this->pElement[i]<<endl;
}
}

classCVector;//前向声明
//矩阵类
/**//*
**author:phinecos
**date:7/17/2008
*/
classCMatrix
{
public:
CMatrix(unsignedintr=0,unsignedintc=0);//由矩阵的行数和列数创建矩阵类对象，并为矩阵元素分配存储空间，将矩阵初始化为单位矩阵；
CMatrix(constchar*pFileName);//由矩阵存储文件名创建矩阵类对象，文件格式可参考附录，但不限于采用此格式；
CMatrix(constCMatrix&m);//复制构造函数，由已有矩阵类对象创建新的矩阵类对象；

unsignedintGetRowsNum()const;//获取矩阵的行数
unsignedintGetColumnsNum()const;//获取矩阵的列数

double&operator()(unsignedintr,unsignedintc);//重载运算符()，用于提取指定行(r)列(c)的元素值
CMatrix&operator=(constCMatrix&m);//重载运算符=,用于矩阵之间相互赋值
CMatrix&operator+()const;//重载一元运算符+,即取矩阵本身
CMatrixoperator-()const;//重载一元运算符-，即矩阵元素取相反数
CMatrixoperator+(constCMatrix&m)const;//重载二元运算符+,即两个矩阵求和
CMatrixoperator-(constCMatrix&m)const;//重载二元运算符-,即两个矩阵求差
CMatrixoperator*(constCMatrix&m)const;//重载二元运算符*,即两个矩阵求积
CMatrixoperator*(constdouble&x)const;//重载二元运算符*,即矩阵与数相乘
CMatrixoperator*(CVector&v)const;//重载二元运算符*，即矩阵与向量相乘
CMatrixoperator/(constdouble&x)const;//重载二元运算符/,即矩阵与数相除
CMatrixoperator^(constint&t)const;//重载二元运算符^,即矩阵求t次幂
voidoperator-=(constCMatrix&m);//重载二元运算符-=，即自减运算
voidoperator+=(constCMatrix&m);//重载二元运算符-=，即自加运算
voidoperator*=(constdoublex);//重载二元运算符*=，即自乘运算（矩阵）

CMatrixTranpose()const;//求矩阵的转置的函数
CMatrixInvert()const;//求矩阵的逆的函数
voidZeros();//矩阵归零化，将当前矩阵的所有元素归零
voidUnit();//矩阵单元化，将当前矩阵转换为单位矩阵
intAddRow(double*pe,unsignedintnr);//在矩阵的nr行位置插入一行，数据存放地址为pe，返回行标
intAddColumn(double*pe,unsignedintnc);//在矩阵的nc列位置插入一列，数据存放地址为pe，返回列标
double*DeleteRow(unsignedintnr);//删除矩阵nr行，返回该行元素值（临时存储地址）
double*DeleteColumn(unsignedintnc);//删除矩阵nc列，返回该列元素值（临时存储地址）

friendclassCLinearEquation;//声明友元类CLinearEquation
friendclassCVector;//声明友元类CVector

//析构函数，释放存储矩阵元素的空间
virtual~CMatrix(void);
//voidprintMatrix()const;
private:
double*pElement;//元素存储区
unsignedintnRow;//行数
unsignedintnColumn;//列数
voiddestroy();
voidalloc(unsignedintn=0);//分配空间
protected:
virtualvoidInitFromFile(constchar*pFileName);//从文件中初始化
public:
virtualintMatOut(constchar*pFileName);//将矩阵以pFileName为文件名进行文件输出
};

矩阵类实现
<!--<br /><br />Code highlighting produced by Actipro CodeHighlighter (freeware)<br />http://www.CodeHighlighter.com/<br /><br />-->#include"stdafx.h"
#include"Matrix.h"
#include"Vector.h"
#include<iostream>
#include<cassert>
#include<fstream>
usingnamespacestd;


CMatrix::CMatrix(constCMatrix&m)
{//复制构造函数，由已有矩阵类对象创建新的矩阵类对象；
unsignedinti;
this->nRow=m.GetRowsNum();
this->nColumn=m.GetColumnsNum();
//this->pElement=newdouble[nRow*nColumn];
this->alloc(nRow*nColumn);
for(i=0;i<nRow*nColumn;++i)
{
this->pElement[i]=m.pElement[i];
}
}
CMatrix::CMatrix(constchar*pFileName)
{//由矩阵存储文件名创建矩阵类对象，文件格式可参考附录，但不限于采用此格式；
if(pFileName!=NULL)
{
this->InitFromFile(pFileName);
}
}

voidCMatrix::InitFromFile(constchar*pFileName)
{//从文件中初始化
ifstreaminFile(pFileName);
if(!inFile)
{
cerr<<"unabletoopeninputfile:"<<pFileName<<"--bailingout!/n";
return;
}

inFile>>this->nRow>>this->nColumn;
this->alloc(nRow*nColumn);//分配空间
unsignedinti;
doublenum;
for(i=0;i<nRow*nColumn;++i)
{
inFile>>num;
this->pElement[i]=num;
}
inFile.close();
}

CMatrix::CMatrix(unsignedintr,unsignedintc):nRow(r),nColumn(c)
{//由矩阵的行数和列数创建矩阵类对象，并为矩阵元素分配存储空间，将矩阵初始化为单位矩阵；
if(nRow>0&&nColumn>0)
{
this->nRow=r;
this->nColumn=c;
//this->pElement=newdouble[nRow*nColumn];
this->alloc(nRow*nColumn);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]=0;
}
}
}

CMatrix&CMatrix::operator=(constCMatrix&m)
{//重载运算符=,用于矩阵之间相互赋值
if(this!=&m)
{
this->destroy();
this->nRow=m.GetRowsNum();
this->nColumn=m.GetColumnsNum();
this->alloc(nRow*nColumn);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]=m.pElement[i];
}
}
return*this;
}

double&CMatrix::operator()(unsignedintr,unsignedintc)
{//重载运算符()，用于提取指定行(r)列(c)的元素值
returnthis->pElement[r*nColumn+c];
}
CMatrix&CMatrix::operator+()const
{//重载一元运算符+,即取矩阵本身
returnconst_cast<CMatrix&>(*this);
}
CMatrixCMatrix::operator-()const
{//重载一元运算符-，即矩阵元素取相反数
CMatrixresult(nRow,nColumn);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
result.pElement[i]=-pElement[i];
}
returnresult;
}
CMatrixCMatrix::operator+(constCMatrix&m)const
{//重载二元运算符+,即两个矩阵求和
assert(this->nRow==m.GetRowsNum()&&this->nColumn==m.GetColumnsNum());
CMatrixresult(*this);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
result.pElement[i]+=m.pElement[i];
}
returnresult;
}
CMatrixCMatrix::operator-(constCMatrix&m)const
{//重载二元运算符-,即两个矩阵求差
assert(this->nRow==m.GetRowsNum()&&this->nColumn==m.GetColumnsNum());
CMatrixresult(*this);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
result.pElement[i]-=m.pElement[i];
}
returnresult;
}
CMatrixCMatrix::operator*(constCMatrix&m)const
{//重载二元运算符*,即两个矩阵求积
CMatrixresult(*this);
intct=0,cm=0,cw=0;//computew(i,j)foralliandj
for(unsignedinti=1;i<=nRow;++i)
{
//computerowiofresult
for(unsignedintj=1;j<=m.GetColumnsNum();++j)
{
//computefirsttermofw(i,j)
doublesum=pElement[ct]*m.pElement[cm];//addinremainingterms
for(unsignedintk=2;k<=nColumn;++k)
{
ct++;//nextterminrowiof*this
cm+=m.GetColumnsNum();//nextincolumnjofm
sum+=pElement[ct]*m.pElement[cm];
};
result.pElement[cw++]=sum;//savew(i,j)
//resettostartofrowandnextcolumn
ct-=nColumn-1;
cm=j;
}//resettostartofnextrowandfirstcolumn
ct+=nColumn;
cm=0;
}
returnresult;
}
CMatrixCMatrix::operator*(constdouble&x)const
{//重载二元运算符*,即矩阵与数相乘
CMatrixresult(*this);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
result.pElement[i]=this->pElement[i]*x;
}
returnresult;
}
CMatrixCMatrix::operator*(CVector&v)const
{//重载二元运算符*，即矩阵与向量相乘
assert(this->nColumn==v.GetDegree());
CMatrixresult(nRow,1);
intct=0,cm=0,cw=0;
for(unsignedinti=0;i<nRow;++i)
{
doublesum=0.0f;
for(unsignedintj=0;j<v.GetDegree();++j)
{
sum+=pElement[ct++]*v[j];
}
result(i,0)=sum;
}
returnresult;

}
CMatrixCMatrix::operator^(constint&t)const
{//重载二元运算符^,即矩阵求t次幂
CMatrixresult(*this);
CMatrixtmp(*this);

for(inti=1;i<t;++i)
{
result=result*(tmp);
}
returnresult;
}
voidCMatrix::operator-=(constCMatrix&m)
{//重载二元运算符-=，即自减运算
assert(nRow==m.GetRowsNum()&&nColumn==m.GetColumnsNum());
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]-=m.pElement[i];
}
}
voidCMatrix::operator+=(constCMatrix&m)
{//重载二元运算符-=，即自加运算
assert(nRow==m.GetRowsNum()&&nColumn==m.GetColumnsNum());
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]+=m.pElement[i];
}
}
voidCMatrix::operator*=(constdoublex)
{//重载二元运算符*=，即自乘运算（矩阵）
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]*=x;
}
}

CMatrixCMatrix::operator/(constdouble&x)const
{//重载二元运算符/,即矩阵与数相除
CMatrixresult(*this);
for(unsignedinti=0;i<nRow*nColumn;++i)
{
result.pElement[i]=this->pElement[i]/x;
}

returnresult;
}
CMatrixCMatrix::Tranpose()const
{//求矩阵的转置的函数
CMatrixresult(*this);
for(unsignedinti=0;i<nRow;++i)
{
for(unsignedintj=0;j<nColumn;++j)
{
result(i,j)=this->pElement[j*nColumn+i];
}
}
returnresult;
}
CMatrixCMatrix::Invert()const
{//求矩阵的逆的函数
CMatrixresult(nRow,nColumn*2);
unsignedinti,j,p,q;
for(i=0;i<nRow;++i)
{
for(j=0;j<nColumn;++j)
{
result(i,j)=this->pElement[i*nColumn+j];
}
for(j=nColumn;j<nColumn*2;++j)
{
if((j-i)==nColumn)
{
result(i,j)=1.0f;
}
else
{
result(i,j)=0.0f;
}
}
}

for(i=0;i<nRow;++i)
{
if(result(i,i)!=1.0f)
{
doublebs=result(i,i);
result(i,i)=1.0f;
for(j=i+1;j<nColumn*2;++j)
{
result(i,j)/=bs;
}
}
for(q=0;q<nRow;++q)
{
if(q!=i)
{
doublebs=result(q,i);
for(p=0;p<nColumn*2;++p)
{
result(q,p)-=bs*result(i,p);
}
}
else
{
continue;
}
}
}
returnresult;
}
voidCMatrix::Zeros()
{//矩阵归零化，将当前矩阵的所有元素归零
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]=0.0f;
}
}

voidCMatrix::Unit()
{//矩阵单位化，将当前矩阵转换为单位矩阵
for(unsignedinti=0;i<nRow;++i)
{
for(unsignedintj=0;j<nColumn;++j)
{
if(i==j)
{
(*this)(i,j)=1.0f;
}
else
{
(*this)(i,j)=0.0f;
}
}
}
}
intCMatrix::AddRow(double*pe,unsignedintnr)
{//在矩阵的nr行位置插入一行，数据存放地址为pe，返回行标
CMatrixtmp(*this);
unsignedintoldRow,oldColumn;
oldRow=this->nRow;
oldColumn=this->nColumn;
this->destroy();
this->nRow=oldRow+1;
this->nColumn=oldColumn;
this->alloc(nRow*nColumn);
unsignedinti,j,k;
for(i=0;i<nr*nColumn;++i)
{
this->pElement[i]=tmp.pElement[i];
}
k=i;
for(j=0;j<nColumn;++j)
{
this->pElement[i++]=pe[j];
}
for(;i<nRow*nColumn;++i)
{
this->pElement[i]=tmp.pElement[k++];
}
returnnr;
}
intCMatrix::AddColumn(double*pe,unsignedintnc)
{//在矩阵的nc列位置插入一列，数据存放地址为pe，返回列标
CMatrixtmp(*this);

unsignedintoldRow,oldColumn;
oldRow=this->nRow;
oldColumn=this->nColumn;
this->destroy();
this->nRow=oldRow;
this->nColumn=oldColumn+1;
this->alloc(nRow*nColumn);

unsignedinti,j,k;
for(i=0;i<nc;++i)
{
for(j=0;j<nRow;++j)
{
this->pElement[j*nColumn+i]=tmp.pElement[j*oldColumn+i];
}
}
k=i;
for(j=0;j<nRow;++j)
{
this->pElement[j*nColumn+i]=pe[j];
}
for(i=k+1;i<nColumn;++i)
{
for(j=0;j<nRow;++j)
{
this->pElement[j*nColumn+i]=tmp.pElement[j*oldColumn+k];
}
k++;
}
returnnc;
}
double*CMatrix::DeleteRow(unsignedintnr)
{//删除矩阵nr行，返回该行元素值（临时存储地址）
double*pTmp=newdouble[nColumn];
CMatrixtmp(*this);

unsignedintoldRow,oldColumn;
oldRow=this->nRow;
oldColumn=this->nColumn;
this->destroy();
this->nRow=oldRow-1;
this->nColumn=oldColumn;
this->alloc(nRow*nColumn);

unsignedinti,j,ct=0;
for(i=0;i<nr;++i)
{
for(j=0;j<nColumn;++j)
{
this->pElement[ct]=tmp(i,j);
ct++;
}
}
for(j=0;j<nColumn;++j)
{
pTmp[j]=tmp(nr,j);
}
for(i=nr+1;i<oldRow;++i)
{
for(j=0;j<nColumn;++j)
{
this->pElement[ct]=tmp(i,j);
ct++;
}
}
returnpTmp;

}
double*CMatrix::DeleteColumn(unsignedintnc)
{//删除矩阵nc列，返回该列元素值（临时存储地址）
double*pTmp=newdouble[nRow];
CMatrixtmp(*this);

unsignedintoldRow,oldColumn;
oldRow=this->nRow;
oldColumn=this->nColumn;
this->destroy();
this->nRow=oldRow;
this->nColumn=oldColumn-1;
this->alloc(nRow*nColumn);

unsignedinti,j,ct=0,cw=0;
for(i=0;i<nRow;++i)
{
for(j=0;j<oldColumn;++j)
{
if(j!=nc)
{
this->pElement[ct]=tmp(i,j);
ct++;
}
else
{
pTmp[cw]=tmp(i,j);
cw++;
}
}
}
returnpTmp;
}


unsignedintCMatrix::GetRowsNum()const
{//获取矩阵的行数
returnthis->nRow;
}
unsignedintCMatrix::GetColumnsNum()const
{//获取矩阵的列数
returnthis->nColumn;
}
voidCMatrix::alloc(unsignedintn/**//*=0*/)
{
this->pElement=newdouble[n];
}
voidCMatrix::destroy()
{
if(this->pElement!=NULL)
{
delete[]this->pElement;
this->pElement=NULL;
this->nColumn=0;
this->nRow=0;
}
}

CMatrix::~CMatrix(void)
{
this->destroy();
}
intCMatrix::MatOut(constchar*pFileName)
{//将矩阵以pFileName为文件名进行文件输出
if(pFileName!=NULL)
{
ofstreamoutFile(pFileName);
if(!outFile)
{
cerr<<"unabletoopenoutputfile:"<<pFileName<<"--bailingout!/n";
return-1;
}
outFile<<this->nRow<<""<<this->nColumn<<endl;
unsignedinti;
for(i=0;i<nRow*nColumn;++i)
{
if((i+1)%nColumn==0)
{
outFile<<this->pElement[i]<<endl;
}
else
outFile<<this->pElement[i]<<"";

}
outFile.close();
}
return0;
}

#include"matrix.h"

classCVector;
//线性方程组类
/**//*
**author:phinecos
**date:7/17/2008
*/
classCLinearEquation:publicCMatrix
{
public:
CLinearEquation(void);
CLinearEquation(CMatrix&coe,CVector&con);//通过系数矩阵和常数向量创建线性方程组
CLinearEquation(constchar*pFile);//通过数据文件创建线性方程组
CLinearEquation(unsignedintne,unsignedintnv);//通过方程个数和未知数创建线性方程组
~CLinearEquation(void);

public:
intAddVariable(double*pcoe,unsignedintpc);//为方程组增加一个变量，其系数矩阵增加一列，pcoe为增加系数地址，pc为增加系数列的序号，为最前面，默认追加在尾部，成功返回未知数个数，否则为-1
intAddEquation(double*pcoe,doublecon,unsignedintpr);//为方程组增加一个方程，其系数矩阵增加一行，pcoe为增加系数地址，con为增加方程的常数项，pr为增加系数行的序号，为最前面，默认追加在尾部，成功返回方程个数，否则为-1
intDeleteVariable(unsignedintpc);//删除方程组中第pc个未知数，并删除其系数列
intDeleteEquation(unsignedintpr);//删除方程组中第pr个方程
CVectorGaussian();//高斯消元法解线性方程组
intCLinearEquation::MatOut(constchar*pFileName);
public:
CVector*Constant;//方程组的常数向量，新增成员
//方程组系数矩阵，继承成员，不可见
//方程组中方程个数，继承成员，矩阵行数，不可见
//方程组中未知数个数，继承成员，矩阵列数，不可见
private:
voidInsertConstant(doublenum,intpos);
voidDeleteConstant(unsignedintpos);
voidInitFromFile(constchar*pFileName);
};

线性方程组类实现
<!--<br /><br />Code highlighting produced by Actipro CodeHighlighter (freeware)<br />http://www.CodeHighlighter.com/<br /><br />-->#include"stdafx.h"
#include"LinearEquation.h"
#include"Vector.h"
#include<cmath>
#include<iostream>
#include<fstream>
usingnamespacestd;

/**//*
**author:phinecos
**date:7/17/2008
*/
CLinearEquation::CLinearEquation(void)
{
}
CLinearEquation::CLinearEquation(constchar*pFile)
{//通过数据文件创建线性方程组
if(pFile!=NULL)
{
this->InitFromFile(pFile);
}

}
voidCLinearEquation::InitFromFile(constchar*pFileName)
{
ifstreaminFile(pFileName);
if(!inFile)
{
cerr<<"unabletoopeninputfile:"<<pFileName<<"--bailingout!/n";
return;
}
inFile>>this->nRow>>this->nColumn;
this->alloc(nRow*nColumn);//分配空间
unsignedinti;
doublenum;
for(i=0;i<nRow*nColumn;++i)
{
inFile>>num;
this->pElement[i]=num;
}
this->Constant=newCVector(nRow);
for(i=0;i<nRow;++i)
{
inFile>>num;
(*this->Constant)[i]=num;
}
inFile.close();
}
CLinearEquation::~CLinearEquation(void)
{
this->destroy();
if(this->Constant!=NULL)
{
deleteConstant;
this->Constant=NULL;
}
}
CLinearEquation::CLinearEquation(unsignedintne,unsignedintnv)
{//通过方程个数和未知数创建线性方程组
this->nRow=ne;
this->nColumn=nv;
this->pElement=newdouble[nRow*nColumn];
this->Constant=newCVector(nRow);
}

CLinearEquation::CLinearEquation(CMatrix&coe,CVector&con)
{//通过系数矩阵和常数向量创建线性方程组
this->nRow=coe.GetRowsNum();
this->nColumn=coe.GetColumnsNum();
this->pElement=newdouble[nRow*nColumn];
for(unsignedinti=0;i<nRow*nColumn;++i)
{
this->pElement[i]=coe.pElement[i];
}
this->Constant=newCVector(con);
}

intCLinearEquation::AddVariable(double*pcoe,unsignedintpc)
{//为方程组增加一个变量，其系数矩阵增加一列，pcoe为增加系数地址，pc为增加系数列的序号，为最前面，默认追加在尾部，成功返回未知数个数，否则为-1
if(pc>=0)
this->AddColumn(pcoe,pc);
else
this->AddColumn(pcoe,this->nColumn);//默认追加在尾部
returnthis->nColumn;
}

intCLinearEquation::AddEquation(double*pcoe,doublecon,unsignedintpr)
{//为方程组增加一个方程，其系数矩阵增加一行，pcoe为增加系数地址，con为增加方程的常数项，pr为增加系数行的序号，为最前面，默认追加在尾部，成功返回方程个数，否则为-1
if(pr>=0)
{
this->AddRow(pcoe,pr);
this->InsertConstant(con,pr);
}
else
{
this->InsertConstant(con,this->nRow);
this->AddRow(pcoe,this->nRow);
}
returnthis->nRow;
}
intCLinearEquation::DeleteVariable(unsignedintpc)
{//删除方程组中第pc个未知数，并删除其系数列
this->DeleteColumn(pc);
return0;
}
intCLinearEquation::DeleteEquation(unsignedintpr)
{//删除方程组中第pr个方程
this->DeleteRow(pr);
this->DeleteConstant(pr);
return0;
}

voidCLinearEquation::DeleteConstant(unsignedintpos)
{
CVectortmp(*this->Constant);
unsignedintoldDegree=this->Constant->GetDegree();
this->Constant->destory();
this->Constant=newCVector(oldDegree-1);
inti,k;
for(i=0;i<(int)pos;++i)
{
(*this->Constant)[i]=tmp[i];
}
for(k=pos+1;k<(int)oldDegree;++k,++i)
{
(*this->Constant)[i]=tmp[k];
}
}


voidCLinearEquation::InsertConstant(doublenum,intpos)
{
CVectortmp(*this->Constant);
unsignedintoldDegree=this->Constant->GetDegree();
this->Constant->destory();
this->Constant=newCVector(oldDegree+1);
inti,k;
for(i=0;i<pos;++i)
{
(*this->Constant)[i]=tmp[i];
}
k=i;
(*this->Constant)[i]=num;
++i;
for(;k<(int)oldDegree;++k,++i)
{
(*this->Constant)[i]=tmp[k];
}
}


CVectorCLinearEquation::Gaussian()
{//高斯消元法解线性方程组
int*js,l,k,i,j,is,p,q;
intn=(int)this->nRow;//方程组阶数
doubled,t;
js=newint[n];
l=1;
for(k=0;k<=n-2;k++)
{
d=0.0f;
for(i=k;i<=n-1;i++)
for(j=k;j<=n-1;j++)
{
t=fabs(this->pElement[i*n+j]);
if(t>d){d=t;js[k]=j;is=i;}
}
if(d+1.0==1.0)
l=0;
else
{
if(js[k]!=k)
{
for(i=0;i<=n-1;i++)
{
p=i*n+k;
q=i*n+js[k];
t=this->pElement[p];
this->pElement[p]=this->pElement[q];
this->pElement[q]=t;
}
}
if(is!=k)
{
for(j=k;j<=n-1;j++)
{
p=k*n+j;
q=is*n+j;
t=this->pElement[p];
this->pElement[p]=this->pElement[q];
this->pElement[q]=t;
}
t=(*this->Constant)[k];
(*this->Constant)[k]=(*this->Constant)[is];
(*this->Constant)[is]=t;
}
}
if(l==0)
{
delete[]js;
return(0);
}
d=this->pElement[k*n+k];
for(j=k+1;j<=n-1;j++)
{
p=k*n+j;
this->pElement[p]=this->pElement[p]/d;
}
(*this->Constant)[k]=(*this->Constant)[k]/d;
for(i=k+1;i<=n-1;i++)
{
for(j=k+1;j<=n-1;j++)
{
p=i*n+j;
this->pElement[p]=this->pElement[p]-this->pElement[i*n+k]*this->pElement[k*n+j];
}
(*this->Constant)[i]=(*this->Constant)[i]-this->pElement[i*n+k]*(*this->Constant)[k];
}
}
d=this->pElement[(n-1)*n+n-1];
if(fabs(d)+1.0==1.0)
{
delete[]js;
return(0);
}
(*this->Constant)[n-1]=(*this->Constant)[n-1]/d;
for(i=n-2;i>=0;i--)
{
t=0.0;
for(j=i+1;j<=n-1;j++)
{
t=t+this->pElement[i*n+j]*(*this->Constant)[j];
}
(*this->Constant)[i]=(*this->Constant)[i]-t;
}
js[n-1]=n-1;
for(k=n-1;k>=0;k--)
{
if(js[k]!=k)
{
t=(*this->Constant)[k];
(*this->Constant)[k]=(*this->Constant)[js[k]];
(*this->Constant)[js[k]]=t;
}
}
delete[]js;
CVectorresult(*this->Constant);
returnresult;
}

intCLinearEquation::MatOut(constchar*pFileName)
{//将矩阵以pFileName为文件名进行文件输出
if(pFileName!=NULL)
{
ofstreamoutFile(pFileName);
if(!outFile)
{
cerr<<"unabletoopenoutputfile:"<<pFileName<<"--bailingout!/n";
return-1;
}
unsignedinti;
for(i=0;i<this->nRow;++i)
{
if(i==nRow-1)
{
outFile<<(*this->Constant)[i]<<endl;
}
else
outFile<<(*this->Constant)[i]<<"";

}
outFile.close();
}
return0;
}

测试程序
<!--<br /><br />Code highlighting produced by Actipro CodeHighlighter (freeware)<br />http://www.CodeHighlighter.com/<br /><br />-->//Test.cpp:Definestheentrypointfortheconsoleapplication.
//
/**//*
**author:phinecos
**date:7/17/2008
*/
#include"stdafx.h"
#include"Vector.h"
#include"Matrix.h"
#include"LinearEquation.h"
#include<iostream>
#include<string>
usingnamespacestd;

intmain(intargc,char*argv[])
{
charch='n',chSave='y';
do
{
stringstrFileName,strSaveFile;
cout<<"请输入数据文件名称:";
cin>>strFileName;

CLinearEquationeq1(strFileName.c_str());
CVectorv21=eq1.Gaussian();
cout<<"高斯法解方程得到的解是:"<<endl;
v21.printElements();
cout<<"是否要保存结果(y/n)?:";
cin>>chSave;
if(chSave!='n'&&ch!='N')
{//保存结果
cout<<"请输入保存文件名称:";
cin>>strSaveFile;
eq1.MatOut(strSaveFile.c_str());
}
cout<<"是否退出程序(y/n)?:";
cin>>ch;
}while(ch!='y'&&ch!='Y');
system("pause");
return0;
}