《数值计算方法精》word版

《《数值计算方法精》word版》由会员上传分享，免费在线阅读，更多相关内容在应用文档-天天文库。

1、数值计算方法5.1引言本章将花较大的篇幅讨论若干常见数值计算问题：线性分析、一元和多元函数分析、微积分、数据分析、以及常微分方程求解等。但与一般数值计算教科书不同，本章的讨论重点是：如何利用现有的世界顶级数值计算资源MATLAB。至于数学描述，本章将遵循“最低限度自封闭”的原则处理，以最简明的方式阐述理论数学、数值数学和MATLAB计算指令之间的内在联系及区别。对于那些熟悉其他高级语言（如FORTRAN，Pascal，C++）的读者来说，通过本章，MATLAB卓越的数组处理能力、浩瀚而灵活的M函数指令、丰富而友善的图形显示指令将使他们体验到解题视野的豁然开

2、朗，感受到摆脱烦琐编程后的眉眼舒展。对于那些经过大学基本数学教程的读者来说，通过本章，MATLAB精良完善的计算指令，自然易读的程序将使他们感悟“教程”数学的基础地位和局限性，看到从“理想化”简单算例通向科学研究和工程设计实际问题的一条途径。对于那些熟悉MATLAB基本指令的读者来说，通过本章，围绕基本数值问题展开的内容将使他们体会到各别指令的运用场合和内在关系，获得综合运用不同指令解决具体问题的思路和借鉴。由于MATLAB的基本运算单元是数组，所以本章内容将从矩阵分析、线性代数的数值计算开始。然后再介绍函数零点、极值的求取，数值微积分，数理统计和分析，拟

3、合和插值，Fourier分析，和一般常微分方程初值问题。本章的最后讨论稀疏矩阵的处理，因为这只有在大型问题中，才须特别处理。从总体上讲，本章各节之间没有依从关系，即读者没有必要从头到尾系统阅读本章内容。读者完全可以根据需要阅读有关节次。除特别说明外，每节中的例题指令是独立完整的，因此读者可以很容易地在自己机器上实践。5.1LU分解和恰定方程组的解5.1.1LU分解、行列式和逆（1）LU分解（2）行列式和逆5.1.2恰定方程组的解【*例5.2.2-1】“求逆”法和“左除”法解恰定方程的性能对比（1）为对比这两种方法的性能，先用以下指令构造一个条件数很大的高阶

4、恰定方程。rand('state',12);%选定随机种子，目的是可重复产生随机阵A。A=rand(100,100)+1.e8;%rand(100,100)生成（100×100）均匀分布随机矩阵。%每个随机阵元素加的目的是使A阵条件数升高。x=ones(100,1);%令解向量x为全1的100元列向量。b=A*x;%为使Ax=b方程一致，用A和x生成b向量。cond(A)%求A阵的条件数。ans=1.4426e+012（2）“求逆”法解恰定方程的误差、残差、运算次数和所用时间flops(0);tic%浮点运算计数器置0；启动计时器StopwatchTime

5、r48xi=inv(A)*b;%xi是用“求逆”法解恰定方程所得的解。ti=toc%关闭计时器，并显示解方程所用的时间。ci=flops%“求逆”法解方程所用的运算次数eri=norm(x-xi)%解向量xi与真解向量x的范-2误差。rei=norm(A*xi-b)/norm(b)%方程的范-2相对残差ti=0.9300ci=2070322eri=3.0708e-004rei=6.6280e-007（3）“左除”法解恰定方程的误差、残差、运算次数和所用时间flops(0);tic;xd=Ab;%是用“左除”法解恰定方程所得的解。td=toc,cd=flo

6、ps,erd=norm(x-xd),red=norm(A*xd-b)/norm(b)td=0.2200cd=741872erd=3.2243e-004red=2.0095e-0165.1.1范数、条件数和方程解的精度【*例5.2.3-1】Hilbert矩阵是著名的病态矩阵。MATLAB中有专门的Hilbert矩阵及其准确逆矩阵的生成函数。本例将对方程近似解和准确解进行比较。所谓n阶Hilbert矩阵的形式是：。N=[68101214];%本例计算的矩阵阶数fork=1:length(N)n=N(k);%矩阵的阶H=hilb(n);%产生n阶Hilbert矩

7、阵Hi=invhilb(n);%产生完全准确的n阶逆Hilbert矩阵b=ones(n,1);%生成n阶全1向量x_approx=Hb;%利用左除H求近似解x_exact=Hi*b;%利用准确逆Hilbert矩阵求准确解ndb=norm(H*x_approx-b);nb=norm(b);ndx=norm(x_approx-x_exact);nx=norm(x_approx);er_actual(k)=ndx/nx;%实际相对误差K=cond(H);%计算Hilbert矩阵的条件数er_approx(k)=K*eps;%最大可能的近似相对误差er_max(

8、k)=K*ndb/nb;%最大可能的相对误差enddisp('Hi