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

<数值计算方法>实验报告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]);}。

《数值计算方法》实验报告专业：年级：学号：姓名：成绩：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 分解法的优点，以及方程组系数矩阵和右端向量的微小变化对解向量的影响。

二、实验内容：解下列两个线性方程组（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编写代码，实现列主元消去法。

首先，我们需要输入一个线性方程组的系数矩阵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```接下来，我们将使用一个具体的例子来说明列主元消去法的求解过程。

数值分析上机实验报告（二）一、问题描述：利用列主元消去法求解下列方程组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课题名称用列主元高斯消去法和列主元三角分解法解线性方程目的和意义高斯消去法是一个古老的求解线性方程组的方法,但由它改进得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效方法;用高斯消去法解线性方程组的基本思想时用矩阵行的初等变换将系数矩阵A 约化为具有简单形式的矩阵上三角矩阵、单位矩阵等,而三角形方程组则可以直接回带求解 用高斯消去法解线性方程组b Ax =其中A ∈Rn ×n 的计算量为：乘除法运算步骤为32(1)(1)(21)(1)(1)262233n n n n n n n n n n nMD n ----+=+++=+-,加减运算步骤为(1)(21)(1)(1)(1)(25)6226n n n n n n n n n n AS -----+=++=;相比之下,传统的克莱姆法则则较为繁琐,如求解20阶线性方程组,克莱姆法则大约要19510⨯次乘法,而用高斯消去法只需要3060次乘除法;在高斯消去法运算的过程中,如果出现absAi,i 等于零或过小的情况,则会导致矩阵元素数量级严重增长和舍入误差的扩散,使得最后的计算结果不可靠,所以目前计算机上常用的解低阶稠密矩阵方程的快速有效的方法时列主元高斯消去法,从而使计算结果更加精确; 2、列主元三角分解法高斯消去法的消去过程,实质上是将A 分解为两个三角矩阵的乘积A=LU,并求解Ly=b 的过程;回带过程就是求解上三角方程组Ux=y;所以在实际的运算中,矩阵L 和U 可以直接计算出,而不需要任何中间步骤,从而在计算过程中将高斯消去法的步骤进行了进一步的简略,大大提高了运算速度,这就是三角分解法采用选主元的方式与列主元高斯消去法一样,也是为了避免除数过小,从而保证了计算的精确度计算公式1、 列主元高斯消去法设有线性方程组Ax=b,其中设A 为非奇异矩阵;方程组的增广矩阵为第1步k=1：首先在A 的第一列中选取绝对值最大的元素1l a ,作为第一步的主元素:111211212222112[,]n n n l n nn n a a a a b a a a b a a a b ⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦a b然后交换A,b 的第1行与第l 行元素,再进行消元计算;设列主元素消去法已经完成第1步到第k -1步的按列选主元,交换两行,消元计算得到与原方程组等价的方程组 Akx=bk第k 步计算如下：对于k=1,2,…,n -11按列选主元：即确定t 使 2如果t ≠k,则交换A,b 第t 行与第k 行元素; 3消元计算消元乘数mik 满足:4回代求解2、 列主元三角分解法 对方程组的增广矩阵 经过k -1步分解后,可变成如下形式：111max 0l i i n a a ≤≤=≠(1)(1)(1)(1)(1)1112111(2)(2)(2)(2)22222()(()1)()()()()()1,1()(,)()[,][,] k k k k nk k nk n k k k k k kk kn k k k k n k k k n nn a a a a b a a a b a a b a b b a a a +++⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥→=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦A b A b ()()max 0k k tk ik k i na a ≤≤=≠,(1,,)ik ik ik kka a m i k n a ←=-=+, (,1,,), (1,,)ij ij ik kji i ik k a a m a i j k n b b m b i k n ←+=+⎧⎨←+=+⎩⎪⎪⎩⎪⎪⎨⎧--=-←←∑+=)1,,2,1(,)(1n n i a x a b x a b x ii n i j j ij i i nnn n [,]A A b =11121,11111222,122221,11,1,1,211,11,2121,112,112,1k k k k k k k j n k k j n k k k i i i k n n kk kj kn k ik ij in i nknjk k k j k n n nnk k n a a a b A a u u u u u u y l l l l l l ll l l l u u u u u y u u u u y a a b a a b l a -------------⎡→⎣⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎦第k 步分解,为了避免用绝对值很小的数kku 作除数,引进量1111 (,1,,;1,2,,) ()/ (1,2,,;1,2,,)k kj kj km mj m k ik ik im mk kkm u a l u j k k n k n l a l u u i k k n k n -=-=⎧=-=+=⎪⎪⎨⎪=-=++=⎪⎩∑∑11(,1,,)k i ik im mk m s a l u i k k n -==-=+∑,于是有kk u =ks ;如果 ,则将矩阵的第t 行与第k 行元素互换,将i,j 位置的新元素仍记为jjl 或jja ,然后再做第k 步分解,这时列主元高斯消去法程序流程图max t ik i n s s ≤≤= ()/ 1,2,,)1 (1,2,,),kk k k t iki k ik u s s s l s s i k k n l i k k n ===++≤=++即交换前的，（且列主元高斯消去法Matlab主程序function x=gauss1A,b,c %列主元法高斯消去法解线性方程Ax=bif lengthA~=lengthb %判断输入的方程组是否有误disp'输入方程有误'return;enddisp'原方程为AX=b：' %显示方程组Abdisp'------------------------'n=lengthA;for k=1:n-1 %找列主元p,q=maxabsAk:n,k; %找出第k列中的最大值,其下标为p,qq=q+k-1; %q在Ak:n,k中的行号转换为在A中的行号if absp<cdisp'列元素太小,detA≈0';break;elseif q>ktemp1=Ak,:; %列主元所在行不是当前行,将当前行与列主Ak,:=Aq,:; 元所在行交换包括bAq,:=temp1;temp2=bk,:;bk,:=bq,:;bq,:=temp2;end%消元for i=k+1:nmi,k=Ai,k/Ak,k; %Ak,k将Ai,k消为0所乘系数Ai,k:n=Ai,k:n-mi,kAk,k:n; %第i行消元处理bi=bi-mi,kbk; %b消元处理endenddisp'消元后所得到的上三角阵是'A %显示消元后的系数矩阵bn=bn/An,n; %回代求解for i=n-1:-1:1bi=bi-sumAi,i+1:nbi+1:n/Ai,i;endclear x;disp'AX=b的解x是' x=b;调用函数解题列主元三角分解法程序流程图列主元三角分解法Matlab主程序①自己编的程序：function x=PLUA,b,eps %定义函数列主元三角分解法函数if lengthA~=lengthb %判断输入的方程组是否有误disp'输入方程有误'return;enddisp'原方程为AX=b：' %显示方程组Abdisp'------------------------'n=lengthA;A=A b; %将A与b合并,得到增广矩阵for r=1:nif r==1for i=1:nc d=maxabsA:,1; %选取最大列向量,并做行交换if c<=eps %最大值小于e,主元太小,程序结束break;elseendd=d+1-1;p=A1,:;A1,:=Ad,:;Ad,:=p;A1,i=A1,i;endA1,2:n=A1,2:n;A2:n,1=A2:n,1/A1,1; %求u1,ielseur,r=Ar,r-Ar,1:r-1A1:r-1,r; %按照方程求取ur,iif absur,r<=eps %如果ur,r小于e,则交换行p=Ar,:;Ar,:=Ar+1,:;Ar+1,:=p;elseendfor i=r:nAr,i=Ar,i-Ar,1:r-1A1:r-1,i; %根据公式求解,并把结果存在矩阵A中endfor i=r+1:nAi,r=Ai,r-Ai,1:r-1A1:r-1,r/Ar,r; %根据公式求解,并把结果存在矩阵A中endendendy1=A1,n+1;for i=2:nh=0;for k=1:i-1h=h+Ai,kyk;endyi=Ai,n+1-h; %根据公式求解yiendxn=yn/An,n;for i=n-1:-1:1h=0;for k=i+1:nh=h+Ai,kxk;endxi=yi-h/Ai,i; %根据公式求解xiendAdisp'AX=b的解x是'x=x'; %输出方程的解②可直接得到P,L,U并解出方程解的的程序查阅资料得子函数PLU1,其作用是将矩阵A分解成L乘以U的形式;PLU2为调用PLU1解题的程序,是自己编的Ⅰ.function l,u,p=PLU1A %定义子函数,其功能为列主元三角分解系数矩阵A m,n=sizeA; %判断系数矩阵是否为方阵if m~=nerror'矩阵不是方阵'returnendif detA==0 %判断系数矩阵能否被三角分解error'矩阵不能被三角分解'endu=A;p=eyem;l=eyem; %将系数矩阵三角分解,分别求出P,L,Ufor i=1:mfor j=i:mtj=uj,i;for k=1:i-1tj=tj-uj,kuk,i;endenda=i;b=absti;for j=i+1:mif b<abstjb=abstj;a=j;endendif a~=ifor j=1:mc=ui,j;ui,j=ua,j;ua,j=c;endfor j=1:mc=pi,j;pi,j=pa,j;pa,j=c;endc=ta;ta=ti;ti=c;endui,i=ti;for j=i+1:muj,i=tj/ti;endfor j=i+1:mfor k=1:i-1ui,j=ui,j-ui,kuk,j;endendendl=trilu,-1+eyem;u=triuu,0Ⅱ.function x=PLU2A,b %定义列主元三角分解法的函数l,u,p=PLU1A %调用PLU分解系数矩阵A m=lengthA; %由于A左乘p,故b也要左乘p v=b;for q=1:mbq=sumpq,1:mv1:m,1;endb1=b1 %求解方程Ly=b for i=2:1:mbi=bi-sumli,1:i-1b1:i-1;endbm=bm/um,m; %求解方程Ux=y for i=m-1:-1:1bi=bi-sumui,i+1:mbi+1:m/ui,i;endclear x;disp'AX=b的解x是' x=b;调用函数解题①②编程疑难这是第一次用matlab编程,对matlab的语句还不是非常熟悉,因此在编程过程中,出现了许多错误提示;并且此次编程的两种方法对矩阵的运算也比较复杂;问题主要集中在循环控制中,循环次数多了一次或者缺少了一次,导致数据错误,一些基本的编程语句在语法上也会由于生疏而产生许多问题,但是语句的错误由于系统会提示,比较容易进行修改,数据计算过程中的一些逻辑错误,比如循环变量的控制,这些系统不会提示错误,需要我们细心去发现错误,不断修正,调试;。

《数值分析》实验报告实验序号：实验二 实验名称: 列主元消去法解方程组 学号: 姓名:任课教师: 专业班级:）1、 实验目的：用列主元Gauss 消元法解n 阶线性代数方程组：⎪⎪⎩⎪⎪⎨⎧=+⋯++⋯⋯⋯⋯⋯=+⋯++=+⋯++nn nn 2n21n12n 2n 2221211n 1n 212111b a a a b a a a b a a a x x x x x x x x x 其基本做法是把上述方程组通过列主元Gauss 消元转化为一个等价的三角形方程组，然后再进行回代就可以求出方程组的解。

列主元消元的基本做法是选取系数矩阵的每一列中绝对值最大的作为主元，然后采取和顺序Gauss 消元法相同的步骤进行 ，求得方程组的解。

要求显示出每一个列主元以及每一大步消元后的系数矩阵),...,2,1(n k =(k )A 和常数项),...,2,1(n k =(k )b ,最后显示出方程组的解),...,2,1(n i x i =。

2、 实验内容：(1)实验分析：1. 列主元Gauss 消元法的算法思想：1． 输入增广矩阵B ；。

2． 对k =1,2,…,n ，循环：（a ） 按列选主元||:max ik ni j a a ≤≤=保存主元所在行的指标k i 。

（b ） 若a=0，则系数矩阵奇异，计算停止；否则，顺序进行。

（c ） 若k i =k 则转向（d ）；否则换行ki kj i b b nj a a k j k ↔=↔,...,2,1 ,（d ） 计算乘子.,...,1,/n k i a a a m ik kk ik ik +=⇒=（e ） 消元： nk i b m b b nk j i a m a a k ik i i kj ij ij ij ,...,1;:,...,1,;:+=-=+=-=3. 回代 1,...,1, ,/:1-=⎪⎪⎭⎫ ⎝⎛-=∑+=n n i a b a b b ii n i j j ij i i 用右端项b 来存放解x 。

贵州师范大学数学与计算机科学学院学生实验报告课程名称： 数值分析 班级： 实验日期： 年 月 日 学 号： 姓名： 指导教师： 实验成绩：一、实验名称实验五:线性方程组的数值解法二、实验目的及要求1. 让学生掌握用列主元gauss 消去法、超松弛迭代法求解线性方程组.2. 培养Matlab 编程与上机调试能力.三、实验环境每人一台计算机，要求安装Windows XP 操作系统，Microsoft office2003、MATLAB6.5(或7.0).四、实验内容1. 编制逐次超松弛迭代(SOR 迭代)函数(子程序),并用于求解方程组⎪⎪⎩⎪⎪⎨⎧=-++=+-+=++-=+++-141414144321432143214321x x x x x x x x x x x x x x x x取初始向量T x )1,1,1,1()0(=,迭代控制条件为 5)1()(1021||||--⨯≤-k k x x 请绘制出迭代次数与松弛因子关系的函数曲线,给出最佳松弛因子.SOR 迭代的收敛速度是否一定比Gauss-Seidel 迭代快?2. 编制列主元 Gauss 消去法函数(子程序),并用于解 ⎪⎩⎪⎨⎧=++-=-+-=+-615318153312321321321x x x x x x x x x要求输出方程组的解和消元后的增广矩阵. 注：题2必须写实验报告五、算法描述及实验步骤Gauss 消去法：功能 解方程组b Ax = .输入 n ,n n ij a A ⨯=)(,T n b b b b ),,,(21 =.输出 方程组的解T n x x x x ),,,(21 =或失败信息.步1 对1,,2,1-=n k 执行步2→步4 .步2 调选列主元模块 .步3 若0=kk a ，则=x “消去法失败”，结束 .步4 对n k k i ,,2,1 ++=执行步5→步6 .步5 对n k k j ,,2,1 ++=执行ij kj kk ik ij a a a a a +⨯-⇐/ .步6 i k kk ik i b b a a b +⨯-⇐/ .步7 nn n n a b x /⇐ .步8 对1,,2,1 --=n n i 执行ii n i j j ij i i a x a b x /)(1∑+=-⇐ .步9 输出T n x x x x ),,,(21 = .选列主元模块：功能 选列主元 .输入 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .输出 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .步1 kk a m ⇐；k l ⇐ .步2 对n k k i ,,2,1 ++=执行若m a ik >则ik a m ⇐；i l ⇐ .步3 若k l ≠，则交换kj a 和lj a ，n k k j ,,1, +=；交换k b 和l b .步4 返回主模块 .六、调试过程及实验结果>> A=[12,-3,3;-18,3,-1;1,1,1];>> b=[15;-15;6];>> x=Gauss1(A,b)Ab =-18.0000 3.0000 -1.0000 -15.00000 1.1667 0.9444 5.16670 0 3.1429 9.4286 index = 1x = 1.0000 2.0000 3.0000七、总结由于数)1(-k kka 在Gauss 消去法中有着突出的作用，第k 步消元时，要用)1(-k kk a 作除数，如果)1(-k kk a =0消元会失败，即使主元)1(-k kk a ≠0，但很小时，舍入误差也会使计算结果面目全非，避免这种缺陷的基本方法就是选主元。