計(jì)算方法實(shí)習(xí)報(bào)告

一．給出一種有效的算法和無(wú)效的算法計(jì)算積分y(n)=∫(x^n)/(4x+1)dx,n=0,1,2,…,10,積分限為（0，1）有效算法運(yùn)用遞推公式y(tǒng)(n)=-y(n-1)/4+1/(4n)，取y0=(log5)/4程序?yàn)椋?include<iostream.h>#include<math.h>voidmain(){doubley0,y1; y0=1/4.0*log(5.0); cout<<"y0="<<y0<<""; for(intn=1;n<=10;n++) { y1=-1.0/4.0*y0+1.0/(4.0*n); cout<<"y"<<n<<"="<<y1<<"";y0=y1; if(n%3==0)cout<<endl; }}其成果為：y0=0.402359y1=0.14941y2=0.0876475y3=0.0614215y4=0.0471446y5=0.0382138y6=0.0321132y7=0.027686y8=0.0243285y9=0.0216957y10=0.0195761Pressanykeytocontinue無(wú)效算法運(yùn)用遞推公式y(tǒng)(n-1)=-4y(n)+1/n,又由廣義積分中值定理可得y(n)=1/((4n+1)(4ζ+1)),ζ∈(0,1),則1/(5(n+1))<y(n)<1/(n+1),因此可取y(n)≈[1/(5(n+1))+1/(n+1)]/2=3/(5(n+1)程序?yàn)椋?include<iostream.h>#include<math.h>voidmain(){floaty9,y10; y10=3.0/55.0; cout<<"y10="<<y10<<""; for(intn=9;n>=0;n--) { y9=-4.0*y10+1.0/(n+1); cout<<"y"<<n<<"="<<y9<<"";y10=y9; if(n%3==0)cout<<endl; }}其成果為：y10=0.0545455y9=-0.118182y8=0.583838y7=-2.21035y6=8.98427y5=-35.7704y4=143.282y3=-572.877y2=2291.84y1=-9166.86y0=36668.4心得體會(huì)由有效算法與無(wú)效算法的成果可以懂得：有效算法的誤差的傳遞是逐漸縮小的，而無(wú)效算法的誤差的傳遞是逐漸擴(kuò)大的。因此有效算法是數(shù)值穩(wěn)定的算法。二．用牛頓法求下列方程lnx+x-2=0的根。1．算法給定初值X0，ε為根的容許誤差，η為|F（X）|的容許誤差，N為迭代次數(shù)的容許值。(1).假如F（X0）的導(dǎo)數(shù)等于零或者迭代次數(shù)不小于N，則算法失敗，結(jié)束；否則執(zhí)行2。(2).計(jì)算X1=X0-F（X0）/F’（X0）。(3).若|X1-X0|<ε,或者|F（X1）|<η,則輸出X1，程序結(jié)束，否則執(zhí)行4。(4).令X0=X1，轉(zhuǎn)向1。2.程序#include<iostream.h>#include<math.h>#defineEPS1e-8#defineETA1e-8#defineN100floatf(floatx){ returnlog(x)+x-2;}floatf1(floatx){ return1+1/x;}voidmain(){ floatx1,x0=1.5,d; cout<<"x0="<<x0<<""; for(inti=0;i<=N;) { if(f1(x0)==0)break; elsex1=x0-f(x0)/f1(x0); d=(x1-x0); if(fabs(f(x1))<ETA&&fabs(d)<EPS) { cout<<x1<<endl; break; } else x0=x1; i++; cout<<"x"<<i<<"="<<x0<<""; while(i%3==0)cout<<endl; } cout<<"Therootoftheequationisx="<<x0<<endl;}3.成果x0=1.5x1=1.55672x2=1.557151.55715Therootoftheequationisx=1.557154．心得體會(huì)通過(guò)牛頓迭代算法的編程實(shí)習(xí)，我理解到牛頓法把非線性方程線性化，加緊了收斂速度。編寫一種用牛頓前插公式計(jì)算函數(shù)值的程序，規(guī)定先輸出差分表，再計(jì)算x點(diǎn)的函數(shù)值，并應(yīng)用于下面的問(wèn)題：xi2021222324yi1.301031.322221.342421.361731.38021求x=21.4時(shí)的插值多項(xiàng)式的值。算法(1).輸入n=4,xi,yi(i=0,1,2,3,4).(2).計(jì)算各階差分f00,f10,f20,f30,f40。(3).計(jì)算函數(shù)值：f00+f10*t+f20*t*(t-1)/2+f30*t*(t-1)*(t-2)/6+f40*t*(t-1)*(t-2)*(t-3)/24;程序#include<iostream.h>floatNewTon(floatf00,floatf10,floatf20,floatf30,floatf40,floatt){ returnf00+f10*t+f20*t*(t-1)/2+f30*t*(t-1)*(t-2)/6+f40*t*(t-1)*(t-2)*(t-3)/24;}voidmain(){ floatf1[4],f2[3],f3[2],f4;floatf0[5]={1.30103,1.32222,1.34242,1.36173,1.38021};floatx[5]={20,21,22,23,24}; cout<<""<<"差分表"<<endl; for(inti=0;i<=4;i++) { cout<<"x"<<i<<"="<<x[i]; } cout<<endl<<"";for(i=0;i<=4;i++) { cout<<""<<f0[i]; } cout<<endl<<"1階差分"<<""; for(i=0;i<=3;i++) { f1[i]=f0[i+1]-f0[i]; cout<<f1[i]<<""; } cout<<endl<<"2階差分"<<"";for(i=0;i<=2;i++) { f2[i]=f1[i+1]-f1[i]; cout<<f2[i]<<""; } cout<<endl<<"3階差分"<<"";for(i=0;i<=1;i++) { f3[i]=f2[i+1]-f2[i]; cout<<f3[i]<<""; } cout<<endl<<"4階差分"<<""; f4=f3[1]-f3[0]; cout<<f4<<endl; floatt=21.4-x[0]; floatN4=NewTon(f0[0],f0[1],f0[2],f0[3],f0[4],t); cout<<"N(21.4)="<<N4<<endl;}成果差分表x0=20x1=21x2=22x3=23x4=241.301031.322221.342421.361731.380211階差分0.02118990.02020.019310.01848012階差分-0.-0.-0.3階差分9.98974e-0056.00815e-0054階差分-3.98159e-005N(21.4)=3.48268心得體會(huì)由于此題目選用的節(jié)點(diǎn)是等距的，因此可用等距節(jié)點(diǎn)的插值公式替代節(jié)點(diǎn)的任意分布的牛頓插值公式，這樣就可以簡(jiǎn)化插值公式，同步可以防止作除法運(yùn)算，并且程序算法的設(shè)計(jì)也相對(duì)比較簡(jiǎn)樸。用改善歐拉措施（取h=0.1）求解下列微分方程初值問(wèn)題：y’=(x^2+y^2),x∈[0,1.0]y(0)=0算法解一階常微分方程初值問(wèn)題y’=f(x,y),y(x0)=y0,x∈[a,b],將區(qū)間[a,b]做n等分，取步長(zhǎng)h=(b-a)/n.歐拉公式為：y(i+1)=y(i)+hf(x(i),y(i))梯形公式為：y(i+1)=y(i)+h[f(x(i),y(i))+f(x(i+1),y(i+1))]/2改善歐拉法，采用公式：?(i+1)=y(i)+hf(x(i),y(i))y(i+1)=y(i)+h[f(x(i),y(i))+f(x(i+1),?(i+1))]/2程序#include<iostream.h>#defineN10floatf(floatx,floaty){ returnx*x+y*y;}voidmain(){ floatx0=0,y0=0,a=0,b=1.0; floath=(b-a)/N; cout<<"x[0]="<<x0<<""<<"y[0]="<<y0<<endl; floatx1=x0,y1=y0; for(inti=0;i<=N;) { y1=y0+h*f(x0,y0);/*歐拉公式*/ x1=x0+h; y1=y0+h*(f(x0,y0)+f(x1,y1))/2;/*梯形公式*/ y0=y1; x0=x1; i++; cout<<"x["<<i<<"]="<<x1<<""<<"y["<<i<<"]="<<y1<<endl; }}成果x[0]=0y[0]=0x[1]=0.1y[1]=0.0005x[2]=0.2y[2]=0.00300013x[3]=0.3y[3]=0.00950303x[4]=0.4y[4]=0.0220247x[5]=0.5y[5]=0.0426214x[6]=0.6y[6]=0.0734421x[7]=0.7y[7]=0.116817x[8]=0.8y[8]=0.175396x[9]=0.9y[9]=0.252374x[10]=1y[10]=0.35183x[11]=1.1y[11]=0.479294心得體會(huì)改善歐拉公式和歐拉公式相比，它們同為單步顯式公式，不過(guò)前者的局部截?cái)嗾`差比后者的局部截?cái)嗾`差高一階。歐拉公式具有一階精度，梯形公式和改善歐拉公式具有二階精度。五．用阿當(dāng)姆斯措施（取h=0.1）求解下列微分方程初值問(wèn)題（用四階龍格-庫(kù)塔公式提供出發(fā)值）：y’=0.1(x^3+y^2),x∈[0,1]y(0)=1算法概要阿當(dāng)姆斯措施是一種線性多步法，其四階顯式公式為y(i+1)=y(i)+h(55f(i)-59f(i-1)+37f(i-2)-9f(i-3))/24四階隱式公式為y(i+1)=y(i)+h(9f(i+1)+19f(i)-5f(i-1)+f(i-2))/24以上兩個(gè)公式均具有四階精度，將它們構(gòu)成阿當(dāng)姆斯預(yù)測(cè)-校正系統(tǒng)：?(i+1)=y(i)+h(55f(i)-59f(i-1)+37f(i-2)-9f(i-3))/24f’(i+1)=f(x(i+1),?(i+1))y(i+1)=y(i)+h(9f’(i+1)+19f(i)-5f(i-1)+f(i-2))/24f(i+1)=f(x(i+1),y(i+1))這是一種四步措施，計(jì)算y(i+1)時(shí)要用到y(tǒng)(i),y(i-1)，y(i-2),y(i-3),因此，它不是自開始的，一般借助于同階的龍格-庫(kù)塔公式為其提供出發(fā)值：y1,y2,y3.程序#include<iostream.h>#include<conio.h>#defineN10floatf(floatx,floaty){ return0.1*(x*x*x+y*y*y);}voidmain(){ floatyy[3],y,x,h,a,b,k1,k2,k3,k4,f1; a=0;b=1.0;y=1;x=0;yy[0]=1; h=(b-a)/N;inti; cout<<"x[0]="<<x<<""<<"y[0]="<<y<<endl; for(i=1;i<=3;i++) {k1=f(x,y); k2=f(x+h/2,y+h*k1/2); k3=f(x+h/2,y+h*k2/2); k4=f(x+h,y+h*k3); y=y+h*(k1+2*k2+2*k3+k4)/6; x=x+h;yy[i]=y; cout<<"x["<<i<<"]="<<x<<""<<"y["<<i<<"]="<<y<<endl; } for(i=3;i<=N;i++) { x=a+i*h; y=yy[3];yy[3]=yy[3]+h*(55*f(x,yy[3])-59*f(x-h,yy[2])+37*f(x-2*h,yy[1])-9*f(x-3*h,yy[0]))/24; f1=f(x+h,yy[3]); yy[3]=y+h*(9*f1+19*f(x,yy[3])-5*f(x-h,yy[2])+f(x-2*h,yy[1]))/24; yy[0]=yy[1]