數(shù)值分析也稱計(jì)算方法

A第一章誤差與范數(shù)數(shù)值分析也稱計(jì)算方法，它不僅是一種研究并解決數(shù)學(xué)問題的數(shù)值近似解的方法，而且是在計(jì)算機(jī)上使用的解數(shù)學(xué)問題的方法.它的計(jì)算對象是那些在理論上有解，而無法用手工計(jì)算的數(shù)學(xué)問題.1.1誤差的來源在運(yùn)用計(jì)算方法解決實(shí)際問題的過程中，會(huì)出現(xiàn)各種各樣的誤差，必須注重誤差分析.否則，一種合理的計(jì)算也可能得出錯(cuò)誤的結(jié)果.例1.1.1用差商f'(a)Mf(a+")—f(a)求f(x)=Inx在x=3處導(dǎo)數(shù)的近似值.h取h=0.1和h=0.0001，用手工計(jì)算，取五位數(shù)字計(jì)算；取h=0.1，h=0.0001，h=0.000000000000001和h=0.0000000000000001分別用MATLAB軟件計(jì)算，取十五位數(shù)字計(jì)算；比較以上的運(yùn)算結(jié)果，說明是否|h|越小則計(jì)算結(jié)果越準(zhǔn)確.解根據(jù)導(dǎo)數(shù)定義，可以用差商求f(x)=lnx在x=3處導(dǎo)數(shù)的近似值f'(3)mf'(3)mf(3+h)-f(3)ln(3+h)—ln3

h從理論上講，Ih從理論上講，Ih|越小則計(jì)算結(jié)果越準(zhǔn)確.如果用手工計(jì)算，取五位數(shù)字計(jì)算.當(dāng)h=0.1時(shí)，得1.131f'(3)m4-1.09861.131f'(3)m4-1.0986=0.3280，0.1與導(dǎo)數(shù)的精確值f'(3)=1=0.33333比較這項(xiàng)計(jì)算還是可取的.但是當(dāng)h=0.0001時(shí)，得f'(3)Mf'(3)M1.0986—1.0986=0，0.0001算出的結(jié)果反而毫無價(jià)值.如果應(yīng)用MATLAB軟件計(jì)算，取十五位數(shù)字計(jì)算，結(jié)果就完全不同了.在MATLAB工作窗口輸入下面程序>>a=3;h=0.1;y=log(a+h)-log(a);yx=y/h運(yùn)行后得yx=0.32789822822991將此程序中h改為0.0001，運(yùn)行后得yx=0.33332777790385后者比前者好.再取h=0.000000000000001，運(yùn)行后得yx=0.44408920985006不如前者好.取h=0.0000000000000001，運(yùn)行后得yx=0算出的結(jié)果反而毫無價(jià)值.

例1.1.2分別求方程組AX=b在下列情況時(shí)的解，其中A=""1 i.01)（1）b="2.01j解（1）首先將方程組AX=b化為同解方程X=A-1b，然后在MATLAB工作窗口輸入程序>>b=[2,2]';A=[1,1;1,1.01];X=A\b運(yùn)行后輸出當(dāng)b=時(shí)，AX運(yùn)行后輸出當(dāng)b=時(shí)，AX=b的解為X=（2）同理可得，當(dāng)bAX=b的解為X例1.1.3計(jì)算e的近似值.解泰勒級(jí)數(shù)取X=1，得這是一個(gè)無限過程，誤差.如果取有限項(xiàng)1_

這是一個(gè)無限過程，誤差.如果取有限項(xiàng)1_

一+3!計(jì)算機(jī)無法求到精確值.只能在（1.2）（1.2）取有限項(xiàng)時(shí)計(jì)算，再估計(jì)s(1)=1+1+—+—+—+...+n 2! 3! 4!(0<9<1).作為e的值必然會(huì)有誤差，根據(jù)泰勒余項(xiàng)定理可知其截?cái)嗾`差為(0<9<1).Sn^—(n+1)!如果取(1.2)的前九項(xiàng)，輸入程序>>n=8;s=1;S=1;fork=1:ns=s*k;S=S+1/s,ends,S,R=3/(s*(n+1))或>>S1=1+1+1/2+1/(1*2*3)+1/(1*2*3*4)+1/(1*2*3*4*5)+1/(1*2*3*4*5*6)+1/(1*2*3*4*5*6*7)+1/(1*2*3*4*5*6*7*8),R1=3/(1*2*3*4*5*6*7*8*9)運(yùn)行后結(jié)果R=8.267195767195768e-006R=8.267195767195768e-0062.71827876984127因?yàn)榻財(cái)嗾`差為所以e的近似值e a3<一a所以e的近似值e a3<一a8.267(8+1)! 9!1S8⑴=1+1+2!+196x10-6(0<9<1),1+3!+6!1-2.71828.8!1.2誤差和有效數(shù)字例1.2.1取2.71828作為e的四舍五入近似值時(shí)，求其絕對誤差和相對誤差.解在MATLAB工作窗口輸入程序>>juewu=exp(1)-2.71828運(yùn)行后輸出結(jié)果為juewu=1.828459045505326e-006例1.2.2計(jì)算』：Sin*dx的近似值，并確定其絕對誤差和相對誤差.解因?yàn)楸环e函數(shù)的原函數(shù)不是初等函數(shù)，故用泰勒級(jí)數(shù)求之.x(一3<x<3), (1.5)TOC\o"1-5"\h\zsinx x2 x4(一3<x<3), (1.5)=1- + - + +x 3! 5! 7! 9!這是一個(gè)無限過程，計(jì)算機(jī)無法求到精確值.可用(這是一個(gè)無限過程，計(jì)算機(jī)無法求到精確值.可用(1.5)的前四項(xiàng)］一M+x4x±代7!"』：土dx提:(1-土+"』：土dx提:(1-土+0x 0 3!(-)5q 55!x6)dx=壬一一+5! 7! 2 33!根據(jù)泰勒余項(xiàng)定理和交錯(cuò)級(jí)數(shù)收斂性的判別定理，得到絕對誤差兀n()9R=Iy-yk-^-WU，

9-9!在MATLAB命令窗口輸入計(jì)算程序如下：symsxf=1-xA2/(1*2*3)+xA4/(1*2*3*4*5)-xA6/(1*2*3*4*5*6*7)y=int(f,x,0,pi/2),y1=double(y)y11=pi/2-(pi/2)A3/(3*3*2)+(pi/2)A5/(5*5*4*3*2)-(pi/2)A7/(7*7*6*5*4*3*2)infd=double(inf),inf=int(sin(x)/x,x,0,pi/2)WU=(pi/2)A9/(9*9*8*7*6*5*4*3*2),R=infd-y11因?yàn)檫\(yùn)行后輸出結(jié)果為：y=1.37076216815449，y-1.37074466418938，R=1.75039651049147e-005，WU-1.782679830970664e-005<10-4.所以，y的絕對誤差為￡=10-4，故y=J：-S!宜dx*1.3707.y的相對誤差為0x￡ 10-4￡^=——=E0—<0.0073%.1.3誤差估計(jì)的基本方法下面的例子揭示了初始誤差的穩(wěn)定性傳播和不穩(wěn)定性傳播.例1.3.4設(shè)計(jì)三種算法求方程2x2+x-15=0在(2,3)的一個(gè)正根x*的近似值，并研究每種算法的誤差傳播情況.解為解已知方程，我們可以設(shè)計(jì)如下三種算法，然后將計(jì)算結(jié)果與此方程的精確解X*=2.5比較，觀察誤差的傳播.算法1將已知方程化為同解方程x=15-2x2.取初值x0=2，按迭代公式\+1=15-2X2依次計(jì)算X],X2,…,Xn,…，結(jié)果列入表1-3中.算法2將已知方程化為同解方程x=;%.取初值x0=2，按迭代公式15x= k+1 2x+1依次計(jì)算X,X,...,X,...，結(jié)果列入表1-3中.1 2 n算法3將已知方程化為同解方程x=x-之X2+X-15.取初值x0=2，按迭代公式為2x2+x-15X+1=Xk- "4x」1依次計(jì)算x,x,...,x,...，結(jié)果列入表1-3中.12 n我們?yōu)檫@種三算法的計(jì)算編寫兩套MATLAB程序如下：(1)MATLAB主程序function[k,juecha,xiangcha,xk]=liti112(x0,x1,limax)%輸入的量--x0是初值,limax是迭代次數(shù)和精確值x;%輸出的量--每次迭代次數(shù)k和迭代值xk,% --每次迭代的絕對誤差juecha和相對誤差xiangcha，x(1)=x0;fori=1:limaxx(i+1)=fl(x(i));%程序中調(diào)用的fl.mjuecha=abs(x(i)-x1);xiangcha=juecha/(abs(x(i))+eps);xk=x(i);k=i-1;[(i-1),juecha,xiangcha,xk]end(2)MATLAB調(diào)用函數(shù)程序及其計(jì)算結(jié)果算法2的MATLAB調(diào)用函數(shù)程序functiony1=fl(x)y1=15/(2*x+1);將MATLAB主程序和調(diào)用函數(shù)程序分別命名liti112.m和fl.m，分別保存為M文件，然后在MATLAB工作窗口輸入命令>>[k,juecha,xiangcha,xk]=liti112(2,2.5,100)運(yùn)行后輸出計(jì)算結(jié)果列入表1-3和表1-4中.將算法2的MATLAB調(diào)用函數(shù)程序的函數(shù)分別用y1=15-2*xA2和y1=x-(2*xA2+x-15)/(4*x+1)代替，得到算法1和算法3的調(diào)用函數(shù)程序，將其保存，運(yùn)行后將三種算法的前8個(gè)迭代值x『x2,…，x8列在一起(見表1-3)，進(jìn)行比較.將三種算法的x1,x2,…，x8對應(yīng)的絕對誤差和相對誤差的值列在一起(見表1—4)，進(jìn)行比較.表1-3例1.3.4中三種算法的計(jì)算結(jié)果-算法迭代次數(shù)算法1的迭代結(jié)果算法2的迭代結(jié)果算法3的迭代結(jié)果022.000000002.00000000

第一章誤差與范數(shù)173.000000002.555555562-832.142857142.500550063-137632.837837842.500000064-3788403232.246963562.500000005-2.8704x10172.246963562.500000006-1.6478x10352.321774842.500000007-5.4307x10702.657901652.50000000:99:-Inf:2.50000001:2.50000000由表1-3和表1-4可以看出，算法1的初始絕對誤差和相對誤差分別是4.50000000和0.64285714,它們是三種算法中最大者，由算法1得到的迭代序列b.｝對應(yīng)的絕對誤差是不穩(wěn)定的，且當(dāng)nrs時(shí)，按指數(shù)無界增長，相對誤差趨近于不定值，x.T-3，即ix.｝不收斂；算法2的初始絕對誤差和相對誤差分別是0.50000000和0.16666667,它們是三種算法中居于中間，由算法2得到的迭代序列｛x.｝對應(yīng)的絕對誤差和相對誤差是穩(wěn)定的，且當(dāng)nr3時(shí)，按指數(shù)緩慢遞減且趨近于0，k｛xk｝緩慢地?cái)[動(dòng)式收斂于精確解2.5（參見算法2的運(yùn)行結(jié)果）；算法3的初始絕對誤差和相對誤差分別是0.05555560和0.02173913,它們是三種算法中最小者，由算法3得到的迭代序列｛x.｝對應(yīng)的絕對誤差和相對誤差是穩(wěn)定的，且當(dāng)nrs時(shí)，按指數(shù)快速遞減且趨近于0，｛x.｝快速地遞減且收斂于精確解2.5. *表1-4例1.3.4中三種算法計(jì)算結(jié)果的誤差算法迭代次數(shù)算法1的i沃差算法2的誤差算法3的誤差絕對誤差相對誤差絕對誤差相對誤差絕對誤差相對誤差01234567:990.50000000000000000013765.500003788403262.87039981x10171.64783901x10355.43074680x1070:Inf0.250000000.642857141.030120480.000100021.00000001111:NaN0.500000000.500000000.357142860.337837840.253036440.230287040.178225160.15790165:0.000000010.250000000.166666670.1666666700.119047620.112612620.084345480.076762470.05940839:0.000000000.500000000.055555600.000550100.000000060.00000000000:00.250000000.021739130.000219970.000000020.00000000000:01.4數(shù)值計(jì)算中應(yīng)注意的問題從例1.3.4可以看出，一個(gè)問題的解決，往往要經(jīng)過多次運(yùn)算.每一步運(yùn)算都可能產(chǎn)生誤差，在反復(fù)多次計(jì)算的過程中，必然產(chǎn)生誤差的傳播和積累.顯然，當(dāng)誤差積累偏大時(shí)，會(huì)使計(jì)算結(jié)果失真.因此，在每一步計(jì)算中，都應(yīng)該防止產(chǎn)生誤差升級(jí)的現(xiàn)象.例1.4.1求數(shù)x=715x（5+8-19-1）的近似值.解（1）直接用MATLAB命令>>x=(7A15)*(sqrt(1+8A(-19))-1)運(yùn)行后輸出結(jié)果x=0問題出現(xiàn)在兩個(gè)相近的數(shù)V'1+8*與1相減時(shí)，計(jì)算機(jī)運(yùn)行程序>>sqrt(1+8A(-19))-1運(yùn)行后輸出結(jié)果ans=0由于計(jì)算機(jī)硬件只支持有限位機(jī)器數(shù)的運(yùn)算，因此在計(jì)算中可能引入和傳播舍入誤差.因?yàn)橛行?shù)字的嚴(yán)重?fù)p失，導(dǎo)致輸出司+8-19-1的結(jié)果為0，計(jì)算機(jī)不能再與數(shù)715繼續(xù)進(jìn)行真實(shí)的計(jì)算，所以，最后輸出的結(jié)果與X的精確值不符.(2)如果化為: 715x8-19x=715x(K1+8-19—1)= —，<1+8-19+1再用MATLAB命令>>x=(7A15)*((8A(-19))/(sqrt(1+8A(-19))+1))運(yùn)行后輸出結(jié)果x=1.6471e-005這是因?yàn)樗?8-19-1化為 一后，計(jì)算機(jī)運(yùn)行程序V1+8-19+1>>x=(8A(-19))/(sqrt(1+8A(-19))+1)運(yùn)行后的結(jié)果為x=3.4694e-018由于有效數(shù)字的損失甚少，所以運(yùn)算的結(jié)果3.4694x10-18再與715繼續(xù)計(jì)算，最后輸出的結(jié)果與X的精確值相差無幾.例1.4.2求數(shù)y=ln(30-2302—1)的近似值.解(1)直接用MATLAB程序>>x=30;x1=sqrt(xA2-1)運(yùn)行后輸出結(jié)果x1=29.9833輸入MATLAB程序>>x=30;x1=29.9833;y=log(x-x1)運(yùn)行后輸出結(jié)果y=-4.0923(2)因?yàn)閘n(30-202-1)中的x=30很大，如果采用倒數(shù)變換法z=X-\X2-1= j，即1 : ln(30-*302—1)=ln ； =-ln(30+V900-1).輸入MATLAB程序>>x=30;y=-log(x+sqrt(xA2-1))運(yùn)行后輸出結(jié)果y=-4.0941(3)輸入MATLAB程序>>x=30;y=log(x-sqrt(xA2-1))運(yùn)行后輸出結(jié)果y=-4.0941可見，（2）計(jì)算的近似值比（1）的誤差小.參加計(jì)算的數(shù)，有時(shí)數(shù)量級(jí)相差很大.如果不注意采取相應(yīng)的措施，在它們的加減法運(yùn)算中，絕對值很小的那個(gè)數(shù)經(jīng)常會(huì)被絕對值較大的那個(gè)數(shù)“吃掉”，不能發(fā)揮其作用，造成計(jì)算結(jié)果失真.例1.4.4請?jiān)?6位十進(jìn)制數(shù)值精度計(jì)算機(jī)上利用軟件MATLAB計(jì)算下面的兩個(gè)數(shù)x*=111111111111111+0.1+0.3和y*=1111111111111111+0.1+0.3將計(jì)算結(jié)果與準(zhǔn)確值比較，解釋計(jì)算結(jié)果.解在MATLAB工作窗口輸入下面程序>>x=11111111111111+0.1+0.3,y=1111111111111111+0.1+0.3運(yùn)行后輸出結(jié)果x=1.11111111111114e+014，y=1.111111111111111e+015從輸出的結(jié)果可以看出，x=x*，而y豐y*.為什么y*僅僅比X*多一位1，而y豐y*呢？這是因?yàn)橛?jì)算機(jī)進(jìn)行運(yùn)算時(shí)，首先要把參加運(yùn)算的數(shù)寫成絕對值小于1而“階碼”相同的數(shù)，這一過程稱為數(shù)的“對階”.在16位十進(jìn)制數(shù)值精度計(jì)算機(jī)上利用軟件MATLAB計(jì)算這兩個(gè)數(shù)，把運(yùn)算的數(shù)x*寫成浮點(diǎn)規(guī)格化形式為X*=0.1111111111111110X1015+0.0000000000000001X1015+0.0000000000000003X1015，在16位十進(jìn)制數(shù)值精度計(jì)算機(jī)上，三項(xiàng)的數(shù)都表示為小數(shù)點(diǎn)后面16位數(shù)字的數(shù)與1015之積，所以，計(jì)算機(jī)沒有對數(shù)進(jìn)行截?cái)啵前丛瓉淼娜齻€(gè)數(shù)進(jìn)行計(jì)算.因此，計(jì)算的結(jié)果x=x*.而y*=0.1111111111111111X1016+0.00000000000000001X1016+0.00000000000000003X1016三項(xiàng)的數(shù)都表示寫成絕對值小于1而“階碼”都為1016的數(shù)以后，第一項(xiàng)的純小數(shù)的小數(shù)點(diǎn)后面有16位數(shù)字.但是，后兩項(xiàng)數(shù)的純小數(shù)的小數(shù)點(diǎn)后面有17位數(shù)字，超過了16位十進(jìn)制數(shù)值精度計(jì)算機(jī)的存儲(chǔ)量，計(jì)算機(jī)對后兩項(xiàng)的數(shù)都進(jìn)行截?cái)嘧詈笠晃?，即后兩?xiàng)的數(shù)都是16位機(jī)上的零，再進(jìn)行計(jì)算，所以計(jì)算結(jié)果與實(shí)際不符.1.5向量和矩陣的范數(shù)在研究方程組近似解的誤差估計(jì)和迭代法收斂性的過程中，需要對向量和矩陣的“大小”加以度量，向量范數(shù)和矩陣范數(shù)正是這種度量指標(biāo)，在數(shù)值分析中起著重要作用.1.5.1向量范數(shù)與MATLAB命令在3維空間中，將向量x=（x,x,x）T的長度記作||x||，通常用||x||=（蕓x2）：來1 2 3 ii=1定義，向量范數(shù)的MATLAB命令及功能列入表1-5：表1-5向量范數(shù)的MATLAB命令及功能MATLAB命令等價(jià)的命令功能Xp=norm(X,P)sum(abs(X).AP)A(1/P)輸入X為向量，輸出為X的p范數(shù)1<p<3

第一章誤差與范數(shù)X2=norm(X)或X2=norm(X,2)sum(abs(X).A2)A(1/2)輸入X為向量，輸出為X的2范數(shù)X1=norm(X,1)sum(abs(X))輸入X為向量，輸出為X的1范數(shù)Xw=norm(X,inf)max(a