高斯消元法_实验报告

合集下载

C++实验报告高斯消元法

高斯肖元法C++上机实验报告学生姓名：学号：专业班级：实验类型：综合一实验项目名称全选主元高斯消去法解线性方程组二实验原理设有n 元线性方程组(考虑便于C++程序数组表示,方程的下标从0开始),0000110,1100000110,111101,111,111n n n n n n n n n n a x a x a x b a x a x a x b a x a x a x b ---------+++=⎧⎪+++=⎪⎨⎪⎪+++=⎩写为矩阵形式为Ax=b,其中A 为线性方程组的系数矩阵,x 为列向量,是方程组的解,b 也是列向量.一般来讲，可以假定矩阵A 是非奇异阵。

（n 阶矩阵A 的行列式不为零，即 |A|≠0，则称A 为非奇异矩阵）00010,10111,1,01,11,1n n n n n n a a a a a a A a a a ----⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦，011n x x x x -⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦ ，011n b b b b -⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦将系数矩阵A 和向量b 放在一起，形成增广矩阵B ：00010,010111,11,01,11,11(,)n n n n n n n a a a b a a a b b A b a a a b -----⎡⎤⎢⎥⎢⎥==⎢⎥⎢⎥⎣⎦全选主元消去就在矩阵B 上进行，整个过程分为如下两个步骤：第一步：消去过程。

对于k 从0开始到n-2结束，进行以下三步。

1. 首先，从系数矩阵A 的k 行k 列开始的子矩阵中选取绝对值最大的元素作为主元素。

例如：11,max 0i j ij k i n k j na a ≤<≤<=≠然后交换B的第k行与第1i行，第k列与第1k列，这样，这个子矩阵中具有最大绝对值的元素被交换到k行k列的位置上.2.其次，进行归一化计算。

计算方法为：/,1,,1/kj kj kkk k kka a a j k nb b a==+-⎧⎪⎨=⎪⎩3.最后进行消去计算：,,1,,1,1,,1 ij ij ik kji i ik ka a a a j i k nb b a b i k n=-=+-⎧⎪⎨=-=+-⎪⎩第二步，回带过程：111,111/,2,,1,0 n n n nni i ij jj ix b ax b a x i n-----=+=⎧⎪⎨=-=-⎪⎩∑三代码的实现整个程序分为5个独立文件，Matrix.h文件中包括矩阵类Matrix的定义，Matrix.cpp文件中包括该类成员函数的实现，LinearEqu.h文件中包括线性方程组类LinearEqu的定义，LinearEqu.cpp文件中包括该类的成员函数实现文件；7-9.cpp文件包括程序的主函数，主函数中定义了一个类LinearEqu的对象，通过这个对象求解一个四元线性方程组。

解线性方程组的列主元素高斯消去法和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五、实验结果分析：实验的数学原理很容易理解，也容易上手。

高斯消去法的数值稳定性实验报告

end
end
%½«³£ÊýÏî´ÓAÖÐ·ÖÀë
X=A(:,length(A(1,:)));
A=A(:,1:length(A(1,:))-1);
%µÃµ½LÓëU¾ØÕó
fori=1:length(A(:,1))
forj=1:length(A(i,:))
if(i==j)
L(i,j)=1;
U(i,j)=A(i,j);
00-2.83542034035461-0.439277018213440
0000.724838760535495
X2=
5.90739256961544
1.88284104193886
-22.2540722833854
14.5809976298923
(3)beselect.m文件与noselect.m文件
%½«³£ÊýÏî¼Ó½øÈ¥
A(:,length(A(1,:))+1)=Y;
%´ÓÉÏÍùÏÂÑ »·
fori=1:length(A(:,1))
%Ô¤¼ÆËãÃ¿ÐÐµÚÒ»¸öÊý×Ö
fork=i:length(A(:,i))
form=1:i-1
A(k,i) = A(k,i)-A(k,m)*A(m,i);
end
end
end
%¼ÆËãXÖµ
fori=length(X):-1:1
forj=i+1:length(U(i,:))
X(i)=X(i)-U(i,j)*X(j);
end
X(i)=X(i)/U(i,i);
end
end
noselect.m文件：
function[ L,U,X ] = noselect( A,Y )

guess消元法

for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
cin>>a[i][j];
}
}
for(k=1;k<=n-1;k++)
{
for(i=k+1;i<=n;i++)
{
for(j=k+1;j<=n+1;j++)
{
a[i][j]=a[i][j]-(a[i][k]/a[k][k])*a[k][j];
《数值分析》实验报告
实验序号：实验一题目名称:Gauss消元法解方程组
学号:XXXXXXX姓名:XXX
任课教师:XXXX专业班级:算机科学与技术（师范）
1、实验目的：
编写用Gauss消元法解线形方程组的程序
2、实验分析：
计算方法分析：
Gauss消元法的基本做法就是把方程组转化成为一个如下图的等价的三角方程组，这个过程叫做消元。得到三角方程组后，就可以逐个求出Xn，Xn-1，…，X1，这个过程叫回代。
程序代码分析：
建立两个数组a和b，通过循环语句将ｎ阶增广矩阵输入进去，通过对列的循环对每一列进行消去未知数，通过ｎ小步ｎ大步把矩阵化简成上三角形矩阵，最后通过迭代法解得方程组得解。
3、函数分析：
程序中只用到一个主函数，求解线形方程组得算法都放在主函数中，
利用以下函数进行求解：
a[i][j]=a[i][j]-a[i][k]/a[k][k]*a[k][j];
b[i]=b[i]-a[i][k]/a[k][k]*b[k];
迭代：
for(i=n-1;i>0;i--)

实验三编程实现列主元高斯消去法

实验三编程实现列主元高斯消去法1.实验目的：实现高斯主消元法，对计算过程加深理解。

2.实验内容：编写c++程序，实现对角占优方程组的编程求解。

3.实验步骤1、设计方程组的存储为二位数组，最大方程组数为100，第i行j列的元素值代表第i个方程的第j个系数，输入时没有的系数项填0。

2 对于每个方程按主消元法从0~n-1依次消元。

3迭代求解，对于第i个未知数的值，依次迭代i+1~n-1已求出的结果4．实验结果分析：对于除数很小的情况程序不能很好的解决，对于不是对角占优的也不能搜索出主元素，以后在进一步解决上述问题。

#include<iostream>using namespace std;const int MAX=100;double str[MAX][MAX];int main(){while(1){cout<<"输入方程组个数0<n<100"<<endl;int n,m,i,j,k,flag=0;cin>>n;if(n<=0||n>=100)break;cout<<"输入方程组的系数"<<endl;for(i=0;i<n;i++){for(j=0;j<=n;j++){cin>>str[i][j];if(i==j&&str[i][j]==0)flag=1;}}if(flag==1){cout<<"方程组不是对角占优的，此程序不能解决"<<endl;system("pause");return 0;}for(i=0;i<n;i++){for(j=n;j>=i;j--)//对角变一str[i][j]/=str[i][i];for(k=i+1;k<n;k++)//方程消元{double tem=str[k][i]/str[i][i];for(j=i;j<=n;j++){str[k][j]-=str[i][j]*tem;}}}double ans[MAX]={0};ans[n-1]=str[n-1][n];for(i=n-2;i>=0;i--){ans[i]=str[i][n];for(j=n-1;j>i;j--)ans[i]-=str[i][j]*ans[j];}for(i=0;i<n;i++)cout<<"X"<<i+1<<" = "<<ans[i]<<endl;system("pause");}}。

计算方法实验报告(1)----Gauss消去法

计算方法实验报告实验名称：实验（一）Gauss 消去法班级：学生姓名：学号：班级序号：课内序号:指导老师：2018-2019学年第2学期一、实验名称：Gauss消去法二、实验学时: 2学时三、实验目的和要求1、掌握高斯消去法基础原理2、掌握高斯消去法解方程组的步骤3、能用程序语言对Gauss消去法进行编程实现四、实验过程代码及结果1、代码：using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace ConsoleApplication_Gauss{class Program{//回带求值的过程static void CalcX(double[,] a, double[] x, int n){for (int i = n - 1; i >= 0; i--){double sum = 0;for (int j = i + 1; j < n; j++){sum += a[i, j] * x[j];}x[i] = (a[i, n] - sum) / a[i, i];}}//消元的过程static void CalcA(double[,] a, int n){for (int k = 0; k < n - 1; k++){for (int i = k + 1; i < n; i++){//double Lik = a[i, k] / a[k, k];// for (int j = k ; j <= n; j++)for (int j = n; j >= k; j--){a[i, j] = a[i, j] - a[i, k] / a[k, k] * a[k, j];}//a[i, k] = 0;}//Output}}//输出未知数x的值static void Output(double[] x, int n){for (int i = 0; i < n; i++){Console.WriteLine("x[{0}]={1}", i, x[i]);}}static void Output(double[,] a, int n){for (int i = 0; i < n; i++){//string s="";for (int j = 0; j <= n; j++){//s += string.Format("{0,-4}", a[i, j]);Console.Write("{0,6}", a[i, j]);}Console.WriteLine();}}//输入函数，表示输入一串值作为方程组的系数static void Input(double[,] a, int n){for (int i = 0; i <= n - 1; i++){string s = Console.ReadLine();string[] ss = s.Split(' ');for (int j = 0; j <= n; j++){a[i, j] = Convert.ToDouble(ss[j]);}}}static void Main(string[] args){Console.WriteLine("请输入矩阵的维数:");int n =Convert.ToInt32( Console.ReadLine());double[,] a = new double[n,n+1];Console.WriteLine("请输入矩阵的各个元;");Input(a, n);Console.WriteLine("------A（i,j)----------");Output(a, n);CalcA(a, n);Console.WriteLine("------消元之后A（i,j)----------");Output(a, n);double[] x = new double[n];CalcX(a, x, n);Output(x, n);Console.ReadLine();}}}2、结果：…。

数值方法高斯列主元消元法试验报告

1）用高斯列主元消元法求解下面的方程组#include <iostream>#include <cmath>#include <iomanip>using namespace std;int main(){//double a[4][5]={1,-1,1,-4,2,5,-4,3,12,4,2, 1,1,11,3,2,-1,7,-1,0};//int i,j,k;int Line,Row;double temp[4][5];//中间量Row=4;Line=5;//Inintial//Result//double X[5];cout<<"最初的矩阵为:"<<endl;for(i=0;i<4;i++){for(j=0;j<5;j++){cout<<a[i][j]<<' ';}cout<<endl;}cout<<endl;/////////////////////////////////for(k=0;k<Row-1;k++){//for(i=k+1;i<Row;i++){//temp[i][k]=a[i][k]/a[k][k];for(j=k;j<Line;j++){a[i][j]-=temp[i][k]*a[k][j];}}}//////////////////////////////////cout<<"经过消元后的增广矩阵为:"<<endl;for(i=0;i<4;i++){for(j=0;j<5;j++){cout<<a[i][j]<<' ';}cout<<endl;}cout<<endl;/////////////////////////////////for(i=Row-1;i>=0;i--){//double temp_new;temp_new=0;for(j=i+1;j<=Row-1;j++){a[i][Row]-=a[i][j]*a[j][Row];}a[i][Row]/=a[i][i];}///////////////////////////////////Print;cout<<"最后的解为:"<<endl;for(i=0;i<Row;i++){//cout<<"X"<<i+1<<"="<<a[i][Row]<<endl;}/////////////////////////////////return 0;}1.2列主元方法#include <iostream>#include <cmath>#include <iomanip>using namespace std;void Print(double a[][5]){int i,j;for(i=0;i<4;i++){for(j=0;j<5;j++){cout<<a[i][j]<<' ';}cout<<endl;}cout<<endl;}int main(){//double a[4][5]={1,-1,1,-4,2,5,-4,3,12,4,2, 1,1,11,3,2,-1,7,-1,0};//double b[4][5]={0.3e-15,int i,j,k,n;int Line,Row;double temp[4][5];//中间量Row=4;Line=5;//Inintial//////////////////////////////////////////////cout<<"最初的矩阵为:"<<endl;Print(a);////////////////////////////////////////////////the main process is underint kk;//flagsdouble max;//bool flag=false;double t;//the temp of change;for(k=0;k<Row-1;k++){/////////////////search the max_num//flag=false;max=a[k][k];kk=k;for(i=k;i<Row;i++){if(abs(a[i][k])>max){max=a[i][k];kk=i;//flag=true;}}//////////////////change the linefor(j=0;j<Line;j++){t=a[k][j];a[k][j]=a[kk][j];a[kk][j]=t;}cout<<"第"<<k+1<<"次换行结果:"<<endl;Print(a);cout<<endl;cout<<"第"<<k+1<<"次消元结果:"<<endl;//////////////////消元的过程for(i=k+1;i<Row;i++){//temp[i][k]=a[i][k]/a[k][k];for(j=k;j<Line;j++){a[i][j]-=temp[i][k]*a[k][j];}}///////////////////Print(a);//}//回带的过程n=Row-1;for(i=n;i>=0;i--){//double temp_new;temp_new=0;for(j=i+1;j<=n;j++){a[i][n+1]-=a[i][j]*a[j][n+1];}a[i][n+1]/=a[i][i];}/////////////////////////////////////////////cout<<"经过消元后的增广矩阵为:"<<endl;Print(a);////////////////////////////////////////////////Print;cout<<"最后的解为:"<<endl;for(i=0;i<Row;i++){//cout<<"X"<<i+1<<"="<<a[i][Row]<<endl;}/////////////////////////////////return 0;}2）分别用列主元消元法与不选主元消元法求解，分析对结果的影响#include <iostream>#include <cmath>#include <iomanip>using namespace std;void Print(double b[][5]){int i,j;for(i=0;i<4;i++){for(j=0;j<5;j++){cout<<b[i][j]<<' ';}cout<<endl;}cout<<endl;}int main(){double b[4][5]={0.3e-15,59.14,3,1,59.17,5.291,-6.130,-1,2,46.78,11.2,9,5,2,1,1,2,1,1,2};int i,j,k,n;int Line,Row;double temp[4][5];//中间量Row=4;Line=5;//Inintial//////////////////////////////////////////////cout<<"最初的矩阵为:"<<endl;Print(b);////////////////////////////////////////////////the main process is underint kk;//flagsdouble max;//bool flag=false;double t;//the temp of change;for(k=0;k<Row-1;k++){/////////////////search the max_num//flag=false;max=b[k][k];kk=k;for(i=k;i<Row;i++){if(abs(b[i][k])>max){max=b[i][k];kk=i;//flag=true;}}//////////////////change the linefor(j=0;j<Line;j++){t=b[k][j];b[k][j]=b[kk][j];b[kk][j]=t;}cout<<"第"<<k+1<<"次换行结果:"<<endl;Print(b);cout<<endl;cout<<"第"<<k+1<<"次消元结果:"<<endl;//////////////////消元的过程for(i=k+1;i<Row;i++){//temp[i][k]=b[i][k]/b[k][k];for(j=k;j<Line;j++){b[i][j]-=temp[i][k]*b[k][j];}}///////////////////Print(b);//}//回带的过程n=Row-1;for(i=n;i>=0;i--){//double temp_new;temp_new=0;for(j=i+1;j<=n;j++){b[i][n+1]-=b[i][j]*b[j][n+1];}b[i][n+1]/=b[i][i];}/////////////////////////////////////////////cout<<"经过消元后的增广矩阵为:"<<endl; Print(b);////////////////////////////////////////////////Print;cout<<"最后的解为:"<<endl;for(i=0;i<Row;i++){//cout<<"X"<<i+1<<"="<<b[i][Row]<<endl;}/////////////////////////////////return 0;}Ax （迭代法收敛速度实验）注意修改不同的A、B的数组2、用迭代法求解;b3、//雅可比迭代法/*@auther luozhengxiao*/#include <iostream>#include <cmath>#include <iomanip>using namespace std;//void Print(double x[]){for(int i=0;i<3;i++){cout<<setprecision(8)<<fixed<<x[i]<<endl;}}int main(){//double B1[3]={-3,2,4};double B2[3]={100,-200,345};double A[3][3]={6,2,-1,1,4,-2,-3,1,4};double x[3],x_old[3],temp;int i,j,k;for(i=0;i<3;i++){cout<<"请输入第"<<i+1<<"个数:";cout<<"\t x["<<i<<"]=";cin>>x[i];//x_old[i]=x[i];}int n;cout<<"请输入迭代次数：";cin>>n;/////////////////////////////////for(k=0;k<n;k++){//for(i=0;i<3;i++){//temp=0;j=0;while(j<3){if(j==i) {j++;continue;}temp+=A[i][j]*x[j];j++;}x_old[i]=B1[i]-temp;x_old[i]/=A[i][i];}for(j=0;j<3;j++){x[j]=x_old[j];}}Print(x);return 0;}。

消元法实验报告12

西京学院数学软件实验任务书课程名称数学软件实验班级*** 学号*** 姓名***实验课题线性方程组高斯消去法，高斯列主元消去法，高斯全主元消去法实验目的熟悉线性代数方程组高斯消去法，高斯列主元消去法，高斯全主元消去法实验要求运用Matlab/C/C++/Java/Maple/Mathematica等其中一种语言完成实验内容线性方程组高斯消去法线性方程组高斯列主元消去法线性方程组高斯全主元消去法成绩教师1. 实验目的掌握高斯消去法的基本思路和迭代步骤；熟悉线性代数方程组高斯消去法，高斯列主元消去法，高斯全主元消去法；培养编程与上机调试能力。

2. 算法描述注：本实验以3行4列的增广矩阵为例1. 高斯消去法基本思路设有方程组A*x=b，设A是可逆矩阵。

高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵，将其中的A 变换成一个上三角矩阵，然后求解这个三角形方程组。

2. 高斯顺序消去法计算步骤将方程组用增广矩阵B={A:b}表示。

1.消元过程(1) a[0][0]!=0.(2)如果a[0][0]=0，则矩阵A奇异，程序结束。

(3)消元每一行都先与第一行消元通式为：a[i][j]=a[i][j]+a[k][j]*（-a[i][k]/a[k][k]）2. 回代过程(1) 若a[k][k]=0,则矩阵奇异，方程组解不唯一，程序结束；(2) 从下往上一步步回代通式为：a[i][3]=a[i][3]-a[i][j]*x[j]x[i]=a[i][3]/a[i][i]3.高斯列主元消去法计算步骤将方程组用增广矩阵B={A:b}表示1.第i次选出i列中最大的行与第i行交换循环同时进行顺序消元过程2.回代过程与顺序法相同4.高斯全主元消去法计算步骤将方程组用增广矩阵B={A:b}表示1.找出所有未知量系数的最大元素记下最大元素所在的行与列2.将最大元素所在的行换到第i行3.将最大元素所在的列换到第i列4.记下列的变换5.回代过程与顺序法相同6.将列变换交换回来7.输出结果3 实验内容解方程组x1+x2+x3=6x2-x3=52x1-2x2+x34 实验步骤C语言代码1.高斯顺序消元法：#include"stdio.h"void main(){int i,j,k,s,x[3],a[3][4];//input matrixprintf("请注意输入的增广矩阵A为3行4列\n");for(i=0;i<3;i++){printf("第%d行\n",i+1);for(j=0;j<4;j++){// printf("%d :",j+1);scanf("%d",&a[i][j]);}// printf("\n");}//outputprintf("\n线性方程组的增广矩阵A为：\n");for(i=0;i<3;i++){for(j=0;j<4;j++)printf("%-5d",a[i][j]);printf("\n");}// gsfor(k=0;k<2;k++){for(i=k+1;i<3;i++){s=-a[i][k]/a[k][k];for(j=k;j<4;j++){a[i][j]=a[i][j]+a[k][j]*s;}}}printf("\n");//outprintf("\n线性方程组的增广矩阵经过高斯消元得到的矩阵为：\n");for(i=0;i<3;i++){for(j=0;j<4;j++)printf("%-5d",a[i][j]);printf("\n");}//solutionprintf("\n线性方程组的解为：\n");x[2]=a[2][3]/a[2][2];for(i=3-2;i>=0;i--){if(a[i][i]!=0){for(j=i+1;j<3;j++){a[i][3]=a[i][3]-a[i][j]*x[j];}x[i]=a[i][3]/a[i][i];}else printf("\n方程组解不唯一");}for(i=1;i<=3;i++){printf("x%d=%-5d",i,x[i-1]);}printf("\n");}2.高斯列主元消去法：#include"stdio.h"void main(){int i,j,k,s,p,l,m,x[3],a[3][4];//input matrixprintf("请注意输入的增广矩阵A为3行4列\n");for(i=0;i<3;i++){printf("第%d行\n",i+1);for(j=0;j<4;j++){// printf("%d :",j+1);scanf("%d",&a[i][j]);}// printf("\n");}//outputprintf("\n线性方程组的增广矩阵A为：\n");for(i=0;i<3;i++){for(j=0;j<4;j++)printf("%-5d",a[i][j]);printf("\n");}// gsfor(k=0;k<2;k++){for(l=k+1;l<3;l++)if(a[k][k]<a[l][k])for(m=k;m<4;m++){p=a[k][m];a[k][m]=a[l][m];a[l][m]=a[k][m];}for(i=k+1;i<3;i++){s=-a[i][k]/a[k][k];for(j=k;j<4;j++){a[i][j]=a[i][j]+a[k][j]*s;}}}printf("\n");//outprintf("\n线性方程组的增广矩阵经过高斯消元得到的矩阵为：\n");for(i=0;i<3;i++){for(j=0;j<4;j++)printf("%-5d",a[i][j]);printf("\n");}//solutionprintf("\n线性方程组的解为：\n");x[2]=a[2][3]/a[2][2];for(i=3-2;i>=0;i--){if(a[i][i]!=0){for(j=i+1;j<3;j++){a[i][3]=a[i][3]-a[i][j]*x[j];}x[i]=a[i][3]/a[i][i];}else printf("\n方程组解不唯一");}for(i=1;i<=3;i++){printf("x%d=%-5d",i,x[i-1]);}printf("\n");}3.高斯全主元消去法：#include"stdio.h"#include"math.h"void main(){int i,j,i1,j1,k,q,l,m,p[3]={0,1,2};float s,x[3],a[3][4],max;printf("请注意输入的增广矩阵A为3行4列\n");for(i=0;i<3;i++){printf("第%d行\n",i+1);for(j=0;j<4;j++){scanf("%f",&a[i][j]);}}printf("\n");for(k=0;k<2;k++){max=abs(a[k][k]);for(i=k;i<3;i++){for(j=k;j<3;j++){if(abs(a[i][j])>max){max=a[i][j];i1=i;j1=j;}}}for(m=k;m<4;m++){q=a[k][m];a[k][m]=a[i1][m];a[i1][m]=q;}for(l=k;l<3;l++){q=a[l][k];a[l][k]=a[l][j1];a[l][j1]=q;}printf("\n");q=p[k];p[k]=p[j1];p[j1]=q;for(i=k+1;i<3;i++){s=-a[i][k]/a[k][k];for(j=k;j<4;j++){a[i][j]=a[i][j]+a[k][j]*s;}}}printf("\n");printf("\n线性方程组的增广矩阵经过高斯全主元消元得到的矩阵为：\n");for(i=0;i<3;i++){for(j=0;j<4;j++)printf("%-10f",a[i][j]);printf("\n");}printf("\n线性方程组的解为：\n");x[2]=a[2][3]/a[2][2];for(i=3-2;i>=0;i--){if(a[i][i]!=0){for(j=i+1;j<3;j++){a[i][3]=a[i][3]-a[i][j]*x[j];}x[i]=a[i][3]/a[i][i];}else printf("\n方程组解不唯一");}for(i=0;i<3;i++){printf("x%d=%-10f",p[i]+1,x[i]);}printf("\n");}2.实验结果：1.高斯顺序消元法：2.高斯列主元消去法：3.高斯全主元消去法：5.实验总结：通过本次实验再次熟悉了高斯主元消元法的思想，加深了对C语言的理解，简洁明了，学会了C语言编写的函数通用实用。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

华中科技大学数值分析实验报告系、年级研究生院2012级学号姓名类别硕士2013年5月6日实验6.1实验要求：根据教材实验6.1做出相应改编：分别使用Gauss 消元、列选主元。

全选主元的方法求解线性方程组，分别比较三种消元方法的结果和算法的区别，并说明主元的选取在Gauss 消元的中的作用。

问题提出：Gauss 消去法是我们在线性代数中已经熟悉的。

但由于计算机的数值运算是在一个有限的浮点数集合上进行的，如何才能确保Gauss 消去法作为数值算法的稳定性呢？Gauss 消去法从理论算法到数值算法，其关键是主元的选择。

主元的选择从数学理论上看起来平凡，它却是数值分析中十分典型的问题。

一般来说书本上采用的列选主元的办法对其线性方程组进行求解的，那么我们是否可以选择一种行列都选取主元消去的办法来减小相应的误差呢？全主元消元法和列主元消元法一样都是由高斯消元法演变而来。

只不过选取主元的范围有所加大。

全选主元相对于列选主元的更加复杂化了，因为在运算的过程中导致了元的位置发生了变化，这样我们就不得不追踪每个元的位置。

本次实验就几个问题进行了matlab 实验分析，比较几种计算方法的优劣性。

实验内容：考虑线性方程组n n n R b R A b Ax ∈∈=⨯,,编制一个程序：分别能进行Gauss 消去、列选主元Gauss 消去、全选主元Gauss 消去法进行解线性方程组。

对三种算法所得到的结果进行比较，分析三种计算方法的准确性。

具体内容：（1）取矩阵⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=1415157,6816816816 b A ，则方程有解T x )1,,1,1(* =。

取n=10、n=20计算矩阵的条件数。

分别编写利用matlab 编写运算程序，实现Gauss 消去、列选主元消去以及全选主元消去的方法。

比较三种计算方法的运算结果。

在列选主元的过程中分别采用每步消去过程总选取按模最小或按模尽可能小的元素作为主元或每步消去过程总选取按模最大的元素作为主元，结果又如何？分析实验的结果。

1 采用普通Gauss 消元法进行计算Gauss 消去法的基本思想是，通过将一个方程乘或除某个数以及两个方程相加减这两种运算手续，逐步减少方程组中变元的数目，最终使某个方程只含有一个变元，从而得出所求的解。

对于Ax =b ，G auss 消去法的求解思路为：（1）若a 11(1)≠0，先让第一个方程组保持不变，利用它消去其余方程组中的x 1,使之变成一个关于变元x 2,x 3……x n 的n -1阶方程组。

（2）按照（1）中的思路继续运算得到更为低阶的方程组。

（3）经过n -1步的消元后，得到一个三角方程。

（4）利用求解公式回代得到线性方程组的解。

根据这个思路编写matlab 程序如下：function x=gauss(n)disp('请输入构造的矩阵的阶数（10/20）')n=input('');disp('构造矩阵为A=')A = diag(6*ones(1,n))+diag(ones(1,n-1),1)+diag(8*ones(1,n-1),-1) disp('由A 构造出矩阵b=')b = A*ones(n,1)[m,n]=size(A);disp('增广矩阵为:')Ab=[A b]for i=1:n-1yuan=Ab(i,i);for k=i+1:nAb(k,i:(n+1))=Ab(k,i:(n+1))-(Ab(k,i)/yuan)*Ab(i,i:(n+1));enddisp(Ab);pauseendx=zeros(n,1);x(n)=Ab(n,(n+1))/Ab(n,n);for i=n-1:-1:1x(i)=(Ab(i,(n+1))-Ab(i,i+1:n)*x(i+1:n))/Ab(i,i);end取n=10运行的结果为：Ab =6 1 0 0 0 0 0 0 0 0 78 6 1 0 0 0 0 0 0 0 150 8 6 1 0 0 0 0 0 0 150 0 8 6 1 0 0 0 0 0 150 0 0 8 6 1 0 0 0 0 150 0 0 0 8 6 1 0 0 0 150 0 0 0 0 8 6 1 0 0 150 0 0 0 0 0 8 6 1 0 150 0 0 0 0 0 0 8 6 1 150 0 0 0 0 0 0 0 8 6 14消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000取n=20运行的结果为：Ab =消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00002 采用列选主元消去法进行计算在采用Gauss消去法的时候，主元绝对值的大小将影响到计算结果，主元的绝对值越大，算法的稳定性越好。

列选主元消去法matlab程序的计算思路为：（1）先构造需要计算的矩阵，得到增广矩阵Ab。

（2）将系数矩阵A的每一列的绝对值最大的元素换至对角线上，矩阵b中的元素也随之改变，先判断主元是否为0。

然后利用此主元逐行消去此主元所在的列中的元素，矩阵b中的元素也随之改变。

（3）经过n-1步运算过后，矩阵A就变换成为一个三角矩阵。

（4）逐次回代，就能计算出方程组的解。

分别每步消去过程总选取按模最大的元素作为主元和每步消去过程总选取按模最小的元素作为主元选取列选主元消去法的matlab程序为：function x=gauss(n)disp('请输入构造的矩阵的阶数（10/20）')n=input('');disp('构造矩阵为A=')A = diag(6*ones(1,n))+diag(ones(1,n-1),1)+diag(8*ones(1,n-1),-1)disp('由A构造出矩阵b=')b = A*ones(n,1)[m,n]=size(A);disp('增广矩阵为:')Ab=[A b]for i=1:n-1if way==1i=idisp('输入每一列的主元所在的行')hang=input('');Ab([i hang],:)=Ab([hang i],:);disp(Ab);pauseendzhuyuan=Ab(i,i);for k=i+1:nAb(k,i:(n+1))=Ab(k,i:(n+1))-(Ab(k,i)/zhuyuan)*Ab(i,i:(n+1));enddisp(Ab);pauseendx=zeros(n,1);x(n)=Ab(n,(n+1))/Ab(n,n);for i=n-1:-1:1x(i)=(Ab(i,(n+1))-Ab(i,i+1:n)*x(i+1:n))/Ab(i,i);end（i（（（（（（（（（（（（（（（（（（（（（（:取n=10手动执行该程序：消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000取n=20运行的结果为：Ab =解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000（ii）每步消去过程总选取按模最小的元素作为主元:取n=10手动执行该程序：消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000取n=20运行的结果为：消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000由上述两种选取主元的方法比较可知，采用每一步消去过程总选取按模最大或按模尽可能大的元素作为主元的选取方法最终得到的方程组的解更加精确，这是因为在计算过程中人为的避免了大数除小数而造成的误差的扩大，而每一步消去过程总选取按模最小或按模尽可能小的元素作为主元的选取方法最终得到的解会有误差，而误差的产生就是由于上诉原因而造成了。

所以在实际计算过程中，必须考察计算模型的可行性。

3 采用全选主元消去进行计算全选主元和列选主元消去法的区别在于全选主元不仅要用到列选主元，而且需要行选绝对值最大的元，这样的麻烦之处就在于不仅需要进行线性方程顺序的变换，而且在解方程的过程中未知数的顺序也发生了相应的改变，这就给程序的编写和运算的时间的上带来了很大的麻烦。

全选主元消去法的基本思路为：（1）先构造需要计算的矩阵，得到增广矩阵Ab（2）将每次更新后系数矩阵的每一列的绝对值最大的元素换至对角线上（绝对值最大的元素必须在对角线的下方），矩阵b中的元素也随之改变。

（3）然后将更新后系数矩阵的每行中的绝对值最大元素换至对角线上（该绝对值最大元素必须要在主元的右侧）。

先判断主元是否为0，然后利用此主元逐行消去此主元所在的列中的元素，矩阵b中的元素也随之改变。

（4）经过n-1步运算过后，矩阵A就变换成为一个三角矩阵。

（5）逐次回代，就能计算出方程组的解。

采用全选主元Gauss消去法编写的matlab程序为：function x=gauss(n)disp('请输入构造的矩阵的阶数（10/20）')n=input('');disp('构造矩阵为A=')A = diag(6*ones(1,n))+diag(ones(1,n-1),1)+diag(8*ones(1,n-1),-1)disp('由A构造出矩阵b=')b = A*ones(n,1)[m,n]=size(A);disp('增广矩阵为:')Ab=[A b]for i=1:n-1zhuyuan = max(max(abs(A(i:n,i:n))));for r = i:nfor t = i:nif zhuyuan == abs(A(r,t))zhuhang = r;zhulie = t;endendendph = A(i,:);A(i,:) = A(zhuhang,:);A(zhuhang,:) = ph;b = b';bb = b(i);b(i) = b(zhuhang);b(zhuhang) = bb;b = b';pl = A(:,i);A(:,i) = A(:,zhulie);A(:,zhulie) = pl;Ab = [A,b];zhuyuan = Ab(i,i);for k=i+1:nAb(k,i:(n+1))=Ab(k,i:(n+1))-(Ab(k,i)/zhuyuan)*Ab(i,i:(n+1));endA = Ab(:,[1:n]);b = Ab(:,[(n+1)]);disp(Ab);pauseendx=zeros(n,1);x(n)=Ab(n,(n+1))/Ab(n,n);for i=n-1:-1:1x(i)=(Ab(i,(n+1))-Ab(i,i+1:n)*x(i+1:n))/Ab(i,i);end取n=10手动执行该程序：消去结束后的矩阵为：解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000取n=20运行的结果为：Ab =解得线性方程组的解为：ans =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000根据消去后的矩阵的可以看出：全选主元消去得到的矩阵十分稳定。