Gauss列主元素消去法实验
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]);}。
高斯列主元消去法实验报告

《数值计算方法》实验报告专业:年级:学号:姓名:成绩: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;}。
实验一-Guass列主元消去法及其应用

上机作业总体要求:1.开发语言可用任一种高级语言2.作业包括1)一份实验报告(书面,内容同电子版实验报告)2)电子版作业的全套(压缩后提交在Webcc上),包括:程序源代码;可执行程序;电子版实验报告(内容包括:一、实验目的二、模型建立3.1 开发环境3.2 程序设计说明(要求设计为通用的)3.3 源代码3.4 程序使用说明3.5 模型的解四、小结(可含个人心得体会))第二章应用题:投入产出(上机:列主元求解)一个城镇有三个主要企业:煤矿、电厂和铁路作为它的经济系统。
●生产价值1元的煤(产品1),需消耗0.3元的电费(中间产品2)和0.2元的运输费(中间产品3);●生产价值1元的电(产品2),需消耗0.4元的煤费(中间产品1)、0.1元的电费和0.1元的运输费;●提供价值1元的铁路运输服务(产品3),需消耗0.3元的煤费和0.2元的电费和0.2元运输费。
在某个星期内,除了三个企业间的彼此需求,煤矿得到6万元的订单(即最终产品1),电厂得到3万元的电量供应要求(即最终产品2),而地方铁路得到价值5万元的运输需求(即最终产品3)。
试问:(1)这三个企业在这星期各应生产多少产值才能满足内外需求?(提示: 中间产品+最终产品(外部需求)=总产品)(2)除了内部需求,试求这星期各企业之间的消耗需求,同时求出各企业新创造的价值。
(提示:内部消耗+新创价值=总产值)(3)若三个企业的外部需求分别增长1个单位(万元),则各企业的总产值分别增长多少?1.列主元算法编为一函数,在main()中调用。
实验报告的要求见“上机作业总要求”。
2.取4位小数计算.。
2-2 Gauss列主元消去法

S2 若ann 0,则输出“ A是奇异矩阵”;停机 . S3 置xn an,n1 / ann ;
对i n 1, n 2,...1,
ai,n1 n aij x j
置xi
j i 1
aii
S4 输出x1, x2,..., xn ;停机.
作业:
P50 习题3
k in
aik
;
S12 若aik ,k 0,则输出“ A是奇异矩阵”;停机 .
S13 若ik k,则
akj aik , j j k,...,n 1;
S14 对i k 1,..., n
置aik aik / akk ; 对j k 1,..., n 1
置aij aij aik akj.
§2-2 Gauss列主元消去法
一、Gauss列主元消去法的引入 例1. 用3位浮点数运算,求解线性方程组
0.0001xx11
x2 x2
1 2
解: 本方程组的精度较高的解为
x* (1.00010001 ,0.99989999 )T
用Gauss消去法求解
A ( A,b)
0.000100 1
1 1
21
0.000100
m2110 000
0
回代后得到
1
1
1.00 104 1.00 104
x1 0.00 , x2 1.00
与精确解相比,该结果显然是错误的 究其原因,在求行乘数时用了很小的数0.0001作除数
如果在求解时将1,2行交换,即
A ( A,b)
1 0.000100
1 1
a(2) i2
,
交换第2行和第i2行,
2in
然后进行消元,得[ A(3) , b(3) ].
Gauss列主元消去法程序设计

《Gauss列主元消去法》实验报告实验名称:Gauss列主元消去法程序设计成绩:___________专业班级:数学与应用数学1202班姓名:王晓阳学号:28实验日期:2014 年11月10日实验报告日期:2014年11月10日一.实验目的1.学习Gauss消去法的基本思路和迭代步骤.2.学会运用matlab编写高斯消去法和列主元消去法程序,求解线性方程组.3.当()ka绝对值较小时,采用高斯列主元消去法.kk4.培养编程与上机调试能力.二、实验内容用消去法解线性方程组的基本思想是用逐次消去未知数的方法把原线性方程组Ax b=化为与其等价的三角形线性方程组,而求解三角形线性方程组可用回代的方法求解.1.求解一般线性方程组的高斯消去法.(1)消元过程:设()0k kk a ≠,第i 个方程减去第k 个方程的()()/k k ik ik kk m a a =倍,(1,,)i k n =+L ,得到()()11k k A x b ++=.()()()()()()()11,,1,,k k k ij ij ik kj k k k ii ik k a a m a i j k n b b m b ++⎧=-=+⎪⎨=-⎪⎩L 经过n-1次消元,可把方程组()()11A x b =化为上三角方程组()()n n A x b =.(2)回代过程:()()()()()1//,1,,1n n n n nn n i i i i i ij j ii j i x b a x b a x a i n =+⎧=⎪⎛⎫⎨=-=- ⎪⎪⎝⎭⎩∑L 以解如下线性方程组为例测试结果.1212312310773264556x x x x x x x x -=⎧⎪-++=⎨⎪-+=⎩2.列主元消去法由高斯消去法可知,在消元过程中可能出现()0k kk a =的情况,这是消去法将无法进行,即使主元素()0k kk a ≠但很小时,用其作除数,会导致其他元素数量级的严重增长和舍入误差的扩散,最后也使得计算解不可靠.这时就需要选取主元素,假定线性方程组的系数矩阵A 是菲奇异的.(1)消元过程:对于1,2,,1k n =-L ,进行如下步骤:1) 按列选主元,记max pk ik k i na a ≤≤=2) 交换增广阵A 的p,k 两行的元素。
列主元素Gauss消去法Jacobi迭代法原理及计算方法

一、 列主元素Gauss 消去法、Jacobi 迭代法原理及计算方法1. 列主元素Gauss 消去法:1.1 Gauss 消去法基本原理设有方程组Ax b =,设A 是可逆矩阵。
高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵[]B A b = ,将其中的A 变换成一个上三角矩阵,然后求解这个三角形方程组。
1.2 列主元Gauss 消去法计算步骤将方程组用增广矩阵[]()(1)ijn n B A b a ⨯+== 表示。
1). 消元过程对1,2,,1k n =-(1) 选主元,找{},1,,k i k k n ∈+ 使得 ,max k i k ik k i na a ≤≤= (2) 如果,0k i k a =,则矩阵A 奇异,程序结束;否则执行(3)。
(3) 如果k i k ≠,则交换第k 行与第k i 行对应元素位置,k kj i j a a ↔,,,1j k n =+ 。
(4) 消元,对,,i k n = ,计算/,ik ik kk l a a =对1,,1j k n =++ ,计算.ij ij ik kj a a l a =-2). 回代过程(1) 若0,nn a =则矩阵奇异,程序结束;否则执行(2)。
(2) ,1/;n n n nn x a a +=对1,,2,1i n =- ,计算,11/n i i n ij j ii j i x a a x a +=+⎛⎫=- ⎪⎝⎭∑2. Jacobi 迭代法2.1 Jacobi 迭代法基本原理Jacobi 迭代法的基本思想是对n 元线性方程组b Ax =,.,n n R b R A ∈∈将其变形为等价方程组f Bx x +=,其中.,,n n n n R x R f R B ∈∈∈⨯B 成为迭代矩阵。
从某一取定的初始向量)0(x 出发,按照一个适当的迭代公式 ,逐次计算出向量f Bx x k k +=+)()1( ( 1,0=k ),使得向量序列}{)(k x 收敛于方程组的精确解.(1)输入1,,,,)0(=k n xb A ε,. (2) )(1,1)0()1(∑≠=-=n j i i j ij i iii x a b a x )1,0(n i = (3)判断 ε≤--≤≤)0()1(10max i i n i x x ,若是,输出1)1(2)1(1,,n x x x ,若否,置1+=k k ,)1()0(i i x x =,)2,1(n i =。
Gauss列主元消去法

贵州师范大学数学与计算机科学学院学生实验报告课程名称: 数值分析 班级: 实验日期: 年 月 日 学 号: 姓名: 指导教师: 实验成绩:一、实验名称实验五:线性方程组的数值解法二、实验目的及要求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,但很小时,舍入误差也会使计算结果面目全非,避免这种缺陷的基本方法就是选主元。
Gauss消去法

一、实验题目用顺序Gauss 消去法和列主元Gauss 消去法解线性方程组:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡--⨯-2178.4617.5911212592.112113.6291.51314.59103.0432115x x x x 二、 实验目的 1.掌握顺序Gauss 消去法和列主元Gauss 消去法的基本思路和迭代步骤;2.培养编程与上机调试能力。
三、实验原理1. 顺序Gauss 消去法的原理:(1).消元计算(k=1,2,….,n-1)()()(1)()()(1)()()/,1,...,,,,1,...,,,1,...,.k k ik ik kk k k k ij ij ik kj k k k i i ik k m a a i k n a a m a i j k n b b m b i k n ++==+=-=+=-=+(2).回带公式()()()()()1/,()/,1,...,2,1.n n n n nn ni i i i iii j ii j i x b a x b a x a i n =+==-=-∑2. 列主元Gauss 消去法的原理:设Ax=b.本算法用A 的具有行交换的列主元素法,消元结果冲掉A,乘数m ij 冲掉a ij ,计算解x 冲掉常数项b ,行列式存放在det 中。
(1)det ←1(2)对于k=1,2,…,n-1按列选主元:,max ik ik k i n a k a ≤≤=如果,0ik a k =,则停止(det (A )=0)如果k i =k 则转(a )换行:1*(*)/*det *det/ikkkij ij ik kjni i ij j ii ik kj i nn nn n n nnm a a a m a b b a b a m b a a b b a =+←←←←-←←=←∑(a ).消元计算对于 i=k+1,…,n① ik a ←ik m =ik a /kk a② 对于j=k+1,…,n*ij ij ik kj a a m a ←←③ *i i ik k b b m b ←←det *det kk a ←(3).如果0nn a =,则计算停止(det (A )=0)(4).回代求解四、实验内容及结果原始数据:A=[0.3e-15,59.14,3,1;5.291,-6.13,-1,2;11.2,9,5,2;1,2,1,1]; b=[59.17;46.78;1;2];1.顺序Gauss 消去法程序源代码%magauss.mfunction x=magauss(A,b,flag)%用途:顺序Gauss 消去法解线性方程组Ax=b%格式:x=magauss(A,b,flag), A 为系数矩阵, b 为右端项, 若flag=0, 则不 % 显示中间过程,否则显示中间过程, 默认为0, x 为解向量if nargin<3,flag=0;endn=length(b);%消元for k=1:(n-1)m=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n);b(k+1:n)=b(k+1:n)-m*b(k);A(k+1:n,k)=zeros(n-k,1);if flag~=0, Ab=[A,b], endend%回代x=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);end输入: A=[0.3e-15,59.14,3,1;5.291,-6.13,-1,2;11.2,9,5,2;1,2,1,1] b=[59.17;46.78;1;2]magauss2(A,b,0)结果:A =0.0000 59.1400 3.0000 1.00005.2910 -6.1300 -1.0000 2.000011.2000 9.0000 5.0000 2.00001.00002.0000 1.0000 1.0000b =59.170046.78001.00002.0000ans =3.84571.6095-15.476110.41132. 列主元Gauss消去法程序源代码%magauss2.mfunction x=magauss2(A,b,flag)%用途:列主元Gauss消去法解线性方程组Ax=b%格式:x=magauss(A,b,flag), A为系数矩阵, b为右端项, 若flag=0, % 则不显示中间过程,否则显示中间过程, 默认为0, x为解向量if nargin<3,flag=0;endn=length(b);for k=1:(n-1)%选主元[ap,p]=max(abs(A(k:n,k))); p=p+k-1;if p>kt=A(k,:); A(k,:)=A(p,:); A(p,:)=t;t=b(k); b(k)=b(p); b(p)=t;end%消元m=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n);b(k+1:n)=b(k+1:n)-m*b(k);A(k+1:n,k)=zeros(n-k,1);if flag~=0, Ab=[A,b], endend%回代x=zeros(n,1);x(n)=b(n)/A(n,n);for k=n-1:-1:1x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k);end输入: A=[0.3e-15,59.14,3,1;5.291,-6.13,-1,2;11.2,9,5,2;1,2,1,1] b=[59.17,46.78,1,2]'magauss(A,b,0)结果:A =0.0000 59.1400 3.0000 1.00005.2910 -6.1300 -1.0000 2.000011.2000 9.0000 5.0000 2.00001.00002.0000 1.0000 1.0000b =59.170046.78001.00002.0000ans =23.68481.0005五、实验结果分析顺序高斯消去法求解线性方程组具有较高的稳定性,但计算过程繁琐。
Gauss消去法和列主元消去法

max=abs(C(s,s));big=0;
if det(C(s:n,s:n))==0
disp('此方程无解');
answer=0;
break;
end
for i=s:n
if max<abs(C(i,s))
max=abs(C(i,s));
k=i;
big=1;
else continue
for i=n:(-1):1
X(i,1)=C(i,n+1);
for j=(i+1):n
X(i,1)=X(i,1)-E(i,j)*X(j,1);
end
X(i,1)=X(i,1)/E(i,i);
ห้องสมุดไป่ตู้end
disp('此方程的解为:')
X
end
5、实验结果
请输入未知数系数矩阵A:
A=[2,-1,3;4,2,5;1,2,0]
3、实验原理
高斯列主元消去法
4、实验内容
clc;clear;format short
disp('请输入未知数系数矩阵A:');
A=input('A=');
disp('请输入常数项列向量B:');
B=input('B=');
C=[A,B];
[m,n]=size(A);
s=1;answer=1;P=zeros(1,n);L=zeros(n);I=eye(n);
y=a(i,k:n+1);a(i,k:n+1)=a(k,k:n+1);a(k,k:n+1)=y;
break;
实验二:Gauss列主元消去法

实验二:Gauss列主元消去法程序1:Gauss列主元消去法A=input('请输入线性方程组的增广矩阵A=');n=length(A)-1;x=zeros(n,1);aa=zeros(n,1);for j=1:nfor i=1:(n+1)AA(j,i)=abs(A(j,i));endendfor k=1:(n-1)for i=k:naa(i-(k-1))=AA(i,k);endfor i=k:nif AA(i,k)==max(aa)breakendendif AA(i,k)==0breakfprintf('方程组系数矩阵奇异\n');elsefor j=k:(n+1)jh=A(i,j);A(i,j)=A(k,j);A(k,j)=jh;endendfenzi=A(k,k);for j=k:(n+1)A(k,j)=A(k,j)/fenzi;endfor p=(k+1):njj=A(p,k);for j=k:(n+1)A(p,j)=A(p,j)-jj*A(k,j);endendendif k==(n-1)x(n)=A(n,(n+1))/A(n,n);for i=(n-1):(-1):1he=0;for j=(i+1):nhe=he+A(i,j)*x(j);endx(i)=A(i,(n+1))-he;endendx用Gauss列主元消去法解方程组:1.请输入线性方程组的增广矩阵A=[1e-008,2,3,1;-1,3.172,4.623,2;-2,1.072,5. 643,3]x =-0.4653-0.07000.38002.请输入线性方程组的增广矩阵A=[4,-2,4,10;-2,17,10,3;-4,10,9,-7];x =2.94640.6071-0.14293.请输入线性方程组的增广矩阵A=[0.3e-020,1,0.7;1,1,0.9]x =0.20000.7000程序2:不选主元的高斯消去法A=input('请输入线性方程组的增广矩阵A=');n=length(A)-1;x=zeros(n,1);for k=1:(n-1)if A(k,k)==0breakfprintf('方程组不能用普通的高斯消去法解\n');elsefenzi=A(k,k);for j=k:(n+1)A(k,j)=A(k,j)/fenzi;endfor p=(k+1):njj=A(p,k);for j=k:(n+1)A(p,j)=A(p,j)-jj*A(k,j);endendx(n)=A(n,(n+1))/A(n,n);for i=(n-1):(-1):1he=0;for j=(i+1):nhe=he+A(i,j)*x(j);endx(i)=A(i,(n+1))-he;endendendx用不选主元的Gauss消去法解方程组:1.请输入线性方程组的增广矩阵A=[4,-2,4,10;-2,17,10,3;-4,10,9,-7];x =2.94640.6071-0.14292.请输入线性方程组的增广矩阵A=[1e-008,2,3,1;-1,3.172,4.623,2;-2,1.072,5. 643,3];x =-0.4653-0.07000.38003.请输入线性方程组的增广矩阵A=[0.3e-020,1,0.7;1,1,0.9]x =0.7000。
高斯消元法与列主元消去法实验报告

实验报告:Gauss消元法小组成员:李岚岚、邱粉珊、缪晓浓、杨水清学号:0917020040、0917010078、0917010073、0917010112一、实验问题编写两个程序,分别利用Gauss消元法和列主元消去法求解方程组二、分析及其计算过程Gauss顺序消元法:源程序:function [x]=gaussl(A,b)[n1,n2]=size(A);n3=size(b);if n1~=n2|n2~=n3|n1~=n3disp('A的行和列的维数不同!');return;endif det(A)==0disp('系数矩阵A奇异');return;end%消元过程L=eye(n1);for j=2:n1for i=j:n1L(i,j-1)=A(i,j-1)/A(j-1,j-1);A(i,:)=A(i,:)-L(i,j-1)*A(j-1,:);b(i)=b(i)-L(i,j-1)*b(j-1);endend%回代过程x(n1)=b(n1)/A(n1,n1);for t=n1-1:-1:1for k=n1:-1:t+1b(t)=b(t)-A(t,k)*x(k);endx(t)=b(t)/A(t,t);end程序的运行以及结果:>>A=[1 2/3 1/3;9/20 1 11/20;2/3 1/3 1];>>b=[2 2 2];>> [x]=gaussl(A,b)x =1 1 1Gauss列主元消去法:源程序:function [x]=gaussll(A,b) [n1,n2]=size(A);n3=size(b);if n1~=n2|n1~=n3|n2~=n3disp('输入的方程错误!');return;endif det(A)==0disp('系数矩阵A奇异');return;endmax=zeros(n1);for m=1:n1%找主元for i=m:n1if abs(A(i,m))>maxmax=A(i,:);A(i,:)=A(m,:);A(m,:)=max;maxb=b(i);b(i)=b(m);b(m)=maxb;endend%消元过程L=eye(n1);for j=2:n1for i=j:n1L(i,j-1)=A(i,j-1)/A(j-1,j-1);A(i,:)=A(i,:)-L(i,j-1)*A(j-1,:);b(i)=b(i)-L(i,j-1)*b(j-1);endendend%回代过程x(n1)=b(n1)/A(n1,n1);for t=n1-1:-1:1for k=n1:-1:t+1b(t)=b(t)-A(t,k)*x(k);endx(t)=b(t)/A(t,t);end程序的运行以及结果:>>A=[-0.002 2 2;1 0.78125 0;3.996 5.5625 4]; >>b=[0.4 1.3816 7.4178];>>[x]= gaussll(A,b)x =1.9273 -0.6985 0.9004。
数学实验题目5 相对Gauss列主元消去法

数学实验题目5 相对Gauss 列主元消去法摘要由一般线性方程组在使用Gauss 消去法求解时,从求解过程中可以清楚地看到,若(1)0k kk a -=,必须施以行交换的手续,才能使消去过程继续下去。
有时既使(1)0k kk a -≠,但其绝对值很小,由于舍入误差的影响,消去过程也会出现不稳定现象。
因此,为使这种不稳定现象发生的可能性减至最小,在施行消去过程时每一步都要选主元素,即要寻找行r ,使(1)(1)||max ||k k rk ik i ka a -->=并将第r 行与第k 行交换,以使(1)k kka -的当前值(即(1)k ika -的数值)远大于0。
这种列主元消去法的主要步骤如下: 1.消元过程 对1,2,,1k n =-,做1º 选主元,记||max ||rk ik i ka a >=若0rk a =,说明方程组系数矩阵奇异,则停止计算,否则进行2º。
2º 交换A (增广矩阵)的,r k 两行元素,,1rj kja a j k n ↔=+3º 计算/ij ij ik kj kk a a a a a =-1,,i k n =+1,,1j k n =++2.回代过程 对,1,,2,1k n n =-,计算,11(/)nk k n kjj kk j k x a ax a +=-=-∑前言利用Gauss 列主元消去法、显式相对Gauss 列主元消去法、隐式相对Gauss 列主元消去法求解线性方程组程序设计流程是否否是开 始输入A (增广矩阵)1k = ||max ||rk ik i ka a >=rk a =1,,1,,1/ij ij ik kj kk i k n j k n a a a a a =+=++=-交换A 中,r k 两行1k n <-,11,1,,2,1()/nk k n kjj kkj k k n n x a ax a +=+=-=-∑输出x结 束 1k =问题1(1)程序运行如下:x = GaussSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000 (2)程序运行如下:x = GaussSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000 (3)程序运行如下:x = GaussSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000 (4)程序运行如下:x = GaussSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000问题2(1)程序运行如下:= GaussSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008x = GaussExpSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008x = GaussIneSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008 (2)程序运行如下:x = GaussSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365x = GaussExpSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365x = GaussIneSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365 (3)程序运行如下:x = GaussSysSolve(Mat2_3,b2_3)x = 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat2_3,b2_3)x = 1 1 1x = GaussIneSysSolve(Mat2_3,b2_3)x = 1.0000 1.0000 1.0000(4)程序运行如下:x = GaussSysSolve(Mat2_4,b2_4)x = 1 1 1x = GaussExpSysSolve(Mat2_4,b2_4)x = 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat2_4,b2_4)x = 1 1 1使用的函数function x = GaussSysSolve(A, b)% GaussSysSolve 用Gauss消去法解线性方程组Ax = b%% Synopsis: x = GaussSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵for i = 1:n[amax, imax] = max(Ab(i:n, i)); %选择主元if amax == 0 %主元为0,矩阵奇异error('Tne Linear System is singular!');elseif i ~= imax+i-1 %主元行数与i不同时,交换这两行Ab([i imax+i-1],:) = Ab([imax+i-1 i], :);endfor j = i+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(i,:) * Ab(j,i)/amax;endendx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算x x(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);endfunction x = GaussExpSysSolve(A, b)% GaussExpSysSolve 用显式Gauss列主元消去法解线性方程组Ax = b%% Synopsis: x = GaussExpSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵for i = 1:n %显式平衡技术s = max(Ab(i,1:n));Ab(i,:) = Ab(i,:)/s;endfor i = 1:n[amax, imax] = max(Ab(i:n, i)); %选择主元if amax == 0 %主元为0,矩阵奇异error('Tne Linear System is singular!');elseif i ~= imax+i-1 %主元行数与i不同时,交换这两行Ab([i imax+i-1],:) = Ab([imax+i-1 i], :);endfor j = i+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(i,:) * Ab(j,i)/amax;endendx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算xx(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);endfunction [x,det] = GaussIneSysSolve(A, b)% GaussIneSysSolve 用隐式Gauss列主元消去法解线性方程组Ax = b%% Synopsis: x = GaussIneSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量% det = 系数矩阵行列式的值[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵det = 1; %初始化系数矩阵行列式为1for i = 1:n %隐式平衡技术s(i) = max(abs(Ab(i,1:n)));if s(i) == 0error('Tne Linear System is singular!'); %系数矩阵某行全为0时,矩阵奇异endends = s(:);for k = 1:n-1[c, kmax] = max(abs(Ab(k:n, k)./s(k:n))); %选择主元if c == 0 %主元为0,矩阵奇异det = 0;error('Tne Linear System is singular! det(A) = 0');elseif k ~= kmax+k-1 %主元行数与k不同时,交换这两行s([k kmax+k-1]) = s([kmax+k-1 k]);Ab([k kmax+k-1],:) = Ab([kmax+k-1 k], :);det = -det;endfor j = k+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(k,:) * Ab(j,k)/Ab(k,k);enddet = Ab(k,k)*det;endif Ab(n,n) == 0det = 0;error('Tne Linear System is singular!'); %最后一行唯一非0元素为0时,矩阵奇异endx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算x x(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);enddet = Ab(n,n)*det;思考题(1)在各主元不是非常小的时候,三种方法结果一致(2)隐式平衡列选主元法最好,应为当主元很小时,普通的Gauss消元法会产生很大误差;显式平衡列选主元法每一行除以其绝对值最大元素时会引入额外的舍入误差。
计算方法实验报告_列主元高斯消去法

row_first=A[i][i]; for(int j=0;j<n+1;j++)
计算方法实验报告
{ A[i][j]=A[i][j]/row_first;
} }
for(int k=n-1;k>0;k--) {
for(int i=0;i<N;i++) {
for(int j=0;j<N;j++) {
A_B[i][j]=A[i][j]; } A_B[i][N]=B[i][0]; } return A_B; }
3
//输出矩阵 A 的 row x col 个元素 void Show_Matrix(double **A,int row,int col) {
for(int i=0;i<N;i++)
{
int row=Choose_Colum_Main_Element(N,A_B,i);
if(Main_Element<=e) goto A_0;
Exchange(A_B,N+1,row,i);
Elimination(N,A_B,i);
cout<<"选取列主元后第"<<i+1<<"次消元:"<<endl;
double factor; for(int i=start+1;i<n;i++) {
factor=A[i][start]/A[start][start]; for(int j=start;j<n+1;j++) {
Gauss列主元消去法程序设计

《Gauss 列主元消去法》实验报告实验名称:Gauss 列主元消去法程序设计专业班级: 信息与计算科学1701姓名:梅礼坤 学号:1508060117 实 验 日 期 : 2019 年10月14日一.实验目的1.学习Gauss 列主元消去法的基本思路和迭代步骤.2.学会运用matlab 编写高斯列主元消去法程序,求解线性方程组.3.培养编程与上机调试能力.二、实验内容用消去法解线性方程组的基本思想是用逐次消去未知数的方法把原线性方程组Ax b =化为与其等价的三角形线性方程组,而求解三角形线性方程组可用回代的方法求解.1.列主元消去法由高斯消去法可知,在消元过程中可能出现()0k kk a =的情况,这是消去法将无法进行,即使主元素()0k kk a ≠但很小时,用其作除数,会导致其他元素数量级的严重增长和舍入误差的扩散,最后也使得计算解不可靠.这时就需要选取主元素,假定线性方程组的系数矩阵A 是菲奇异的.(1)消元过程:对于1,2,,1k n =-,进行如下步骤:1) 按列选主元,记max pk ik k i na a ≤≤=2) 交换增广阵A 的p,k 两行的元素。
A(k,j)=A(p,j) (j=k,…,n +1)3)交换常数项b 的p,k 两行的元素。
b(k)=b(p)4)计算消元()()()()()()()11,,1,,k k k ij ij ik kj k k k ii ik k a a m a i j k n b b m b ++⎧=-=+⎪⎨=-⎪⎩(2)回代过程 ()()()()()1//,1,,1n n n n nn n i i i i i ij j ii j i x b a x b a x a i n =+⎧=⎪⎛⎫⎨=-=- ⎪⎪⎝⎭⎩∑三、实验环境MATLAB R2010b四、实验步骤1.高斯列主元消去法流程图:2.程序设计:(一)高斯列主元消去法:a=input('请输入系数阵:'); b=input('请输入常数项:'); n=length(b);A=[a,b];x=zeros(n,1); %初始值for k=1:n-1if abs(A(k,k))<10^(-4);%判断是否选主元 y=1elsey=0;endif y; %选主元for i=k+1:n;if abs(A(i,k))>abs(A(k,k))p=i;else p=k;endendif p~=k;for j=k:n+1;s=A(k,j);A(k,j)=A(p,j);%交换系数A(p,j)=s;endt=b(k);b(k)=b(p);%交换常数项b(p)=t;endendfor i=k+1:nm(i,k)=A(i,k)/A(k,k); %第k次消元for j=k+1:nA(i,j)=A(i,j)-A(k,j)*m(i,k);endb(i)=b(i)-m(i,k)*b(k);endendx(n)=b(n)/A(n,n); %回代for i=n-1:-1:1;s=0;for j=i+1:n;s=s+A(i,j)*x(j);endx(i)=(b(i)-s)/A(i,i)end五、实验结果六、实验讨论、结论本实验通过matlab程序编程实现了高斯列主元消去法的求解,能加深对高斯消去法基本思路与计算步骤的理解。
数值分析实验二(列主元Gauss消去法)

《数值分析》实验报告实验编号:实验二课题名称:列主元Gauss消去法一、算法介绍1、输入矩阵的阶数n,方程组的增广矩阵A;2、对k=0,1,…,n-2,循环:选取列中绝对值最大的元素,将主元所在的行的元素保存在数组temp[n+1]中。
若主元为零,则系数矩阵奇异,计算停止;否则,顺序进行。
如果绝对值最大的元素就在矩阵的对角线上,则进行普通高斯消元法的第一大步,否则将方程组系数换行之后再进行普通高斯消元法的第一大步;3、然后利用回代法求解线性方程组。
二、程序代码#include<iostream>#include<cmath>#include<iomanip>using namespace std;int main(){int n=0,k=0,i=0,j=0,h=0,g=0,flag=0,i1,j1;double max=0,m=0;cout<<"***利用列主元Gauss消元法求解线性方程组***"<<endl;cout<<"请输入矩阵的阶数:"<<endl;cin>>n;double a[n][n+1];double t[n+1];double x[n];memset(a,0,sizeof(a));memset(x,0,sizeof(x));cout<<"请输入方程组的增广矩阵:"<<endl;for(i=0;i<n;i++){for(j=0;j<n+1;j++){cin>>a[i][j];}}for(k=0;k<n-1;k++){max=0;j1=0;for(i=k;i<n;i++){if(fabs(a[i][k])>max){max=fabs(a[i][k]);i1=i;j1=k;}}if(max==0){cout<<"该系数矩阵为奇异矩阵,计算停止"<<endl;flag=1;break;}else{cout<<"第"<<j1+1<<"列中绝对值最大的元素是"<<a[i1][j1]<<",在线性方程组的第"<<i1+1<<"行"<<endl;if(i1!=k){for(j=0;j<=n;j++){t[j]=a[i1][j];a[i1][j]=a[k][j];a[k][j]=t[j];}}for(i=k+1;i<=n-1;i++){m=a[i][k]/a[k][k];for(j=k;j<=n;j++)a[i][j]=a[i][j]-m*a[k][j];for(g=0;g<n;g++){for(h=0;h<n+1;h++)cout<<setiosflags(ios::fixed)<<setprecision(2)<<a[g][h]<<" ";cout<<endl;}cout<<endl;}}}if(flag==0){x[n-1]=a[n-1][n]/a[n-1][n-1] ;double sum=0;for(k=n-2;k>=0;k--){sum=0;for(i=n-1;i>=k;i--)sum+=a[k][i]*x[i];x[k]=(a[k][n]-sum)/a[k][k];}cout<<"该线性方程组的解为:"<<endl;for(i=0;i<n;i++)cout<<"x"<<i+1<<"="<<setiosflags(ios::fixed)<<setprecision(2)<<x[i]<<endl;}system("pause");return 0;}三、运算结果截屏四、算法分析列主元Gauss消元法避免了普通高斯消元法中出现的问题:遇到某个主元为零或者当主元绝对值很小时,计算将会停止或求出的结果将与其实际结果相差很远。
高斯列主消元数值分析实验报告

数值分析实验报告之高斯列主消元法一、实验目的:清楚高斯列主元消去法与高斯主元素消去法的区别,以及它提出的必要性;掌握高斯列主元消去法的原理及推导过程,会用其解简单的线性方程组。
二、实验内容:用高斯列主元消去法解方程组⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--643.5072.1000.2623.4712.3000.1000.3000.2001.0⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡321x x x =⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡000.3000.2000.1 三、实验原理:在采用高斯消去法解方程组时,小主元可能产生麻烦,即用其做除数,会导致其他元素数量级的严重增长和舍入误差的扩散,最终使得计算的解不可靠。
故应避免采用绝对值小的主元素。
在消元之前,选择一个绝对值最大的元素作为主元,用其做除数来进行消元,这样就具有较好的数值稳定性。
这就是选主元消去法。
下面详细说明列主元素消去法。
第一步:在Ax=b 即)1()1(b x A =的系数矩阵)1(A 的第一列元素中选择一个绝对值最大的元素,不妨设为)1(1l a 。
对调)1(1j a 和)1(lj a 及)1(1b 和)1(l b (j=1,2,……,n ,1≤l ≤n)。
以)1(1l a 作为新的)1(11a 进行消元(消去对调后的第2~n 个方程中的1x )。
第k 步:(1≤k ≤n-1)设第k-1步消元过程完成,得到)()(k k b x A =,检查)(k A 中第k 列的后n-k+1个元素)(k kk a ,)(1k k k a +,…,)(k nk a ,从中选出绝对值最大者,不妨设是)(k pk a ,称它为第k 列主元素。
若p=k ,则取)(k kk a 做除数直接进行消元。
若p ≠k,则将第p 个方程与第k 个方程对调,使)(k pk a 成为新的)(k kk a ,然后以其作为除数进行消元,继续这一过程,直至得到等价的三角形方程组)()(n n b x A =,下一阶段的回代过程不变。
列主元素高斯消去法

实验报告题目:列主元素高斯消去法学生姓名陈玉霞学号20111325015学院信息与控制学院专业信息工程(系统工程方向)指导教师殷传洋2013 年11 月18日一、实验目的(1)掌握高斯消去法的基本思路和迭代步骤;(2)培养编程与上机调试能力。
二、实现功能本程序采用GAUSS列主元消去法求解线性方程组。
AX=b.其中A为N阶矩阵,X,b均为N围列向量。
三、算法描述高斯消去法基本思路设有方程组Ax=b,设A是可逆矩阵。
高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵B [A b],将其中的A变换成一个上三角矩阵,然后求解这个三角形方程组。
四、实验内容X1-0.5*X2+1.5*X3=0.5X2-0.25*X3=0.5X3=-6五、实验步骤(1) 写出增广矩阵A;(2) 选主元;(3) 判断是否奇异;(4) 交换对应行元素;(5) 消元;(6) 回代,计算出X[1], X[2], X[3].六、代码# include<stdio.h># include<math.h># define delta 1e-6#define N 100void main(){int i,j,t,r,n,u,c=0;float p,L,max,s;float X[N];float a[N][N+1];printf("请输入方程的阶数\n"); scanf("%d",&n);printf("输入的原方程系数,中间用空格隔开\n"); for(i=0;i<n;i++) for(j=0;j<n+1;j++)scanf("%f",&a[i][j]);printf("方程系数为\n");for(i=0;i<n;i++){for(j=0;j<n+1;j++){printf("%.2f ",a[i][j]);if(j==n) {printf("\n");}}}for(j=0;j<n-1;j++){{max=fabs(a[j][j]);r=j;}for(i=j+1;i<n;i++)if(fabs(a[i][j])>max){max=a[i][j];r=i;}if(fabs(a[i][j])<delta)printf("矩阵奇异");for(t=j;t<n+1;t++){} p=a[j][t]; a[j][t]=a[r][t]; a[r][t]=p; } for(i=j+1;i<N;i++) { L=a[i][j]/a[j][j]; for(t=j;t<n+1;t++) a[i][t]=a[i][t]-L*a[j][t]; }printf("输出原方程的解为:\n");X[n-1]=a[n-1][n]/a[n-1][n-1];for(i=n-2;i>=0;i--){s=a[i][n];for(j=i+1;j<n;j++)s=s-a[i][j]*X[j];X[i]=s/a[i][i];}for(u=0;u<n;u++)for(j=0;j<n+1;j++){printf("%12f",a[u][j]); c++;if(c%(n+1)==0) printf("\n");}for(i=0;i<n;i++){printf("x(%d)=%.4f\n",i+1,X[i]); if(i==n-1)printf("\n");}}七、实验结果八、实验体会从实验中,我学会了很多,首先要弄懂高斯消去法的原理和计算步骤,使我更加了解Gauss列主元素消去法的实质,掌握了新知识的同时,也为以后的学习打下一个好的基础。
Gauss完全主元素消去法解方程组完全

计算方法实验报告(三)班级:地信10801 序号:姓名:一、实验题目:Gauss完全主元素消去法解方程组二、实验学时: 2学时三、实验目的和要求1、掌握高斯完全主元素消去法基础原理2、掌握高斯完全主元素消去法解方程组的步骤3、能用程序语言对高斯完全主元素消去法进行编程实现四、实验过程代码及结果1.代码#include<stdio.h>#include<iostream.h>#include"math.h"float a[100][101];float x[10];int N; //阶数void shuchu(){for(int i=1;i<=N;i++){for(int j=1;j<=N+1;j++){cout<<a[i][j];cout<<" "<<" ";}cout<<endl;}}void initdata(){cout<<"请输入阶数N:";cin>>N;cout<<endl;cout<<"请输入N*(N+1)个数"<<endl; //输入矩阵中的数for(int i=1;i<=N;i++)for(int j=1;j<=N+1;j++){cin>>a[i][j];}cout<<endl;cout<<"建立的矩阵为:"<<endl; //打印出矩阵shuchu();}void main(){int z[10];int maxi,maxj;initdata();for(int i=1;i<=N;i++)z[i]=i;for(int k=1;k<N;k++){maxi=k;maxj=k;float maxv=abs(a[k][k]);for(i=k;i<=N;i++)for(int j=k;j<=N;j++)if(abs(a[i][j])>maxv){maxv=abs(a[i][j]);maxi=i;maxj=j;}if(maxi!=k) //换行{for(int j=1;j<=N+1;j++){float t=a[k][j];a[k][j]=a[maxi][j];a[maxi][j]=t;}}if(maxj!=k) //换列{for(i=1;i<=N;i++){float t=a[i][k];a[i][k]=a[i][maxj];a[i][maxj]=t;}int t=z[k]; z[k]=z[maxj];z[maxj]=t;}for(int i=k+1;i<=N;i++) //消元{float l=a[i][k]/a[k][k];for(int j=k;j<=N+1;j++){a[i][j]+=-l*a[k][j];}}}//回代for(i=N;i>0;i--){float s=0;for(int j=i+1;j<=N;j++){s+=a[i][j]*x[z[j]];}x[z[i]]=(a[i][N+1]-s)/a[i][i];}cout<<"用完全主元素消去法后的矩阵为:" <<endl; shuchu();for(i=1;i<=N;i++) // 打印出x[i]cout<<"x["<<i<<"]="<<x[i]<<endl;}2.结果。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
Lab06.Gauss 列主元素消去法实验【实验目的和要求】1.使学生深入理解并掌握Gauss 消去法和Gauss 列主元素消去法步骤; 2.通过对Gauss 消去法和Gauss 列主元素消去法的程序设计,以提高学生程序设计的能力;3.对具体问题,分别用Gauss 消去法和Gauss 列主元素消去法求解。
通过对结果的分析比较,使学生感受Gauss 列主元素消去法优点。
【实验内容】1.根据Matlab 语言特点,描述Gauss 消去法和Gauss 列主元素消去法步骤。
2.编写用不选主元的直接三角分解法解线性方程组Ax=b 的M 文件。
要求输出Ax=b 中矩阵A 及向量b ,A=LU 分解的L 与U ,det A 及解向量x 。
3.编写用Gauss 列主元素消去法解线性方程组Ax=b 的M 文件。
要求输出Ax=b 中矩阵A 及向量b 、PA=LU 分解的L 与U 、det A 及解向量x ,交换顺序。
4.给定方程组(1) ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--11134.981.4987.023.116.427.199.103.601.3321x x x(2) ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡----15900001.582012151********.23107104321x x x x 先用编写的程序计算,再将(1)中的系数3.01改为3.00,0.987改为0.990;将(2)中的系数2.099999改为2.1,5.900001改为9.5,再用Gauss 列主元素消去法解,并将两次计算的结果进行比较。
【实验仪器与软件】1.CPU 主频在1GHz 以上,内存在128Mb 以上的PC ;2.Matlab 6.0及以上版本。
实验讲评:实验成绩:评阅教师:200 年 月 日Lab06.Gauss 列主元素消去法实验第一题:1、算法描述:Ⅰ、Gauss 消去法由书上定理5可知 设Ax=b ,其中A ∈R^(n(1)如果()0(1,2,....,1)k kka k n ≠=-,则可通过高斯消去法将Ax=b 约化为等价的 角形线性方程组,且计算公式为:① 消元计算(k=1,2,….,n-1)()()(1)()()(1)()()/,1,...,,,,1,...,,,1,...,.k k ik ik kk k k k ij ij ik kj k k k iiik k m a a i k n a a m a i j k n b b m b i k n ++==+=-=+=-=+② 回带公式()()()()()1/,()/,1,...,2,1.n n n n nn ni i i i iii j ii j i x b a x ba x a i n =+==-=-∑(2)如果A 为非奇异矩阵,则可通过高斯消去法将方程组Ax=b 约化方程组为上三角矩阵以上消元和回代过程总的乘除法次数为332333nn nn +-≈,加减法次数为32353263nnn n+-≈以上过程就叫高斯消去法。
Ⅱ、Gauss 列主元素消去法 设Ax=b.本算法用A 的具有行交换的列主元素消去法,消元结果冲掉A ,乘法ik m 冲掉ij a ,计算解x 冲掉常数项b ,行列式存放在det 中。
1、 det ←12、对于k=1,2,…,n-1 (1) 按列主元,max ik ikk i na k a ≤≤=(2) 如果,0ik a k =,则停止(det (A )=0)(3) 如果k i =k 则转(4)换行:1*(*)/*det *det 0/ik kkij ij ik kjni i ij j ii ik kj i nn nn n n nnm a a a m a b b a b a m b a a b b a =+←←←←-←←=←∑(4) 消元计算对于 i=k+1,…,n ① ik a ←ik m =ik a /kk a② 对于j=k+1,…,n*ij ij ik kj a a m a ←←③ *i i ik k b b m b ←←(5) d e t *d e tkka ← 3、如果0nn a =,则计算停止(det (A )=0) 4、回代求解 (1)/n n nnb b a ← (2)对于i=n-1,…,2,11(*)/ni i ij j ii j i b b a b a =+←-∑5、det *det nn a ←第二题:编写用不选主元的直接三角分解法解线性方程组Ax=b 的M 文件。
要求输出Ax=b 中矩阵A 及向量b ,A=LU 分解的L 与U ,det A 及解向量x 。
编写M 文件如下:function [x,L,U]=Doolittle(A,b) N = size(A); n = N(1);L = eye(n,n); U = zeros(n,n);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i)endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy = SolveDownTriangle(L,b);x = SolveUpTriangle(U,y);function y=SolveDownTriangle(L,b)N=size(L);n=N(1);for i=1:nif(i>1)s=L(i,1:(i-1))*y(1:(i-1),1);elses=0;endy(i,1)=(b(i)-s)/L(i,i);endfunction x=SolveUpTriangle(U,y)N=size(U);n=N(1);for i=n:-1:1if(i<n)s=U(i,(i+1):n)*x((1+i):n,1);elses=0;endx(i,1)=(y(i)-s)/U(i,i);end第三题:编写用Gauss列主元素消去法解线性方程组Ax=b的M文件。
要求输出Ax=b中矩阵A及向量b、PA=LU分解的L与U、det A及解向量x,交换顺序。
编写M文件如下:function [x,det,index]=gauss(A,b)[n,m]=size(A);nb=length(b);if n~=merror('The rows and columns of natrix A must be equal !'); return;endif m~=nberror('The columns of A must be equal the length b !'); return;endindex=1;det=1;x=zeros(n,1);for k=1:n-1maxa=0;for i=k:nif abs(A(i,k))>maxamaxa=abs(A(i,k));r=i;endendif maxa<1e-10index =0;return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n);if abs(A(n,n))<1e-10index=0;return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);end第四题:Ⅰ、运用Gauss列主元素消去法解上述程序计算矩阵:⑴、计算第一个矩阵①、计算源程序:计算程序如下:function [x,det,index]=gauss(A,b)A=[3.01,6.03,1.99;1.27,4.16,-1.23;0.987,-4.81,9.34];b=[1;1;1];[n,m]=size(A);nb=length(b);if n~=merror('The rows and columns of natrix A must be equal !');return;endif m~=nberror('The columns of A must be equal the length b !');return;endindex=1;det=1;x=zeros(n,1);for k=1:n-1maxa=0;for i=k:nif abs(A(i,k))>maxamaxa=abs(A(i,k));r=i;endendif maxa<1e-10index =0;return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n);if abs(A(n,n))<1e-10index=0;return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);end运行结果如下:ans =1.0e+003 *1.5926-0.6319-0.4936>>②、计算修改后的程序:计算程序如下:function [x,det,index]=gauss(A,b)A=[3.00,6.03,1.99;1.27,4.16,-1.23;0.990,-4.81,9.34];b=[1;1;1]; [n,m]=size(A);nb=length(b);if n~=merror('The rows and columns of natrix A must be equal !'); return;endif m~=nberror('The columns of A must be equal the length b !'); return;endindex=1;det=1;x=zeros(n,1);for k=1:n-1maxa=0;for i=k:nif abs(A(i,k))>maxamaxa=abs(A(i,k));r=i;endendif maxa<1e-10index =0;return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n);if abs(A(n,n))<1e-10index=0;return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);end计算结果如下:ans =119.5273-47.1426-36.8403>>⑵:计算第二个矩阵:①:计算源程序:计算程序如下:function [x,det,index]=gauss(A,b)A=[10,-7,0,1;-3,2.09999,6,2;5,-1,5,-1;2,1,0,2];b=[8;5.900001;5;1]; [n,m]=size(A);nb=length(b);if n~=merror('The rows and columns of natrix A must be equal !'); return;endif m~=nberror('The columns of A must be equal the length b !');return;endindex=1;det=1;x=zeros(n,1);for k=1:n-1maxa=0;for i=k:nif abs(A(i,k))>maxamaxa=abs(A(i,k));r=i;endendif maxa<1e-10index =0;return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z; endz=b(k);b(k)=b(r);b(r)=z;det=-det; endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n);if abs(A(n,n))<1e-10index=0;return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);end计算结果如下:ans =0.0000-1.00001.00001.0000>>②、计算修改后的程序:计算程序如下:function [x,det,index]=gauss(A,b)A=[10,-7,0,1;-3,2.1,6,2;5,-1,5,-1;2,1,0,2];b=[8;9.5;5;1]; [n,m]=size(A);nb=length(b);if n~=merror('The rows and columns of natrix A must be equal !'); return;endif m~=nberror('The columns of A must be equal the length b !'); return;endindex=1;det=1;x=zeros(n,1);for k=1:n-1maxa=0;for i=k:nif abs(A(i,k))>maxamaxa=abs(A(i,k));r=i;endendif maxa<1e-10index =0;return;endif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:nm=A(i,k)/A(k,k);for j=k+1:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);enddet=det*A(k,k);enddet=det*A(n,n);if abs(A(n,n))<1e-10index=0;return;endfor k=n:-1:1for j=k+1:nb(k)=b(k)-A(k,j)*x(j);endx(k)=b(k)/A(k,k);end计算结果如下:ans =-0.3543-1.42521.38271.5669>>Ⅱ、运用不选主元的直接三角分解法解上述程:⑴、计算第一个矩阵:①、计算源程序:计算程序如下:function [x,L,U]=Doolittle(A,b)A=[3.01,6.03,1.99;1.27,4.16,-1.23;0.987,-4.81,9.34];b=[1,1,1];N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i) endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy = SolveDownTriangle(L,b);x = SolveUpTriangle(U,y);function y=SolveDownTriangle(L,b)N=size(L);n=N(1);for i=1:nif(i>1)s=L(i,1:(i-1))*y(1:(i-1),1);elses=0;endy(i,1)=(b(i)-s)/L(i,i);endfunction x=SolveUpTriangle(U,y)N=size(U);n=N(1);for i=n:-1:1if(i<n)s=U(i,(i+1):n)*x((1+i):n,1);elses=0;endx(i,1)=(y(i)-s)/U(i,i);end计算结果如下:U =3.0100 6.0300 1.99000 1.6158 00 0 0U =3.0100 6.0300 1.99000 1.6158 -2.06960 0 0L =1.0000 0 00.4219 1.0000 00.3279 -4.2006 1.0000U =3.0100 6.0300 1.99000 1.6158 -2.06960 0 -0.0063ans =1.0e+003 *1.5926-0.6319-0.4936>>②、计算修改后的程序:计算程序如下:function [x,L,U]=Doolittle(A,b)A=[3.00,6.03,1.99;1.27,4.16,-1.23;0.990,-4.81,9.34];b=[1,1,1]; N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i) endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy = SolveDownTriangle(L,b);x = SolveUpTriangle(U,y);function y=SolveDownTriangle(L,b)N=size(L);n=N(1);for i=1:nif(i>1)s=L(i,1:(i-1))*y(1:(i-1),1);elses=0;endy(i,1)=(b(i)-s)/L(i,i);endfunction x=SolveUpTriangle(U,y)N=size(U);n=N(1);for i=n:-1:1if(i<n)s=U(i,(i+1):n)*x((1+i):n,1);elses=0;endx(i,1)=(y(i)-s)/U(i,i);end计算结果:U =3.0000 6.0300 1.99000 1.6073 00 0 0U =3.0000 6.0300 1.99000 1.6073 -2.07240 0 0L =1.0000 0 00.4233 1.0000 00.3300 -4.2306 1.0000U =3.0000 6.0300 1.99000 1.6073 -2.07240 0 -0.0844ans =119.5273-47.1426-36.8403>>⑵、计算第二个矩阵的程序:①、计算源程序:计算程序如下:function [x,L,U]=Doolittle(A,b)A=[10,-7,0,1;-3,2.09999,6,2;5,-1,5,-1;2,1,0,2];b=[8,5.900001,5,1]; N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i)endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy = SolveDownTriangle(L,b);x = SolveUpTriangle(U,y);function y=SolveDownTriangle(L,b)N=size(L);n=N(1);for i=1:nif(i>1)s=L(i,1:(i-1))*y(1:(i-1),1);elses=0;endy(i,1)=(b(i)-s)/L(i,i);endfunction x=SolveUpTriangle(U,y)N=size(U);n=N(1);for i=n:-1:1if(i<n)s=U(i,(i+1):n)*x((1+i):n,1);elses=0;endx(i,1)=(y(i)-s)/U(i,i);end计算结果如下:U =10.0000 -7.0000 0 1.00000 -0.0000 0 00 0 0 00 0 0 0U =10.0000 -7.0000 0 1.00000 -0.0000 6.0000 00 0 0 00 0 0 0U =10.0000 -7.0000 0 1.0000 0 -0.0000 6.0000 2.3000 0 0 0 00 0 0 0L =1.0e+005 *0.0000 0 0 0 -0.0000 0.0000 0 0 0.0000 -2.5000 0.0000 0 0.0000 0 0 0.0000L =1.0e+005 *0.0000 0 0 0-0.0000 0.0000 0 0 0.0000 -2.5000 0.0000 0 0.0000 -2.4000 0 0.0000U =1.0e+006 *0.0000 -0.0000 0 0.0000 0 -0.0000 0.0000 0.0000 0 0 1.5000 00 0 0 0U =1.0e+006 *0.0000 -0.0000 0 0.00000 -0.0000 0.0000 0.00000 0 1.5000 0.57500 0 0 0L =1.0e+005 *0.0000 0 0 0-0.0000 0.0000 0 00.0000 -2.5000 0.0000 00.0000 -2.4000 0.0000 0.0000U =1.0e+006 *0.0000 -0.0000 0 0.00000 -0.0000 0.0000 0.00000 0 1.5000 0.57500 0 0 0.0000ans =0.0000-1.00001.00001.0000>>②、计算修改后的程序:计算程序如下:function [x,L,U]=Doolittle(A,b)A=[10,-7,0,1;-3,2.1,6,2;5,-1,5,-1;2,1,0,2];b=[8,9.5,5,1]; N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i) endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k) endendy = SolveDownTriangle(L,b);x = SolveUpTriangle(U,y);function y=SolveDownTriangle(L,b)N=size(L);n=N(1);for i=1:nif(i>1)s=L(i,1:(i-1))*y(1:(i-1),1);elses=0;endy(i,1)=(b(i)-s)/L(i,i);endfunction x=SolveUpTriangle(U,y)N=size(U);n=N(1);for i=n:-1:1if(i<n)s=U(i,(i+1):n)*x((1+i):n,1);elses=0;endx(i,1)=(y(i)-s)/U(i,i);end计算结果如下:U =10 -7 0 10 0 0 00 0 0 00 0 0 0U =10 -7 0 10 0 6 00 0 0 00 0 0 0U =10.0000 -7.0000 0 1.00000 0 6.0000 2.30000 0 0 00 0 0 0Warning: Divide by zero.> In hghj465 at 14L =1.0000 0 0 0-0.3000 1.0000 0 00.5000 Inf 1.0000 00.2000 0 0 1.0000Warning: Divide by zero.> In hghj465 at 14L =1.0000 0 0 0-0.3000 1.0000 0 00.5000 Inf 1.0000 00.2000 Inf 0 1.0000U =10.0000 -7.0000 0 1.00000 0 6.0000 2.30000 0 -Inf 00 0 0 0 U =10.0000 -7.0000 0 1.00000 0 6.0000 2.30000 0 -Inf -Inf0 0 0 0 L =1.0000 0 0 0-0.3000 1.0000 0 00.5000 Inf 1.0000 00.2000 Inf NaN 1.0000 U =10.0000 -7.0000 0 1.00000 0 6.0000 2.30000 0 -Inf -Inf0 0 0 NaNWarning: Divide by zero.> In hghj465>SolveUpTriangle at 39In hghj465 at 18ans =NaNNaNNaNNaN >>。