线性方程组AXB的数值计算方法实验Word格式文档下载.docx
- 文档编号:21820568
- 上传时间:2023-02-01
- 格式:DOCX
- 页数:40
- 大小:397.98KB
线性方程组AXB的数值计算方法实验Word格式文档下载.docx
《线性方程组AXB的数值计算方法实验Word格式文档下载.docx》由会员分享,可在线阅读,更多相关《线性方程组AXB的数值计算方法实验Word格式文档下载.docx(40页珍藏版)》请在冰豆网上搜索。
由于新计算出的分量比旧分量准确些,因此设想一旦新分量
求出,马上就用新分量
代替雅可比迭代法中
来求
这就是高斯-赛德尔(Gauss-Seidel)迭代法。
把矩阵A分解成
(6)其中
分别为
的主对角元除外的下三角和上三角部分,于是,方程组
(1)便可以写成
即
其中
(7)
以
为迭代矩阵构成的迭代法(公式)
(8)
称为高斯—塞德尔迭代法(公式),用变量表示的形式为
(9)
4.稀疏矩阵
矩阵中非零元素的个数远远小于矩阵元素的总数,并且非零元素的分布没有规律,则称该矩阵为稀疏矩阵(sparsematrix);
与之相区别的是,如果非零元素的分布存在规律(如上三角矩阵、下三角矩阵、对角矩阵),则称该矩阵为特殊矩阵。
常见于进行大量数据计算。
3、实验内容
1.P1081
许多科学应用包含的矩阵带有很多的零。
在实际情况中很重要的三角形线性方程组有如下形式:
d1x1+c1x2=b1
a1x1+d2x2+c2x3=b2
a2x2+d3x3+c3x4=b3
·
aN-2xN-2+dN-1xN-1+cN-1xN=bN-1
aN-1xN-1+dNxN=bN
构造一个程序求解三角形线性方程组。
可假定不需要行变换,而且可用第k行消去第k+1行的xk。
2.P1201
求解线性方程组AX=B,其中
A=
B=
使用三角分解法求解X。
3.P1202
求解线性方程组AX=B,其中A=[aij]N×
N,,aij=ij-1;
而且B=[bij]N×
1,b11=N,当i≥时,bij=(iN-1)/(i-1)。
对N=3,7,11的情况分别求解。
精确解为X=[11…11]’。
对得到的结果与精确解的差异进行解释。
4.P1203
通过重复求解N各线性方程组
ACJ=EJ,其中J=1,2,…,N
来得到A-1,
则
A[C1C2…CN]=[E1E2…EN]
而且
A-1=[C1C2…CN]
保证对LU分解只计算一次!
5.P1293
设有如下三角线性方程组,而且系数矩阵具有严格对角优势:
d1x1+c1x2=b1
(i)根据方程组
(1),式
(2)和式(3),设计一个算法来求解上述方程组。
算法必须有效地利用系数矩阵的稀疏性。
a11x1+a12x2+·
·
a1jxj+·
+a1NxN=b1
a21x1+a22x2+·
a2jxj+·
+a2NxN=b2
方程组
(1)
aj1x1+aj2x2+·
ajjxj+·
+ajNxN=bj
aN1x1+aN2x2+·
aNjxj+·
+aNNxN=bN
=
,j=1,2,·
,N·
式
(2)
式(3)
(ii)根据(i)中设计的算法构建一个C++程序,并求解下列上三角线性方程组
(a)4m1+m2=3(b)4m1+m2=1
m1+4m2+m3=3m1+4m2+m3=2
m2+4m3+m4=3m2+4m3+m4=1
m3+4m4+m5=3m3+4m4+m5=2
m48+4m49+m50=3m48+4m49+m50=1
m49+m50=3m49+m50=2
6.P1294
利用高斯—赛德尔迭代法求解下列带状方程。
12x1-2x2+x3=5
-2x1+12x2-2x3+x4=5
x1-2x2+12x3-2x4+x5=5
x2-2x3+12x4-2x5+x6=5
x46-2x47+12x48-2x49+x50=5
x47-2x48+12x49-2x50=5
x48-2x49+12x50=5
4、实验结果及分析
实验描述:
本次实验使用系数矩阵的第k行消去第k+1行的xk,消除方法为第k行减去第k-1行乘上系数ak-1/bk-1,待消至第N行时,求解出xN,并依次会带求出各xN-1至x1,为了检验结果的正确性使用上面的方程组组
(1)及方程组
(2)进行验证。
方程组
(1)的结果为
,方程组
(2)的结果为
算法流程图:
Y
实验结果:
结果截图
(1)
图1
输入矩阵
,输出结果为X=
,与预期结果一致
(2)
图2
实验结论:
通过对系数矩阵的增广矩阵进行高斯消元和回带容易得到线性方程组的解,同时,利用这种方法可以求得矩阵的逆。
本次实验的解法为使用LU矩阵求解X,该解法的内容为将系数矩阵A分解为一上三角矩阵及一下三角矩阵,且有A=LU,之后由LY=B,UX=Y分别求解出Y,X。
图3
输入矩阵A=
,B=
将X代入AX后结果与矩阵B一致,运行结果正确无误,该程序正确,且有X=
(1)本次实验仍使用三角矩阵求解矩阵X的值,求解方法与实验二大致相同;
(2)矩阵A编写函数buildA完成,buildA的输入为矩阵A的阶数N,输出结果为矩阵A的地址,对于矩阵A有A=[aij]N×
(3)矩阵B编写函数buildB完成,buildB的输入为矩阵B的阶数N,输出结果为矩阵B的地址,对于矩阵B有B=[bij]N×
(1)N=3时
图4
(2)N=7时
图5
(3)N=11时
图6
实验结论:
由图4,图5,图6与结果对比可知,程序运行结果与预期结果相一致,程序正确无误。
(1)本次实验的目的为求解A的逆矩阵A-1,求解方法为利用A[C1C2…CN]=[E1E2…EN],A-1=[C1C2…CN]求解,故可将之分解为对于ACk=Ek,k=1,2,3,·
,N中对于Ck的求解,之后另A-1=[C1C2…CN];
(2)由于A-1=[C1C2…CN],A[C1C2…CN]=[E1E2…EN],固有[E1E2…EN]=I,故Ek=[a1j],j=1,2,·
,N,其中a1k=1,other=0;
(3)对于ACk=Ek,k=1,2,3,·
,N的求解方法使用LU三角分解法求解,用此方法求解出各个Ek对应的Ck,最后以此构成A-1;
(4)求出LU后,应判断矩阵对角线上是否存在为0的元素,若存在,则A不存在逆矩阵;
若不存在,则可求解逆矩阵A-1;
(5)上述方法中的LU分解只需要进行一次;
(6)对于程序的正确性使用矩阵
及矩阵
进行验证,其中,矩阵
不存在逆矩阵,矩阵
的逆矩阵为
实验结果:
(1)输入为矩阵
图7
(2)输入矩阵为
图8
如果系数矩阵能分解为LU的形式,其中L为下三角矩阵,U为上三角矩阵,通过对系数矩阵的分解再求解可应用简单的迭代进行求解x。
实验描述:
(1)本实验的目的为使用高斯—赛德尔迭代求解带状方程组;
(2)由于实验中的带状方程有极强的稀疏性和相似性,故编程时应考虑该矩阵的以上特点以减少运算量及运行时占用的内存;
(3)为了减少程序的运行次数,故选择式(3)作为运行的程序;
(4)应无明确的结束标志,故选择delta=
<
1e-5作为结束迭代的标志,此时再进行一轮迭代后输出结果。
a)4m1+m2=3
m1+4m2+m3=3
m2+4m3+m4=3
m3+4m4+m5=3
m48+4m49+m50=3
m49+m50=3
图9
(b)4m1+m2=1
m1+4m2+m3=2
m2+4m3+m4=1
m3+4m4+m5=2
m48+4m49+m50=1
m49+m50=2
图10
求解带状线性方程组的解可使用高斯-赛德尔迭代法。
表1.通过高斯赛德尔迭代得到的x的解
1~10
0.46379552381655
0.5372846051999656
0.5090229246013306
0.4982216344361741
0.4989418602397619
0.4999853512481308
0.5000887238901357
0.500015318846052
0.4999947932669753
0.4999978569134675
11~20
0.5000001084251992
0.5000002015766873
0.500000022610945
0.4999999862385722
0.499999995873979
0.500000*********8
0.5000000004390388
0.5000000000239392
0.499999999966016
0.4999999999925573
21~30
0.5000000000017429
0.5000000000009169
0.5
0.4999999999999205
0.4999999999999212
0.5
0.5000000000009176
0.5000000000017445
31~40
0.499999999992555
0.499999999966017
0.500000000023939
0.5000000004390391
0.500000*********6
0.4999999862385723
0.5000000226109451
0.500000108425199
41~50
0.5372846051999657
0.4637955238165501
图11
X的值如上图所示。
附件(代码):
#include<
iostream>
stdlib.h>
usingnamespacestd;
//此函数用于计算矩阵的解
float*uptrbk(float**A,intN)
{
floatc;
float*x=newfloat[N-1];
//生成一维数组x[N]
intn;
for(n=1;
n<
N;
n++)//进行高斯消元法计算x[k],k=1,2,3·
,N
{
c=A[n][n-1]/A[n-1][n-1];
A[n][n-1]=0;
A[n][n]=A[n][n]-c*A[n-1][n];
A[n][N]=A[n][N]-c*A[n-1][N];
}
x[N-1]=A[N-1][N]/A[N-1][N-1];
for(n=N-2;
n>
=0;
n--)
A[n][N]=A[n][N]-x[n+1]*A[n][n+1];
x[n]=A[n][N]/A[n][n];
returnx;
//带回计算结果数组x[N]的地址x
}
//实验的main函数
intmain()
intN;
cout<
"
请输入矩阵的阶数:
;
//输入矩阵的阶数用于生成动态矩阵
cin>
>
inti,k;
float*uptrbk(float**,int);
float*x;
float**A=newfloat*[N-1];
//生成动态增广矩阵
for(i=0;
i<
i++)
A[i]=newfloat[N];
请输入增广矩阵的值"
endl;
//输入增广矩阵的值
for(k=0;
k<
=N;
k++)
cin>
A[i][k];
x=uptrbk(A,N);
//计算矩阵的解列向量X
x的值为:
i++)//输出矩阵的解的列向量X
cout<
x["
i+1<
]="
x[i]<
system("
pause"
);
return0;
intN,i,k;
//输入矩阵的阶数,用于生成动态矩阵
N=N-1;
float*B=newfloat[N];
float*lufact(float**,float*,int);
float**A=newfloat*[N];
//生成动态系数矩阵A
请输入矩阵A的值:
//输入系数矩?
阵的值
请输入矩阵B的值:
//生成动态矩阵B
cin>
B[i];
//输入矩阵B的值
x=lufact(A,B,N);
x的值为a:
i++)//输出AX=B解的列向量X
//使用LU法求解X的函数,A为系数矩阵,B为计算结果,N为矩阵阶数
float*lufact(float**A,float*B,intN)
inti,j,k;
float**U;
//生成二维数组U
float*x=newfloat[N];
//生产保存结果的列矩阵X
float*y=newfloat[N];
//生产用于保存中间值的列矩阵Y
U=A;
i++)//将A转换为LU
for(j=i+1;
j<
j++)
{
c=U[j][i]/U[i][i];
for(k=i;
U[j][k]=U[j][k]-c*U[i][k];
U[j][i]=c;
}
i++)//计算中间矩阵Y的值
for(j=0;
i;
B[i]=B[i]-y[j]*U[i][j];
y[i]=B[i];
for(i=N;
i>
i--)//计算解X的值
for(j=N;
j>
j--)
y[i]=y[i]-x[j]*U[i][j];
x[i]=y[i]/U[i][i];
cmath>
intN,i;
double*x;
double**A;
double*B;
double**buildA(int);
double*buildB(int);
double*lufact(double**,double*,int);
A=buildA(N);
B=buildB(N);
i++)//输出AX=B解的列矩阵X
double*lufact(double**A,double*B,intN)
doublec;
double**U;
//生成二维数组U,用于储存L&
U矩阵
double*x=newdouble[N];
//生成保存结果的列矩阵X
double*y=newdouble[N];
//生成用于保存中间值的列矩阵Y
i++)//将A转换为LU矩阵
//用于产生阶数为N的矩阵A的函数
double**buildA(intN)
inti,j;
double**A=newdouble*[N];
//生成二维动态数组A[N-1][N-1]
A[i]=newdouble[N];
i++)//产生矩阵A的元素并存储到A[N-1][N-1]中
A[i][j]=pow(double(i+1),j);
returnA;
//用于产生1*N的列矩阵B
double*buildB(intN)
inti;
double*B=
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 线性方程组 AXB 数值 计算方法 实验