线性方程组数值解法.docx
《线性方程组数值解法.docx》由会员分享,可在线阅读,更多相关《线性方程组数值解法.docx(19页珍藏版)》请在冰豆网上搜索。
线性方程组数值解法
计算方法实验
题目:
班级:
学号:
姓名:
1实验目的
1.通过编程加深对列主元高斯消去法、LU三角分解法和雅克比迭代法等求解多项式方程方法的理解
2.观察上述三种方法的计算稳定性和求解精度并比较各种方法利弊
2实验步骤
2.1环境配置:
VS2013,C++控制台程序
2.2添加头文件
#include"stdio.h"
#include"stdlib.h"
#include"stdafx.h"
#include
2.3主要模块
程序一共分成三层,最底层是数据结构部分,负责存储数据,第二层是交互部分,即多项式方程部分,负责输入输出获得数据,最上层是核心的算法部分,负责处理已获得的数据。
具体功能如下:
●数据结构部分
数据结构部分是整个程序的最底层,负责存储部分。
因数组作为数据元素插入和删除操作较少,而顺序表空间利用率大且查看方便,故此程序选用二维顺序表保存系数。
数据结构文件中写的是有关其的所有基本操作以供其他文件调用。
●多项式方程部分
多项式方程部分是程序的第二层,内容是有关方程组的所有函数、构建方程、输出方程等等,同时在此文件中获得方程系数并储存,同时此文件还负责显示菜单部分。
●算法部分
此文件负责核心算法,处于整个程序最上层部分,负责列主元高斯消去法、LU三角分解法和雅克比迭代法的具体实现过程。
通过调用方程文件的函数,将获得的数据进行处理运算,可以得到结果返回给方程主函数和输出的第二层。
总结:
主函数负责获取方程系数并显示,算法和方程作为后台程序,顺序表作为存储手段。
3代码
3.1主程序部分
//Solutionoflinearquations.cpp:
定义控制台应用程序的入口点。
//
#include"stdio.h"
#include"stdlib.h"
#include"stdafx.h"
#include"squencelist.h"
#include"equation.h"
#include"algorithm.h"
#include
int_tmain(intargc,_TCHAR*argv[])
{
while(Exflag)
{
GetEquation();
ShowMenu();
}
return0;
}
3.2多项式方程部分
●方程部分头文件
#ifndef_EQUATION_H
#define_EQUATION_H
#include"stdio.h"
#include"stdlib.h"
#include"squencelist.h"
externintXnumbers;
externintFnumber;
externintExflag;
externdatacoa*A;
voidGetEquation(void);
voidShowMenu(void);
voidprintfunction(datacoa*A);
voidprintresx(datacoa*A);
voidTip(void);
#endif
●方程部分CPP文件
#include"stdafx.h"
#include"equation.h"
#include"math.h"
#include"algorithm.h"
#include"stdio.h"
#include"stdlib.h"
#include
#include
usingnamespacestd;
//全局变量
intXnumbers=0;
intFnumber=0;
intExflag=1;
datacoa*A;
////////////////////////多项式函数系数/////////////////////////
voidGetEquation(void)
{
inti,j,flag=1;
floatx;
A=InitStruct();
while(flag)
{
cout<<"方程未知量和解总个数:
"<cin>>Xnumbers;
cout<<"方程个数:
"<cin>>Fnumber;
cout<<"输入方程系数,输入00结束(如输入215或者21534100"<cout<<"34100):
"<cin>>x;
while(x!
=00)
{
for(i=1;i<=Fnumber;i++)
{
for(j=1;j<=Xnumbers;j++)
{
if(!
Insert(A,x,i,j))exit(0);
cin>>x;
}
}
}
j=1;
printfunction(A);
if(Xnumbers==Fnumber+1)flag=0;
else
{
cout<<"方程可能无解"<A->m=0;
A->n=0;
}
}
}
//////////////////////////显示交互/////////////////////////////
voidShowMenu(void)
{
intx;
cout<<"选择求解方程根的方法:
"<cout<<"1.列主元高斯消去法"<cout<<"2.三角分解法"<cout<<"3.迭代法"<cin>>x;
switch(x)
{
case1:
ColumnGaussmethod(A);
Tip();
break;
case2:
LUmethod(A);
Tip();
break;
case3:
ISODATAmethod(A);
Tip();
break;
default:
break;
}
}
////////////////////////打印输出函数///////////////////////////
voidprintfunction(datacoa*A)
{
inti,j;
cout<<"矩阵="<for(i=1;i<=A->m;i++)
{
for(j=1;j<=A->n;j++)
{
cout<data[i][j];
}
cout<}
}
////////////////////////打印结果///////////////////////////
voidprintresx(datacoa*A)
{
inti;
for(i=1;i<=A->m;i++)
{
cout<<"X"<data[i][Xnumbers+1];
cout<}
}
////////////////////////返回提示///////////////////////////
voidTip(void)
{
intflag;
cout<<"输入000退出,其余返回:
"<cin>>flag;
if(flag==000)Exflag=0;
}
3.3核心算法部分
●算法部分头文件
#ifndef_ALGORITHM_H
#define_ALGORITHM_H
#include"stdio.h"
#include"stdlib.h"
#include"stdafx.h"
intJudge(doublex1,doublex0,doublee);
voidColumnGaussmethod(datacoa*A);
voidLUmethod(datacoa*A);
voidISODATAmethod(datacoa*A);
#endif
●算法部分CPP文件
#include"algorithm.h"
#include"stdafx.h"
#include"squencelist.h"
#include"equation.h"
#definemax100
//////////////////////////误差判别////////////////////////////
inlineintJudge(doublex1,doublex0,doublee)
{
if(e==0)
return0;
if(x1==0)
{
if(fabs(x0)elsereturn0;
}
if(x0==0)
{
if(fabs(x1)elsereturn0;
}
if(fabs(x1-x0)return1;
elsereturn0;
}
//////////////////////////列主元高斯消元法////////////////////////////
voidColumnGaussmethod(datacoa*A)
{
cout<<"///////////////////列主元高斯消元法///////////////////"<inti,j,i2,flagc,k,j2;
floattemp,res;
for(i=1;i{
flagc=i;
for(i2=i+1;i2<=Fnumber;i2++)
if((fabs(A->data[i2][i]))>(fabs(A->data[flagc][i])))
flagc=i2;
if(flagc!
=i)
for(k=i;k<=Xnumbers;k++)
{
temp=A->data[i][k];
A->data[i][k]=A->data[flagc][k];
A->data[flagc][k]=temp;
}
for(i2=i+1;i2<=Fnumber;i2++)
{
temp=A->data[i2][i]/A->data[i][i];
for(j2=i;j2<=Xnumbers;j2++)
A->data[i2][j2]=A->data[i2][j2]-temp*A->data[i][j2];
}
}
for(i=Fnumber;i>=1;i--)
{
for(j=Fnumber;j>=i+1;j--)
A->data[i][Xnumbers]=A->data[i][Xnumbers]-A->data[i][j]*A->data[j][Xnumbers+1];
res=A->data[i][Xnumbers]/A->data[i][i];
Insert(A,res,i,Xnumbers+1);
}
printresx(A);
}
//////////////////////////直接三角分解法////////////////////////////
voidLUmethod(datacoa*A)
{
datacoa*L,*U;
inti,j,k,q;
floattemp=0.0;
L=InitStruct();
U=InitStruct();
for(i=1;i<=Fnumber;i++)
{
for(j=1;j<=Xnumbers;j++)
{
Insert(L,0,i,j);
Insert(U,0,i,j);
}
}
for(j=1;j<=Xnumbers;j++)
Insert(U,A->data[1][j],1,j);
for(i=2;i<=Fnumber;i++)
Insert(L,A->data[i][1]/A->data[1][1],i,1);
for(k=2;k<=Fnumber;k++)
{
for(j=k;j<=Xnumbers;j++)
{
for(q=1;qtemp=temp+L->data[k][q]*U->data[q][j];
Insert(U,A->data[k][j]-temp,k,j);
temp=0;
}
for(i=k+1;i<=Fnumber;i++)
{
for(q=1;qtemp=temp+L->data[i][q]*U->data[q][k];
temp=A->data[i][k]-temp;
Insert(L,temp/U->data[k][k],i,k);
temp=0;
}
}
Insert(U,U->data[Fnumber][Xnumbers]/U->data[Fnumber][Xnumbers-1],Fnumber,Xnumbers+1);
for(k=Fnumber-1;k>=1;k--)
{
for(q=k+1;qtemp=temp+U->data[k][q]*U->data[q][Xnumbers+1];
temp=U->data[k][Xnumbers]-temp;
Insert(U,temp/U->data[k][k],k,Xnumbers+1);
temp=0;
}
printresx(U);
}
//////////////////////////迭代法////////////////////////////
voidISODATAmethod(datacoa*A)
{
inti=1,j,k=0;
floatx0=0,x1=1,temp=0;
for(i=1;i<=Fnumber;i++)
Insert(A,0,i,Xnumbers+1);
while
(1)
{
for(i=1;i<=Fnumber;i++)
{
for(j=1;j<=Fnumber;j++)
{
if(j==i)continue;
temp=temp+A->data[i][j]*A->data[j][Xnumbers+1];
}
temp=A->data[i][Xnumbers]-temp;
temp=temp/A->data[i][i];
Insert(A,A->data[i][Xnumbers+1],i,Xnumbers+2);
Insert(A,temp,i,Xnumbers+1);
temp=0;
}
k++;
for(i=1;i<=Fnumber;i++)
temp=temp+fabs(A->data[i][Xnumbers+1]-A->data[i][Xnumbers+2]);
if((temp<0.000001)||k>max)
break;
temp=0;
}
DeleteLie(A,Xnumbers+2);
printfunction(A);
cout<printresx(A);
cout<}
3.4数据结构部分
●数据结构头文件
#ifndef_SQUENCELIST_H
#define_SQUENCELIST_H
#include"stdio.h"
#include"stdlib.h"
#include"stdafx.h"
#include
usingnamespacestd;
#definemaxsize1024
/**
*sequenlist
*/
typedeffloatdatatype;
typedefstruct
{
datatypedata[maxsize][maxsize];
intm,n;
}datacoa;
datacoa*InitStruct();
intLength(datacoa*);
intInsert(datacoa*,datatype,int,int);
voidDeleteLie(datacoa*L,intj);
voidDeleteLine(datacoa*L,inti);
#endif
●数据结构CPP文件
#include"stdafx.h"
#include"squencelist.h"
///////////////////////////////////数据结构部分///////////////////////////////////////////
///////////////////////////////////sequenlist///////////////////////////////////////////
datacoa*InitStruct()
{
datacoa*L=(datacoa*)malloc(sizeof(datacoa));
L->m=0;
L->n=0;
returnL;
//datacoa*L=newdatacoa;
}
intLength(datacoa*L)
{
returnL->m*L->n;
}
intInsert(datacoa*L,datatypex,inti,intj)
{
intk;
if((L->m>=maxsize-1)||(L->n>=maxsize-1))
{
cout<<"表已满"<return0;
}
for(k=L->n;k>=j;k--)
L->data[i][k+1]=L->data[i][k];
L->data[i][j]=x;
if(i>L->m)L->m++;
if(j>L->n)L->n++;
return1;
}
voidDeleteLie(datacoa*L,intj)
{
intk,i;
if((j<1)||(j>L->n))
{
cout<<"非法删除位置"<}
for(i=1;i<=L->m;i++)
{
for(k=j;k<=L->n;k++)
L->data[i][j]=L->data[i][j+1];
}
L->n--;
}
voidDeleteLine(datacoa*L,inti)
{
intk,j;
if((i<1)||(i>L->m))
{
cout<<"非法删除位置"<}
for(j=1;j<=L->n;j++)
{
for(k=i;k<=L->m;k++)
L->data[i][j]=L->data[i+1][j];
}
L->m--;
}
4运行结果
4.1列主元高斯消去法运行结果
4.2LU三角分解法运行结果
4.3雅克比迭代法运行结果
边界情况调试
●超定方程等可能无解的情况
●迭代法迭代次数超出100次的情况
5总结
这次的程序设计题看似简单实则临界代码区很多,调试时很容易出错。
输入输出
使用二维顺序表存储数据对于此题更加方便,也解决了C++/C语言数组下标从0开始的不便,将下标转换成习惯下标。
在迭代法的时候。
将两次的x值分别存储在二维表后面,不需要另开空间,节省变量。
列主元高斯消元法
1.循环的顺序很重要
2.边界条件单独处理
雅克比迭代法
1.设置最大循环次数不能进入死循环。
6参考资料
《计算方法与实习》孙志忠等著第五版---东南大学出版社
《软件技术基础》周大为等著第一版---西安电子科技大学出版社