数值分析上机题matlab线性方程组,数值分析上机实验报告 - 线性方程组部分实验题1...

s=A(i，(i+1):n)*x((i+1):n，1); else s=0; end

x(i，1)=(b(i)-s)/A(i，i);end

%Cholosky分解方法********************************************************************

----------------------------------------------------------------------------------------------------------

%利用对称正定矩阵之Cholesky分解求解线性方程组Ax=b

function x=Chol_Solve(A，b)

clear;clc

% A=input('输入A=');

% b=input('输入b=');n=input('请输入系数矩阵的维数n=');%定义hilbert矩阵

A=hilb(n);x_exact=ones(n，1);

b=A*x_exact;n=length(b);

l=Cholesky(A);x=ones(n，1);

y=ones(n，1);

y(1)=b(1)/l(1，1);%求Ly=b

for i=2:n z=0;

for k=1:i-1

z=z+l(i，k)*y(k); end

y(i)=(b(i)-z)/l(i，i);endx(n)=y(n)/l(n，n);%求L'x=y

for i=n-1:-1:1 z=0;

for k=i+1:n

z=z+l(k，i)*x(k); end

x(i)=(y(i)-z)/l(i，i);endx

wc=x-x_exact %误差计算

a=input('是否继续，继续请按1，结束请按0:');

if a==1

Chol_Solve(A，b);

else end

%对称正定矩阵之Cholesky分解法：

function L=Cholesky(A)

% n=input('请输入系数矩阵的维数n=');%定义hilbert矩阵

% A=hilb(n);n=length(A);

L=zeros(n);

for k=1:n

delta=A(k，k);

for j=1:k-1

delta=delta-L(k，j)^2; end

if delta<1e-10

return; end

L(k，k)=sqrt(delta);

for i=k+1:n

L(i，k)=A(i，k);

for j=1:k-1

L(i，k)=L(i，k)-L(i，j)*L(k，j); end

L(i，k)=L(i，k)/L(k，k); endend3.2.条件数的计算

根据2.3叙述原理，写出条件数的计算程序如下：

%tiaojianshu.m

m=input ('input m:=') ;

N=[1:m];

for i=1:mn=N(i);

H=hilb(n);

k=cond(H);

disp('矩阵的阶数')

disp(n)

disp('矩阵')

disp(H)

disp('矩阵的条件数')

disp(k)end4.结果与误差分析4.1.线性方程组计算结果

分别用Gauss法和Cholesky分解法计算线性方程组的结果，并进行误差分析结果如下：表 1

n=2、5、8、10、15 时线性方程组的结果nx分量

Gauss法

Cholesky分解法求解误差(1.0e-15 *)求解误差(1.0e-15 *)211.0000 0.44411.0000 0.444121.0000

-0.66611.0000

-0.7772nx分量

Gauss法

Cholesky分解法求解误差(1.0e-10 *)求解误差(1.0e-10 *)511.0000

-0.00211.0000

-0.000921.0000 0.0371.0000 0.016731.0000

-0.15481.0000

-0.073241.0000 0.22891.0000 0.111651.0000

-0.11021.0000

-0.055 nx分量

Gauss法

Cholesky分解法求解误差( 1.0e-06 *)求解误差(1.0e-06 *)811.0000 0.0000 1.0000 021.0000 0.0023 1.0000 0.000531.0000

-0.0304 1.0000

-0.006341.0000 0.1644 1.0000 0.035351.0000

-0.4418 1.0000

-0.097961.0000 0.6242 1.0000 0.14271.0000

-0.4437 1.0000

-0.103281.0000 0.1250 1.0000 0.0296nx分量

Gauss法

Cholesky分解法求解误差求解误差1011.0000 0.0000 NaNNaN21.0000 0.0000 NaN NaN31.0000 0.0000 NaN NaN41.0001 0.0001 NaN NaN50.9997

-0.0003 NaN NaN61.0008 0.0008 NaN NaN70.9986

-0.0014 NaN NaN81.0013 0.0013 NaN NaN90.9993

-0.0007 -Inf -Inf101.0002 0.0002 Inf Infnx分量

Gauss法

Cholesky分解法求解误差求解误差1511.0000 0.0000 NaN NaN20.9999

-0.0001 NaN NaN31.0023 0.0023 NaN NaN40.9601

-0.0399 NaN NaN51.3780 0.3780 NaN NaN6

-1.1655

-2.1655 NaN NaN78.9425 7.9425 NaN NaN8

-18.0672

-19.0672 NaN NaN930.4852 29.4852 NaN NaN10

-25.6924

-26.6924 NaN NaN119.0097 8.0097 NaN NaN1211.3971 10.3971 NaN NaN13

-12.5904

-13.5904 NaN NaN147.5529 6.5529 NaN NaN15

-0.2123

-1.2123 NaN NaN4.2.条件数计算结果

数值分析实验报告姓名：学号：实验1：1.实验项目的性质和任务通过上机实验，对病态问题、线性方程组求解和函数的数值逼近方法有一个初步理解。2．教学内容和要求1)对高阶多多项式编程求下面方程的解并绘图演示方程的解与扰动量 的关系。(实验2.6)2)对，生成对应的Hilbert矩阵，计算矩阵的条件数；通过先确定解获得常向量b的方法,确定方程组最后，用矩阵分解方法求解方程组，并分析计算结果。(第三章，实验题4)3)对函数的Chebyshev点编程进行Lagrange插值，并分析插值结果。(第四章 实验1)项目涉及核心知识点病态方程求解、矩阵分解和方程组求解、Lagrange插值。重点与难点算法设计和matlab编程。1)a．实验方案：先创建一个20*50的零矩阵X，然后利用Matlab中的roots()和poly()函数将50个不同的ess扰动值所产生的50个解向量分别存入X矩阵中。然后再将ess向量分别和X的20个行向量绘图。即可直观的看出充分小的扰动值会产生非常大的偏差。即证明了这个问题的病态性。b．编写程序：X=zeros(20,50);ve=zeros(1,21);ess=linspace(0,0.00001,50);k=1;while kve(2)=ess(k);X(1:20,k)=roots(poly(1:20)+ve);k=k+1;endm=1;while mfigure(m),plot(ess,X(m,:));m=m+1;endC．实验结果分析和拓展由上面的实验结果可以看出一个充分小的扰动值可以让方程的解产生非常大的偏差，而且这个偏差随着ess的变大偏差也随即变大。但可以看出在相对小的根处根比较稳定，也就是说这些根关于ess并不敏感，而在较大根处时，根很不稳定，即这些解关于ess的变化是敏感的。这就说明了这个问题本身就是一个病态问题，与算法好坏无关。若扰动在x^18处，只要把程序中的ve(2)改为ve(3)即可，其图形和此类似。d．实验结论：高次多项式扰动求方程解问题是一个病态问题。2)a．实验方案：先创建一个20*20的零矩阵A，再通过给定解x和Hilbert矩阵求出列向量b，然后通过LU分解法求出方程HX=b的解X，然后将x-X’这一行向量存入A矩阵中，形成一循环，最后，如果Hilbert矩阵非病态的话，则可输出一个20*20的对角矩阵。b.编写程序：n=2;A=zeros(20,20);while nx=1:n;H=hilb(n);b=H*x';[L U]=lu(H);y=L\\b;X=U\\y;A(n,1:n)=x-X';n=n+1;endWarning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=4.455948e-017.Warning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=7.948463e-017.Warning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=1.798429e-016.Warning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=7.626119e-018.Warning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=6.040620e-017.Warning:Matrix is close to singular or badly scaled.Results may be inaccurate.RCOND=5.444860e-017.AA=1.0e+003*Columns 1 through 100 0 0 0 0 0 0 0 0 00.0000 0.0000 0 0 0 0 0 0 0 00.0000 0.0000-0.0000 0 0 0 0 0 0 00.0000 0.0000-0.0000 0.0000 0 0 0 0 0 00.0000-0.0000 0.0000-0.0000 0.0000 0 0 0 0 00.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0 0 0 00.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.0000 0 0 00.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.0000 0 00.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 00.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.00000.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.0000-0.0000 0.00000.0000 0.0000-0.0000 0.0000-0.0000 0.0001-0.0003 0.0006-0.0007 0.00050.0000-0.0000 0.0000-0.0001 0.0005-0.0027 0.0096-0.0223 0.0348-0.03610.0000-0.0000 0.0000-0.0004 0.0030-0.0098 0.0080 0.0593-0.2570 0.51540.0000-0.0000 0.0000-0.0001 0.0005-0.0029 0.0095-0.0171 0.0086 0.03470.0000-0.0000 0.0000-0.0000 0.0003-0.0016 0.0059-0.0133 0.0145 0.00940.0000-0.0000 0.0000-0.0001 0.0009-0.0042 0.0118-0.0182 0.0082 0.01850.0000 0.0000-0.0000 0.0002-0.0027 0.0187-0.0762 0.1806-0.2249 0.08130.0000 0.0000-0.0000 0.0001-0.0017 0.0120-0.0497 0.1224-0.1699 0.10640.0000-0.0000 0.0000-0.0003 0.0028-0.0137 0.0371-0.0464-0.0164 0.1243Columns 11 through 200 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 00.0000 0 0 0 0 0 0 0 0 00.0002 0.0000 0 0 0 0 0 0 0 00.0238-0.0091 0.0015 0 0 0 0 0 0 00.6091 0.4336-0.1727 0.0296 0 0 0 0 0 00.0944 0.1170-0.0824 0.0318-0.0053 0 0 0 0 00.0624 0.1107-0.1110 0.0674-0.0232 0.0035 0 0 0 00.0289 0.0059 0.0103 0.0082-0.0263 0.0181-0.0042 0 0 00.0524 0.1690-0.3743-0.1862 1.0944-1.2171 0.6004-0.1156 0 00.0327 0.1652-0.3051-0.0485 0.7195-0.9387 0.5714-0.1699 0.0191 00.1120-0.0421 0.0883 0.0222-0.0628 0.1013-0.2902 0.3783-0.2173 0.0469C．实验结果分析和拓展：当Hilbert矩阵的阶数比较小时，其解X和给定解x偏差不大；但当Hilbert矩阵的阶数变大时，偏差就会变大。这就说明了Hilbert矩阵是一组病态矩阵，从Matlab运行中的Warning可以看出，其条件数相当大。d．实验结论：Hilbert矩阵是一组病态矩阵，用它来做线性方程的系数矩阵时，往往会得出与精确解相差较大的解。3)a．实验方案：在区间【-1，1】上取点，先按Chebyshev取点，即xk=cos((2k-1)pi/2/(n+1))取点，然后再进行拉格朗日插值，绘出图和插值点。而后再进行均匀取点再拉格朗日插值。将两种插值结果进行比较。b．编程实现：for a=1:10b=a+1;for c=1:bX(c)=cos((2*c-1)*pi/2/(a+1));Y(c)=1/(1+25*X(c)^2);x=-1:0.05:1;endm=length(x);for i=1:mz=x(i);s=0;for k=1:bL=1;for j=1:bif j~=kL=L*(z-X(j))/(X(k)-X(j));endends=s+L*Y(k);endy(i)=s;endfigure(1)plot(x,y,'r');hold on;figure(2)plot(X,Y,'b*')hold onendfor a=2:2:10b=a+1;X=linspace(-1,1,b);Y=1./(1+25*X.^2);x=-1:0.05:1;m=length(x);for i=1:mz=x(i);s=0;for k=1:bL=1;for j=1:bif j~=kL=L*(z-X(j))/(X(k)-X(j));endends=s+L*Y(k);endy(i)=s;endfigure(1)plot(x,y,'r');hold on;figure(2)plot(X,Y,'b*')hold onendC．实验结果分析及拓展：均匀插值时，当n比较大时，就会出现多项式插值的Runge现象，即当插值节点的个数n增加时，Lagrange插值多项式对原来函数的近似并非越来越好。当进行非等距节点插值时，其近似效果明显要比均匀插值是要好。原因是非均匀插值时，在远离原点处的插值节点比较密集，所以其插值近似效果要比均匀插值时的效果要好。d．实验结论：利用Chebyshev点进行非等距节点插值的对原函数的近似效果要比均匀节点插值的好内容来自www.book6789.com请勿采集。

声明：本网内容收集自互联网，旨在传播知识仅供参考，不代表本网赞同其观点，文字及图片版权归原网站所有。