列主元消去法解方程组实验报告
Gauss列主元消去法实验

<数值计算方法>实验报告1.实验名称实验2 Gauss 列主元消去法2.实验题目用Gauss 列主元消去法求解线性方程组。
0.0011 2.0002 3.0003 1.0001.0001 3.7122 4.6233 2.0002.0001 1.0722 5.6433 3.000x x x x x x x x x ++=⎧⎪-++=⎨⎪-++=⎩3.实验目的加深自己对Gauss 列主元消去法的理解和认识,并且通过做实验或做练习来加强自己Gauss 列主元消去法的掌握,学会并灵活运用Gauss 列主元消去法来求解方程组。
4.基础理论-------Gauss 列主元消去法1.Gauss 列主元消去法的基本思想是:在进行第k (k=1,2,...,n-1)步消元时,从第k 列的kk a 及以下的各元素中选取绝对值最大的元素,然后通过行变换将它交换到主元素kk a 的位置上,再进行消元。
2.Gauss 列主元消去法的优点:当kk a (k=1,2,...,n-1)的绝对值很小时,用Gauss 列主元消去法来求解方程组时,可以避免所的数值结果产生较大误差或失真。
5.实验环境实验系统:Win 7实验平台:VisualC++语言6.实验过程写出算法→编写程序→计算结果Gauss 列元消去法的算法Input:方程组未知量的个数n;增广矩阵()()1,2,...,T ij A a A A An ==,其中i=1,2,…,n; j=1,2,…,n+1Output:方程组的解x1,x2,…,xn,或失败信息。
1. for i ←1ton-1 do;2. temp ←|ii a |;3. p ←I;4. for j ←i+1 to n do5. if ||ji a >temp then6. p ←j;8. end9. end10. if temp=0 then11. |return False;12. end13. if p ≠I then14. p A ⇔i A ;//i,p 两行交换15. end//列选主元16. for j ←i+1 to n do17.*j ji i A m A -ji m ←/ji ii a a ;18. j A ←*j ji i A m A -;//消元19. end7.实验结果原方程组的解为:X1=-0.490396 , x2=-0.051035 ,x3=0.3675208.附录程序清单#include<iostream.h> #include"stdio.h"#include"math.h"void main ( ){ int n=3,i,j,k,p;doubleA[10][10]={{0.001,2.000,3.000,1.000},{-1.000,3.712,4.623,2.000},{-2.0 00,1.072,5.643,3.000}},temp,m,x[100];for(i=0;i<n;i++){ //选主元temp=fabs(A[i][i]); p=i;for(k=i+1;k<n;k++)if(fabs(A[k][i])>temp){temp=fabs(A[k][i]); p=k;}if(temp==0){ printf("\n无法求解:");return;}if(p!=i)for(j=0;j<n+1;j++){ temp=A[i][j];A[i][j]=A[p][j];A[p][j]=temp;}//消元for(k=i+1;k<n;k++){ m=A[k][i]/A[i][i];for(j=i+1;j<=n;j++)A[k][j]=A[k][j]-m*A[i][j];}}//回代for(i=n-1;i>=0;i--){x[i]=A[i][n];for(j=i+1;j<n;j++)x[i]=x[i]-A[i][j]*x[j];x[i]=x[i]/A[i][i];}printf("\nx=\n");for(i=0;i<n;i++)printf("%f \n",x[i]);}。
高斯列主元消元法解方程组

format long
n=50;
for a=1:n;
for c=1:n;
A(a,c)=1/(a+c-1);
end
end
b=[1:n]';
x=gauss(A,b);
fprintf('%12.8f ',x)
实验结果为:
23964.74834020-4122991.37806186176597502.11690632-3260288682.1972771031898227709.80876200-181566275672.69705000 619956768219.91504000-1230092760464.934800001193508577123.06670000 -71714652666.11485300-587044042244.59863000-454093769462.97290000 1744941129911.71780000-2715153157131.341300002706193814657.00730000 758555481068.28345000-2729174782832.09330000-650678553386.91565000 770544482197.08948000661071297466.389400003114329102228.31640000 -2419453655729.64110000-2498342518894.922900001991167370388.64180000 806469536595.88586000-3357159715176.254400003434745294165.19480000 -1115456797533.65920000835226816241.179930001614362560020.15380000 -2300534244207.34180000781205861296.42358000-3724002810480.48780000 -101978267755.295900003815588397071.54350000-269576755745.41956000 -2354847317920.947800005144417865217.228500001355869718979.75980000 1213622862842.97950000-11349162850658.63300000367430197567.42474000 2362009479539.30520000712955508308.938350005794833430082.82320000 -2649289143449.391600002040834706833.95900000-2941986228491.25980000 -2402163300368.774400002234821162570.22750000
高斯列主元消去法实验报告

《数值计算方法》实验报告专业:年级:学号:姓名:成绩:1.实验名称实验2高斯列主元消去法2. :用Gauss列主消去法求解线性方程组0.001*X1+2.000*X2+3.000*X3=1.000-1.000*X1+3.217*X2+4.623*X3=2.000-2.000*X1+1.072*X2+5.643*X3=3.0003.实验目的a.熟悉运用已学的数值运算方法求解线性方程—Gauss列主消去法;b.加深对计算方法技巧的认识,正确使用计算方法来求解方程;c.培养用计算机来实现科学计算和解决问题的能力。
4.基础理论列主元消去法:a.构造增广矩阵b.找到每列绝对值的最大数;c.行变换;d.消去;e.回代5.实验环境Visual C++语言6.实验过程实现算法的流程图:7.结果分析a.实验结果与理论一致;b.由于数值设置成双精度浮点型,所以初值对计算结果影响不大;c.运用程序能更好的实现计算机与科学计算的统一和协调。
8. 附录程序清单#include<stdio.h>#include<math.h>int main(){int n=3,i,j,k,p;double a[4][4];double b[4];double x[4];double m[4][4];double temp;a[1][1]=0.001; a[1][2]=2.000; a[1][3]=3.000; b[1]=1.000;a[2][1]=-1.000; a[2][2]=3.1712; a[2][3]=2.000; b[2]=2.000;a[3][1]=-2.000; a[3][2]=1.072; a[3][3]=5.643; b[3]=3.000;for(i=1;i<=n-1;i++){temp=a[i][i];p=i;for(j=i+1;j<=n;j++)if(fabs(a[j][i])>temp){temp=a[j][i];p=j;}if(temp==0)return 0;if(p!=i) //换行{for(j=1;j<=n;j++)a[0][j]=a[i][j];for(j=1;j<=n;j++)a[i][j]=a[p][j];for(j=1;j<=n;j++)a[p][j]=a[0][j];b[0]=b[i];b[i]=b[p];b[p]=b[0];}for(j=i+1;j<=n;j++){m[j][i]=a[j][i]/a[i][i];for(k=i;k<=n;k++)a[j][k]=a[j][k]-m[j][i]*a[i][k];}}if(a[n][n]==0)return 0;x[n]=b[n]/a[n][n];for(i=n-1;i>=1;i--)//回代{temp=0;for(j=i+1;j<=n;j++)temp=temp+a[i][j]*x[j];temp=b[i]-temp;x[i]=temp/a[i][i];}for(i=1;i<=n;i++)//输出结果{printf("输出结果为:x[%d]=%lf ",i,x[i]);}printf("\n");return 0;}。
列主元消去法解方程组实验报告

实验名称:列主元消去法解方程组1 引言我们知道,高斯消去法是一个古老的解线性方程组的方法。
而在用高斯消去法解Ax=b时,其中设A为非奇异矩阵,可能出现的情况,这时必须进行带行交换的高斯消去法。
但在实际计算中即使但其绝对值很小时,用作除数,会导致中间结果矩阵元素数量级严重增长和舍入误差的扩散,使得最后的结果不可靠。
因此,小主元可能导致计算的失败,我们应该避免采用绝对值很小的主元素。
为此,我们在高斯消去法的每一步应该在系数矩阵或消元后的低阶矩阵中选取绝对值最大的元素作为主元素,保持乘数,以便减少计算过程中舍入误差对计算解的影响。
一种方式是完全主元消去法,这种消去法是在每次选主元时,选择为主元素。
这种方法是解低阶稠密矩阵方程组的有效方法,但这种方法在选取主元时要花费一定的计算机时间。
实际计算中我们常采用部分选主元的的消去法。
列主元消去法即在每次选主元时,仅依次按列选取绝对值最大的元素作为主元素,且仅交换两行,再进行消元计算。
2 实验目的和要求运用matlab编写一个.m文件,要求用列主元消去法求解方程组(实现PA=LU):要求输出以下内容:(1)计算解x;(2) L,U;(3)整形数组IP(i)(i=1,2,…,n-1)(记录主行信息)3 算法原理与流程图(1)算法原理设有线性方程组Ax=b,其中设A为非奇异矩阵。
方程组的增广矩阵为第1步(k=1):首先在A的第一列中选取绝对值最大的元素,作为第一步的主元素:,然后交换(A,b)的第1行与第i1行元素,再进行消元计算。
设列主元素消去法已经完成第1步到第k-1步的按列选主元,交换两行,消元计算得到与原方程组等价的方程组第k步计算如下:对于k=1,2,…,n-1(1)按列选主元:即确定ik使(2)如果,则A为非奇异矩阵,停止计算。
(3)如果ik≠k,则交换[A,b]第ik行与第k行元素。
(4)消元计算消元乘数满足:(5)回代求解计算解在常数项b(n)内得到。
解线性方程组的列主元素高斯消去法和LU分解法实验报告

解线性方程组的列主元素高斯消去法和LU 分解法一、实验目的:通过数值实验,从中体会解线性方程组选主元的必要性和LU 分解法的优点,以及方程组系数矩阵和右端向量的微小变化对解向量的影响。
二、实验内容:解下列两个线性方程组(1)⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4987.023.116.427.199.103.601.3321x x x(2)⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----15900001.582012151526099999.23107104321x x x x 三、实验要求:(1) 用你熟悉的算法语言编写程序用列主元高斯消去法和LU 分解求解上述两个方程组,输出Ax=b 中矩阵A 及向量b, A=LU 分解的L 及U ,detA 及解向量x.(2) 将方程组(1)中系数3.01改为3.00,0.987改为0.990,用列主元高斯消去法求解变换后的方程组,输出列主元行交换次序,解向量x 及detA ,并与(1)中结果比较。
(3) 将方程组(2)中的2.099999改为2.1,5.900001改为5.9,用列主元高斯消去法求解变换后的方程组,输出解向量x 及detA ,并与(1)中的结果比较。
(4)用MATLAB的内部函数inv求出系数矩阵的逆矩阵,再输入命令x=inv(A)*b,即可求出上述各个方程组的解,并与列主元高斯消去法和LU分解法求出的解进行比较,体会选主元的方法具有良好的数值稳定性。
用MATLAB的内部函数det求出系数行列式的值,并与(1)、(2)、(3)中输出的系数行列式的值进行比较。
四、实验过程:(1)列主元高斯消去法的主程序为function [RA,RB,n,X]=liezhuY(A,b)B=[A b]; n=length(b); RA=rank(A);RB=rank(B);zhica=RB-RA;D=det(A)if zhica>0,disp('请注意:因为RA~=RB,所以此方程组无解.')returnendif RA==RBif RA==ndisp('请注意:因为RA=RB=n,所以此方程组有唯一解.')X=zeros(n,1); C=zeros(1,n+1);for p= 1:n-1[Y,j]=max(abs(B(p:n,p))); C=B(p,:);B(p,:)= B(j+p-1,:); B(j+p-1,:)=C;for k=p+1:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1);endendb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);endelsedisp endend解方程组(1)在MATLAB工作窗口输入>>A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.819.34];b=[1;1;1];[RA,RB,n,X]=liezhuY(A,b)运行后输出结果为请注意:因为RA=RB=n,所以此方程组有唯一解. D=-0.1225RA =3 RB =3 n =3X = 397.8654-157.6242-123.1120解方程组(2)在MATLAB工作窗口输入>>A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];b=[8;5.900001;5;1];[RA,RB,n,X]=liezhu(A,b)运行后输出结果为请注意:因为RA=RB=n,所以此方程组有唯一解. D=-762.0000RA =4 RB =4 n =4X =0.0000-1.00001.00001.0000LU分解法及MATLAB主程序为function hl=zhjLU(A)[n n] =size(A); RA=rank(A);D=det(A)if RA~=ndisp('请注意:因为A的n阶行列式hl等于零,所以A不能进行LU分解.A的秩RA如下:'), RA,hl=det(A);returnendif RA==nfor p=1:nh(p)=det(A(1:p, 1:p));endhl=h(1:n);for i=1:nif h(1,i)==0disp('请注意:因为A的r阶主子式等于零,所以A不能进行LU分解.A 的秩RA和各阶顺序主子式值hl依次如下:'), hl;RAreturnendendif h(1,i)~=0disp('请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:')for j=1:nU(1,j)=A(1,j);endfor k=2:nfor i=2:nfor j=2:nL(1,1)=1;L(i,i)=1;if i>jL(1,1)=1;L(2,1)=A(2,1)/U(1,1); L(i,1)=A(i,1)/U(1,1); L(i,k)=(A(i,k)- L(i,1:k-1)*U(1:k-1,k))/U(k,k);elseU(k,j)=A(k,j)-L(k,1:k-1)*U(1:k-1,j);endendendendhl;RA,U,Lendend解方程组(1)在MATLAB工作窗口输入>>A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.819.34];h1=zhjLU(A)运行输出结果为请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:D=9.8547RA =3U =3.0100 6.0300 1.99900 4.1600 -2.07340 0 5.3016L =1.0000 0 00.4219 1.0000 00.3279 -1.6316 1.0000h1 =3.0100 4.8635 -0.1225解方程组(2)在MATLAB工作窗口输入>>A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 02];h1=zhjLU(A)运行后输出结果为请注意:因为A的各阶主子式都不等于零,所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下:D=-762.0000RA =4U =10.0000 -7.0000 0 1.00000 2.1000 6.0000 2.30000 0 -2.1429 -4.23810 -0.0000 0 12.7333L =1.0000 0 0 0-0.3000 1.0000 0 00.5000 1.1905 1.0000 -0.00000.2000 1.1429 3.2000 1.0000h1 =10.0000 -0.0000 -150.0001 -762.0001(2)在MATLAB工作窗口输入>>A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.819.34];b=[1;1;1];A(1,1)=3;A(1,3)=0.990;[RA,RB,n,X]=liezhu(A,b)请注意:因为RA=RB=n,所以此方程组有唯一解.RA =3 RB =3 n =3X = -4.02641.91931.5210hi = 3.0000 4.8219 9.8547在MATLAB工作窗口输入x=[397.8654;-157.6242;-123.1120]';x1=[-4.0264;1.9193;1.5210]';wucha=x1-x运行后输出结果为wucha =-401.8918 159.5435 124.6330(3)在MATLAB工作窗口输入>>A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];A(2,2)=2.1;b(2,1)=5.9;b=[8;5.900001;5;1];[RA,RB,n,X]=lie zhu(A,b)运行后输出结果为请注意:因为RA=RB=n,所以此方程组有唯一解.RA =4 RB =4 n =4X =0.0000-1.00001.00001.0000h1 =10.0000 -0.0000 -150.0000 -762.0000在MATLAB工作窗口输入>>x=[0;-1;1;1]';x1=[0;-1;1;1]';wucha=x1-x运行后输出结果为wucha = 0 0 0 0(4)解方程组(1)在MATLAB工作窗口输入>>A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.81 9.34];B=inv(A)运行后结果为B =-268.9293 538.3418 128.4529106.7599 -213.4281 -50.956183.3992 -166.8022 -39.7090在MATLAB工作窗口输入>>b=[1;1;1];x=inv(A)*b运行后结果为x =397.8654-157.6242-123.1120在MATLAB工作窗口输入>>A=[3.01 6.03 1.999;1.27 4.16 -1.23;0.987 -4.81 9.34];A(1,1)=3;A(1,3)=0.990;B=inv(A)运行输出结果为B = 3.3424 -6.1983 -1.1705-1.3269 2.7442 0.5020-1.0365 2.0682 0.4893在MATLAB工作窗口输入>>b=[1;1;1];x=inv(A)*b运行后输出结果为x =-4.02641.91931.5210解方程组(2)在MATLAB工作窗口输入>>A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];B=inv(A) 运行后结果为B =-0.0223 -0.0984 0.1181 0.1686-0.1601 -0.1181 0.1417 0.26900.0108 0.1063 0.0724 -0.07550.1024 0.1575 -0.1890 0.1969在MATLAB工作窗口输入>>b=[8;5.900001;5;1];x=inv(A)*b运行后输出结果为x = 0-1.00001.00001.0000在MATLAB工作窗口输入>>A=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];A(2,2)=2.1;B=inv(A)运行后输出结果为B =-0.0223 -0.0984 0.1181 0.1686-0.1601 -0.1181 0.1417 0.26900.0108 0.1063 0.0724 -0.07550.1024 0.1575 -0.1890 0.1969在MATLAB工作窗口输入>>b=[8;5.900001;5;1];b(2,1)=5.9;x=inv(A)*b运行后输出结果为x =-0.0000-1.00001.00001.0000五、实验结果分析:实验的数学原理很容易理解,也容易上手。
列主元消去法matlab实验报告

列主元消去法matlab实验报告列主元消去法是一种常用的线性方程组求解方法,它通过选取主元元素来消去其他元素,从而简化方程组的求解过程。
本文将以Matlab为工具,对列主元消去法进行实验研究,并给出相应的实验报告。
我们需要明确列主元消去法的基本原理。
列主元消去法的核心思想是选取每一列的主元素,通过消去其他元素,从而将方程组转化为上三角形或下三角形的形式。
具体来说,通过选取第一列的主元素,将第一列下方的元素消去;然后选取第二列的主元素,将第二列下方的元素消去;依此类推,直到最后一列。
这样,我们就得到了一个上(下)三角形的方程组,可以通过回代(代入)的方法求解。
接下来,我们使用Matlab编写代码,实现列主元消去法。
首先,我们需要输入一个线性方程组的系数矩阵A和常数向量b,其中A 是一个n×n的矩阵,b是一个n×1的向量。
然后,我们通过选取主元素的方式进行消去操作,得到一个上三角形的方程组。
最后,我们通过回代(代入)的方法求解方程组的解。
具体实现的代码如下所示:```matlabfunction x = gauss_elimination(A, b)n = size(A, 1); %方程组的个数% 消元过程for k = 1:n-1[~, p] = max(abs(A(k:n, k))); %选取主元素 p = p + k - 1;% 交换第k行和第p行temp = A(k, :);A(k, :) = A(p, :);A(p, :) = temp;temp = b(k);b(k) = b(p);b(p) = temp;% 消去操作for i = k+1:nfactor = A(i, k) / A(k, k);A(i, :) = A(i, :) - factor * A(k, :);b(i) = b(i) - factor * b(k);endend% 回代(代入)过程x = zeros(n, 1);x(n) = b(n) / A(n, n);for i = n-1:-1:1x(i) = (b(i) - A(i, i+1:n) * x(i+1:n)) / A(i, i);endend```接下来,我们将使用一个具体的例子来说明列主元消去法的求解过程。
试验四 Gauuss列主元消去法解线性方程组实验报告

for(j=0;j<m+1;j++) scanf("%lf",&ab[i][j]); printf("您输入的系数增广矩阵为:\n"); for(i=0;i<m;i++) { for(j=0;j<=m;j++) printf(" %10.9f",ab[i][j]); printf("\n"); } for(i=0;i<m-1;i++)//i 是行 { Change(ab,m,i); for(k=i+1;k<m;k++) { f=ab[i][i]; g=ab[k][i]; e=-(g/f); for(t=0;t<m+1;t++) { ab[k][t]+=ab[i][t]*e; } } } printf("经行处理后原矩阵变为:\n"); for(k=0;k<m;k++) { for(t=0;t<=m;t++) printf(" %10.9f",ab[k][t]); printf("\n"); } printf("方程组的解为:\n"); ab[2][2]=ab[2][3]/ab[2][2]; ab[1][1]=(ab[1][3]-ab[1][2]*ab[2][2])/ab[1][1]; ab[0][0]=(ab[0][3]-ab[0][2]*ab[2][2]-ab[0][1]*ab[1][1])/ab[0][0]; for(i=0;i<m;i++) printf("x%d=%10.9f\n",i+1,ab[i][i]); }
列主元消去法解方程组

数值分析上机实验报告(二)一、问题描述:利用列主元消去法求解下列方程组2X1+5X2+3X3 - 2X4=72X1- 2X2+3X3+5X4=-1X1+3X2+2X3+3X4=0X1+2X2+ X3 - 2X4=4二、算法原理:由高斯消去法知道,在消去过程中可能出现a kk(k)=0的情况,这时候消去法将无法进行,所以最好选取系数矩阵(或消元后的低阶矩阵)中绝对值最大的元素作为主元,以使高斯消去法具有较好的数值稳定性。
三、实验步骤:1、det 1;2、对于k=1,2,···,n-1(1)按列选主元|a ik.k|=max|a ik|(2)如果a i.k=0,则计算停止(det(A)=0)(3)如果i k=k则转(4)换行:a kj a ik·j(j=k,k+1,···,n)b k b ikdet -det(4)消元计算对于i=k+1,···,n○1am ik=a ik/a kk○2对于j=k+1,···,n a ij a ij—m ik*a kj○3b i b i-m ik*b ik(5)det a kk*det3、如果则计算停止(det(A)=0)4、回代求解(1)b n b n/a nn(2)对于i=n-1···,2,1bi(bi-∑aij*bj)/aii5.det ann*det四、实验框图五、源程序# include <stdio.h># include<math.h># define n 4main(){int i,j,k,l;float A[n][n],b[n],x[n],max;//输入系数矩阵及右端项for(i=0;i<n;i++)for(j=0;j<n;j++){printf("A[%d][%d]=",i,j);scanf("%f;",&A[i][j]);}for(i=0;i<n;i++){printf("b[%d]=",i);scanf("%f;",&b[i]);}//列主元消去过程for(k=0;k<n-1;k++){max=abs(A[k][k]);l=k;for(i=k+1;i<n;i++)if(abs(A[i][k])>max){max=abs(A[i][k]);l=i;} if(l>k){for(j=k;j<n;j++){max=A[k][j];A[k][j]=A[l][j];A[l][j]=max;}max=b[k];b[k]=b[l];b[l]=max;}for(i=k+1;i<n;i++){max=A[i][k]/A[k][k];for(j=k+1;j<n;j++)A[i][j]=A[i][j]-max*A[k][j];b[i]=b[i]-max*b[k];}}//回代过程x[n-1]=b[n-1]/A[n-1][n-1];for(k=1;k<n;k++){i=n-k-1;x[i]=b[i];for(j=i+1;j<n;j++)x[i]=x[i]-A[i][j]*x[j];x[i]=x[i]/A[i][i];}//输出解for(i=0;i<n;i++)printf("x[%d]=%f;",i,x[i]);getchar();}六、运行结果。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
实验名称:列主元消去法解方程组
1引言
我们知道,高斯消去法是一个古老的解线性方程组的方法。
而在用高斯消去法解Ax=b 时,其中设A为非奇异矩阵,可能出现a kk)0的情况,这时必须进行带行交换的高斯消去法。
但在实际计算中即使a kk)0但其绝对值很小时,用a;:〉作除数,会导致中间结果矩阵A(k)元素数量级严重增长和舍入误差的扩散,使得最后的结果不可靠。
因此,小主元可能导致计算的失败,我们应该避免采用绝对值很小的主元素。
为此,我们在高斯消去法的每一步应该在系数矩阵或消元后的低阶矩阵中选取绝对值最大的元素作为主元素,保持乘数1,以便减少计算过程中舍入误差对计算解的影响。
一种方式是完全主元消去法,这种消去法是在每次选主元时,选择a(二k max a j k)0
kJ k k in」
k j n
为主元素。
这种方法是解低阶稠密矩阵方程组的有效方法,但这种方法在选取主元时要花费一定的计算机时间。
实际计算中我们常采用部分选主元的的消去法。
列主元消去法即在每次
选主元时,仅依次按列选取绝对值最大的元素作为主元素,且仅交换两行,再进行消元计算
2实验目的和要求
运用matlab编写一个.m文件,要求用列主元消去法求解方程组(实现PA=L)
111111
2 11111
3 2 1111
4 3 2 1 1 1
5 4 3 2 1 1
6 5 4 3 2 1
7 6 5 4 3 2 要求输出以下内容:
(1) 计算解x;
x.
117 1x28 1x310 1x413 1x517
1
x6
22 128 x 7
(2) L,U ;
(3) 整形数组IP (i ) (i=1,2,…,n-1 )(记录主行信息)
3算法原理与流程图
(1) 算法原理
设有线性方程组Ax=b,其中设A为非奇异矩阵。
方程组的增广矩阵为
an a12L L a1n M bi
a21a22L L a2n M b2
[代口
M
M M M
a i1
M M M M
a n1a n2L L a nn M
b n
第1步(k=1):首先在A的第一列中选取绝对值最大的元素a^,作为第一步的主元素3in max
耳!0,然后交换(A,b)的第1行与第i1行元素,再进行消元计算。
1 1 i n
设列主元素消去法已经完成第1步到第k-1步的按列选主元,交换两行,消元计算得到与原方程组等价的方程组A(k)x b(k)
J1) c ⑴
an ai2
L a(k)a1(n)M
⑵
L (2)⑵
M
a22a2k a2n b2()
O M M M [A,b] [A(k),b(k)]a kk)L a kk)M.(k)bk
a k k)u a k k)1,n M(k)
d 1
M M M M
a nk)L a(k)叫
n
M
第k步计算如下:
对于k=1,2,…,n-1
(1)按列选主元:即确定ik使
时max无0
(2) 如果a ik,k 0,则A为非奇异矩阵,停止计算
(3) 如果ik 丰k ,则交换[A , b]第ik 行与第k 行元素。
(4) 消元计算
(5) 回代求解
(3)
4程序代码及注释
for k=1:n-1
a ik
E ik
电,(i
a kk
1L ,n)
a j
b E k a
kj ,
mb.
(i ,j (i
k 1,L ,n)
k 1,L ,n)
消元乘数m ik 满足: E ik
a ik a kk
1,(i k 1,L ,n)
0 / a nn
(当 a
nn
b
(b i
j b n
n a ij b j ) / a ii ,
I
i 1
0) n 1,K ,2,1
计算解(X 1,L ,X n )在常数项 b(n)内得至叽
(2)
流程图见图1
IP(k)=k;
if i<j U(i,j)=A(i,j);
elseif i==j
L(i,j)=1; U(i,j)=A(i,j);
else
L(i,j)=A(i,j);
end
end
end
end
输入n, A , b,£
det=1.0
k=1,2,…,n-1
输出det( A) 0
按列选主元a k ,k max k i n a ik0
c0= a ik ,k
血
停机
(I
(j
c0|<
Y
ik=k
换行
a kj a ik ,j (j
k, L , n)
b k b ik,det det
消元计算
k+1,L
a ij
k+1,L
det
回代求解
b n
det
(b
a nn
a nn(当a nn0)
nn nn
n
a^ b j)/ a i ,I 1,2, L n
j i 1
* det
输岀计算解及行列式值b
(I)( I=1,2,…,n)及det
,n)
,n)
停机
,a ik
* a kj
,b i
a kk det
a ik / a kk
6
附图1
4.0000
5.0000
6.0000 0
5算例分析
1、测试示例 >> A=[ 1 2 3 4 5 6];
>> b=[3 7 11];
>> [x L U IP P]=gauss(A,b) Error! Please in put aga in! >> A=[ 1 2 3 0 0 0 4 5 6]; >> b=[3 7 11];
>> [x L U IP P]=gauss(A,b)
The equati ons have no unique soluti on! x = NaN -Inf Inf
1.0000
0.2500 1.0000
0 1.0000
0 0.7500 1.5000
IP =
3 3
P =
1 0 0
0 0 1
0 1 0
2、计算过程
(1)首先输入系数矩阵A和矩阵b
>> A=[
1 1 1 1 1 1 1
2 1 1 1 1 1 1
3 2 1 1 1 1 1
4 3 2 1 1 1 1
5 4 3 2 1 1 1
6 5 4 3 2 1 1
7 6 5 4 3 2 1];
>> b=[7 8 10 13 17 22 28];
(2)输出结果
>> [x L U IP P]=gauss(A,b)
x =
1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000
1.000000 0000
0.2857 1.000000 000
0.42860.8000 1.000000 00
0.57140.60000.7500 1.0000000
0.71430.40000.50000.6667 1.000000
0.1429-0.2000-0.2500-0.3333-0.5000 1.00000
0.85710.20000.25000.33330.5000 --1.0000 1.0000
7.0000 6.0000 5.0000 4.0000 3.0000 2.0000 1.0000
0 -0.7143 -0.4286 -0.1429 0.1429 0.4286 0.7143
0 -0.8000 -0.6000 -0.4000 -0.2000 0.0000
0 -0.7500 -0.5000 -0.2500 0.0000
0 -0.6667 -0.3333 -0.0000
0 0.5000 1.0000
0 1.0000
IP =
6讨论与结论
1、时间复杂度:
>> tic;[x L U IP P]=gauss(A,b);toc
Elapsed time is 0.000856 sec on ds.
2、程序优化
初次编程时,没有考虑到给一个变量赋初值的情况。
虽然在MATLAB^变量不赋初值是完全允许的,但是由于一个变量中含有多个元素时,每次改变该数组的长度,便会增加计算机时间。
另外,给程序加上一定的判断条件及报错信息,一定程度上有程序优化的作用。
因此,本程序中的以下程序段都起到了程序优化的作用。
参考文献
[1] 易大义,沈云宝,李有法•计算方法(第2版),浙江大学出版社.p.29-53.
[2] 张琨高思超毕靖编著MATLAB2010从入门到精通电子工业出版社
0 0 0 1 0 0 0
0 0 0 0 0 1 0。