线性方程组的数值解法实验报告
- 格式:doc
- 大小:326.00 KB
- 文档页数:9
迭代法解线性方程组数值分析实验报告一、实验目的本次实验旨在深入研究和掌握迭代法求解线性方程组的基本原理和方法,并通过数值实验分析其性能和特点。
具体目标包括:1、理解迭代法的基本思想和迭代公式的推导过程。
2、掌握雅克比(Jacobi)迭代法、高斯赛德尔(GaussSeidel)迭代法和超松弛(SOR)迭代法的算法实现。
3、通过实验比较不同迭代法在求解不同类型线性方程组时的收敛速度和精度。
4、分析迭代法的收敛性条件和影响收敛速度的因素。
二、实验原理1、线性方程组的一般形式对于线性方程组$Ax = b$,其中$A$ 是$n×n$ 的系数矩阵,$x$ 是$n$ 维未知向量,$b$ 是$n$ 维常向量。
2、迭代法的基本思想迭代法是从一个初始向量$x^{(0)}$出发,按照某种迭代公式逐步生成近似解序列$\{x^{(k)}\}$,当迭代次数$k$ 足够大时,$x^{(k)}$逼近方程组的精确解。
3、雅克比迭代法将系数矩阵$A$ 分解为$A = D L U$,其中$D$ 是对角矩阵,$L$ 和$U$ 分别是下三角矩阵和上三角矩阵。
雅克比迭代公式为:$x^{(k+1)}= D^{-1}(b +(L + U)x^{(k)})$。
4、高斯赛德尔迭代法在雅克比迭代法的基础上,每次计算新的分量时立即使用刚得到的最新值,迭代公式为:$x_i^{(k+1)}=(b_i \sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)}\sum_{j=i+1}^{n}a_{ij}x_j^{(k)})/a_{ii}$。
5、超松弛迭代法在高斯赛德尔迭代法的基础上引入松弛因子$\omega$,迭代公式为:$x_i^{(k+1)}= x_i^{(k)}+\omega((b_i \sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)}\sum_{j=i}^{n}a_{ij}x_j^{(k)})/ a_{ii} x_i^{(k)})$。
线性方程组求解的数值实验一、课题名称:双对角线线性方程组的数值实验二、班级和姓名:09101900 张争(学号 0910190161) 三、实验问题摘要:考虑一种特殊的对角线元素不为零的双对角线线性方程组(以n=7为例)⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7665544332211d a d a d a d a d a da d ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7654321x x x x x x x =⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡7654321b b b b b b b 写出解一般的n (奇数)阶方程组的程序(不要用消元过程,因为不用它可以十分方便的解出这个方程组)。
四、过程中的问题及处理方法:解方程组时应当先解出x1和xn,r 然后根据公式求出其他的解,由于矩阵的特殊性,应该从头尾开始解,最后解出21+n x 。
编程序时使用循环语句要分情况讨论。
五、计算公式: 1x =11d b kk k k k d x a b x 11---=(k<n/2)n x =nnd b kk k k k d x a b x 1+-=(k>n/2)21212123212121+--++++--=n n n n n n n d x a x a b x六、结构程序设计输入矩阵阶数→输入各系数→解出前(n-1)/2个解→解出后(n-1)/2个解→解出第(n+1)/2个解七、程序运行结果:这里令所有的系数等于1,实验结果如下:八、主要程序段介绍:解方程函数:由于矩阵具有好的对称性,解的形式也一样,用循环语句。
但是第一个和最后一个解以及中间一个解的形式不同,需要单独解出,程序如下:void jiefangcheng(){int i;x[1]=b[1]/d[1];for(i=2;i<(i+1)/2;i++)x[i]=(b[i]-a[i-1]*x[i-1])/d[i];x[n]=b[n]/d[n];for(i=n-1;i>(i+1)/2;i--)x[i]=(b[i]-a[i]*x[i+1])/d[i];x[(n+1)/2]=(a[(n-1)/2]*x[(n-1)/2]-a[(n+1)/2]*x[(n+3)/2])/d[(n+1)/2];}九、源程序如下://*****************张争******学号0910190161**********#include<iostream.h>#define N 10000int n;int a[N],b[N],d[N],x[N];void input(){int i;cout<<"请输入方程组的阶数:";cin>>n;cout<<"请输入方程组的系数:"<<endl;for(i=1;i<=n;i++){cout<<"d"<<i<<"=";cin>>d[i];cout<<"b"<<i<<"=";cin>>b[i];}for(i=1;i<n;i++){cout<<"a"<<i<<"=";cin>>a[i];}}void jiefangcheng(){int i;x[1]=b[1]/d[1];for(i=2;i<(i+1)/2;i++)x[i]=(b[i]-a[i-1]*x[i-1])/d[i];x[n]=b[n]/d[n];for(i=n-1;i>(i+1)/2;i--)x[i]=(b[i]-a[i]*x[i+1])/d[i];x[(n+1)/2]=(a[(n-1)/2]*x[(n-1)/2]-a[(n+1)/2]*x[(n+3)/2])/d[(n+1)/2];}void output(){ int i;cout<<"方程组的解为:"<<endl;for(i=1;i<=n;i++)cout<<"x"<<i<<"="<<x[i]<<endl;}void main(){input();jiefangcheng();output();}。
实验报告——线性方程组的数值解法姓名:班级:学号:日期:一 实践目的1. 熟悉求解线性方程组的有关理论和方法。
2. 会编列主元消去法,全主元消去法,雅克比迭代法及高斯-赛德尔迭代法的程序。
3.通过实际计算,进一步了解各种方法的优缺点,选择合适的数值方法。
4. 进一步应用数学知识,扩展数学思维,提高编程能力。
二 问题定义及题目分析1.求解线性方程是实际中常遇到的问题。
而一般求解方法有两种,一个是直接求解法,一个是迭代法。
2.直接法是就是通过有限步四则运算求的方程准确解的方法。
但实际计算中必然存在舍入误差,因此这种方法只能得到近似解。
3.迭代法是先给一个解的初始近似值,然后按一定的法则求出更准确的解,即是用某种极限过程逐步逼近准确解的方法。
4.这次实践,将采用直接法:列主元高斯消去法,全主元高斯消去法;迭代法:雅克比迭代法和高斯-赛德尔迭代法。
三 详细设计 设有线性方程组11112211n n a x a x a x b +++=21122222n n a x a x a x b +++= 1122n n nn n n a x a x a x b +++=令A= 111212122212n n n n nn a a a aa a a a a ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦ x=12n x x x ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦ b= 12n b b b ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦一. 上三角方程组的求解很简单。
所以若能将方程组化为上三角方程组,那就很容易得到方程的解,高斯消去法则是一种简洁高效的方法,可以将方程组化为上三角方程组,但是这种方法必须满足每一步时0kk a =才能求解,还有当kk a 绝对值很小时,作分母会引起较大的舍入误差。
因此在消元过程中应尽量选择绝对值较大的系数作为主元素。
1.列主元消去法很好的解决这一问题。
将1x 的系数1(1)k a k n ≤≤中绝对值最大者作为主元素,交换第一行和此元素所在的行,然后主元消去非主元项1x 的系数,。
1大学数学实验 实验报告 | 2014/4/5一、 实验目的1、学习用Matlab 软件数值求解线性代数方程组,对迭代法的收敛性和解的稳定性作初步分析;2、通过实例学习用线性代数方程组解决简化问题。
二、 实验内容项目一:种群的繁殖与稳定收获:种群的数量因繁殖而增加,因自然死亡而减少,对于人工饲养的种群(比如家畜)而言,为了保证稳定的收获,各个年龄的种群数量应维持不变。
种群因雌性个体的繁殖而改变,为方便起见以下种群数量均指其中的雌性。
种群年龄记作k=1,2,…,n ,当年年龄k 的种群数量记作x k ,繁殖率记作b k (每个雌性个体1年的繁殖的数量),自然存活率记作s k (s k =1−d k ,d k 为1年的死亡率),收获量记作ℎk ,则来年年龄k 的种群数量x ̌k 应该为x ̌k =∑b k n k=1x k , x ̌k+1=s k x k −ℎk , (k=1,2,…,n -1)。
要求各个年龄的种群数量每年维持不变就是要求使得x ̌k =x k , (k=1,2,…,n -1).(1) 如果b k , s k 已知,给定收获量ℎk ,建立求各个年龄的稳定种群数量x k 的模型(用矩阵、向量表示).(2) 设n =5,b 1=b 2=b 5=0,b 3=5,b 4=3,s 1=s 4=0.4,s 2=s 3=0.6,如要求ℎ1~ℎ5为500,400,200,100,100,求x 1~x 5.(3) 要使ℎ1~ℎ5均为500,如何达到?问题分析:该问题属于简单的种群数量增长模型,在一定的条件(存活率,繁殖率等)下为使各年龄阶段的种群数量保持不变,各个年龄段的种群数量将会满足一定的要求,只要找到种群数量与各个参量之间的关系,建立起种群数量恒定的方程就可以求解出各年龄阶段的种群数量。
模型建立:根据题目中的信息,令x ̌k =x k ,得到方程组如下:{x ̌1=∑b k nk=1x k =x 1x ̌k+1=s k x k −ℎk =x k+1整理得到:{−x 1∑b k nk=1x k =0−x k+1+s k x k =ℎk2 大学数学实验 实验报告 | 2014/4/52写成系数矩阵的形式如下:A =[b 1−1b 2b 3s 1−100s 2−1…b n−1b n0000⋮⋱⋮000000000⋯00−10s n−1−1]令h =[0, ℎ1,ℎ2,ℎ3,…,ℎn−2,ℎn−1]Tx =[x n , x n−1,…,x 1]T则方程组化为矩阵形式:Ax =h ,即为所求模型。
第12 章 实验十一线性方程组的数值解法实验目的:理解线性方程组基本定理,掌握求解线性方程组的各种方法,逆矩阵求法、高斯消去法、带选主元的分解法、迭代法。
12.1 线性方程组(基本定理) 已知m 个方程n 个未知量的方程组:⎪⎪⎩⎪⎪⎨⎧=+++=+++=+++m n mn m m n n n n b x a x a x a b x a x a x a b x a x a x a 22112222212********* (12.1)其中ij a 为已知系数,j x 未知,j b 为已知项,当),,2,1(n j b j =全为零时,称(12.1)方程组为线性齐次方程组,当),,2,1(n j b j =不全为零时,称(12.1)方程组为线性非齐次方程组。
上述方程组也可表示如下矩阵形式:b Ax = 其中A 、x 和y 分别定义为:⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=mn m n a a a a A LL L L L 1111⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=n x x x 1 ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=m b b b 1。
方程组(12.1)中的m 与n 有三种情况:m=n 或m<n 或m>n ,m=n 是最常见的。
当m=n 时,有系数矩阵的行列式0≠A ,方程组(12.1)有唯一解。
用矩阵的秩可以得到更满意的结果:(1)当n B R A R ==)()(时,方程组(12.1)有唯一解, (2)当n B R A R <=)()(时,方程组(12.1)有无穷多组解, (3)当)()(B R A R ≠时,方程组(12.1)无解。
12.2举例例12.1 用MATLAB 求解线性方程组:⎪⎪⎩⎪⎪⎨⎧=+++-=+++=+++=+++62332022428340213424321432143214321x x x x x x x x x x x x x x x x解法一:逆矩阵求法设 A=[1 2 1 4 ; 2 0 4 3 ;4 2 2 1 ;-3 1 3 2 ]; det(A)=-180, B=[13 28 20 6]’; 则:X=A\B % A\表示矩阵A 的逆,或用1-A 表示得到X =3.0000 -1.00004.0000 2.0000也可以用如下命令: X=B ’/A ’ 得到X =3.0000 -1.00004.0000 2.0000解法二:在MA TLAB 上,用高斯消去法求解方程组:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-62028132313122434024214321x x x x (12.2) 首先写出方程组对应的增广矩阵:a=[ 1 2 1 4 13; ... 2 0 4 3 28 ; ... 4 2 2 1 20 ; ...-3 1 3 2 6 ]其中前四列是系数矩阵,最后一列是方程组的常数列。
线性代数方程组的数值解法【实验目的】1. 学会用MATLAB 软件数值求解线性代数方程组,对迭代法的收敛性和解的稳定性作初步分析;2. 通过实例学习用线性代数方程组解决简化的实际问题。
【实验内容】【题目1】通过求解线性方程组A1x=b1和A2x=b2,理解条件数的意义和方程组的性态对解的影响。
其中A1是n阶范德蒙矩阵,即⎡1x0⎢1x1⎢A1=⎢⎢⎢⎣1xn-12x0x12 2xn-1n-1⎤ x0⎥ x1n-1⎥1,...,n-1 ,xk=1+0.1k,k=0,⎥ n-1⎥ xn-1⎥⎦A2是n阶希尔伯特矩阵,b1,b2分别是A1,A2的行和。
(1)编程构造A1(A2可直接用命令产生)和b1,b2;你能预先知道方程组A1x=和A2x=。
b2的解吗?令n=5,用左除命令求解(用预先知道的解可检验程序)b1(2)令n=5,7,9,…,计算A1,A2的条件数。
为观察它们是否病态,做以下试验:b1,b2不变,A1和A2的元素A1(n,n),A2(n,n)分别加扰动ε后求解;A1和A2不变,b1,b2的分量b1(n),分析A和b的微小扰动对解的影响。
b2(n)分别加扰动ε求解。
ε取10-1010,-8,10-6。
(3)经扰动得到的解记做x~,计算误差-x~x,与用条件数估计的误差相比较。
1.1构造A1,A2和b1,b2首先令n=5,构造出A1,A2和b1,b2。
首先运行以下程序,输出A1。
运行以下程序对A1,A2求行和:由于b1,b2分别是A1,A2的行和,所以可以预知x1=运行下列程序,用左除命令对b1,b2进行求解:得到以下结果: T。
x2=(1,1, ,1)1.2 计算条件数并观察是否为病态1.不加扰动,计算条件数。
运行以下程序:由此可知,A1,A2的条件数分别是3.574∗10, 4,766∗10。
2.b1,b2不变,A1(n,n),A2(n,n)分别加扰动(1)n=5时设x11,x12,x13分别为A1添加扰动10−10,10−8,10−6后的解。
线性方程组的数值解法一 实验目的1 学会用MATLAB 软件求解线性代数方程组,对迭代法的收敛性和解的稳定性做初步分析;2 通过实例学习线性代数方程组解决简化的实际问题。
二 实验内容1 已知方程组Ax =b ,其中A ∈ℝ20×20,定义为:31/21/41/231/21/41/41/231/21/41/41/231/21/41/23--⎡⎤⎢⎥---⎢⎥⎢⎥---⎢⎥-⎢⎥⎢⎥---⎢⎥--⎣⎦试通过迭代法求解此方程组,认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。
实验要求:(1) 选取不同的初始向量x (0)和不同的方程组右端项向量b ,给定迭代误差要求,用雅可比迭代法和高斯-赛德尔迭代法计算,观测得到的迭代向量序列是否均收敛?若收敛,记录迭代次数,分析计算结果并得出你的结论。
(2) 取定右端项向量b 和初始向量x (0),将A 的主对角线元素成倍增长若干次,非对角线元素不变,每次用雅可比迭代法计算,要求迭代误差满足‖x (k+1)−x (k)‖∞<10−5,比较收敛速度,分析现象并得出你的结论。
初步解决:首先建立利用雅克比迭代法和高斯-赛德尔迭代法计算的函数M 文件。
雅克比迭代法:高斯-赛德尔迭代法:关于迭代法的收敛性:由原式Ax=b化简可得x=Bx+f,也就是说,研究此处迭代法是否收敛实际上就是研究有原矩阵A经变换之后所得矩阵B的性质。
而由于该题第一问中,A保持不变,所以对于收敛性只需研究一次。
在命令栏中输入以下命令:(1)构造题目中的矩阵A:(2)提取对角矩阵、上三角矩阵、下三角矩阵:(3)对于两种迭代方法分别求B:(4)分别计算两个矩阵的三种范数:n1和m1、n2和m2、n3和m3分别是矩阵B1、B2的1-范数、2-范数和∞范数。
由书上定理,矩阵的谱半径不超过任何一个范数,即ρ(B)≤‖B‖,而由图中可以看出,两个矩阵的六个范数没有一个大于1,所以两个矩阵的谱半径一定都小于1,所以此时两种迭代方法均收敛。
数值分析实验报告(一)一、实验名称:求解线性方程组二、实验问题:从2到9的每一个n 值,解n 阶方程组Ax=b 。
这里A 和b 如下定义:a ij =(i+j-1)-1,b i =p(n+i-1)-p(i-1)式中p(x)=(x 2/24)(2+x 2(-7+n 2(14+n(12+3n))))解释发生的现象。
三、计算公式(数学模型):Gauss 消去法是将方程组Ax=b ,用逐次消去未知数的方法,化为等价的三角方程组,再回代求解的方法。
设Ax=b ,A n n R ⨯∈,如果),...,2,1(,0)(n k a k kk=≠,则可通过Gauss 消去法将Ax=b 约化为等价的三角形方程组⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡)()2(2)1(121)()2(2)2(22)1(1)1(12)1(11n n n n nn n n b b b x x x a a a a a a且计算公式为:(a)消元计算(k=1,2,…,n)()()()⎪⎪⎪⎩⎪⎪⎪⎨⎧+=-=+=-=+==++n k i a m b b n k j i a m a a n k i a a m k k ik k i k i k kj ik k ij k ijk kn k ik ik ,,1,,1,,1)()()1()()()1()()( (b)回代计算:()⎪⎪⎩⎪⎪⎨⎧-=-==∑+=n i j i ii ji i ij i i in nn n n n n ia x ab x a b x 1)()()()()(1,2,,1四、结构程序设计:1、数据流向图2五、计算程序C++程序//改进的平方根法#include<iostream.h>double p(double x,double n){return x*x/24*(2+x*x*(-7+n*n*(14+n*(12+3*n))));}void main(void){cout<<"*************************线性方程组求解的数值实验题*****************************"<<'\n';double A[9][9]={0};double L[9][9]={0};double B[9]={0};double d[9]={0};double x[9]={0};double y[9]={0};double dd=0;int n,i,j,k;for(n=2;n<10;n++){cout<<n<<"阶方程组Ax=b的解为"<<'\n';for(i=0;i<n;i++){for(j=0;j<n;j++){A[i][j]=1.0/(i+j+1); //赋Ax=b中系数A的初值(A是对称阵)}B[i]=p((double)n+i,(double)n)-p((double)i,(double)n); //赋Ax=b中系数b的初值(调用函数p)}d[0]=A[0][0]; //d[]用来记录A[k][k]for(i=0;i<n;i++) L[i][0]=A[i][0]/d[0]; //LDL'分解中L[i][j]的第一列for(i=1;i<n;i++){for(k=0;k<i;k++) dd+=L[i][k]*L[i][k]*d[k];d[i]=A[i][i]-dd; //LDL'分解中的d[i]dd=0;for(j=i;j<n;j++){for(k=0;k<j;k++) dd+=L[i][k]*d[k]*L[j][k]; //LDL'分解中除第一列外的L[i][j]L[i][j]=(A[i][j]-dd)/d[j];dd=0;}}y[0]=B[0]; //对y的回代过程for(i=1;i<n;i++){for(k=0;k<i;k++) dd+=L[i][k]*y[k];y[i]=B[i]-dd;dd=0;}x[n-1]=y[n-1]/d[n-1]; //对x的回代过程for(i=n-2;i>=0;i--){for(k=i+1;k<n;k++) dd+=L[k][i]*x[k];x[i]=y[i]/d[i]-dd;dd=0;}for(i=0;i<n;i++){cout<<'x'<<i+1<<"="<<x[i]<<'\n';}cout<<'\n';}} //在C++数组中是以开始的,所以i,j,k取值相应变到从0开始运行结果如下:六、结果和讨论:对于2到9阶的线性方程组的解来看,Gauss 消去法有其自己的弊端,使有些结果有误差。
实验5 线性代数方程组的数值解法化工系 分0班 2010011805 张亚清【实验目的】1、 学会用MATLAB 软件数值求解线性代数方程组,对迭代法的收敛性和解的稳定性作初步分析;2、 通过实例学习用线性代数方程组解决简化的实际问题。
【实验内容】 1、 题目3已知方程组Ax=b ,其中2020⨯ℜ∈A ,定义为试通过迭代法求解此方程组,认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。
实验要求:(1)选取不同的初始向量)0(x和不同的方程组右端向量b ,给出迭代误差要求,用雅克比迭代法和高斯-赛德尔迭代法计算,观测得到的迭代向量序列是否均收敛?若收敛,记录迭代次数,分析计算结果并得出你的结论;(2)取定右端向量b 和初始向量)0(x,将A 的主对角线元素成倍增长若干次,非主对角线元素不变,每次用雅克比迭代法计算,要求迭代误差满足5)()1(10-∞+<-k k x x ,比较收敛速度,分析现象并得出你的结论。
【问题分析】对于线性方程组Ax=b ,满足0≠ii a (i=1,2,...,n ),将A 分解为A=D-L-U ,则雅克比迭代公式等价于如下的矩阵形式,...2,1,0,)(1)(1)1(=++=--+k b D x D L D x k k或b D f U L D B f x B xJ J J k J k 11)()1(),(,--+=+=+=。
类似的,Ax=b 的高斯-赛德尔迭代公式等价于如下矩阵形式b L D f U L D B f x B x J S G S G k S G k 11)()1()(,)(,-----+-=-=+=。
【问题解答】(1)选取初始向量T x )1,...,1,1()0(=,T b )1,...,1,1(=,迭代要求为4)()1(10-∞+<-k k x x 。
将A 按A=D-L-U 分解为如下三个矩阵:①对方程组进行雅克比迭代,利用MATLAB 编程计算。
实验报告
——线性方程组的数值解法
姓名:
班级:
学号:
日期:
一 实践目的
1. 熟悉求解线性方程组的有关理论和方法。
2. 会编列主元消去法,全主元消去法,雅克比迭代法及高斯-赛德尔迭代法的程序。
3.通过实际计算,进一步了解各种方法的优缺点,选择合适的数值方法。
4. 进一步应用数学知识,扩展数学思维,提高编程能力。
二 问题定义及题目分析
1.求解线性方程是实际中常遇到的问题。
而一般求解方法有两种,一个是直接求解法,一个是迭代法。
2.直接法是就是通过有限步四则运算求的方程准确解的方法。
但实际计算中必然存在舍入误差,因此这种方法只能得到近似解。
3.迭代法是先给一个解的初始近似值,然后按一定的法则求出更准确的解,即是用某种极限过程逐步逼近准确解的方法。
4.这次实践,将采用直接法:列主元高斯消去法,全主元高斯消去法;迭代法:雅克比迭代法和高斯-赛德尔迭代法。
三 详细设计 设有线性方程组
11112211n n a x a x a x b +++=
21122222n n a x a x a x b +++= 1122n n nn n n a x a x a x b +++=
令A= 1112121
2221
2
n n n n nn a a a a
a a a a a ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥
⎣⎦ x=
12n x x x ⎡⎤⎢⎥
⎢⎥⎢⎥⎢⎥⎣⎦ b= 12n b b b ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦
一. 上三角方程组的求解很简单。
所以若能将方程组化为上三角方程组,那就很容易得到方程的解,高斯消去法则是一种简洁高效的方法,可以将方程组化为上三角方程组,但是这种方法必须满足每一步时0kk a =才能求解,还有当kk a 绝对值很小时,作分母会引起较大的舍入误差。
因此在消元过程中应尽量选择绝对值较大的系数作为主元素。
1.列主元消去法很好的解决这一问题。
将1x 的系数1(1)k a k n ≤≤中绝对值最大者作为主元素,交换第一行和此元素所在的行,然后主元消去非主元项1x 的系数,。
一般的,在k x 的一些相关系数之前,先从,1,,,k k k k n k a a a + 中选出绝对值最大者作为主元素,交换第k 行和次主元素所在的行,再消去k x 的一些相关系数,只要方程组A 非奇异,则每个主元定不为0,故消元过程一定能将方程组化为上三角方程组。
将方程用增广矩阵[](1)()ij n n A b a ⨯+= 表示。
(1) 消元过程
对01,2,1k n =- ,;
1. 选主元,找{},1,k i k k n ∈+ ,使得,max ()k
i k ik a a k i n =≤≤
2. ,=0k
i k a ,则矩阵A 奇异,程序结束,否则执行下一步。
3. 如果k i k ≠,则交换第k 行与第k i 行对应的元素位置,即 (,1,,1)k
kj i j a a j k k n ↔=++
4. 消元,对1,2,,i k k n =++ 计算/ik ik kk l a a = 对1,2,,1j k k n =+++ 计算ij ij ik kj a a l a =- (2) 回代过程
,1/n n n nn x a a +=,对1,,2,1i n =- 计算 ,11()/n
i i n ij i ii j i x a a x a +=+=-∑
2. 全主元消去法(对列主元消去法的改进)
首先,从方程组系数ij a 中选出绝对值最大者作为主元,然后通过行交换及列交换,是主元素移到第一行第一列,然后消去第一列的非主元项的系数。
一般的,在第k 次消元之前,先从n-k+1行和n-k+1列中选取绝对值最大者作为主元素,交换第k 行和此主元素所在之行,再交换第k 列和此主元素所在之列。
然后消元,消去第k 列的一些系数,但每一步都必须跟踪记录12,n x x x 之间顺序的变化。
for(l=i;l<n;l++) 从剩余的n-i 列和n-i 列找主元 {
for(j=i;j<n;j++) 并将最大者记为A[k1][k2] if(fabs(A[l][j])>fabs(A[k1][k2]))k1=l,k2=j; }
if(k1!=i) 如果k1不等于i ,则第i 行与第k1行互换 for(j=0;j<=n;j++) {
p=A[i][j];
A[i][j]=A[k1][j]; A[k1][j]=p; }
for(i=0;i<n;i++)Y[i]=i; 初始化 Y[n]
if(k2!=i) if 成立则第i 列与第k2列互换 {
p=Y[i]; 相应的将Y[n]的第i 列与第k2列互换
Y[i]=Y[k2]; 则Y[n]起到了跟踪记录列交换的作用
Y[k2]=p;
for(j=0;j<n;j++) 第i 列与第k2列互换 {
p=A[j][i];
A[j][i]=A[j][k2]; A[j][k2]=p; } }
二 迭代法经常用于求解阶数高且矩阵稀疏的方程组。
1. 雅克比迭代法
设方程组系数矩阵A 的对角元素0ii a ≠(1,2,,)i n = ,max 为迭代的最大次数,e 为误差线。
(1)区初始向量000
12
(,,)n x x x x = ,令k=0; (2) 对1,2,,i n = ,计算
1
1,()/n
k k i
i ij j ii j j i
x
b a x a +=≠=-
∑
(3) 11
n
k k i i i x x e +=-<∑,则输出1k i x +,程序结束,否则进行下一步。
(4) 如果max k ≥,则并不收敛,程序结束,否则执行第2步。
2. 高斯-赛德尔迭代法(对雅克比迭代法的改进)
这种迭代法在于充分利用每一小步所求出的新分量。
因此在新分量一旦求出之时,马上、就用新分量代替。
即在雅克比迭代法法中求1
k i x +时用111121,k k k i x x x +++- 分别代替121,k k k
i x x x - 。
使得程序效率更高。
for(i=0;i<n;i++) for 循环进行迭代
{ 并记录每次迭代的误差 p=0;
for(j=0;j<n;j++) {
if(i==j) continue; 如果i 等于j 则不执行下一步 p=p+A[i][j]*X[j];
}
Y[i]=(A[i][n]-p)/A[i][i]; 新分量的计算
error=error+fabs(Y[i]-X[i]); 误差计算
X[i]=Y[i]; 新分量代替旧分量
}
四调试与测试分析
由于for循环嵌套语之多句及数组下标控制之难,在编程过程中困难重重。
因此就先将程序的大致算法写出来,再按照算法编程。
而C 语言数组下标从0 开始,因此按算法中所有下标减一编程。
在列主元消去法中,每次交换第i行和主元所在的第k行(按下标从1开始),没有交换第1,2…i-1列,但不影响结果。
而在编写全主元消去法时,每次能正确跟踪,结果总不对。
经过几次详细检查发现在交换行或列时,都是从第i行或第i列开始交换,发现列交换不能类比列主元消去法中的交换方法,必须考虑全局。
1.直接法测试用例:1 1 0 3 4
2 1 -1 1 1
3 -1 -1 3 -3
-1 2 3 -1 4
列主元高斯消去法输出全主元高斯消去法输出
X[1]=-1.33333 X[4]=1
X[2]=2.33333 X[3]=-0.333333
X[3]=-0.333333 X[1]=-1.33333
X[4]= 1 X[2]=2.33333
真实值为:X[1]=- 4/3,X[2]= 7/3,X[3]= -1/3, X[4]= 1。
2.迭代发测试用例:max=100
e=0.0001
10 -2 -1 3
-2 10 -1 15
-1 -2 5 10
雅克比迭代法输出高斯-赛德尔迭代法输出
迭代次数为12迭代次数为7
X[1]=0.999991 X[1]=0.999994
X[2]=1.99999 X[2]=2
X[3]=2.99999 X[3]=3
真实值为:X[1]=1,X[2]=2,X[3]=3。
两种直接法并无很大明显差别。
高斯-赛德尔迭代法明显优于雅克比迭代法。
五结果分析及总结
实验结果与预期一致,程序设计基本合理,能正确运行出结果。
在设计过程中能够及时发现问题并解决问题。
通过这次计算实践,对方程组求解的方法有了深入了解和掌握。
可视化实验结果为:
图1.列主元高斯消去法
图2.全主元高斯消去法
图3.雅克比迭代法
图4.高斯-赛德尔迭代法。