專題3-4 插值與擬合+數(shù)值微分積分

1、專題3插值與擬合說明：專題 3作業(yè)插值與擬合-數(shù)學(xué)實驗P81-實驗內(nèi)容 1（d）,3,5,6 1、p81-1-d（數(shù)值、畫圖都要體現(xiàn)出來）第一題d函數(shù)lagr文件function y=lagr1(x0,y0,x)n=length(x0);m=length(x);for i=1:m z=x(i); s=0; for k=1:n p=1; for j=1:n; if j=k p=p*(z-x0(j)/(x0(k)-x0(j); end end s=p*y0(k)+s; end y(i)=s;end具體程序：x0=-2:2;%取n=5y0= exp(-x0.2);x=-2:0.1:2;%取m=41y

2、=exp(-x.2);y1=lagr1(x0,y0,x);y2=interp1(x0,y0,x);y3=spline(x0,y0,x);for k=1:11 xx(k)=x(19+2*k);%x>=0且間隔0.2的插值點 yy(k)=y(19+2*k);%x>=0且間隔0.2的函數(shù)值 yy1(k)=y1(19+2*k);%x>=0且間隔0.2的拉格朗日插值 yy2(k)=y2(19+2*k);%x>=0且間隔0.2的分段線性插值 yy3(k)=y3(19+2*k);%x>=0且間隔0.2的三次樣條插值endxx;yy;yy1;yy2;yy3'z=0*x;p

3、lot(x,z,x,y,'b-', x,y1,'r')pauseplot(x,z,x,y,'k-', x,y2,'r')pauseplot(x,z,x,y,'k-',x,y3,'r')y=exp(-x2),-2<=x<=2,三種插值結(jié)果的計較ans = 0 1.0000 1.0000 1.0000 1.0000 0.2000 0.9608 0.9698 0.8736 0.9623 0.4000 0.8521 0.8815 0.7472 0.8617 0.6000 0.6977 0.742

4、7 0.6207 0.7168 0.8000 0.5273 0.5657 0.4943 0.5459 1.0000 0.3679 0.3679 0.3679 0.3679 1.2000 0.2369 0.1714 0.2980 0.2011 1.4000 0.1409 0.0036 0.2281 0.0642 1.6000 0.0773 -0.1035 0.1581 -0.0243 1.8000 0.0392 -0.1126 0.0882 -0.04572.0000 0.0183 0.0183 0.0183 0.0183一、題目：2000年世界人口預(yù)測60年代末世界人口增長情況：年： 1960

5、 1961 1962 1963 1964 1965 1966 1967 1968 人口：29.72 30.61 31.51 32.13 32.34 32.85 33.56 34.20 34.83記人口數(shù)為N(t)，采用指數(shù)函數(shù)N=exp(a+bt)對數(shù)據(jù)進行擬合，畫出擬合曲線及散點圖。(1)使用Matlab中的polyfit； (2) 直接使用最小二乘法求解（超定）方程組。2、 p81-3題程序：clc,clearx=1:0.5:10; y=x.3-6*x.2+5*x-3; y0=y+rand(1); figure(1)f1=polyfit(x,y0,3) %輸出多項式系數(shù)y1=polyval

6、(f1,x); %計算各x點的擬合值plot(x,y,'+',x,y1) grid on title('三次擬合曲線'); c1=1 -6 5 -3'%原三次多項式系數(shù)xishubijiao3=c1,f1',c1-f1'%三次擬合系數(shù)比較wucha3=y',y0',y1',y1'-y',y1'-y0'%三次擬合誤差figure(2); f2=polyfit(x,y0,2) %2次多項式擬合y2=polyval(f2,x); plot(x,y,'+',x,y2); gr

7、id on title('二次擬合曲線'); wucha2=y',y0',y2',y2'-y',y2'-y0'%二次擬合誤差figure(3); f3=polyfit(x,y0,4) %4次多項式擬合y3=polyval(f3,x); plot(x,y,'+',x,y3) grid on title('四次擬合曲線');wucha4=y',y0',y3',y3'-y',y3'-y0'%四次擬合誤差bjiao3_4=(wucha3(:,3

8、)'*(wucha3(:,3)-(wucha4(:,3)'*(wucha4(:,3)%比較三次四次的誤差平方和（擬合數(shù)據(jù)和原來沒有干擾的數(shù)據(jù)進行比較）結(jié)果：f1 = Columns 1 through 3 0.999999999999998 -5.999999999999961 4.999999999999779 Column 4 -2.450276391708465xishubijiao3 = 1.000000000000000 0.999999999999998 0.000000000000002 -6.000000000000000 -5.999999999999961

9、-0.000000000000039 5.000000000000000 4.999999999999779 0.000000000000221 -3.000000000000000 -2.450276391708465 -0.549723608291535wucha3 = 1.0e+002 * Columns 1 through 3 -0.030000000000000 -0.024502763917089 -0.024502763917087 -0.056250000000000 -0.050752763917089 -0.050752763917087 -0.09000000000000

10、0 -0.084502763917089 -0.084502763917088 -0.123750000000000 -0.118252763917089 -0.118252763917088 -0.150000000000000 -0.144502763917089 -0.144502763917088 -0.161250000000000 -0.155752763917089 -0.155752763917089 -0.150000000000000 -0.144502763917089 -0.144502763917089 -0.108750000000000 -0.1032527639

11、17089 -0.103252763917089 -0.030000000000000 -0.024502763917089 -0.024502763917089 0.093750000000000 0.099247236082911 0.099247236082911 0.270000000000000 0.275497236082911 0.275497236082911 0.506250000000000 0.511747236082911 0.511747236082911 0.810000000000000 0.815497236082911 0.815497236082911 1.

12、188750000000000 1.194247236082912 1.194247236082911 1.650000000000000 1.655497236082911 1.655497236082911 2.201250000000000 2.206747236082912 2.206747236082911 2.850000000000000 2.855497236082911 2.855497236082909 3.603750000000000 3.609247236082911 3.609247236082910 4.470000000000000 4.475497236082

13、911 4.475497236082909 Columns 4 through 5 0.005497236082913 0.000000000000002 0.005497236082913 0.000000000000001 0.005497236082912 0.000000000000001 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082911 0 0.005497236082911 -0.000000000000000 0.005497236082911 -0.

14、000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 -0.000000000000001 0.005497236082911 -0.000000000000001 0.005497236082911 -0.0000

15、00000000001 0.005497236082909 -0.000000000000002 0.005497236082911 -0.000000000000001 0.005497236082909 -0.000000000000002f2 = 10.499999999999993 -72.299999999999940 89.949723608291023wucha2 = 1.0e+002 * Columns 1 through 3 -0.030000000000000 -0.024502763917089 0.281497236082911 -0.056250000000000 -

16、0.050752763917089 0.051247236082911 -0.090000000000000 -0.084502763917089 -0.126502763917089 -0.123750000000000 -0.118252763917089 -0.251752763917089 -0.150000000000000 -0.144502763917089 -0.324502763917089 -0.161250000000000 -0.155752763917089 -0.344752763917089 -0.150000000000000 -0.14450276391708

17、9 -0.312502763917089 -0.108750000000000 -0.103252763917089 -0.227752763917088 -0.030000000000000 -0.024502763917089 -0.090502763917089 0.093750000000000 0.099247236082911 0.099247236082911 0.270000000000000 0.275497236082911 0.341497236082911 0.506250000000000 0.511747236082911 0.636247236082911 0.8

18、10000000000000 0.815497236082911 0.983497236082910 1.188750000000000 1.194247236082912 1.383247236082910 1.650000000000000 1.655497236082911 1.835497236082911 2.201250000000000 2.206747236082912 2.340247236082910 2.850000000000000 2.855497236082911 2.897497236082911 3.603750000000000 3.6092472360829

19、11 3.507247236082909 4.470000000000000 4.475497236082911 4.169497236082910 Columns 4 through 5 0.311497236082911 0.305999999999999 0.107497236082911 0.102000000000000 -0.036502763917089 -0.042000000000000 -0.128002763917089 -0.133500000000000 -0.174502763917089 -0.180000000000000 -0.183502763917089

20、-0.189000000000000 -0.162502763917089 -0.168000000000000 -0.119002763917088 -0.124500000000000 -0.060502763917089 -0.066000000000000 0.005497236082911 -0.000000000000000 0.071497236082911 0.066000000000000 0.129997236082911 0.124500000000000 0.173497236082910 0.167999999999999 0.194497236082910 0.18

21、8999999999999 0.185497236082911 0.179999999999999 0.138997236082910 0.133499999999999 0.047497236082910 0.041999999999999 -0.096502763917091 -0.102000000000002 -0.300502763917091 -0.306000000000002f3 = Columns 1 through 3 0.000000000000000 0.999999999999998 -5.999999999999990 Columns 4 through 5 4.9

22、99999999999997 -2.450276391708882wucha4 = 1.0e+002 * Columns 1 through 3 -0.030000000000000 -0.024502763917089 -0.024502763917089 -0.056250000000000 -0.050752763917089 -0.050752763917089 -0.090000000000000 -0.084502763917089 -0.084502763917089 -0.123750000000000 -0.118252763917089 -0.118252763917089

23、 -0.150000000000000 -0.144502763917089 -0.144502763917088 -0.161250000000000 -0.155752763917089 -0.155752763917088 -0.150000000000000 -0.144502763917089 -0.144502763917088 -0.108750000000000 -0.103252763917089 -0.103252763917088 -0.030000000000000 -0.024502763917089 -0.024502763917088 0.093750000000

24、000 0.099247236082911 0.099247236082912 0.270000000000000 0.275497236082911 0.275497236082912 0.506250000000000 0.511747236082911 0.511747236082912 0.810000000000000 0.815497236082911 0.815497236082911 1.188750000000000 1.194247236082912 1.194247236082912 1.650000000000000 1.655497236082911 1.655497

25、236082911 2.201250000000000 2.206747236082912 2.206747236082912 2.850000000000000 2.855497236082911 2.855497236082912 3.603750000000000 3.609247236082911 3.609247236082912 4.470000000000000 4.475497236082911 4.475497236082912 Columns 4 through 5 0.005497236082911 -0.000000000000000 0.005497236082911

26、 -0.000000000000000 0.005497236082911 -0.000000000000000 0.005497236082911 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082912 0.0000000

27、00000000 0.005497236082912 0.000000000000000 0.005497236082912 0.000000000000000 0.005497236082911 0 0.005497236082911 0 0.005497236082911 0 0.005497236082911 0 0.005497236082912 0.000000000000001 0.005497236082912 0.000000000000001 0.005497236082912 0.000000000000001bjiao3_4 = -5.238689482212067e-0

28、10圖像：分析：在本題中利用了二次擬合，三次擬合，四次擬合分別擬合出函數(shù)y=x.3-6*x.2+5*x-3的多項式，并畫出相應(yīng)的擬合曲線。運行之后，對二次擬合曲線，三次擬合曲線，四次擬合曲線的圖像進行比較分析，發(fā)現(xiàn)：相比較而言，三次擬合的吻合性更高，即精度更好。二次效果最差。3、p81 第五題：（要求程序、計算結(jié)果和圖形）n=0000 1153 2045 2800 3466 4068 4621 5135 5619 6152;t=0 20 40 60 80 100 120 140 160 183.5;R=n'.2,n'aa=Rt'y=aa(1)*n'.2+aa(2

29、)*n'plot(n,t,'b+',n,y,'r')xlabel('n'),ylabel('T')t=0 20 40 60 80 100 120 140 160 184;n=0 1141 2019 2760 3413 4004 4545 5051 5525 6061;A=(n.2)' n'X= At' 注意：（1）、如果直接用polyfit做多項式擬合，會得出一般二次多項式；（2）、另外，可以做變換，兩邊同時除以n，然后用一次函數(shù)擬合。3、p81-6題（本答案在閭宇提供的基礎(chǔ)上修改）分析：本題已經(jīng)建

30、立好了模型！只需對其進行求解： 由題目所給的方程 （已知）可見，v(t)與（T）成指數(shù)變化關(guān)系，所以在通過曲線擬合的時候，使用指數(shù)曲線。（注：這是一種非線性曲線擬合，首先要進行變量代換，為方便計算做以下變量代換） 用代替，是擬合后的曲線方程 對 變形后取對數(shù)，有 。 令 則： 程序：t=0.5 1 2 3 4 5 7 9 ;v1=6.36 6.48 7.26 8.22 8.66 8.99 9.43 9.63;y=log(10-v1);a=polyfit(t,y,1);a1=a(1);a2=a(2);T=-1/a1v0=10-exp(a2)v2=10-(10-v0)*exp(-t/T);擬合后的

31、函數(shù)值plot(t,v1,'k+',t,v2,'r')結(jié)果：t0 = 3.5269v0 =5.6221專題4 數(shù)值積分微分專題4作業(yè)：數(shù)值積分與數(shù)值微分-數(shù)學(xué)實驗P97-實驗內(nèi)容：2(d),3,9+補充完成p96人口增長率二、數(shù)值積分與微分部分0（未布置）97-1題利用三種公式計算積分，另外利用三點公式計算每點處的導(dǎo)數(shù)，數(shù)值計算結(jié)果都與精確值進行比較。(數(shù)值比較需列表、另外還要求畫圖比較)程序：syms xy=(x+sin(x/3);z=int(y,0.3,1.5)z=zxx=0.3:0.2:1.5;y=xx+sin(xx./3);z1=sum(y(1:6)*0

32、.2z2=sum(y(2:7)*0.2z3=quad('x+sin(x./3)',0.3,1.5)z4=trapz(xx,y)u1=z-z1,u2=z-z2,u3=z-z3,u4=z-z4h=0.2;x=0.3:h:1.5;y=0.3895 0.6598 0.9147 1.1611 1.3971 1.6212 1.8325;format longz1=sum(y(1:6)*h,z2=sum(y(2:7)*h,z3=trapz(y)*h,z4=quad('y',0.3,1.5),z5=quadl('y',0.3,1.5),format short1

33、、 p97-2題選擇一些函數(shù)用梯形、辛普森和隨機模擬三種方法計算積分。改變步長（對梯形公式），改變精度要求（對辛普森公式），改變隨機點數(shù)目（對隨機模擬），進行比較、分析。如下函數(shù)供選擇參考：d.（1組同學(xué)提供答案）程序：編輯函數(shù)M文件：function y=fun(x)y=(1./sqrt(2*pi)*exp(-x.2./2)；程序：clc,clearformat longh=4/50; %步長x=-2:h:2;y=fun(x);z1=trapz(y)*h %梯形公式z2=quad('fun',-2,2) %辛普森公式n=1000; %隨機模擬方法x1=rand(1,n);y1

34、=fun(x1.*2);z3=sum(y1)*4/n當對梯形公式取h=4/50,4/100,4/10000;對辛普森公式分別取精度為 ，對隨機點數(shù)目取n=1000,10000,100000,得到的結(jié)果如下：梯形公式辛普森公式蒙特卡羅方法0.954384567677890.954499438241540.939992110595860.954470941689640.954499736107350.957943385948660.954499733224120.954491037397360.95427317381756分析：對于梯形公式，步長越小，結(jié)果越精確；對于辛普森公式，在一般情況精度下的

35、結(jié)果都很精確，辛普森有很好的優(yōu)越性，但是需要函數(shù)表達式；蒙特卡羅方法，雖然結(jié)果具有隨機性，但是隨著n的增大，其計算結(jié)果越來越接近準確值。2.p98-3題（要求15位小數(shù)顯示）x=0:0.1:1;y=(1-x.2).(0.5);z1=4*trapz(x,y)z2=4*quad('(1-x.2).(0.5)',0,1)z3=4*quad8('(1-x.2).(0.5)',0,1)n=10000000;x3=rand(1,n);yy=(1-x3.2).(0.5);z4=4*sum(yy)/n以下是閭宇同學(xué)1組（1）蒙特卡羅方法創(chuàng)建M文件function y=f(x)m

36、=0;for n=1:x if rand(1)2+rand(1)2<=1 m=m+1; end;end;4*m/x輸入Format long>> f(10000);ans = 3.132400000000000>> f(50000);ans = 3.139920000000000>> f(100000);ans = 3.147320000000000（2）數(shù)值積分法 設(shè) ，將區(qū)間0,1n等分，取 梯形法： 辛普森法： 創(chuàng)建M文件fun.mfunction y=fun(x)y=4./(1+x.2);創(chuàng)建M文件f(x)function y=f(x)for

37、n=1:x-1 a(n)=2*fun(n/x);end;vpa(1/(2*k)*(sum(a)+fun(0)+fun(1)輸入：digits(30) %保留小數(shù)點后30位>> f(100)ans =3.14157598692312900467982217378>> f(500)ans =3.14159198692313035294887413329>> f(10000)ans =3.1415926519231400781961838223. p99-9題編程并最終回答問題。 題目：圖7是一個國家的地圖，為了算出它的國土面積，首先對地圖作如下測量：以由西向東方

38、向為x軸，由南到北方向為y軸，選擇方便的原點，并將從最西邊界點到最東邊界點在X軸上的區(qū)間適當?shù)貏澐譃槿舾啥?，在每個分點的y軸方向測出南邊界點和北邊界點的y坐標y1和y2，這樣就得到了表中的數(shù)據(jù)（單位mm）。  根據(jù)地圖的比例我們知道18mm相當于40km，試由測量數(shù)據(jù)計算該國土的近似面積，與它的精確值41288km2比較。（2012本-代雄，雷雨錦，李穎藝提供）程序1：>> x=7.0 10.5 13.0 17.5 34.0 40.5 44.5 48.0 56.0 61.0 68.5 76.5 80.5 91.0 96.0 101.0 104.0 106.5 111.5

39、118.0 123.5 136.5 142.0 146.0 150.0 157.0 158.0;y1=44 45 47 50 50 38 30 30 34 36 34 41 45 46 43 37 33 28 32 65 55 54 52 50 66 66 68;y2=44 59 70 72 93 100 110 110 110 117 118 116 118 118 121 124 121 121 121 122 116 83 81 82 86 85 68;figure(1)plot(x,y1,'o')hold on plot(x,y2,'rp')hold o

40、nxi=10:10:100yi=interp1(x,y1,xi,'spline')yi=interp1(x,y2,xi,'spline')xx=min(x):0.1:max(x);yy1=interp1(x,y1,xx,'spline');plot(xx,yy1)hold onyy2=interp1(x,y2,xx,'spline');plot(xx,yy2,'r')hold onxi = 10 20 30 40 50 60 70 80 90 100yi = 44.6613 51.4824 53.2436 39.1

41、722 30.9541 35.8817 34.5265 44.6028 46.3324 38.2206yi = 56.3424 72.5187 87.6339 98.9170 108.8743 115.7368 117.3082 117.7779 117.8652 124.2328另外：兩種方法，可以參考。結(jié)果：方法一本題利用蒙特卡羅方法計算國土面積，通過for循環(huán)產(chǎn)生隨機點，if語句限定曲線所圍區(qū)域內(nèi)。 通過for循環(huán)多次實驗，求的平均值。>> x=7.0 10.5 13.0 17.5 34.0 40.5 44.5 48.0 56.0 61.0 68.5 76.5 80

42、.5 91.0 . 96.0 101.0 104.0 106.5 111.5 118.0 123.5 136.5 142.0 146.0 150.0 . 157.0 158.0; y1=44.0 45.0 47.0 50.0 50.0 38.0 30.0 30.0 34.0 36.0 34.0 41.0 45.0 46.0 . 43.0 37.0 33.0 28.0 32.0 65.0 55.0 54.0 52.0 50.0 66.0 66.0 68.0; y2=44.0 59.0 70.0 72.0 93.0 100.0 110.0 110.0 110.0 117.0 118.0 116.0

43、 . 118.0 118.0 121.0 124.0 121.0 121.0 121.0 122.0 116.0 83.0 81.0 . 82.0 86.0 85.0 68.0;L=max(x)-min(x);H=max(y2)-min(y1);s=0;for k=1:10 n=10000;m=0;for i=1:nu(i)=unifrnd(min(x),max(x);v(i)=unifrnd(min(y1),max(y2);f1=interp1(x,y1,u(i),'cubic');f2=interp1(x,y2,u(i),'cubic');if(v(i)&g

44、t;=f1&&v(i)<=f2)m=m+1;endendS(k)=1600*L*H*m/(n*182);s=s+S(k);end s=s/10s = 4.2552e+004方法二：x=7.0 10.5 13.0 17.5 34.0 40.5 44.5 48.0 56.0 61.0 68.5 76.5 80.5 91.0 . 96.0 101.0 104.0 106.5 111.5 118.0 123.5 136.5 142.0 146.0 150.0 . 157.0 158.0; y1=44.0 45.0 47.0 50.0 50.0 38.0 30.0 30.0 34.

45、0 36.0 34.0 41.0 45.0 46.0 . 43.0 37.0 33.0 28.0 32.0 65.0 55.0 54.0 52.0 50.0 66.0 66.0 68.0; y2=44.0 59.0 70.0 72.0 93.0 100.0 110.0 110.0 110.0 117.0 118.0 116.0 . 118.0 118.0 121.0 124.0 121.0 121.0 121.0 122.0 116.0 83.0 81.0 . 82.0 86.0 85.0 68.0;newx=7:0.1:158;newy1=interp1(x,y1,newx,'lin

46、ear');newy2=interp1(x,y2,newx,'linear');area=sum(newy2-newy1)*0.1/182*1600)area = 4.2414e+004結(jié)果分析：對于方法一（插值更光滑），采用隨機投點法結(jié)合分段三次插值模擬上下邊界，求國土面積；方法二（簡單，實用，不光滑），同樣采取隨機投點法，但結(jié)合分段線性插值模擬邊界求解?？梢钥闯龇椒ǘ贸龅慕Y(jié)果要比方法一得出的結(jié)果相對精確值更接近些，說明國土邊界是不規(guī)則，不光滑的。方法一、二都各有所長。補充4.完成p96人口增長率 要求：1）編程計算p96表格顯示的增長率 2）編程計算1980年人口

47、x=1900:10:2000y=76.0 92.0 106.5 123.2 131.7 150.7 179.3 204.0 226.5 251.4 281.4r(1)=(-3*y(1)+4*y(2)-y(3)/20/y(1)r(11)=(y(9)-4*y(10)+3*y(11)/20/y(11)for i=2:10 r(i)=(y(i+1)-y(i-1)/20/y(i)endr年份 1970 1972 1974 1976 1978 1980年增長率（%） 0.87 0.85 0.89 0.91 0.95 1.10為了從1.3表2估計1980年的人口，仍用上面的記號，人口增長率滿足微分方程它在初

48、始條件下的解為 因為在題目中增長率以離散型數(shù)據(jù)給出（表3），所以上式要用數(shù)值積分計算.用梯形公式計算的結(jié)果是，1980年該地區(qū)人口為230.2（百萬）.%p96 renkou(2)clc,cleary(1)=210r=0.87 0.85 0.89 0.91 0.95 1.10*0.01；t=0:2:10；for i=2:6 y(i)=y(1)*exp(2*trapz(r(1:i);endy'或者%p96 renkou(2)clc,cleary(1)=210r=0.87 0.85 0.89 0.91 0.95 1.10*0.01t=0:2:10for i=2:6 y(i)=y(1)*ex

49、p(trapz(t(1:i),r(1:i);endy'第七次 第四章Lingo練習(xí)程序1： max 72x1+64x2st2) x1+x2<503) 12x1+8x2<4804) 3x1<100Endmax 72x1+64x2st2) x1+x2<503) 12x1+8x2<4804) 3x1<100end LP OPTIMUM FOUND AT STEP 2 OBJECTIVE FUNCTION VALUE 1) 3360.000 VARIABLE VALUE REDUCED COST X1 20.000000 0.000000 X2 30.000000 0.000000 ROW SLACK OR SURPLUS DUAL PRICES 2) 0.000000 48.000000 3) 0.000000 2.0