《计算方法矩阵直接法与迭代法》由会员分享，可在线阅读，更多相关《计算方法矩阵直接法与迭代法（20页珍藏版）》请在金锄头文库上搜索。
1、需要两份东西1实验程序：输入输出注释2实验报告：问题描述方法描述方案设计结果分析结论谢谢，麻烦写的详细些实验五 解线性方程组的直接方法实验 5.1 （主元的选取与算法的稳定性）问题提出：Gauss 消去法是我们在线性代数中已经熟悉的。但由于计算机的数值运算是在一个有限的浮点数集合上进行的，如何才能确保 Gauss 消去法作为数值算法的稳定性呢？Gauss 消去法从理论算法到数值算法，其关键是主元的选择。主元的选择从数学理论上看起来平凡，它却是数值分析中十分典型的问题。实验内容：考虑线性方程组nnRbAbx,编制一个能自动选取主元，又能手动选取主元的求解线性方程组的 Gauss消去过程。实验要求
2、：（1）取矩阵 ，则方程有解 。1457,6816MObA Tx)1,(*L取 n=10 计算矩阵的条件数。让程序自动选取主元，结果如何？（2）现选择程序中手动选取主元的功能。每步消去过程总选取按模最小或按模尽可能小的元素作为主元，观察并记录计算结果。若每步消去过程总选取按模最大的元素作为主元，结果又如何？分析实验的结果。（3）取矩阵阶数 n=20 或者更大，重复上述实验过程，观察记录并分析不同的问题及消去过程中选择不同的主元时计算结果的差异，说明主元素的选取在消去过程中的作用。（4）选取其他你感兴趣的问题或者随机生成矩阵，计算其条件数。重复上述实验，观察记录并分析实验结果。1） 首先编写 g
3、auss 消元法代码；输入为矩阵 A，向量 b，精度 ptol，输出为方程组的解 x。程序如下function out = gausle( A,b,flag,ptol )%UNTITLED 高斯消元法 可选择的主元消去法% A n x n 矩阵% b n x 1 % flag 标志主元为自动选取还是手动选取，% 0，自动（对角为主元） 1，选取最小模为主元，% 2，选取模最大， （误差最小）3 ，次最小的模为主元% ptol 精度% out 输出值，为 nx1 的解 Ax=b 的解% if nargin0,disp(请注意：因为 RA=RB，所以此方程组无解.)returnend% 消元过程f
4、or i=1:n-1 ;i;if flag=0;pivot=Ab(i,i);ri=i;elseif flag=1;% 按照每列模最小的选取pivot,min_index=min(abs(A(i:n,i);ri=i+min_index-1;elseif flag=2;% 按照模最大的选取 pivot,max_index=max(abs(A(i:n,i);ri=i+max_index-1;elseif flag=3 %方程最小非 0 数tA=A;pivot,min_index=min(abs(tA(i:n,i); while pivot=0;tA(min_index+i-1,i)=inf;pivo
5、t,min_index=min(abs(tA(i:n,i);end ri=i+min_index-1;endif (pivot=0)|(pivotptol) ;warning(系数矩阵奇异！ );return;endif ri=i; % 交换行tmp=A(i,:);A(i,:)=A(ri,:);A(ri,:)=tmp;t=b(i);b(i)=b(ri);b(ri)=t;endfor kk=i+1:nL(kk,i)=A(kk,i)/A(i,i);A(kk,i+1:n)=A(kk,i+1:n)-L(kk,i)*A(i,i+1:n);b(kk)=b(kk)-L(kk,i)*b(i);endendif
6、 A(n,n)=0warning(系数矩阵奇异！ );return;end% % 回代求解x=zeros(n,1);for k=n:-1:1 % k% A(k,k+1:n)% x(k+1:n)if k=nx(n)=b(n)/A(n,n);elsex(k)=(b(k)-sum( A(k,k+1:n).*x(k+1:n) ) )/(A(k,k);endendout=x; end主程序为N=10 时，方程解为1.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000A 的条件数如下：1-条件数 2.5575e+0032-条件数 1.727
7、6e+003无穷 条件数 2.5575e+003代码在 里，运行计算即可。Flag=0；，时为自动选取2） Flag=1，为选取最小模的值为主元，flag=2，为选取最大模为主元，代码是一样的，修改 flag 值和维数 n 即可。3） N=20 ,只需要更改 n，重复上述步骤即可。4）思考题一：（Vadermonde 矩阵）设，niniininnnxxbxxxxA002102212100MLML，其中， ，kxk ,10,.1（1）对 n=2,5,8，计算 A 的条件数；随 n 增大，矩阵性态如何变化？（2）对 n=20，解方程组 Ax=b；设 A 的最后一个元素有扰动 10-4，再求解Ax=
8、b（3）计算（2）扰动相对误差与解的相对偏差，分析它们与条件数的关系。（4）你能由此解释为什么不用插值函数存在定理直接求插值函数而要用拉格朗日或牛顿插值法的原因吗？（5）尝试做出范德蒙矩阵的预优。答 1） n=2，5，8 ；A 的条件数为如下所示，随 n 增大，条件数越大，矩阵的病态性越严重。1277 71010 101 22,()6;()=4;()6;. .59;().8;8;5;345ncodcondcodAn 代码为 sikaoti_111.m2,. e+6 * 0.7- .12536 0.9-8 1.432 0.69-87 .53 -0.14 .Txx） 无 扰 动 时 为 的 向 量
9、 ；有 扰 动 时 ， 0.- . 0.,sikaot_12.代 码 文 件 为3) 由此看来，当系数矩阵的条件数越大，则病态性越严重，因而系数矩阵很小相对误差也会让解产生很大的相对误差。4） 插值函数存在定理直接求插值函数而要用拉格朗日或牛顿插值法的原因是为了补偿计算中所产生偏差，而且迭代次数越多，偏差就越大。5） 预优，相关 MATLAB 函数提示：zeros(m,n) 生成 m 行，n 列的零矩阵ones(m,n) 生成 m 行，n 列的元素全为 1 的矩阵eye(n) 生成 n 阶单位矩阵rand(m,n) 生成 m 行,n 列(0,1)上均匀分布的随机矩阵diag(x) 返回由向量
10、x 的元素构成的对角矩阵tril(A) 提取矩阵 A 的下三角部分生成下三角矩阵triu(A) 提取矩阵 A 的上三角部分生成上三角矩阵rank(A) 返回矩阵 A 的秩det(A) 返回方阵 A 的行列式inv(A) 返回可逆方阵 A 的逆矩阵V， D=eig(A) 返回方阵 A 的特征值和特征向量norm(A,p) 矩阵或向量的 p 范数cond(A,p) 矩阵的条件数L，U，P=lu(A) 选列主元 LU 分解R=chol(X) 平方根分解Hi=hilb(n) 生成 n 阶 Hilbert 矩阵实验六 解线性方程组的迭代法实验 6.1（病态的线性方程组的求解）问题提出：理论的分析表明，求
11、解病态的线性方程组是困难的。实际情况是否如此，会出现怎样的现象呢？实验内容：考虑方程组 Hx=b 的求解，其中系数矩阵 H 为 Hilbert 矩阵，njijihHjinji ,21,1,)(, L这是一个著名的病态问题。通过首先给定解（例如取为各个分量均为 1）再计算出右端 b 的办法给出确定的问题。实验要求：（1）选择问题的维数为 6，分别用 Gauss 消去法、列主元 Gauss 消去法、J 迭代法、 GS 迭代法和 SOR 迭代法求解方程组，其各自的结果如何？将计算结果与问题的解比较，结论如何？（2）逐步增大问题的维数（至少到 100） ，仍然用上述的方法来解它们，计算的结果如何？计算
12、的结果说明了什么？（3）讨论病态问题求解的算法答：1）Guass 消元法的程序在 gausle.m 程序里，如下所示。 Flag=0，是自动消元，Guass 消去法的结果为：x=1，1，1，1，1，1flag=2 为列主元 Gauss 消去法法，Gauss 列主元消去法的结果为：x=1，1，1，1，1，1Jacobi 迭代法程序在 jacobi.m，结果为：Inf ,Inf,NaN,NaN,NaN,NaN 。也就是迭代发散，因为迭代阵 B 的谱半径大于 1，为 4.3085,所以此方程不适合用 Jacobi 方法迭GS 迭代法在 gauseidel.m，其结果为：0.99921.01250.9
13、5901.03061.02670.9715 迭代次数为 2016SOR 超级松弛法程序在 SOR.m 中，w 选择为 1.89 其结果为：0.99911.01800.91361.15860.88001.0304迭代次数为 517 次。主程序如下。% 实验 6.1% 实验 6.1（病态的线性方程组的求解）clc;clear;close all;n=6; % 矩阵维数 A=hilb(n);eps=1e-5;R=max(abs(eig(A) ;%disp(A 的谱半径为 ,num2str(R);b=sum(A,2);x_gauss=gausle(A,b,0)x_gauss_2=gausle(A,b,
14、2)% 编写 jacobi.m 和 gauseidel.mx_jaco_3,n3=jacobi(A,b,zeros(size(b),eps ) %( A,b,omega,x0,eps,N)x_gs_4,n4=gauseidel(A,b,zeros(size(b),eps )% SOR 法omega=1.89;x_SOR_5,n5=SOR(A,b,omega,zeros(size(b),eps)2) 将维数增加为 100，高斯消元法，结果为1.00001.00000.99931.02280.67543.5051-10.564434.3970-60.724877.6571-77.081786.19
15、87-72.494316.549923.9588-7.319232.2195-71.5173-9.4683102.1479-4.3990-49.2259-65.428477.29858.937810.7393-2.8585-94.704818.5097112.925253.2932-77.7178-45.60030.8528-20.4096-56.1519 106.6452-3.433045.76888.1846-62.144310.476826.8102-75.727076.4337-55.6616-22.277789.8267-63.1822111.9219-73.864928.5538-22.5285-44.4054-20.734745.2838-39.7352-57.6491101.373957.695647.3377105.3975-47.4092-159.0873-38.1213-73.969618.9134-69.1624118.195726.6726141.973669.7617-144.1574-13.454667.8139-44.4608-126.1802-39.
