數(shù)學(xué)建模MATLAB教案

1、數(shù)學(xué)建模 MATLAB 教案1 / 61三種插值方法 拉格朗日多項(xiàng)式插值 構(gòu)造基函數(shù)l(X)(x Xo)L (x XiJ(x XiJL (x Xn)i(X Xo)L (XiXi i)(XiXiJL (XiXn)插值多項(xiàng)式nLn(x)yi(x)i 0分段線性插值將每?jī)蓚€(gè)相鄰的節(jié)點(diǎn)用直線連起來，即在每個(gè)小區(qū)間上是線性函數(shù)。有現(xiàn)成命令。 三次樣條插值一根有彈性的細(xì)長(zhǎng)木條固定在節(jié)點(diǎn)上，其他地方自然彎曲，如此稱為樣條曲線。 普遍使用的樣條函數(shù)是分段三次多項(xiàng)式：在每個(gè)小區(qū)間上是三次多項(xiàng)式，在大區(qū)間上二階導(dǎo)數(shù)連續(xù)，通過全部節(jié)點(diǎn)。有現(xiàn)成命令。例子1對(duì) y 亠，5x5，用 11 個(gè)等分節(jié)點(diǎn)作上述三種插值，用 2

2、1 個(gè)等分插值點(diǎn)1 x作圖。比較結(jié)果，spli ne 插值最好。2.數(shù)據(jù)擬合2.1 多項(xiàng)式擬合指令方法x=1 2 3 4 5 6 7 8 9;y=9 7 6 3 -1 2 5 7 20;P=polyfit(x,y,3);xi=0:.2:10;yi=polyval(P,xi);plot(xi,yi, x,y,r*);圖形窗口方法x=1 2 3 4 5 6 7 8 9; y=9 7 6 3 -1 2 5 7 20; plot(x,y,r*);2.2 指定函數(shù)類型擬合例子某次阻尼振蕩實(shí)驗(yàn)中測(cè)得數(shù)據(jù)點(diǎn)x=0;0.4;1.2;2;2.8;3.6;4.4;5.2;6;7.2;8;9.2;10.4;11.6

3、;12.4;13.6;14.4;15;y=1;0.85;0.29;-0.27;-0.53;-0.4;-0.12;0.17;0.28;0.15;-0.03;-0.15;-0.071;0.059;0.08;0.032;-0.015;-0.02;f=fittype(a*cos(k*t)*exp(w*t),i ndepe nden t,t,coefficie nts,a,k,w);cfun=fit(x,y,f)xi=0:.1:20;yi=cfu n( xi);plot(x,y,r*,xi,yi,b-);補(bǔ)充：三維圖形例子數(shù)學(xué)建模 MATLAB 教案2 / 6作曲面z2 2f (x, y)的圖形，zsi

4、y,x2, 7.5 x 7.5, 7.5 y 7.5Jx2y2作螺旋線x sint,y cost,z t練習(xí)矩陣運(yùn)算：+加法；-減法；轉(zhuǎn)置；*乘法；八乘幕；左除；/右除。.*乘法；A乘幕；左除；./右除。矩陣函數(shù)：zeros；ones eye；rand; randn。圖形：線型：-實(shí)線；：點(diǎn)線；-.虛點(diǎn)線；-波折線；圓點(diǎn)；+加號(hào)；*星號(hào)；x x 形; 小圓。顏色：y 黃；r 紅；g 綠；b 藍(lán)；w 白；k 黑；m 紫；c 青。例子車燈光源投影區(qū)域的繪制(CUMCM2002A)3.數(shù)值積分矩形公式Lnn 1hfk,hk 0fk,h梯形公式(相當(dāng)于將兩矩陣公式平均，也相當(dāng)于用分段線性插值函數(shù)作為近

5、似)n 1hTnhfk-(fo+fn),hk 12辛普森公式(拋物線法)m 1m 1( f0+f2m4f2k 12 f2k), h3k 0k 1b a2m例子02si nxdx例子衛(wèi)星軌道長(zhǎng)度人造地球衛(wèi)星軌道可視為平面上的橢圓。我國(guó)第一顆人造地球衛(wèi)星近地點(diǎn)距地球 表面 439km,遠(yuǎn)地點(diǎn)距地球表面 2384km,地球半徑為 6371km,求該衛(wèi)星的軌道 長(zhǎng)度。4.數(shù)值微分f(a)f(a)f(a)f(af(a h)h) f f前差公式 hf ff(af(a h)h)后差公式 hf(af(a h)h) f(af(a h)h)中心差商2h數(shù)學(xué)建模 MATLAB 教案3 / 6數(shù)學(xué)建模 MATLAB

6、教案4 / 6前差形式數(shù)值微分，現(xiàn)成命令：diff(x)。輸入 x 是 n 維數(shù)組，輸出為 n-1 維數(shù)組5常微分方程數(shù)值解設(shè)yf (X, y), y(Xo)乂，其中 f 適當(dāng)光滑，對(duì) y 滿足 Lipschitz 條件，以保證解存在且唯一。在一系列離散點(diǎn)x0 x-ix2L xnL上求y(xn)的近似值yn(n 1,2,L )，通常向前歐拉公式在小區(qū)間Xn,Xm上用差商y(Xy(Xnl) ) y(Xy(Xn) )代替方程左端的導(dǎo)數(shù)，方程右端f(X,y)h中的X取小區(qū)間召，召i的左端點(diǎn)Xn，可得y(Xn i) y(Xn) hf(Xn,y(Xn)從Xo出發(fā)，由初值y(Xo) yo，得到y(tǒng)(xj的近

7、似值y-為(以下用 代替)yiyohf (Xo, yo)再以y-作為y(xi)的近似值，得到y(tǒng)(X2)的近似值y為y2yihf (Xi, yi)繼續(xù)下去，一般地yn iynhf (Xn,yn),n0,i,L龍格庫(kù)塔方法按照微分中值定理有yxyx ww y(xnh),o ih注意到方程y f (x, y)就有y(Xn i) y(Xn) hf (Xnh, y(Xnh)記K f (Xnh, y(Xnh)，稱為區(qū)間Xn,Xni上的平均斜率。向前歐拉公式簡(jiǎn)單地取f (xn, yn)為 K，精度很低。龍格庫(kù)塔方法的基本思想：在區(qū)間Xn,Xni內(nèi)多取幾個(gè)點(diǎn)，將它們的斜率加權(quán)平 均作為 K 取 2 個(gè)點(diǎn)，2

8、階龍格庫(kù)塔公式。取 4 個(gè)點(diǎn)，4 階龍格庫(kù)塔公式。有現(xiàn)成的命令：X2Xi,X3X2,L ,xnXn 1。取等步長(zhǎng) h，即Xnx0nh(n 1,2,L )。數(shù)學(xué)建模 MATLAB 教案5 / 6t,x=ode23(f,ts,x0,optio ns)t,x=ode45(f,ts,x0,optio ns)其中，f 是由待解方程寫成的 m 文件名；ts 為自變量的取值；x0 為函數(shù)的初值； options 用于設(shè)定誤差限(可以缺省，缺省時(shí)設(shè)定為相對(duì)誤差103，絕對(duì)誤差106)， 設(shè)定命令：options=odeset(reltol,rt,abstol,at) rt， at 分別為設(shè)定的相對(duì)誤差和絕 對(duì)

9、誤差。t，x 為輸出的自變量和函數(shù)值。例子單擺運(yùn)動(dòng)方程ml mg sin，初始條件(0)0,牌0) 0在0不大的條件下，可將方程中的 sin 近似為，于是得到線性常系數(shù)微分方程 皺g0，容易算出其在初始條件下的解為(t)0cos t,、g。當(dāng)較大時(shí)，若仍用 近似 sin ，誤差太大了。試用數(shù)值方法在等于10和30o兩種情況下求解(設(shè) l 25cm)，畫出(t)的圖形，并與近似解的結(jié)果比較。先將它化為方程組。令Xi,X2&，則方程化為& x2,x2gsin，初始條件為x1(0)為0,x2(0)0，其中g(shù) 9.8， l 25， 細(xì) 為10o0.1745弧度和30o0.5236弧度兩

10、種情況。對(duì)于近似解，周期 T 2 ,學(xué) 2 ,10(s)可以看出，初始角度為10o時(shí)精確做值)解與近似解相差不大，而初始角度為30o時(shí)，隨著時(shí)間的增加二者差別就很大了。食餌-捕食者模型 食餌甲和捕食者乙在時(shí)刻t的數(shù)量分別記作x(t)和y(t)當(dāng)甲獨(dú)立生存時(shí)它的(相對(duì))增長(zhǎng)率為r，即& rx，而乙的存在使甲的增長(zhǎng)率減小， 設(shè)減小的程度與乙的數(shù)量成正比，于是x(t)滿足2 x(r ay) rx axy，比例系數(shù)a反映捕食者掠取食餌的能力。捕食者離開食餌無法生存，設(shè)乙獨(dú)自存在時(shí)死亡率為d，即& dy，而甲的存在使乙的死亡率降低，設(shè)這個(gè)作用與甲的數(shù)量成正比，于是y(t)滿足&

11、y( d bx) dy bxy，比例系數(shù) b 反映食餌對(duì)捕食者的供養(yǎng)能力。設(shè)初始數(shù)量分別為x(0) Xo, y(0)yo。數(shù)學(xué)建模 MATLAB 教案6 / 6設(shè)r 1,d0.5,a 0.1,b 0.02, Xo25, y。2，求數(shù)值解，畫出函數(shù)x(t)和y(t)的圖形以及相圖(x,y)6.線性規(guī)劃LINPROGLin ear program ming.X=LINPROG(f,A,b) solves the linear programming problem:min f*xxsubject to: A*x = b例子min z2為3x2x34x22x38st3x12x26人叢， 0c=2;3

12、;1;a=1,4,2；3,2,0;b=8;6;x,y=li nprog(c,-a,-b,zeros(3,1)c=2;3;-5;a=-2,5,-1;b=-10;aeq=1,1,1;beq=7;x=li nprog(-c,a,b,aeq,beq,zeros(3,1)value=c*x7.優(yōu)化工具箱help optim 查看優(yōu)化工具箱的幫助信息基本程序名：fmi n, fmi nu,fmi ns,leastsq,curvefit,co nstr,lp,qp,mi ni max, nnl s,co nls,fzero,fsolve等。例子2 2求解min?a(本題精確結(jié)果為 x=y=0)用下面一組數(shù)據(jù)擬合c(t) rekt中的系數(shù) r 和 k。(也可以用指定函數(shù)類型擬合)max z2%3x25x3% x?X37st 25x2X310 xn%2壓0數(shù)學(xué)建模 MATLAB 教案7 / 6t=0.25 0.5 1 1.5 2 3 4 6 8;c=19.21 18.15 15.36 14