数值线性代数大作业报告

合集下载

数值线性代数第二版徐树方高立张平文上机习题第一章实验报告(供参考)

数值线性代数第二版徐树方高立张平文上机习题第一章实验报告(供参考)

上机习题1.先用你所熟悉的的计算机语言将不选主元和列主元Gauss 消去法编写成通用的子程序;然后用你编写的程序求解84阶方程组;最后将你的计算结果与方程的精确解进行比较,并就此谈谈你对Gauss 消去法的看法。

Sol :(1)先用matlab 将不选主元和列主元Gauss 消去法编写成通用的子程序,得到P U L ,,: 不选主元Gauss 消去法:[])(,A GaussLA U L =得到U L ,满足LU A =列主元Gauss 消去法:[])(,,A GaussCol P U L =得到P U L ,,满足LU PA =(2)用前代法解()Pb or b Ly =,得y用回代法解y Ux =,得x求解程序为()P U L b A Gauss x ,,,,=(P 可缺省,缺省时默认为单位矩阵)(3)计算脚本为ex1_1代码%算法(计算三角分解:Gauss 消去法)function [L,U]=GaussLA(A)n=length(A);for k=1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);endU=triu(A);L=tril(A);L=L-diag(diag(L))+diag(ones(1,n));end%算法计算列主元三角分解:列主元Gauss消去法)function [L,U,P]=GaussCol(A)n=length(A);for k=1:n-1[s,t]=max(abs(A(k:n,k)));p=t+k-1;temp=A(k,1:n);A(k,1:n)=A(p,1:n);A(p,1:n)=temp;u(k)=p;if A(k,k)~=0A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n); elsebreak;endendL=tril(A);U=triu(A);L=L-diag(diag(L))+diag(ones(1,n));P=eye(n);for i=1:n-1temp=P(i,:);P(i,:)=P(u(i),:);P(u(i),:)=temp;endend%高斯消去法解线性方程组function x=Gauss(A,b,L,U,P)if nargin<5P=eye(length(A));endn=length(A);b=P*b;for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j); endb(n)=b(n)/L(n,n);y=b;for j=n:-1:2y(j)=y(j)/U(j,j);y(1:j-1)=y(1:j-1)-y(j)*U(1:j-1,j);endy(1)=y(1)/U(1,1);x=y;endex1_1clc;clear;%第一题A=6*eye(84)+diag(8*ones(1,83),-1)+diag(ones(1,83),1);b=[7;15*ones(82,1);14];%不选主元Gauss消去法[L,U]=GaussLA(A);x1_1=Gauss(A,b,L,U);%列主元Gauss消去法[L,U,P]=GaussCol(A);x1_2=Gauss(A,b,L,U,P);%解的比较subplot(1,3,1);plot(1:84,x1_1,'o-');title('Gauss');subplot(1,3,2);plot(1:84,x1_2,'.-');title('PGauss');subplot(1,3,3);plot(1:84,ones(1,84),'*-');title('精确解');结果为(其中Gauss表示不选主元的Gauss消去法,PGauss表示列主元Gauss 消去法,精确解为[]'⨯8411,,1 ):-6-4-202468Gauss050100PGauss 00.20.40.60.811.21.41.61.82精确解由图,显然列主元消去法与精确解更为接近,不选主元的Gauss 消去法误差比列主元消去法大,且不如列主元消去法稳定。

解线性方程组的迭代法数值计算上机实习报告

解线性方程组的迭代法数值计算上机实习报告

解线性方程组的迭代法数值计算上机实习报告一.综述:考虑用迭代法求解线性方程组,取真解为,初始向量取为零,以范数为度量工具,取误差指标为.其中。

分别完成下面各小题:第六题:编制程序实现Jacobi迭代方法和Gauss-Seidel 方法。

对应不同的停机标准(例如残量,相邻误差,后验误差停机标准),比较迭代次数以及算法停止时的真实误差。

其中残量准则:、相邻误差准则:后验误差停机准则:解:为了结果的可靠性,这里我分别对矩阵阶数为400、2500、10000进行试验,得到对应不同的方法、取不同的停机标准,迭代次数和真实误差的数据如下:分析上面数据可知,对应不同的停机标准,GS方法的迭代次数都近似为J方法的一半,这与理论分析一致。

而且从迭代次数可以看出,在这个例子中,作为停机标准,最强的依次为后验误差,再到残量,再到相邻误差。

第七题:编写程序实现SOR 迭代方法。

以真实误差作为停机标准,数值观测SOR 迭代方法中松弛因子对迭代次数的影响,找到最佳迭代因子的取值。

解:本题中考虑n=50,即对2500阶的矩阵A。

由于我们已经知道要使SOR方法收敛,松弛因子需要位于。

下面来求SOR方法中对应的最佳松弛因子。

应用筛选法的思想,第一次我们取松弛因子,间距为0.05,得到的对应的图像如下,从图中可以看出迭代次数随着的增大,先减小后变大,这与理论相符。

同时可以看出最佳松弛因子.第二次将区间细分为10份,即取,可得下面第二幅图像,从图像中可以看出最佳松弛因子第八题:对于J 方法,GS方法和(带有最佳松弛因子的)SOR 方法,分别绘制误差下降曲线以及残量的下降曲线(采用对数坐标系),绘制(按真实误差)迭代次数与矩阵阶数倒数的关系;解:对于J方法,考虑n=50时,采用相邻误差为迭代的终止条件,误差下降曲线及残量的下降曲线如下:对于GS方法,考虑n=50的时候,采用相邻误差作为迭代的终止条件,所得到的残量和误差的下降曲线如下图:从中可以看出,当相邻误差满足误差指标时,真实误差却并不小于误差指标,而为2.6281e-04。

线性代数实验报告[1].doc

线性代数实验报告[1].doc

线性代数实验报告
专业班级姓名学号
实验日期年月日星期
成绩评定教师签名批改日期
题目1:交通流量问题:
下图给出某城市部分街道的交通流量(单位:辆/小时):
假设:(1)全部流入网络的流量等于全部流出网络的流量;
(2)全部流入一个节点的流量等于全部流出此节点的流量. 试建立数学模型,以确定该交通网络未知部分的具体流量.
(要求:1. 模型建立(即:列出线性方程组),2. 求解,3. 输出结果,4. 结果综述.)
题目2:求一个正交变换,将二次型:434241312
1242322211262421993x x x x x x x x x x x x x x f --++-+++=
化为标准型 ,判断此二次型的正定性。

(完整word版)数值线性代数第二版徐树方高立张平文上机习题第四章实验报告

(完整word版)数值线性代数第二版徐树方高立张平文上机习题第四章实验报告

第四章上机习题1考虑两点边值问题⎪⎩⎪⎨⎧==<<=+.1)1(,0)0(10 ,22y y a a dx dy dx y d ε 容易知道它的精确解为ax e e ay x +---=--)1(111εε为了把微分方程离散化,把[0,1]区间n 等分,令h=1/n ,1,,1,-==n i ih x i得到差分方程,21211a hy y h y y y i i i i i =-++-++-ε简化为 ,)2()(211ah y y h y h i i i =++-+-+εεε从而离散化后得到的线性方程组的系数矩阵为⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡+-++-++-++-=)2()2()2()2(h h h h h h h A εεεεεεεεεε 对,100,2/1,1===n a ε分别用Jacobi 迭代法,G-S 迭代法和SOR 迭代法求线性方程组的解,要求有4位有效数字,然后比较与精确解得误差。

对,0001.0,01.0,1.0===εεε考虑同样的问题。

解 (1)给出算法:为解b Ax =,令U L D A --=,其中][ij a A =,),,,(2211nn a a a diag D = ,⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡------=-00001,21323121n n n n a a a a a a L,⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡------=-0000,122311312 n n n n a a a a a a U 利用Jacobi 迭代法,G-S 迭代法,SOR 迭代法解线性方程组,均可以下步骤求解: step1给定初始向量x0=(0,0,...,0),最大迭代次数N ,精度要求c ,令k=1 step2令x=B*x0+gstep3若||x-x0||2<c ,算法停止,输出解和迭代次数k ,否则,转step4 step4若k>=N,算法停止,迭代失败,否则,令x0=x ,转step2在Jacobi 迭代法中,B=D -1*(L+U),g=D -1*b在G-S 迭代法中,B=D -1*(L+U),g=D -1*b在SOR 迭代法中,B=(D-w*L)-1*[(1-w)*D+w*U],g=w*(D-w*L)-1*b另外,在SOR 迭代法中,上面算法step1中要给定松弛因子w ,其中0<w<2 为计算结果,规定w=0.5。

线性代数报告

线性代数报告

建造环境与能源应用工程 1 班陈嘉威 3013214105杜澎磊 3013214106宋子旭 3013214127近几十年来,随着科学技术的发展,特殊是计算机技术的发展,数学的应用领域已由传统的物理领域(包括力学、电子等学科以及土木、机电等工程技术) 迅速扩展到非物理领域(人口、经济、金融、生物、医学等)。

数学在发展高科技、提高生产力水平和实现现代化管理等方面的作用越来越明显。

这就要求我们如何将实际问题经过分析、简化,转化为一个数学问题,然后用一个适当的数学方法去解决。

线性代数是一个数学分支,是代数的一个重要学科它对于培养学生严谨的逻辑推理和抽象思维能力起着不可或者缺的作用。

线性代数研究最多的是矩阵。

矩阵是一个数表,而这个数表可以进行变换,以形成新的数表。

也就是说如果抽象出某种变化规律,就可以用代数的理论对研究的数表进行变换,并得出想要的一些结论。

所以,矩阵是一种方便的计算工具,可以以简单的形式表示复杂的公式,比如数字图象处理、计算机图形学、计算几何学、人工智能、网络通信以及普通的算法设计和分析等。

因此,矩阵的应用日益广泛,不少领域都要用到矩阵的知识。

本文将要探讨的,就是矩阵在实际生活中的一些应用形式。

经过分析和筛选,本文将从以下三个方面展开论述:可逆矩阵在保密通信中的应用,矩阵与成本利润的计算以及矩阵与数字图象。

随着计算机与网络技术的迅猛发展,通信技术中的保密工作显得尤其重要, 怎样确保通信过程中信息的安全变得至关重要, 因此大量各具特色的密码体系不 断涌现。

矩阵作为线性代数的重要组成部份, 其应用领域也从传统的物理领域迅 速扩展到非物理领域,特别是在保密通信中发挥着重要作用。

矩阵的定义: m 行n 列的矩形数表称为m 行n 列矩阵,简称m × n 矩阵,矩阵用大写黑体字母A ,B ,C ,…表示。

如: A=[a 1a 2…a a 1 1a 12……a22… ] 这m × n 个数称为矩阵m1 m2 mnA 的元素, a ij 称为矩阵A 的第i 行第j 列元素,一个m × n 矩阵A 也可简记为A = (a ij ) m×n 或者 A m×n 。

数值代数实验报告(1)绝对经典

数值代数实验报告(1)绝对经典

数值代数实验报告Numerical Linear Algebra And ItsApplications学生所在学院:理学院学生所在班级:计算数学10-1学生姓名:戈东潮指导教师:于春肖教务处2012年12月实验一实验名称: Poisson 方程边值问题的五点差分格式 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件实现Poisson 方程的边值问题的五点差分格式的线性代数方程组来认真解读Poisson 方程边值问题的具体思想与方法,使我们掌握得更加深刻,并且做到学习与使用并存,增加学习的实际动手性,不再让学习局限于书本和纸上,而是利用计算机学习来增加我们的学习兴趣。

二、实验内容:利用Poisson 方程来解下列问题:⎪⎩⎪⎨⎧∂∈=∈=∂∂+∂∂-ΩΩ),(,),(),(,sin sin )(y x y x u y x y x y u x u 0222222πππ 其中}{的边界是ΩΩ∂<<∈Ω,1,0),(y x y x 。

边值问题的解释是y x y x u ππsin sin ),(=(1)用例题中模型问题的方法取N=5,10(也可以取N=20)列出五点差分格式的线性代数方程组三、实验过程:系数矩阵A 的Matlab 实现函数如下:%文件名:poisson.m function T=poisson(N) for i=1:N a(i)=4; end b=diag(a); for j=1:N-1b(j,j+1)=-1;endfor j=2:Nb(j,j-1)=-1;endfor i=1:NI(i)=-1;endI=diag(I);for i=1:N:N*Nfor j=1:NT(i+j-1,i:i+N-1)=b(j,:);endendfor i=1:N:N*N-Nfor j=1:NT(i+j-1,i+N:i+N-1+N)=I(j,:);endendfor i=1+N:N:N*Nfor j=1:NT(i+j-1,i-N:i-1)=I(j,:);endend求解常数项b的Matlab实现函数如下:%函数名:constant.mfunction b=constant(N)n=N+1;h=1/n;i=1;for y=1/n:1/n:N/nfor x=1/n:1/n:N/nf(i)=2*pi*pi*sin(pi*x)*sin(pi*y);i=i+1;endendb=h*h*f;b=b';四、实验结果(总结/方案)在Matlab运行窗口输入>>A=poisson(5) Enter输出结果A =Columns 1 through 114 -1 0 0 0 -1 0 0 0 0 0-1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 -1 0 0 0 -1 0 00 0 0 -1 4 0 0 0 0 -1 0-1 0 0 0 0 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 -1 00 0 0 0 -1 0 0 0 -1 4 00 0 0 0 0 -1 0 0 0 0 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 00 0 0 0 0 0 0 0 0 -1 00 0 0 0 0 0 0 0 0 0 -10 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 Columns 12 through 220 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 00 -1 0 0 0 0 0 0 0 0 00 0 -1 0 0 0 0 0 0 0 00 0 0 -1 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 0 0 04 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 0 0 0 0 -1 0 00 0 0 0 4 -1 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 0 00 0 0 0 -1 0 0 0 0 4 -10 0 0 0 0 -1 0 0 0 -1 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 0 Columns 23 through 250 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 4 -1 0 -1 4 -10 -1 4 在运行窗口输入 >>b=constant(5) Enter 常数项的输出结果b =[ 0.1371 0.2374 0.2742 0.2374 0.1371 0.2374 0.4112 0.4749 0.4112 0.2374 0.2742 0.4749 0.5483 0.4749 0.2742 0.2374 0.4112 0.4749 0.4112 0.2374 0.1371 0.2374 0.2742 0.2374 0.1371]T所以A*u=b 的五点差分格式为j i j i j i j i j i j i f h u u u u u ,,,,,,211114=-----+-+ (i,j=1,2……)实验二实验名称: 用Jacobi 迭代法和SOR 迭代法求解方程组 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件采用Jacobi 迭代法和SOR 迭代法来求解方程组,求解过程中了解两种迭代方法的基本思想与迭代过程,并分析两种迭代方法的收敛性和实用用以及他们的异同点。

线性代数学习报告doc

线性代数学习报告doc

线性代数学习报告篇一:浅谈学习线性代数的心得体会沈阳药科大学选修课结课论文沈阳药科大学浅谈学习线性代数的心得体会学校:沈阳药科大学姓名:郑亚娟学号:10106331 专业:药物制剂年级:XX级班级:03班一、内容摘要线性代数是一门较抽象的数学课程,但是线性代数除了其抽象之外还具有另外一个重要的特点:“实用性”,由于计算机的飞速发展和广泛应用,线性代数已成为越来越多的科技工作者必不可少的数学工具。

掌握线性代数的基本概念、基本理论与基本方法,为解决工科各专业的实际问题,为进一步学习相关课程及扩大数学知识都将奠定必要的数学基础。

在初步学习了高等数学这门课程后,里面涉及了一些线性代数的求解方法,听老师说,某些题目用线性代数的方法求解更容易,但是由于我们还未系统的学习这门课程,老师也是一带而过,并未深讲。

致使我对线性代数这门学科有了浓厚的兴趣,在首先简单了解了这门学科的背景后,发现线性代数是一门丰富多彩充满未知的科学,在看到学校开设了这门课程的选修课后,我义无反顾的叫我们全寝室的人都选修了这门奇妙的课程。

学习线性代数的初步感受就是它的概念多,推理论证多,基本理论与结论多,线性代数在内容上,思想方法上及论证方法上都与“高等数学”有所区别。

它具有较强的逻辑性和抽象性,一开始就要高度重视。

它又与中学所学的代数有一定的联系,所以有些内容并不是完全陌生的。

我相信只要我每节每章地,一步一个脚印的弄懂、弄通,记住有关的概念和结论,并通过反复的应用(练习)来掌握它,循序渐进掌握这门课程是容易的。

关键词:数学线性代数背景应用计算方法感受二、绪论2.1 线性代数的发展史由于费马和笛卡儿的工作,线性代数基本上出现于十七世纪。

直到十八世纪末,线性代数的领域还只限于平面与空间。

十九世纪上半叶才完成了到n维向量空间的过渡,矩阵论始于凯莱,在十九世纪下半叶,因若当的工(本文来自: 小草范文网:线性代数学习报告)作而达到了它的顶点。

1888年,皮亚诺以公理的方式定义了有限维或无限维向量空间。

数值线性代数第二版徐树方高立张平文上机习题第一章实验报告

数值线性代数第二版徐树方高立张平文上机习题第一章实验报告

@上机习题1.先用你所熟悉的的计算机语言将不选主元和列主元Gauss 消去法编写成通用的子程序;然后用你编写的程序求解84阶方程组;最后将你的计算结果与方程的精确解进行比较,并就此谈谈你对Gauss 消去法的看法。

Sol :(1)先用matlab 将不选主元和列主元Gauss 消去法编写成通用的子程序,得到P U L ,,: 不选主元Gauss 消去法:[])(,A GaussLA U L =得到U L ,满足LU A = 列主元Gauss 消去法:[])(,,A GaussCol P U L =得到P U L ,,满足LU PA = (2)用前代法解()Pb or b Ly =,得y用回代法解y Ux =,得x]求解程序为()P U L b A Gauss x ,,,,=(P 可缺省,缺省时默认为单位矩阵)(3)计算脚本为ex1_1 代码%算法(计算三角分解:Gauss 消去法) function [L,U]=GaussLA(A) n=length(A);—for k=1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n); endU=triu(A);L=tril(A);L=L-diag(diag(L))+diag(ones(1,n));end!%算法计算列主元三角分解:列主元Gauss消去法)function [L,U,P]=GaussCol(A)n=length(A);for k=1:n-1[s,t]=max(abs(A(k:n,k)));p=t+k-1;temp=A(k,1:n);¥A(k,1:n)=A(p,1:n);A(p,1:n)=temp;u(k)=p;if A(k,k)~=0A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n); elsebreak;^endendL=tril(A);U=triu(A);L=L-diag(diag(L))+diag(ones(1,n)); P=eye(n);for i=1:n-1temp=P(i,:);P(i,:)=P(u(i),:);{P(u(i),:)=temp;endend%高斯消去法解线性方程组function x=Gauss(A,b,L,U,P)if nargin<5P=eye(length(A));¥endn=length(A);b=P*b;for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j);endb(n)=b(n)/L(n,n);<y=b;for j=n:-1:2y(j)=y(j)/U(j,j);y(1:j-1)=y(1:j-1)-y(j)*U(1:j-1,j);endy(1)=y(1)/U(1,1);x=y;|endex1_1clc;clear;%第一题A=6*eye(84)+diag(8*ones(1,83),-1)+diag(ones(1,83),1); b=[7;15*ones(82,1);14];%不选主元Gauss消去法)[L,U]=GaussLA(A);x1_1=Gauss(A,b,L,U);%列主元Gauss消去法[L,U,P]=GaussCol(A);x1_2=Gauss(A,b,L,U,P);%解的比较subplot(1,3,1);plot(1:84,x1_1,'o-');title('Gauss'); subplot(1,3,2);plot(1:84,x1_2,'.-');title('PGauss');(subplot(1,3,3);plot(1:84,ones(1,84),'*-');title('精确解');结果为(其中Gauss 表示不选主元的Gauss 消去法,PGauss 表示列主元Gauss消去法,精确解为[]'⨯8411,,1 ):8Gauss50100PGauss精确解由图,显然列主元消去法与精确解更为接近,不选主元的Gauss 消去法误差比列主元消去法大,且不如列主元消去法稳定。

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

数值线性代数实验大报告指导老师:赵国忠姓名:1108300001 刘帅1108300004 王敏1108300032 郭蒙一、实验名称:16题P75上机习题二、实验目的:编制通用的子程序,完成习题的计算任务三、实验内容与要求:P75上机习题先用熟悉的计算机语言将算法2.5.1编制成通用的子程序,然后再用所编制的子程序完成下面两个计算任务:(1) 估计5到20阶Hilbert 矩阵的无穷范数条件数。

(2) 设A n = 11...111................1-1 (01)-- 先随机地选取x ∈R n ,并计算出b=An x;然后再用列主元Gauss 消去法求解该方程组,假定计算解为∧x .试对n 从5到30估计计算解∧x 的精度,并且与真实的相对误差作比较。

四、 实验原理:(1)矩阵范数(martix norm )是数学上向量范数对矩阵的一个自然推广。

利用for循环和cond (a )Hilbert 求解Hilbert 矩阵的无穷范数,再利用norm(a,inf)求矩阵的无穷范数条件数。

(2)本题分为4步来求解。

先运用rand 随机选取x ∈R n,输入A n 矩阵,编制一个M 文件计算出b 。

第二步用列主元高斯消去法求解出方程的解X2。

第三步建立M 文件: soluerr.m 估计计算解∧x 的精度。

第四步, 建立M 文件: bijiao.m ,与真实相对误差作比较。

五、 实验过程:(1)程序:clearfor n=5:20for i=1:nfor j=1:na(i,j)=1/(i+j-1);endendc=cond(a);f=norm(c,inf);fprintf('n=%3.0f\nnorm(c,inf)%e\n',n,f) end运行结果:n= 5norm(c,inf)4.766073e+005n= 6norm(c,inf)1.495106e+007n= 7norm(c,inf)4.753674e+008n= 8norm(c,inf)1.525758e+010n= 9norm(c,inf)4.931542e+011n= 10norm(c,inf)1.602467e+013n= 11norm(c,inf)5.224376e+014n= 12norm(c,inf)1.698855e+016n= 13norm(c,inf)3.459404e+017n= 14norm(c,inf)4.696757e+017n= 15norm(c,inf)2.569881e+017n= 16norm(c,inf)7.356249e+017n= 17norm(c,inf)4.362844e+017n= 18norm(c,inf)1.229633e+018n= 19norm(c,inf)9.759023e+017n= 20norm(c,inf)1.644051e+018(2)程序:M文件:matrix1.mfunction [a,b,x1]=matrix1(n) format longA1=-1*ones(n,n)A2=tril(A1)for i=1:nA2(i,i)=1endA2(:,n)=1a=A2x1=rand(n,1)b=A2*x1end运行结果:>> A1 =-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1A2 =-1 0 0 0 0 -1 -1 0 0 0 -1 -1 -1 0 0 -1 -1 -1 -1 0 -1 -1 -1 -1 -1A2 =1 0 0 0 0 -1 -1 0 0 0 -1 -1 -1 0 0 -1 -1 -1 -1 0 -1 -1 -1 -1 -1A2 =1 0 0 0 0-1 -1 -1 0 0 -1 -1 -1 -1 0 -1 -1 -1 -1 -1A2 =1 0 0 0 0 -1 1 0 0 0 -1 -1 1 0 0 -1 -1 -1 -1 0 -1 -1 -1 -1 -1A2 =1 0 0 0 0 -1 1 0 0 0 -1 -1 1 0 0 -1 -1 -1 1 0 -1 -1 -1 -1 -1A2 =1 0 0 0 0 -1 1 0 0 0 -1 -1 1 0 0 -1 -1 -1 1 0 -1 -1 -1 -1 1A2 =1 0 0 0 1 -1 1 0 0 1 -1 -1 1 0 1 -1 -1 -1 1 1 -1 -1 -1 -1 1a =-1 1 0 0 1 -1 -1 1 0 1 -1 -1 -1 1 1 -1 -1 -1 -1 1x1 =0.8147236863931790.9057919370756190.1269868162935060.9133758561390190.632359246225410b =1.4470829326185890.723427496907850-0.961169560949882-0.301767337397875-2.128519049675914a =1 0 0 0 1 -1 1 0 0 1 -1 -1 1 0 1 -1 -1 -1 1 1 -1 -1 -1 -1 1b =1.4470829326185890.723427496907850-0.961169560949882-0.301767337397875-2.128519049675914x1 =0.8147236863931790.9057919370756190.1269868162935060.9133758561390190.632359246225410M文件:LZYgauss.mfunction[x2]=LZYgauss(a,b)format longn=length(a);x2=zeros(n,1);a=[a b];for k=1:n-1max=k;for i=k+1:nif a(i,k)>a(max,k)max=i;endendtemp=a(k,k:n+1);a(k,k:n+1)=a(max,k:n+1);a(max,k:n+1)=temp;for i=k+1:na(i,k)=-a(i,k)/a(k,k);a(i,k+1:n+1)=a(i,k+1:n+1)+a(i,k)*a(k,k+1:n+1);endendx2(n,1)=a(n,n+1)/a(n,n);for i=n-1:-1:1sum=0;for j=i+1:nsum=sum+x2(j,1)*a(i,j);endx2(i,1)=(a(i,n+1)-sum)/a(i,i);end运行结果:>> LZYgauss(a,b)ans =0.8147236863931790.9057919370756190.1269868162935060.913375856139020 0.632359246225410估计计算解x的精度:M文件: soluerr.mfunction [x,error]=soluerr(a,b)format long%估计计算解的精度% 算法:列主元Gauss消去法,其中% A --- 系数矩阵% b-右端项% index --- index=0表示计算成败;index=1表示计算成功%输出结果:error--本算法给出的计算解的估计% normA--逆矩阵无穷范数估计% rnorm--计算解的残量[n,m]=size(a); nb=length(b);if n~=merror('The rows and columns of matrix a must be equal!');return;endif m~=nberror('The columns of a must be equal the dimension of b!'); return;endindex=1;%列主元矩阵三角分解[L,U,u,index_col]=Gauss_col(a);%解下三角方程组Ly=Pb[y,index_low]=Gauss_low(L,b(u));%解上三角方程组Ux=y[x,index_upp]=Gauss_upp(U,y);%输出数值解xpause(0.3)%估计矩阵逆的无穷大范数normA=normAinv(L,U,u);%估计计算解的残量rnorm=norm(b-a*x,inf);%计算右端项bnorm=norm(b,inf);%计算矩阵A的范数Anorm=norm(a,inf);%计算解的精度error=normA*Anorm*rnorm/bnorm;运行结果:>> [x,error]=soluerr(a,b)x =0.8147236863931790.9057919370756190.1269868162935060.9133758561390200.632359246225410x =0.8147236863931790.9057919370756190.1269868162935060.9133758561390200.632359246225410error =5.215941218821592e-016X1与真实相对误差做比较M文件:bijiao.mfunction [x,x1,error,error1]=bijiao(a,b,n) [a,x1]=matrix1(n)[x,error]=soluerr(a,b)error1=abs((x-x1)/x1)运行结果:>>x =0.8147236863931790.9057919370756190.1269868162935060.9133758561390200.632359246225410x1 =0.9578935505663481.6132601689718680.629241553693582-0.799716694069468-1.770467991515150error =5.215941218821592e-016error1 =Columns 1 through 40 0 0 00 0 0 00 0 0 00 0 0 00 0 0 0Column 50.0808655478999350.3995939126190040.2836847318376260.9675930648949141.357170674226221六、结果分析:1、矩阵范数(martix norm)是数学上向量范数对矩阵的一个自然推广。

而条件数在一定程度上刻画了扰动对方方程组解的影响程度。

相关文档
最新文档