文档库 最新最全的文档下载
当前位置:文档库 › 病态希尔伯特方程组解法-高等数值分析大作业

病态希尔伯特方程组解法-高等数值分析大作业

病态希尔伯特方程组解法-高等数值分析大作业
病态希尔伯特方程组解法-高等数值分析大作业

病态方程组的求解

理论分析表明,数值求解病态线性方程组很困难,考虑求解如下的线性方程组的求解:

Hx=b ,其中H 是Hilbert 矩阵,H=(hij)n×n ,h ij =1

i+j?1,ii ,jj =1,2…..nn

1.估计矩阵的2-条件数和矩阵阶数的关系;

2.对不同的n ,取x =(1,1,1…,1)∈R n ,利用Hx=b 得出b.然后分别用Guass 消去,Jacobi 迭代,Gauss-seidel 迭代,SOR 迭代求解方程: Ax=b ;

把近似解x (k )与准确解x 对比,计算结果;

3.综合计算结果,讨论病态方程组的求解。

解:1.用MATLAB 计算1-9阶Hilbert 矩阵的2-条件数得到的结果如下表所示

(cond2(H

nn )表示n 阶Hilbert 矩阵的2-条件数)

对表中数据进行分析可发现,随矩阵的阶数的增加,2-条件数越来越大,后项与前项之比约为30,2-条件数呈近似指数式增长。对cond2(H nn )取以10为底的对数,可发现在前几项log10(cond2(Hn))随n 的增加约呈线性增加关系。 n 为20时,log10(cond2(Hn))与n 的关系与下图所示:

此时可发现当n ≤13时,log10(cond2(Hn))与n 接近于线性关系,n >13后,log10(cond2(Hn))变化十分缓慢。为了对图像的线性进行比较,对n 为1-13部分

1-9阶hilbert 矩阵的2-条件数

n 1

2 3 4 5 6 7 8 9 cond2(H nn ) 1

19.2815 524.0568 1.55E+04 4.77E+05 1.49E+07 4.75E+08 1.53E+10 4.93E+11

进行线性拟合,结果如下图所示:

从上图可看出,n较小时,log10(cond2(Hn))与n成近似线性关系。

继续增大n,n=200时,log10(cond2(Hn))与n如下图所示:

由图可看出,n为200时,log10(cond2(Hn))与n的关系与n=20时相近,在

n较小时log10(cond2(Hn))与n的关系接近线性关系,当n大于一定的值后log10(cond2(Hn))变化十分缓慢,但有上升的趋势。

因此hilbert矩阵的2-条件数与阶数的关系是:hilbert矩阵的2-条件数本身随矩阵的阶数的增加会以近似于指数形式增长,随阶数的增加,其2-条件数会变得很大,但对其做一定处理,如取hilbert矩阵的2-条件数以10为底的对数后,可发现,2-条件数以10为底的对数整体很小,在阶数较小时其值与阶数的关系接近线性关系,当阶数大于一定的值后,其值还是增加得十分缓慢,变化较小。

2.使用Guass消去,Jacobi迭代,Gauss-seidel迭代,SOR迭代求解近似解:

定迭代初值为(0.99,0.99,0.99,0.99)T,迭代次数为20000,判断迭代法是否收敛的标准都是||xx(kk)?xx(kk?1)||2

||xx(kk?1)||2+eeeeee≤1ee?06,即前后两项二范数的相对误差小于

1ee?06,对Guass消去法,当||x?xx?||2||xx?||2+eeeeee≤1ee?06(x为近似解,xx?为准确解)时,认为近似解还是准确的。阶数不同(以n=5,8,9,13为例)时,四种方法得到的近似解如下:

n=5时四种方法得到的近似解

x分量x_Guass x_Jacobi x_gauseid

x_sor

x精确ω=0.67 ω=1.5

x(1) 1.000000 NaN 1.000040 1.000050 1.000016 1 x(2) 1.000000 NaN 0.999675 0.999536 0.999919 1 x(3) 1.000000 NaN 1.000396 1.000652 0.999854 1 x(4) 1.000000 NaN 1.000369 1.000440 1.000694 1 x(5) 1.000000 NaN 0.999491 0.999271 0.999501 1 是否准确/

迭代次数

准确发散496 133 686

n=8时四种方法得到的近似解

x分量x_Guass x_Jacobi x_gauseid

x_sor

x精确ω=1.85 ω=1.5

x(1) 1.000000 NaN 0.999988 0.999937 0.999955 1 x(2) 1.000000 NaN 1.000321 1.000889 1.000729 1 x(3) 1.000000 NaN 0.998996 0.997092 0.997763 1 x(4) 1.000000 NaN 1.000248 1.003317 1.001549 1 x(5) 1.000000 NaN 1.000717 0.998044 1.000227 1 x(6) 1.000000 NaN 1.000512 1.002650 1.000796 1 x(7) 1.000000 NaN 0.999946 0.997242 0.999699 1 x(8) 1.000000 NaN 0.999225 1.000815 0.999233 1 是否准确/

迭代次数

准确发散315 134 284

n=9时四种方法得到的近似解

x分量x_Guass x_Jacobi x_gauseid

x_sor

x精确ω=1.31 ω=1.5

x(1) 1.000000 NaN 0.999964 0.999934 0.999951 1 x(2) 1.000000 NaN 1.000553 1.000887 1.000552 1 x(3) 1.000000 NaN 0.998692 0.997798 0.998791 1 x(4) 1.000002 NaN 1.000007 1.000707 1.000310 1 x(5) 0.999992 NaN 1.000717 1.000677 1.000202 1 x(6) 1.000017 NaN 1.000751 1.000817 1.000602 1 x(7) 0.999980 NaN 1.000371 1.000322 1.000235 1 x(8) 1.000012 NaN 0.999781 0.999739 0.999938 1 x(9) 0.999997 NaN 0.999103 0.999052 0.999389

是否准确/

迭代次数

不准确发散272 257 822

n=13时四种方法得到的近似解

x分量x_Guass x_Jacobi x_gauseid

x_sor

x精确ω=0.19 ω=1.5

x(1) 1.000000 NaN 0.999946 0.999968 0.999986 1 x(2) 0.999989 NaN 1.000365 1.000288 0.999946 1 x(3) 1.000414 NaN 0.999869 0.999642 1.000863 1 x(4) 0.993147 NaN 0.999375 0.999618 0.998665 1 x(5) 1.061246 NaN 0.999667 0.999889 0.999967 1 x(6) 0.669261 NaN 1.000104 1.000198 0.999814 1 x(7) 2.149074 NaN 1.000412 1.000416 1.000440 1 x(8) -1.654171 NaN 1.000530 1.000497 1.000406 1 x(9) 5.118548 NaN 1.000475 1.000435 1.000474 1 x(10) -3.243050 NaN 1.000285 1.000245 1.000271 1 x(11) 3.783005 NaN 0.999999 0.999951 1.000046 1 x(12) -0.051793 NaN 0.999651 0.999578 0.999717 1 x(13) 1.174329 NaN 0.999264 0.999148 0.999366 1 是否准确/

迭代次数

不准确发散585 112 738 n=13时,二范数相对精度为1e-8时GS迭代法和SOR迭代得到的近似解

x分量x_gauseid

x_sor

x精确ω=0.43 ω=1.5

x(1) 1.000002 1.000003 1.000000 1

x(2) 0.999983 0.999954 1.000014 1

x(3) 1.000002 1.000118 0.999857 1

x(4) 1.000126 1.000023 1.000369 1

x(5) 0.999939 0.999892 0.999777 1

x(6) 0.999868 0.999869 0.999931 1

x(7) 0.999922 0.999943 0.999876 1

x(8) 1.000019 1.000051 1.000038 1

x(9) 1.000097 1.000133 1.000082 1

x(10) 1.000124 1.000155 1.000127 1

x(11) 1.000088 1.000101 1.000083 1

x(12) 0.999990 0.999973 0.999996 1

x(13) 0.999839 0.999779 0.999848 1

迭代次数11890 4008 12961

实际计算时对n为2-20的近似解均进行了计算,以上列出了阶数为5,8,9,13时四种方法得到的近似解,以及当相对精度提高时GS迭代和SOR迭代的结果。对结果进行分析可以发现:

(1)对Guass消去法,在n从1变到8时,计算得到的近似解还是准确的,n大于8时,近似解与准确解开始有偏差,且n越大,偏差越大;

(2)对Jacobi迭代,无论阶数n为多少,Jacobi迭代都不收敛,使用Jacobi迭代法无法得到近似解;

(3)使用Gauss-seidel迭代法时,对所有的n都收敛,不同阶数收敛速度不一样,但结果与精确解一直存在误差,在精确解附近波动,当对解的精度要求提高时(如要求精度为1e-08)GS迭代的收敛速度会大大变慢;

(4)使用SOR迭代时,对所有的n也都收敛,ω不同时,收敛速度不同,对不同的阶数,总能找到一个ω,使SOR迭代法的收敛速度比Gauss-seidel迭代法快,当对解的精度的要求提高时,选用合适的迭代因子的效果更明显。

从上面的分析可以看出当病态问题的阶数升高时作为直接法的Gauss消去法由能求解变得不能求解。用Jacobi迭代方法迭代发散,无法求解;而GS和SOR 迭代法在阶数升高时仍能求解,选择合适的ω,可使SOR迭代法收敛速度加快,但在阶数较低时直接法能求得精确解而迭代法却总存在一定的误差,可见直接法与迭代法各有优势与不足。

3.综合计算结果,讨论病态方程组的求解。

由2中的结果可知对系数矩阵为Hilbert矩阵的的病态方程组,在阶数较低时,能使用高斯消去法直接解,阶数高时使用高斯消去会带来很大误差,不能继续使用,这时可换用迭代法求解;迭代法中Jacobi迭代法由于不收敛,不能用于病态方程组的求解;GS法和SOR法都能用来求解,对迭代速度要求高的话,推荐使用最优迭代因子的SOR迭代法。

迭代法之所以在阶数高时还能继续用来求解病态方程组,是因为迭代法本质上相当于对系数矩阵做了预处理,降低了系数矩阵的条件数。

数值分析实验报告1

实验一误差分析 实验1.1(病态问题) 实验目的:算法有“优”与“劣”之分,问题也有“好”与“坏”之别。对数值方法的研究而言,所谓坏问题就是问题本身对扰动敏感者,反之属于好问题。通过本实验可获得一个初步体会。 数值分析的大部分研究课题中,如线性代数方程组、矩阵特征值问题、非线性方程及方程组等都存在病态的问题。病态问题要通过研究和构造特殊的算法来解决,当然一般要付出一些代价(如耗用更多的机器时间、占用更多的存储空间等)。 问题提出:考虑一个高次的代数多项式 显然该多项式的全部根为1,2,…,20共计20个,且每个根都是单重的。现考虑该多项式的一个扰动 其中ε(1.1)和(1.221,,,a a 的输出b ”和“poly ε。 (1(2 (3)写成展 关于α solve 来提高解的精确度,这需要用到将多项式转换为符号多项式的函数poly2sym,函数的具体使用方法可参考Matlab 的帮助。 实验过程: 程序: a=poly(1:20); rr=roots(a); forn=2:21 n form=1:9 ess=10^(-6-m);

ve=zeros(1,21); ve(n)=ess; r=roots(a+ve); -6-m s=max(abs(r-rr)) end end 利用符号函数:(思考题一)a=poly(1:20); y=poly2sym(a); rr=solve(y) n

很容易的得出对一个多次的代数多项式的其中某一项进行很小的扰动,对其多项式的根会有一定的扰动的,所以对于这类病态问题可以借助于MATLAB来进行问题的分析。 学号:06450210 姓名:万轩 实验二插值法

数值分析大作业-三、四、五、六、七

大作业 三 1. 给定初值 0x 及容许误差 ,编制牛顿法解方程f (x )=0的通用 程序. 解:Matlab 程序如下: 函数m 文件:fu.m function Fu=fu(x) Fu=x^3/3-x; end 函数m 文件:dfu.m function Fu=dfu(x) Fu=x^2-1; end 用Newton 法求根的通用程序Newton.m clear; x0=input('请输入初值x0:'); ep=input('请输入容许误差:'); flag=1; while flag==1 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)

while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)=ep flag=0; end end fprintf('最大的sigma 值为:%f\n',sigma); 2.求下列方程的非零根 5130.6651()ln 05130.665114000.0918 x x f x x +??=-= ?-???解:Matlab 程序为: (1)主程序 clear clc format long x0=765; N=100; errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1; while nerrorlim n=n+1; else break ; end x0=x; end disp(['迭代次数: n=',num2str(n)]) disp(['所求非零根: 正根x1=',num2str(x),' 负根x2=',num2str(-x)]) (2)子函数 非线性函数f function y=f(x) y=log((513+0.6651*x)/(513-0.6651*x))-x/(1400*0.0918); end

数值分析实验报告2

实验报告 实验项目名称函数逼近与快速傅里叶变换实验室数学实验室 所属课程名称数值逼近 实验类型算法设计 实验日期 班级 学号 姓名 成绩

512*x^10 - 1280*x^8 + 1120*x^6 - 400*x^4 + 50*x^2 - 1 并得到Figure,图像如下: 实验二:编写程序实现[-1,1]上n阶勒让德多项式,并作画(n=0,1,…,10 在一个figure中)。要求:输入Legendre(-1,1,n),输出如a n x n+a n-1x n-1+…多项式。 在MATLAB的Editor中建立一个M-文件,输入程序代码,实现勒让德多项式的程序代码如下: function Pn=Legendre(n,x) syms x; if n==0 Pn=1; else if n==1 Pn=x; else Pn=expand((2*n-1)*x*Legendre(n-1)-(n-1)*Legendre(n-2))/(n); end x=[-1:0.1:1]; A=sym2poly(Pn); yn=polyval(A,x); plot (x,yn,'-o'); hold on

end 在command Windows中输入命令:Legendre(10),得出的结果为: Legendre(10) ans = (46189*x^10)/256 - (109395*x^8)/256 + (45045*x^6)/128 - (15015*x^4)/128 + (3465*x^2)/256 - 63/256 并得到Figure,图像如下: 实验三:利用切比雪夫零点做拉格朗日插值,并与以前拉格朗日插值结果比较。 在MATLAB的Editor中建立一个M-文件,输入程序代码,实现拉格朗日插值多项式的程序代码如下: function [C,D]=lagr1(X,Y) n=length(X); D=zeros(n,n); D(:,1)=Y'; for j=2:n for k=j:n D(k,j)=(D(k,j-1)- D(k-1,j-1))/(X(k)-X(k-j+1)); end end C=D(n,n); for k=(n-1):-1:1

数值分析作业答案

数值分析作业答案 插值法 1、当x=1,-1,2时,f(x)=0,-3,4,求f(x)的二次插值多项式。 (1)用单项式基底。 (2)用Lagrange插值基底。 (3)用Newton基底。 证明三种方法得到的多项式是相同的。 解:(1)用单项式基底 设多项式为: , 所以: 所以f(x)的二次插值多项式为: (2)用Lagrange插值基底 Lagrange插值多项式为: 所以f(x)的二次插值多项式为: (3) 用Newton基底: 均差表如下: xk f(xk) 一阶均差二阶均差 1 0 -1 -3 3/2 2 4 7/ 3 5/6 Newton插值多项式为: 所以f(x)的二次插值多项式为: 由以上计算可知,三种方法得到的多项式是相同的。 6、在上给出的等距节点函数表,若用二次插值求ex的近似值,要使截断误差不超过10-6,问使用函数表的步长h应取多少? 解:以xi-1,xi,xi+1为插值节点多项式的截断误差,则有 式中 令得 插值点个数

是奇数,故实际可采用的函数值表步长 8、,求及。 解:由均差的性质可知,均差与导数有如下关系: 所以有: 15、证明两点三次Hermite插值余项是 并由此求出分段三次Hermite插值的误差限。 证明:利用[xk,xk+1]上两点三次Hermite插值条件 知有二重零点xk和k+1。设 确定函数k(x): 当或xk+1时k(x)取任何有限值均可; 当时,,构造关于变量t的函数 显然有 在[xk,x][x,xk+1]上对g(x)使用Rolle定理,存在及使得 在,,上对使用Rolle定理,存在,和使得 再依次对和使用Rolle定理,知至少存在使得 而,将代入,得到 推导过程表明依赖于及x 综合以上过程有: 确定误差限: 记为f(x)在[a,b]上基于等距节点的分段三次Hermite插值函数。在区间[xk,xk+1]上有 而最值 进而得误差估计: 16、求一个次数不高于4次的多项式,使它满足,,。

数值分析实验报告

数值分析实验报告 姓名:周茹 学号: 912113850115 专业:数学与应用数学 指导老师:李建良

线性方程组的数值实验 一、课题名字:求解双对角线性方程组 二、问题描述 考虑一种特殊的对角线元素不为零的双对角线性方程组(以n=7为例) ?????????? ?????? ? ???? ?d a d a d a d a d a d a d 766 55 44 3 32 211??????????????????????x x x x x x x 7654321=?????????? ? ???????????b b b b b b b 7654321 写出一般的n (奇数)阶方程组程序(不要用消元法,因为不用它可以十分方便的解出这个方程组) 。 三、摘要 本文提出解三对角矩阵的一种十分简便的方法——追赶法,该算法适用于任意三对角方程组的求解。 四、引言 对于一般给定的d Ax =,我们可以用高斯消去法求解。但是高斯消去法过程复杂繁琐。对于特殊的三对角矩阵,如果A 是不可约的弱对角占优矩阵,可以将A 分解为UL ,再运用追赶法求解。

五、计算公式(数学模型) 对于形如????? ?? ????? ??? ?---b a c b a c b a c b n n n n n 111 2 2 2 11... ... ...的三对角矩阵UL A =,容易验证U 、L 具有如下形式: ??????? ????? ??? ?=u a u a u a u n n U ...... 3 3 22 1 , ?? ????? ? ?? ??????=1 (1) 1132 1l l l L 比较UL A =两边元素,可以得到 ? ?? ??-== = l a b u u c l b u i i i i i i 111 i=2, 3, ... ,n 考虑三对角线系数矩阵的线性方程组 f Ax = 这里()T n x x x x ... 2 1 = ,()T n f f f f ... 2 1 = 令y Lx =,则有 f Uy = 于是有 ()?????-== --u y a f y u f y i i i i i 1 1 11 1 * i=2, 3, ... ,n 再根据y Lx =可得到

数值分析实验报告2

实验名称 插值法 实验目的 (1)学习并熟练掌握MA TLAB 语言的编程; (2)通过课程实习能够应用MATLAB 软件来计算函数的插值,了解函数插值方法。 实验原理 牛顿差商形式多项式 P(x)=f(x0)+f[x0,x1](x-x0)+f[x0,x1,x2](x-x0)(x-x1)+…+f[x0,x1,x2…xn](x-x0)…(x-xn-1) 牛顿插值多项式的余项 Rn(x)=f[x0,x1,x2…xn]wn+1(x) 实验题目 {1}已知函数在下列各点的值为 i x 0.2 0.4 0.6 0.8 1.0 ()i f x 0.98 0.92 0.81 0.64 0.38 试用4次牛顿插值多项式()4P x 及三次样条函数()Q x (自然边界条件)对数据进行插 值。用图给出{(,i i x y ),i x =0.2+0.08i ,i=0,1,11,10},()4P x 及()Q x 。 ①实验过程 x1=[0.2 0.4 0.6 0.8 1.0]; y1=[0.98 0.92 0.81 0.64 0.38]; n=length(y1); c=y1(:); for j=2:n %求差商 for i=n:-1:j c(i)=(c(i)-c(i-1))/(x1(i)-x1(i-j+1)); end end syms x df d ; df(1)=1;d(1)=y1(1); for i=2:n %求牛顿差值多项式 df(i)=df(i-1)*(x-x1(i-1)); d(i)=c(i-1)*df(i); end P4=vpa(sum(d),5) %P4即为4次牛顿插值多项式,并保留小数点后5位数 pp=csape(x1,y1, 'variational');%调用三次样条函数 q=pp.coefs; q1=q(1,:)*[(x-.2)^3;(x-.2)^2;(x-.2);1]; q1=vpa(collect(q1),5) q2=q(1,:)*[(x-.4)^3;(x-.4)^2;(x-.4);1]; q2=vpa(collect(q2),5) q3=q(1,:)*[(x-.6)^3;(x-.6)^2;(x-.6);1];

北航数值分析大作业一

《数值分析B》大作业一 SY1103120 朱舜杰 一.算法设计方案: 1.矩阵A的存储与检索 将带状线性矩阵A[501][501]转存为一个矩阵MatrixC[5][501] . 由于C语言中数组角标都是从0开始的,所以在数组MatrixC[5][501]中检索A的带内元素a ij的方法是: A的带内元素a ij=C中的元素c i-j+2,j 2.求解λ1,λ501,λs ①首先分别使用幂法和反幂法迭代求出矩阵按摸最大和最小的特征值λmax和λmin。λmin即为λs; 如果λmax>0,则λ501=λmax;如果λmax<0,则λ1=λmax。 ②使用带原点平移的幂法(mifa()函数),令平移量p=λmax,求 出对应的按摸最大的特征值λ,max, 如果λmax>0,则λ1=λ,max+p;如果λmax<0,则λ501=λ,max+p。 3.求解A的与数μk=λ1+k(λ501-λ1)/40的最接近的特征值λik (k=1,2,…,39)。 使用带原点平移的反幂法,令平移量p=μk,即可求出与μk最接近的特征值λik。 4.求解A的(谱范数)条件数cond(A)2和行列式d etA。 ①cond(A)2=|λ1/λn|,其中λ1和λn分别是矩阵A的模最大和 最小特征值。

②矩阵A的行列式可先对矩阵A进行LU分解后,detA等于U所有对角线上元素的乘积。 二.源程序 #include #include #include #include #include #include #include #define E 1.0e-12 /*定义全局变量相对误差限*/ int max2(int a,int b) /*求两个整型数最大值的子程序*/ { if(a>b) return a; else return b; } int min2(int a,int b) /*求两个整型数最小值的子程序*/ { if(a>b) return b; else return a; } int max3(int a,int b,int c) /*求三整型数最大值的子程序*/ { int t; if(a>b) t=a; else t=b; if(t

数值分析第一次作业及参考答案

数值计算方法第一次作业及参考答案 1. 已测得函数()y f x =的三对数据:(0,1),(-1,5),(2,-1), (1)用Lagrange 插值求二次插值多项式。(2)构造差商表。(3)用Newton 插值求二次插值多项式。 解:(1)Lagrange 插值基函数为 0(1)(2)1 ()(1)(2)(01)(02)2 x x l x x x +-= =-+-+- 同理 1211 ()(2),()(1)36 l x x x l x x x = -=+ 故 2 20 2151 ()()(1)(2)(2)(1) 23631 i i i p x y l x x x x x x x x x =-==-+-+-++=-+∑ (2)令0120,1,2x x x ==-=,则一阶差商、二阶差商为 011215 5(1) [,]4, [,]20(1) 12 f x x f x x ---= =-= =----- 0124(2) [,,]102 f x x x ---= =- 实际演算中可列一张差商表: (3)用对角线上的数据写出插值多项式 2 2()1(4)(0)1*(0)(1)31P x x x x x x =+--+-+=-+ 2. 在44x -≤≤上给出()x f x e =的等距节点函数表,若用二次插值求x e 的近似值,要使 截断误差不超过6 10-,问使用函数表的步长h 应取多少 解: ()40000(), (),[4,4],,,, 1.x k x f x e f x e e x x h x x h x x th t ==≤∈--+=+≤考察点及

(3) 2000 4 43 4 3 () ()[(()]()[()] 3! (1)(1) (1)(1) 3!3! .(4,4). 6 f R x x x h x x x x h t t t e t h th t h e h e ξ ξ =----+ -+ ≤+??-= ≤∈- 则 4 36 ((1)(1) 100.006. t t t h - -+± << Q在点 得 3.求2 () f x x =在[a,b]上的分段线性插值函数() h I x,并估计误差。 解: 22 22 11 1 111 22 11 11 1 () () k k k k h k k k k k k k k k k k k k k k k k k x x x x x x I x x x x x x x x x x x x x x x x x x x x x ++ + +++ ++ ++ + --- =+= --- ?-? -=+- - [] 2 11 22 11 ()()()[()] 11 ()() 44 h h k k k k k k k k R x f x I x x x x x x x x x x x x x h ++ ++ =-=-+- =--≤-= 4.已知单调连续函数() y f x =的如下数据 用插值法计算x约为多少时() 1. f x=(小数点后至少保留4位) 解:作辅助函数()()1, g x f x =-则问题转化为x为多少时,()0. g x=此时可作新 的关于() i g x的函数表。由() f x单调连续知() g x也单调连续,因此可对() g x的数值进行反插。的牛顿型插值多项式为 1()0.110.097345( 2.23)0.451565( 2.23)( 1.10) 0.255894( 2.23)( 1.10)(0.17) x g y y y y y y y - ==-+++++ -++-

数值分析实验报告模板

数值分析实验报告模板 篇一:数值分析实验报告(一)(完整) 数值分析实验报告 1 2 3 4 5 篇二:数值分析实验报告 实验报告一 题目:非线性方程求解 摘要:非线性方程的解析解通常很难给出,因此线性方程的数值解法就尤为重要。本实验采用两种常见的求解方法二分法和Newton法及改进的Newton法。利用二分法求解给定非线性方程的根,在给定的范围内,假设f(x,y)在[a,b]上连续,f(a)xf(b) 直接影响迭代的次数甚至迭代的收敛与发散。即若x0 偏离所求根较远,Newton法可能发散的结论。并且本实验中还利用利用改进的Newton法求解同样的方程,且将结果与Newton法的结果比较分析。 前言:(目的和意义) 掌握二分法与Newton法的基本原理和应用。掌握二分法的原理,验证二分法,在选对有根区间的前提下,必是收

敛,但精度不够。熟悉Matlab语言编程,学习编程要点。体会Newton使用时的优点,和局部收敛性,而在初值选取不当时,会发散。 数学原理: 对于一个非线性方程的数值解法很多。在此介绍两种最常见的方法:二分法和Newton法。 对于二分法,其数学实质就是说对于给定的待求解的方程f(x),其在[a,b]上连续,f(a)f(b) Newton法通常预先要给出一个猜测初值x0,然后根据其迭代公式xk?1?xk?f(xk) f'(xk) 产生逼近解x*的迭代数列{xk},这就是Newton法的思想。当x0接近x*时收敛很快,但是当x0选择不好时,可能会发散,因此初值的选取很重要。另外,若将该迭代公式改进为 xk?1?xk?rf(xk) 'f(xk) 其中r为要求的方程的根的重数,这就是改进的Newton 法,当求解已知重数的方程的根时,在同种条件下其收敛速度要比Newton法快的多。 程序设计: 本实验采用Matlab的M文件编写。其中待求解的方程写成function的方式,如下 function y=f(x);

数值分析实验报告

实验一、误差分析 一、实验目的 1.通过上机编程,复习巩固以前所学程序设计语言及上机操作指令; 2.通过上机计算,了解误差、绝对误差、误差界、相对误差界的有关概念; 3.通过上机计算,了解舍入误差所引起的数值不稳定性。 二.实验原理 误差问题是数值分析的基础,又是数值分析中一个困难的课题。在实际计算中,如果选用了不同的算法,由于舍入误差的影响,将会得到截然不同的结果。因此,选取算法时注重分析舍入误差的影响,在实际计算中是十分重要的。同时,由于在数值求解过程中用有限的过程代替无限的过程会产生截断误差,因此算法的好坏会影响到数值结果的精度。 三.实验内容 对20,,2,1,0 =n ,计算定积分 ?+=10 5dx x x y n n . 算法1:利用递推公式 151--=n n y n y , 20,,2,1 =n , 取 ?≈-=+=1 00182322.05ln 6ln 51dx x y . 算法2:利用递推公式 n n y n y 51511-= - 1,,19,20 =n . 注意到 ???=≤+≤=10 10202010201051515611261dx x dx x x dx x , 取 008730.0)12611051(20120≈+≈y .: 四.实验程序及运行结果 程序一: t=log(6)-log(5);

n=1; y(1)=t; for k=2:1:20 y(k)=1/k-5*y(k-1); n=n+1; end y y =0.0884 y =0.0581 y =0.0431 y =0.0346 y =0.0271 y =0.0313 y =-0.0134 y =0.1920 y =-0.8487 y =4.3436 y =-21.6268 y =108.2176 y =-541.0110 y =2.7051e+003 y =-1.3526e+004 y =6.7628e+004 y =-3.3814e+005 y =1.6907e+006 y =-8.4535e+006 y =4.2267e+007 程序2: y=zeros(20,1); n=1; y1=(1/105+1/126)/2;y(20)=y1; for k=20:-1:2 y(k-1)=1/(5*k)-(1/5)*y(k); n=n+1; end 运行结果:y = 0.0884 0.0580 0.0431 0.0343 0.0285 0.0212 0.0188 0.0169

数值分析大作业三 四 五 六 七

大作业 三 1. 给定初值 0x 及容许误差 ,编制牛顿法解方程f (x )=0的通用程序. 解:Matlab 程序如下: 函数m 文件:fu.m function Fu=fu(x) Fu=x^3/3-x; end 函数m 文件:dfu.m function Fu=dfu(x) Fu=x^2-1; end 用Newton 法求根的通用程序Newton.m clear; x0=input('请输入初值x0:'); ep=input('请输入容许误差:');

flag=1; while flag==1 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)

while flag==1 sigma=k*eps; x0=sigma; k=k+1; m=0; flag1=1; while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)=ep flag=0;

end end fprintf('最大的sigma 值为:%f\n',sigma); 2.求下列方程的非零根 5130.6651()ln 05130.665114000.0918 x x f x x +?? =-= ?-???解: Matlab 程序为: (1)主程序 clear clc format long x0=765; N=100; errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1;

数值分析作业答案part

6.4.设??? ? ? ??=5010010a b b a A ,0det ≠A ,用a ,b 表示解线性方程组f Ax =的雅可比迭代与 高斯—塞德尔迭代收敛的充分必要条件。 解 雅可比迭代法的迭代矩阵 ? ??? ??? ? ??----=???? ? ??----????? ??=-050100100100000001010101 a b b a a b b a B J , ?? ? ?? -=-1003||2ab B I J λλλ,10||3)(ab B J = ρ。 雅可比迭代法收敛的充分必要条件是3 100 ||

(完整版)哈工大-数值分析上机实验报告

实验报告一 题目:非线性方程求解 摘要:非线性方程的解析解通常很难给出,因此线性方程的数值解法就尤为重要。本实验采用两种常见的求解方法二分法和Newton法及改进的Newton法。 前言:(目的和意义) 掌握二分法与Newton法的基本原理和应用。 数学原理: 对于一个非线性方程的数值解法很多。在此介绍两种最常见的方法:二分法和Newton法。 对于二分法,其数学实质就是说对于给定的待求解的方程f(x),其在[a,b]上连续,f(a)f(b)<0,且f(x)在[a,b]内仅有一个实根x*,取区间中点c,若,则c恰为其根,否则根据f(a)f(c)<0是否成立判断根在区间[a,c]和[c,b]中的哪一个,从而得出新区间,仍称为[a,b]。重复运行计算,直至满足精度为止。这就是二分法的计算思想。

Newton法通常预先要给出一个猜测初值x0,然后根据其迭代公式 产生逼近解x*的迭代数列{x k},这就是Newton法的思想。当x0接近x*时收敛很快,但是当x0选择不好时,可能会发散,因此初值的选取很重要。另外,若将该迭代公式改进为 其中r为要求的方程的根的重数,这就是改进的Newton法,当求解已知重数的方程的根时,在同种条件下其收敛速度要比Newton法快的多。 程序设计: 本实验采用Matlab的M文件编写。其中待求解的方程写成function的方式,如下 function y=f(x); y=-x*x-sin(x); 写成如上形式即可,下面给出主程序。 二分法源程序: clear %%%给定求解区间 b=1.5; a=0;

%%%误差 R=1; k=0;%迭代次数初值 while (R>5e-6) ; c=(a+b)/2; if f12(a)*f12(c)>0; a=c; else b=c; end R=b-a;%求出误差 k=k+1; end x=c%给出解 Newton法及改进的Newton法源程序:clear %%%% 输入函数 f=input('请输入需要求解函数>>','s') %%%求解f(x)的导数 df=diff(f);

数值分析实验报告总结

数值分析实验报告总结 随着电子计算机的普及与发展,科学计算已成为现代科 学的重要组成部分,因而数值计算方法的内容也愈来愈广泛和丰富。通过本学期的学习,主要掌握了一些数值方法的基本原理、具体算法,并通过编程在计算机上来实现这些算法。 算法算法是指由基本算术运算及运算顺序的规定构成的完 整的解题步骤。算法可以使用框图、算法语言、数学语言、自然语言来进行描述。具有的特征:正确性、有穷性、适用范围广、运算工作量少、使用资源少、逻辑结构简单、便于实现、计算结果可靠。 误差 计算机的计算结果通常是近似的,因此算法必有误差, 并且应能估计误差。误差是指近似值与真正值之差。绝对误差是指近似值与真正值之差或差的绝对值;相对误差:是指近似值与真正值之比或比的绝对值。误差来源见表 第三章泛函分析泛函分析概要 泛函分析是研究“函数的函数”、函数空间和它们之间 变换的一门较新的数学分支,隶属分析数学。它以各种学科

如果 a 是相容范数,且任何满足 为具体背景,在集合的基础上,把客观世界中的研究对象抽 范数 范数,是具有“长度”概念的函数。在线性代数、泛函 分析及相关的数学领域,泛函是一个函数,其为矢量空间内 的所有矢量赋予非零的正长度或大小。这里以 Cn 空间为例, Rn 空间类似。最常用的范数就是 P-范数。那么 当P 取1, 2 ,s 的时候分别是以下几种最简单的情形: 其中2-范数就是通常意义下的距离。 对于这些范数有以下不等式: 1 < n1/2 另外,若p 和q 是赫德尔共轭指标,即 1/p+1/q=1 么有赫德尔不等式: II = ||xH*y| 当p=q=2时就是柯西-许瓦兹不等式 般来讲矩阵范数除了正定性,齐次性和三角不等式之 矩阵范数通常也称为相容范数。 象为元素和空间。女口:距离空间,赋范线性空间, 内积空间。 1-范数: 1= x1 + x2 +?+ xn 2-范数: x 2=1/2 8 -范数: 8 =max oo ,那 外,还规定其必须满足相容性: 所以

数值分析大作业

数值分析报大作业 班级:铁道2班 专业:道路与铁道工程 姓名:蔡敦锦 学号:13011260

一、序言 该数值分析大作业是通过C语言程序编程在Microsoft Visual C++ 6.0编程软件上运行实现的。本来是打算用Matlab软间来计算非线性方程的根的。学习Matlab也差不多有一个多月了,感觉自己编程做题应该没什么问题了;但是当自己真心的去编程、运行时才发现有很多错误,花了一天时间修改、调试程序都没能得到自己满意的结果。所以,我选择了自己比较熟悉的C程序语言来编程解决非线性的求值问题,由于本作业是为了比较几种方法求值问题的收敛速度和精度的差异,选择了一个相对常见的非线性函数来反映其差异,程序运行所得结果我个人比较满意。编写C语言,感觉比较上手,程序出现问题也能比较熟练的解决。最终就决定上交一份C程序语言编程的求值程序了!

二、选题 本作业的目的是为了加深对非线性方程求根方法的二分法、简单迭代法、、牛顿迭代法弦截法等的构造过程的理解;能将各种方法的算法描述正确并且能够改编为程序并在计算机上实现程序的正确合理的运行,能得到自己满意的结果,并且能调试修改程序中可能出现的问题和程序功能的增减修改。本次程序是为了比较各种方法在求解同一非线性方程根时,在收敛情况上的差异。 为了达到上面的条件我选择自己比较熟悉的语言—C语言来编程,所选题目为计算方程f(x)=x3-2x-5=0在区间[2,3]内其最后两近似值的差的绝对值小于等于5 ?的根的几种方法的比较。 110- 本文将二分法、牛顿法、简单迭代法、弦截法及加速收敛法这五种方法在同一个程序中以函数调用的方式来实现,比较简洁明了,所得结果能很好的比较,便于分析;发现问题和得出结论。

数值分析实验报告

学生实验报告实验课程名称 开课实验室 学院年级专业班 学生姓名学号 开课时间至学年学期

if(A(m,k)~=0) if(m~=k) A([k m],:)=A([m k],:); %换行 end A(k+1:n, k:c)=A(k+1:n, k:c)-(A(k+1:n,k)/ A(k,k))*A(k, k:c); %消去end end x=zeros(length(b),1); %回代求解 x(n)=A(n,c)/A(n,n); for k=n-1:-1:1 x(k)=(A(k,c)-A(k,k+1:n)*x(k+1:n))/A(k,k); end y=x; format short;%设置为默认格式显示,显示5位 (2)建立MATLAB界面 利用MA TLAB的GUI建立如下界面求解线性方程组: 详见程序。 五、计算实例、数据、结果、分析 下面我们对以上的结果进行测试,求解:

? ? ? ? ? ? ? ? ? ? ? ? - = ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? - - - - - - 7 2 5 10 13 9 14 4 4 3 2 1 13 12 4 3 3 10 2 4 3 2 1 x x x x 输入数据后点击和,得到如下结果: 更改以上数据进行测试,求解如下方程组: 1 2 3 4 43211 34321 23431 12341 x x x x ?? ???? ?? ???? ?? ???? = ?? ???? - ?? ???? - ???? ?? 得到如下结果:

数值分析2016上机实验报告

序言 数值分析是计算数学的范畴,有时也称它为计算数学、计算方法、数值方法等,其研究对象是各种数学问题的数值方法的设计、分析及其有关的数学理论和具体实现的一门学科,它是一个数学分支。是科学与工程计算(科学计算)的理论支持。许多科学与工程实际问题(核武器的研制、导弹的发射、气象预报)的解决都离不开科学计算。目前,试验、理论、计算已成为人类进行科学活动的三大方法。 数值分析是计算数学的一个主要部分,计算数学是数学科学的一个分支,它研究用计算机求解各种数学问题的数值计算方法及其理论与软件实现。现在面向数值分析问题的计算机软件有:C,C++,MATLAB,Python,Fortran等。 MATLAB是matrix laboratory的英文缩写,它是由美国Mathwork公司于1967年推出的适合用于不同规格计算机和各种操纵系统的数学软件包,现已发展成为一种功能强大的计算机语言,特别适合用于科学和工程计算。目前,MATLAB应用非常广泛,主要用于算法开发、数据可视化、数值计算和数据分析等,除具备卓越的数值计算能力外,它还提供了专业水平的符号计算,文字处理,可视化建模仿真和实时控制等功能。 本实验报告使用了MATLAB软件。对不动点迭代,函数逼近(lagrange插值,三次样条插值,最小二乘拟合),追赶法求解矩阵的解,4RungeKutta方法求解,欧拉法及改进欧拉法等算法做了简单的计算模拟实践。并比较了各种算法的优劣性,得到了对数值分析这们学科良好的理解,对以后的科研数值分析能力有了极大的提高。

目录 序言 (1) 问题一非线性方程数值解法 (3) 1.1 计算题目 (3) 1.2 迭代法分析 (3) 1.3计算结果分析及结论 (4) 问题二追赶法解三对角矩阵 (5) 2.1 问题 (5) 2.2 问题分析(追赶法) (6) 2.3 计算结果 (7) 问题三函数拟合 (7) 3.1 计算题目 (7) 3.2 题目分析 (7) 3.3 结果比较 (12) 问题四欧拉法解微分方程 (14) 4.1 计算题目 (14) 4.2.1 方程的准确解 (14) 4.2.2 Euler方法求解 (14) 4.2.3改进欧拉方法 (16) 问题五四阶龙格-库塔计算常微分方程初值问题 (17) 5.1 计算题目 (17) 5.2 四阶龙格-库塔方法分析 (18) 5.3 程序流程图 (18) 5.4 标准四阶Runge-Kutta法Matlab实现 (19) 5.5 计算结果及比较 (20) 问题六舍入误差观察 (22) 6.1 计算题目 (22) 6.2 计算结果 (22) 6.3 结论 (23) 7 总结 (24) 附录

数值计算实验报告2

贵州师范大学数学与计算机科学学院学生实验报告 课程名称: 数值分析 班级: 实验日期:2013年9月27日 学 号: 姓名: 指导教师: 实验成绩: 一、实验名称 实验二: Lagrange 插值与曲线拟合的最小二乘法 二、实验目的及要求 1.让学生掌握Lagrange 插值与曲线拟合的最小二乘法 2.让学生能够用这些方法解决一些实际问题 三、实验环境 每人一台计算机,要求安装Windows XP 操作系统,Microsoft office2003、MATLAB6.5(或7.0). 四、实验内容 题1: 对函数2 11 )(x x f += ,取n+1个等距分布的插值节点,取不同的n ,作n 次Lagrange 插值,把)(x f 和插值多项式的图象绘制在同一张图上进行比较. 题2: 给定数据点 分别用一次,二次,和三次多项式曲线,以及最小二乘法拟合这些数据点,哪一种曲线拟合较好?为什么?你能找出更好的拟合曲线吗? 提示:用残差平方的大小来判断拟合的优劣,越小越好. 五、算法描述及实验步骤 针对实验1: (1)运用Matlab 创建M 文件 (2)在命令窗口调用文件 针对实验2: (1)运用Matlab 作出上面表中的数据的散点图 (2)分别作出一次二次三次多项式拟合曲线 (3)比较三种曲线拟合的精度即比较残差平方的大小

六、调试过程及实验结果 https://www.wendangku.net/doc/c31262768.html,grange插值: (1)命令窗口输入: >> f=shuru(-5,5,10) (回车) f = Columns 1 through 6 0.0385 0.0588 0.1000 0.2000 0.5000 1.0000 Columns 7 through 11 0.5000 0.2000 0.1000 0.0588 0.0385 (2)所得结果为图形所示: 2. 曲线拟合的最小二乘法 (1)先作出该散点图的一次拟合曲线(程序如下):>> x=[-3,-1,0,1,3,5]; >> t=-4:0.01:6; >> y=[-6,-3,-1,0,1,3]; >> subplot(1,3,1) >> scatter(x,y,'filled','r'); >> hold on >> p1=polyfit(x,y,1) p1 = 1.0776 -1.8980 >> y1=polyval(p1,x);

北航数值分析大作业第一题幂法与反幂法

《数值分析》计算实习题目 第一题: 1. 算法设计方案 (1)1λ,501λ和s λ的值。 1)首先通过幂法求出按模最大的特征值λt1,然后根据λt1进行原点平移求出另一特征值λt2,比较两值大小,数值小的为所求最小特征值λ1,数值大的为是所求最大特征值λ501。 2)使用反幂法求λs ,其中需要解线性方程组。因为A 为带状线性方程组,此处采用LU 分解法解带状方程组。 (2)与140k λλμλ-5011=+k 最接近的特征值λik 。 通过带有原点平移的反幂法求出与数k μ最接近的特征值 λik 。 (3)2cond(A)和det A 。 1)1=n λλ2cond(A),其中1λ和n λ分别是按模最大和最小特征值。 2)利用步骤(1)中分解矩阵A 得出的LU 矩阵,L 为单位下三角阵,U 为上三角阵,其中U 矩阵的主对角线元素之积即为det A 。 由于A 的元素零元素较多,为节省储存量,将A 的元素存为6×501的数组中,程序中采用get_an_element()函数来从小数组中取出A 中的元素。 2.全部源程序 #include #include void init_a();//初始化A double get_an_element(int,int);//取A 中的元素函数 double powermethod(double);//原点平移的幂法 double inversepowermethod(double);//原点平移的反幂法 int presolve(double);//三角LU 分解 int solve(double [],double []);//解方程组 int max(int,int); int min(int,int); double (*u)[502]=new double[502][502];//上三角U 数组 double (*l)[502]=new double[502][502];//单位下三角L 数组 double a[6][502];//矩阵A int main() { int i,k; double lambdat1,lambdat2,lambda1,lambda501,lambdas,mu[40],det;

数值分析实验报告2

实验报告 一、实验名称 复合梯形求积公式、复合辛普森求积公式、龙贝格求积公式及自适应辛普森积分。 二、实验目的及要求 1. 掌握复合梯形求积计算积分、复合辛普森求积计算积分、龙贝格求积计算积分和自适应辛普森积分的基本思路和步骤. 2. 培养Matlab 编程与上机调试能力. 三、实验环境 计算机,MATLAB 软件 四、实验内容 1.用不同数值方法计算积分9 4 ln 1 0-=? xdx x 。 (1)取不同的步长h 。分别用复合梯形及复合辛普森求积计算积分,给出误差中关于h 的函数,并与积分精确指比较两个公式的精度,是否存在一个最小的h ,使得精度不能再被改善。 (2)用龙贝格求积计算完成问题(1)。 (3)用自适应辛普森积分,使其精度达到10-4。 五、算法描述及实验步骤 1.复合梯形公式 将区间[a,b]划分为n 等份,分点x k =a+ah,h=(b-a)/h,k=0,1,...,n ,在每个子区间[x k ,x k +1](k=0,1,...,n-1)上采用梯形公式(),得 )]()([2 )(b f a f a b dx x f b a +-≈? () )]()(2)([2)]()([21 1 110b f x f b f h x f x f h T n k k k n k k n ++=+=∑∑-=+-= () ),(),(12 )(' '2b a f h a b f R n ∈-- =ηη () 其中Tn 称为复合梯形公式,Rn 为复合梯形公式的余项。 2.复合辛普森求积公式 将区间[a,b]划分为n 等份,在每个子区间[x k ,x k +1](k=0,1,...,n-1)上采用辛普森公式(),得 )]()2 (4)([6b f b a f a f a b S +++-= ()

相关文档
相关文档 最新文档