計(jì)算方法實(shí)例

第1講計(jì)算方法（數(shù)值分析）計(jì)算方法研究者的鼻祖應(yīng)該是泰勒等，計(jì)算方法的形成與發(fā)展應(yīng)該是與計(jì)算機(jī)的產(chǎn)生與發(fā)展結(jié)伴而形的。記得有部電影《橫空出世》中描寫科技人員為了得到制造原子彈的數(shù)據(jù)，有成百上千的人在用古老的計(jì)算工具一一算盤，晝夜計(jì)算。由于計(jì)算機(jī)只能表示離散的數(shù)（且只能表示有限位、有限個數(shù)），因此，使用計(jì)算機(jī)進(jìn)行計(jì)算，通常會產(chǎn)生誤差，得到的結(jié)果有些是“可接受”的，有些是“不可接受”的。于是，通常要研究以下兩個內(nèi)容：所求結(jié)果的近似值，計(jì)算過程產(chǎn)生的誤差；得到結(jié)果所用的計(jì)算方法是否可行、快速，節(jié)省計(jì)算機(jī)資源.其內(nèi)容大致可分為三個部分：數(shù)值逼近；數(shù)值代數(shù)；微分方程的數(shù)值求解第一章計(jì)算方法的一般概念§1.1引言通過研究下面這個例題給出所要討論的幾個方面的內(nèi)容。TOC\o"1-5"\h\zEx.1.1.1建立/廣』1：5必的遞推公式并求當(dāng)n=0,1,2,…時(shí),/“的值.解I=|1xn+5對T—5對-1dx=j1"1dx—5』1X"Tdx=1-5I,n0x+5 o ox+5n n-1而10=j1—dx=ln6-ln5.I=ln6-ln5(1.1)0(1.1)< 1 _I=—-5I,n=1,2，…nn n-1由補(bǔ)Ex1_1_1(迭代次數(shù))的結(jié)果知：(1)用Maple，以系統(tǒng)默認(rèn)的精度計(jì)算10-120時(shí)：evalf(Ex1_1_1(0)); ^182321557 \# 系統(tǒng)默認(rèn)取十位evalf(Ex1_1_1(10)); .02evalf(Ex1_1_1(20)); 100000.(2)用Maple，以較高精度(取20位)，計(jì)算10-120有：evalf(Ex1_1_1(0),20); .1823215567939546262evalf(Ex1_1_1(10),20); .015367550448evalf(Ex1_1_1(20),20); .00800TOC\o"1-5"\h\zevalf(Ex1_1_1(26),20); .1evalf(Ex1_1_1(27),20); 0.evalf(Ex1_1_1(28),20); 0.(3)用Maple,以更高精度(取30位),有：evalf(Ex1_1_1(26),30); .00621001212evalf(Ex1_1_1(26),50); .0062100121272123878307412230028即不同精度要求，所得結(jié)果相差很大。為什么會這樣呢？哪個結(jié)果是120的近似值呢？這可由式(1.1)本身的理論研究給出答案。事實(shí)上，按理論推導(dǎo)有/1=-5/0+1一一. .. /V .. . . ... ..設(shè)計(jì)算機(jī)計(jì)算/°時(shí)的結(jié)果為/°,則代入上式經(jīng)計(jì)算機(jī)計(jì)算有/=-5/+1

10記/0-1記/0-10=S，則有同理/-/=-5/-/=-5￡

11 00/廣。=一5(、1 -1)=(FWn-2-/n-2)= =(一5"(/0F=(-5)n￡.可見,誤差/-7的大小既與/-I的誤差￡有關(guān),也與迭代次數(shù)n有關(guān).盡管￡可以取得非常nn 0 01小，但誤差的傳播逐步擴(kuò)大，結(jié)果并不可靠.當(dāng)￡<<—時(shí),才足夠精確.5n與書上的思路類似，也可以用下面的關(guān)系研究積分的值的近似值：0</二』1—dx〈L』1xndx= T0(nT3).n0x+5 50 5(n+1)有關(guān)計(jì)算程序參后面的補(bǔ)1.1.1。Ex.1.1.2求方程工2+(a+P)x+109=0的根，其中a=-109,。=-1.解由求根公式得fx1,2-b±\b2-4x1,22a其中a=1,-b=-a-P=109+1=0.1x1010+0.000,000,000,1x1010,c=109=0.1x1010

b2-4ac=0.01x1020-4x1x0.1x1010=0.01x1020-0.000,000,000,04x1020.在只能取到小數(shù)點(diǎn)后面第8位的計(jì)算機(jī)上對階計(jì)算后,將得到-b=0.1x1010+0x1010=0.1x1010，--a吃掉了。-—加減運(yùn)算時(shí)大吃小b2-4ac=0.01x1020-0=0.01x1020--b2吃掉了4ac.這樣得到的根將是0.1x1010+0.1x101。 “ 0.1x1010-0.1x1010八=109, x2= 21 =0,這樣得到的根將是0.1x1010+0.1x1010小 0.1x1010-0.1x1010八x1= 2_1 =109,x2= 2~1 =0,而這相當(dāng)于P,c在方程中沒有出現(xiàn),即所求的是方程x2+以x=0的根.這雖然不正確.為了避免以上情況并考慮到分子部分可能出現(xiàn)三個相近的數(shù)相減而導(dǎo)致有效數(shù)字的嚴(yán)重?fù)p失的不利情況,在計(jì)算機(jī)上求ax2+bx+c=0的根時(shí)，是按下列步驟進(jìn)行的(退化情況a=0或b=0另行考慮)-b-sgn(b)\b2-4ac2ax2a1cx= 2ax1就算b2吃掉了4ac,第一式分子仍是兩個同號數(shù)的和,在求x］時(shí)仍用了c;或4ac吃掉了b2,第一式分子第一項(xiàng)還含有b.這樣一般就不再出現(xiàn)本例直接用求根公式計(jì)算過程中出現(xiàn)的情況.實(shí)驗(yàn)：法1——用求根函數(shù)；法2——按上式編程Ex.1.1.3表示p(x)=axn+axn-1+ +a,并求〃=10,a=^—,x=1時(shí)的值.0 1 n kk2解通常計(jì)算勤⑴要用n+(n-1)+...+2+1=也口次乘法,n次加法而采用我國13世2紀(jì)的數(shù)學(xué)家秦九韶所首創(chuàng)的所謂秦九韶算法p(x)=(…((ax+a)x+a)x+ +a)x+an 0 1 2 n-1 n則只需要n次乘法與n次加法.程序算法：b=a0b^=b^x+a^,k=1,2,—,np(x)=bn n實(shí)驗(yàn)程序見后面：實(shí)驗(yàn)程序一Ex.1.13__Maple:Ex.1.1.4計(jì)算方的近似值，精確到10-5.R(二)= ——，日位于0與1之間n2 (n+1)!2n+1 2解,/ex=2土=X七+R(x),R(x)=——―xn+1,&位于0與工之間.R(二)= ——，日位于0與1之間n2 (n+1)!2n+1 2...=e2=2^^+Rk=0k!.2k n1欲使R(1)<10-5,只要 一一—<——=—1—<-^,故取〃=6.n2 (n+1)!2n+1 (n+1)!-2n+1 (n+1)!-2n 100000(7!*2人6=322560>100000)為保證前7項(xiàng)計(jì)算結(jié)果滿足精度0.00001的要求，每項(xiàng)都要精確到第7位，經(jīng)四舍五入累加到第六位，再四舍五入也能保證最終精度.(因?yàn)樗愕叫?shù)點(diǎn)后的6位再進(jìn)行四舍五入可能會產(chǎn)生大于10-5的誤差.).?一v6 1「. \evk廣kk=0練習(xí)題：一、給出前21個結(jié)果(以任何方式算出)補(bǔ)1.1.1：用Maple計(jì)算：(1)對于浮點(diǎn)數(shù)，Maple將始終按照輸入時(shí)的精度對它進(jìn)行處理，而且會將計(jì)算結(jié)果同算式中小數(shù)點(diǎn)后有效位數(shù)最多的一個因子保持精度一致；(2)對于非浮點(diǎn)運(yùn)算，包括含e,n,sin(x),ln(x)等,Maple會作精確計(jì)算(因?yàn)樗鼤鞣栠\(yùn)算)，必要時(shí)用evalf(%,n)給出輸出結(jié)果；(3)不加位數(shù)n給出15位。restart;Ex1_1_1:=proc(n::numeric)ifn<1thenf(x):=xAn/(x+5);int(f(x),x=0..1);else1/n-5*Ex1_1_1(n-1);endif;end;Ex1_1_1:=proc(n::numeric)ifn<1thenf(x):=xAn/(x+5);int(f(x),x=0..1)else1/n-5xEx1_1_1(n-1)endifendprocevalf(Ex1_1_1(0));.182321558evalf(Ex1_1_1(1));.088392212evalf(Ex1_1_1(2));.05803894evalf(Ex1_1_1(3));.0431386evalf(Ex1_1_1(4));.034307evalf(Ex1_1_1(5));.028466evalf(Ex1_1_1(6));.02434evalf(Ex1_1_1(7));.0212evalf(Ex1_1_1(8));.0191evalf(Ex1_1_1(9));.016evalf(Ex1_1_1(10));.02evalf(Ex1_1_1(11));-.02evalf(Ex1_1_1(12));.3evalf(Ex1_1_1(13));0.evalf(Ex1_1_1(14));5.evalf(Ex1_1_1(15));-20.evalf(Ex1_1_1(16));200.evalf(Ex1_1_1(17));-1000.evalf(Ex1_1_1(18));3000.evalf(Ex1_1_1(19));-20000.evalf(Ex1_1_1(20));100000.結(jié)果分析：1）計(jì)算I20，用系統(tǒng)默認(rèn)位數(shù)給出結(jié)果為134922987458706126014977597520>Ex1_1_1(20); 134922987458706126014977597520卜95367431640625ln(2)+95367431640625ln(3)-95367431640625ln(5)134922987458706126014977597520取精確結(jié)果134922987458706126014977597520>1349229874587061260149/77597520;>evalf(%); .17387538601014取精確結(jié)果后三項(xiàng)單獨(dú)計(jì)算得：>95367431640625*ln(2)+95367431640625*ln(3)-95367431640625*ln(5);95367431640625ln(2)+95367431640625ln(3)-95367431640625ln(5)>evalf(%);>evalf(%);.1738753871014用后一個結(jié)果減去前一個結(jié)果時(shí)，計(jì)算機(jī)是先對階計(jì)算小數(shù)的差，得0.000,000,001，再乘1014因?yàn)檫@個數(shù)有9位小數(shù)，所以乘后的結(jié)果就為10000.2)如果是取20位，分別得.173875386041597883561014.17387538604159796361014后者小數(shù)減去前者小數(shù)后再將小數(shù)點(diǎn)右移14位得0.79636-0.788356=.008004雖然這一結(jié)果與結(jié)果evalf(Ex1_1_1(20),20);.00800有所差異，但已屬正常又：forifrom0to20doevalf(Ex1_1_1(i),20);od;.1823215567939546262.0883922160302268689.058038919848865656.04313873408900506.0343063295549747.0284683522251264.024324905541035.02123261515197.01883692424015.0169264899103.015367550448.014071338668

.01297664000.0120398770.0112291867.010520734.00989634.0093419.0088462.008400.00800實(shí)驗(yàn)程序——Ex.1.13__Maple:法1：restart;>with(linalg):Ex1_1_3:=proc(a,x)>localn,b,k;n:=linalg[vectdim](a);>b:=a[1];forkfrom2tondo>b:=b*x+a[k];od;print(p(x)=b);>print(p(x)=simplify(%));end：c:=array(1..5);c:=array(1..5,[])Ex1_1_3(c,x);p(x)=(((cx+c)x+c)x+c)x+c

1 2 3 4 5p(x)=cx4+x3c+x2c+xc+c1 2 3 4 5當(dāng)輸入多項(xiàng)式系數(shù)序列及X后就可具體計(jì)算出值來.>c:=array(1..5);c:=array(1..5,[])Ex1_1_3(c,x);TOC\o"1-5"\h\zp(X)=cx4+x3c+x2c+xc+c

1 2 3 4 5>c:={seq(1/k人2,k=1??10)};-11 1111111、c:={1 },,, ,,, ,,64811009162536494d:=convert(c