实 验 报 告 课程名称 数值分析 实验项目 解线性方程组的直接法 专业班级 姓 名 学 号 指导教师 成 绩 日 期 月 日 一. 实验目的 1、掌握程序的录入和matlab的使用和操作； 2、了解影响线性方程组解的精度的因素——方法与问题的性态。 3、学会Matlab提供的“\\”的求解线性方程组。 二. 实验要求 1、按照题目要求完成实验内容； 2、写出相应的Matlab 程序； 3、给出实验结果(可以用表格展示实验结果)； 4、分析和讨论实验结果并提出可能的优化实验。 5、写出实验报告。 三. 实验步骤 1、用LU分解及列主元高斯消去法解线性方程组 ?7?10???32.099999a)?5?1??21?1??x1??8??????62??x2??5.900001???， 5?1??x3??5??????????02??x4??1?0输出Ax?b中系数A?LU分解的矩阵L和U，解向量x和det(A);用列主元法的行交换次序解向量x和求det(A)；比较两种方法所得结果。 2、用列主高斯消元法解线性方程组Ax?b。
?3.016.031.99??x1??1???????4.16?1.23??x2???1? （1）、?1.27?0.987?4.819.34??x??1????3????3.006.031.99??x1??1???????4.16?1.23??x2???1? （2）、?1.27?0.990?4.819.34??x??1????3???分别输出A,b,det(A)，解向量x，（1）中A的条件数。分析比较（1）、（2）的计算结果 3、线性方程组Ax?b的A和b分别为 ?10??7A??8??7?77??32????565??23?, b????610933????31?5910????8则解x?(1,1,1,1,)T. 用MATLAB内部函数求det(A)和A的所有特征值和cond(A)2. 若令 78.17.2??10??5??7.085.046, A??A???85.989.899???6.99599.98???求解(A??A)(x??x)?b，输出向量?x和?x2，从理论结果和实际计算两方面分析线性方程组Ax?b解的相对误差?x4、 希尔伯特矩阵Hn?(hij)?Rn?n2/x2以及A的相对误差?A2/A2的关系。 1， i?j?1，其中hij?（1）分别对n?2,3,?,6计算cond(Hn)?，分析条件数作为n的函数如何变化。（2）令x?(1,1,?,1,)T?Rn，计算bn?Hnx，然后用高斯消去法解线性方程组Hnx?bn求出x，计算剩余向量rn?bn?Hnx以及?x?x?x。分析当n增加时解x分量的有效位数如何随n变化，它与条件数有何关系？当n多大时x连一位有效数字也
没有了？ 将每种情形的两个结果进行表格对比，如： n=6时： GAUSS列主消去法求得的x 四、实验结果 五、讨论分析 （对上述算例的计算结果进行比较分析，主要说清matlab的算符与消去法的适用范围不同，自己补充） 六、改进实验建议 （自己补充） x的有效数字 1.列主元的高斯消去法 利用列主元的高斯消去法matlab程序源代码： 首先建立一个gaussMethod.m的文件，用来实现列主元的消去方法。 function x=gaussMethod(A,b) %高斯列主元消去法，要求系数矩阵非奇异的， % n = size(A,1); if abs(det(A))<= 1e-8 error('系数矩阵是奇异的'); return; end % for k=1:n ak = max(abs(A(k:n,k))); index = find(A(:,k)==ak); if length(index) == 0 index = find(A(:,k)==-ak); end %交换列主元 temp = A(index,:);
A(index,:) = A(k,:); A(k,:) = temp; temp = b(index);b(index) = b(k); b(k) = temp; %消元过程 for i=k+1:n m=A(i,k)/A(k,k); %消除列元素 A(i,k+1:n)=A(i,k+1:n)-m*A(k,k+1:n); b(i)=b(i)-m*b(k); end end %回代过程 x(n)=b(n)/A(n,n); for k=n-1:-1:1; x(k)=(b(k)-A(k,k+1:n)*x(k+1:n)')/A(k,k); end x=x'; end 然后调用gaussMethod函数，来实现列主元的高斯消去法。在命令框中输入下列命令： 输出结果如下： 利用LU分解法及matlab程序源代码： function [L,U]=myLU(A) %实现对矩阵A的LU分解，L为下三角矩阵 A[n,n]=size(A);
L=zeros(n,n); U=zeros(n,n); for i=1:n L(i,i)=1; end for k=1:n for j=k:n U(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j)'); end for i=k+1:n L(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)'))/U(k,k); end end 在命令框中输入下列命令： >> a=[10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2] a = 10.0000 -7.0000 0 1.0000 -3.0000 2.1000 6.0000 2.0000 5.0000 -1.0000 5.0000 -1.0000 2.0000 1.0000 0 2.0000 >> [l,u]=lu(a) l = 1.0000 0 0 0 -0.3000 -0.0000 1.0000 0 0.5000 1.0000 0 0 0.2000 0.9600 -0.8000 1.0000 u = 10.0000 -7.0000 0 1.0000 0 2.5000 5.0000 -1.5000 0 0 6.0000 2.3000 0 0 0 5.0800 >> b=[8 5.900001 5 1]' b = 8.0000 5.9000 5.0000 1.0000 >> y=l\\b y = 8.0000 1.0000 8.3000