數(shù)學(xué)實驗與Matlab市公開課獲獎?wù)n件省名師優(yōu)質(zhì)課賽課一等獎?wù)n件

數(shù)學(xué)試驗與Matlab1/229

實驗一

矩陣運算與Matlab命令

2/229日常矩陣及其運算矩陣應(yīng)用實例:欖球防護用具生產(chǎn)管理3/229

應(yīng)用問題一個工廠生產(chǎn)三種橄欖球用具：防護帽、墊肩、臀墊。需要不一樣數(shù)量：硬塑料、泡沫塑料尼龍線、勞動力。為監(jiān)控生產(chǎn)，管理者對它們之間關(guān)系十分關(guān)心。為把握這些量關(guān)系，他列出下面表4/229

原料產(chǎn)品關(guān)系表5/229訂單管理者接到四份訂單如上表所表示。問應(yīng)該怎樣計算每份訂單所需原材料，方便組織生產(chǎn)？6/229將表格寫成矩陣形式7/229計算輸入下面Matlab指令A(yù)=[423;132;133;322],B=[35206045;10155040;20124520]C=A*B請自行計算觀看結(jié)果8/229Matlab基本指令向量創(chuàng)建和運算9/2291.直接輸入向量》x1=[124],x2=[1,2,1],x3=x1’運行結(jié)果

x1=124x2=121x3=124

10/2292.冒號創(chuàng)建向量

x1=3.4:6.7,x2=3.4:2:6.7x3=2.6:-0.8:0運算結(jié)果x1=3.40004.40005.40006.4000x2=3.40005.4000x3=2.60001.80001.00000.11/2293.生成線性等分向量指令x=linspace(a,b,n)在[a,b]區(qū)間產(chǎn)生n個等分點(包含端點)x=linspace(0,1,5)結(jié)果x=00.25000.50000.75001.000012/229工作空間在Matlab窗口創(chuàng)建向量后并運行后，向量就存在于工作空間，能夠被調(diào)用。13/229向量運算設(shè)x=[x1x2x3];y=[y1y2y3];為兩個三維向量，a,b為標量。向量數(shù)乘：a*x=[a*x1a*x2a*x3]向量平移：x+b=[x1+bx2+bx3+b]向量和：x+y=[x1+y1x2+y2x3+y3]向量差：x-y=[x1-y1x2-y2x3-y3]數(shù)乘冪：如a^214/229元素群運算(四則運算）x.*y=[x1*y1x2*y2x3*y3](元素群乘積)x./y=[x1/y1x2/y2x3/y3](元素群右除，右邊y做分母)x.\y=[y1/x1y2/x2y3/x3](元素群左除，左邊x做分母)x.^5=[x1^5x2^5x3^5](元素群乘冪)2.^x=[2^x12^x22^x3](元素群乘冪)x.^y=[x1^y1x2^y2x3^y3](元素群乘冪)15/229元素群運算（函數(shù)計算）Matlab有許多內(nèi)部函數(shù)，可直接作用于向量產(chǎn)生一個同維函數(shù)向量。x=linspace(0,4*pi,100);（產(chǎn)生100維向量x）y=sin(x);(y也自動為100維向量)y1=sin(x).^2;y2=exp(-x).*sin(x);觀察結(jié)果16/229創(chuàng)建矩陣（數(shù)值矩陣創(chuàng)建）直接輸入法創(chuàng)建簡單矩陣。

A=[1234;5678;9101112]B=[-1.3,sqrt(3);(1+2)*4/5,sin(5);exp(2),6]

觀察運行結(jié)果17/229創(chuàng)建矩陣（符號矩陣創(chuàng)建）用指令“syms”說明符號變量。

symsa11a12a13a14a21a22a23a24a31a32a33a34b11b12b13b14b21b22b23b24b31b32b33b34A1=[a11a12a13a14;a21a22a23a24;a31a32a33a34],B1=[b11b12b13b14;b21b22b23b24;b31b32b33b34]運行18/229矩陣運算(矩陣加減、數(shù)乘、乘積)C=A1+B1D=A1-B1symsc,cA=c*A1A2=A1(:,1:3),B1G=A2*B119/229矩陣運算(矩陣加減、數(shù)乘、乘積)A,A_trans=A’H=[123;210;123],K=[123;210;231]h_det=det(H),k_det=det(K),H_inv=inv(H),K_inv=K^-120/229矩陣運算(左除和右除)左除“\”：

求矩陣方程AX=B解；（

A、B行要保持一致）

解為X=A\B；當A為方陣且可逆時有X=A\B=inv(A)*B；右除“/”：

求矩陣方程XA=B解

（A、B列要保持一致）

解為X=B/A，當A為方陣且可逆時有X=B/A=B*inv(A)21/229矩陣運算(左除和右除)求矩陣方程：設(shè)A、B滿足關(guān)系式：AB＝2B+A,求B。其中A=[301;110;014]。解：有(A-2I)B＝A程序：A=[301;110;014];B=inv(A-2*eye(3))*A,B=(A-2*eye(3))\A觀察結(jié)果：22/229分塊矩陣（矩陣標識）1.矩陣元素標識：A(i,j)表示矩陣A第i行j列元素；2.向量標識方式A(vr,vc)：vr=[i1,i2,…,ik]、vc=[j1,j2,…,ju]分別是含有矩陣A行號和列號單調(diào)向量。A(vr,vc)是取出矩陣A第i1,i2,…,ik行與j1,j2,…,ju列交叉處元素所組成新矩陣。23/229分塊矩陣（矩陣標識）取出A1、3行和1、3列交叉處元素組成新矩陣A1程序A=[10112;01-123;30510;23121],vr=[1,3];vc=[1,3];A1=A(vr,vc)觀察結(jié)果24/229分塊矩陣（矩陣標識）將A分為四塊，并把它們賦值到矩陣B中，觀察運行后結(jié)果。程序A11=A(1:2,1:2),A12=A(1:2,3:5),A21=A(3:4,1:2),A22=A(3:4,3:5)B=[A11A12;A21A22]結(jié)果25/229分塊矩陣（矩陣修改和提取）修改矩陣A，將它第1行變?yōu)?。程序：A=[10112;01-123;30510;23121]A(1,:)=[00000];A刪除上面矩陣A第1、3行。程序：A([1,3],:)=[]結(jié)果26/229生成特殊矩陣全1陣ones(n),ones(m,n),ones(size(A))全零陣：zeros(n)，zeros(m,n),zeros(size(A))經(jīng)常用于對某個矩陣或向量賦0初值單位陣：eye(n)，eye(m,n)隨機陣：rand(m,n)，rand(n)=rand(n,n)用于隨機模擬，常和rand('seed',k)配合使用。27/229生成特殊矩陣將rand指令運行屢次，觀察結(jié)果。程序：y1=rand(1,5),y2=rand(1,5),rand('seed',3),x1=rand(1,5),rand('seed',3),x2=rand(1,5)結(jié)果28/229慣用矩陣函數(shù)det(A)：方陣行列式；rank(A)：矩陣秩；eig(A)：方陣特征值和特征向量；trace(A)：矩陣跡；rref(A)：初等變換階梯化矩陣Asvd(A)：矩陣奇異值分解。cond(A)：矩陣條件數(shù)；29/229數(shù)據(jù)簡單分析1.當數(shù)據(jù)為行向量或列向量時，函數(shù)對整個向量進行計算.2.當數(shù)據(jù)為矩陣時，命令對列進行計算，即把每一列數(shù)據(jù)當成同一變量不一樣觀察值。max(求最大)、min(求最小)、mean(求平均值)、sum(求和)、std(求標準差)、cumsum(求累積和)、median(求中值)、diff(差分)、sort(升序排列)、sortrows(行升序排列)等等。30/229數(shù)據(jù)簡單分析觀察：生成一個3×6隨機數(shù)矩陣，并將其各列排序、求各列最大值與各列元素之和。程序rand('seed',1);A=rand(3,6),Asort=sort(A),Amax=max(A),Asum=sum(A)結(jié)果31/229試驗二函數(shù)可視化與Matlab作圖

32/229函數(shù)可視化

f(x),g(x)是周期函數(shù)嗎？觀察它們圖象。程序clf,x=linspace(0,8*pi,100);F=inline('sin(x+cos(x+sin(x)))');y1=sin(x+cos(x+sin(x)));y2=0.2*x+sin(x+cos(x+sin(x)));plot(x,y1,'k:',x,y2,'k-')legend('sin(x+cos(x+sin(x))','0.2x+sin(x+cos(x+sin(x)))',2)令33/229繪制平面曲線(plot指令)plot(x,y)：以x為橫坐標、y為縱坐標繪制二維圖形x,y是同維數(shù)向量；plot(y)：相當于x=[1,2,…,length(y)]時情形。34/229繪制平面曲線（繪制多個圖形)1.plot(x,[y1;y2;…]),x是橫坐標向量，[y1;y2;…]是由若干函數(shù)縱坐標拼成矩陣2.plot(x,y1),holdon,plot(x,y2),holdoff3.plot(x,y1,x,y2,…)4.plotyy兩個坐標系，用于繪制不一樣尺度函數(shù)。35/229繪制平面曲線（線型、點形和顏色控制）plot(x,y,‘顏色＋線型＋點形’)plot(x,y,‘顏色＋線型＋點形’,x,y,‘顏色＋線型＋點形’,…)句柄圖形和set命令改變屬性值，可套用：h=plot(x,y),set(h,‘屬性’,屬性值,‘屬性’,屬性值,…)也可用plot(x,y,'屬性','屬性值')設(shè)置圖形對象屬性。36/229繪制平面曲線（屬性變量和屬性值）線寬：LineWidth點大小：MarkerSize線型：LineStyle顏色：color37/229繪制平面曲線（例）觀察：改變繪圖線型和顏色。用gridon指令為圖形窗口加上網(wǎng)格線，并改變網(wǎng)格線型和字體大小。程序h=plot([0:0.1:2*pi],sin([0:0.1:2*pi]));set(h,'LineWidth',5,'color','red');gridonset(gca,'GridLineStyle','-','fontsize',16)

觀察結(jié)果38/229繪制平面曲線（坐標軸控制）axis指令axis([xminxmaxyminymax])：

設(shè)定二維圖形x和y坐標范圍；

axis([xminxmaxyminymaxzminymax])：

設(shè)定三維圖形坐標范圍；其中xmin<x<xmax，ymin<y<ymax，zmin<z<zmax。39/229繪制平面曲線（gca屬性控制）改變當前軸對象句柄gca屬性用set(gca,‘屬性’,屬性值,…)可改變字體大小、坐標刻度等軸對象內(nèi)容。比如：set(gca,'ytick',[-1-0.500.51])將y坐標按向量[-1-0.500.51]將刻度分成4格；set(gca,'yticklabel','a|b|c|d|e')改變y坐標刻度說明。40/229繪制平面曲線（gca屬性控制,例）設(shè)置y坐標刻度并加以說明，并改變字體大小。程序

plot([0:0.1:2*pi],sin([0:0.1:2*pi]),'k.-',);gridon,axis([06.3-1.11.1])，set(gca,'ytick',[-1-0.500.51]),set(gca,'yticklabel','a|b|c|d|e'),set(gca,'fontsize',20)get(gca)運行結(jié)果41/229繪制平面曲線（文字標注）

title(‘圖形標題’)；xlabel（‘x軸名稱’）；ylabel（‘y軸名稱’）；zlabel（‘z軸名稱’）；text(‘說明文字’)：創(chuàng)建說明文字；gtext（'說明文字'）：用鼠標在特定位置輸入文字。文字標注慣用符號：\pi（π）；\alpha（α）；\beta（β）；\leftarrow（左箭頭）\rightarrow（右箭頭）；\bullet（點號）42/229繪制平面曲線（程序講解，exp2_1.m）clf,t=0:0.1:3*pi;alpha=0:0.1:3*pi;plot(t,sin(t),'r-');holdon;plot(alpha,3*exp(-0.5*alpha),'k:');set(gca,'fontsize',15,'fontname','timesNewRoman'),xlabel('\it{t(deg)}');ylabel('\it{magnitude}');title('\it{sinewaveand{\it{Ae}}^{-\alpha{\itt}}wave}');

43/229繪制平面曲線（程序講解，exp2_1.m）text(6,sin(6),'\fontsize{15}TheValue\it{sin(t)}at{\itt}=6\rightarrow\bullet','HorizontalAlignment','right'),text(2,3*exp(-0.5*2),['\fontsize{15}\bullet\leftarrowTheValueof\it{3e}^{-0.5\it{t}}=',num2str(3*exp(-0.5*2)),'at\it{t}=2']);legend('\itsin(t)','{\itAe}^{-\alphat}')注1：num2str：['string1',num2str,'string2']，用方括號注2：legend請結(jié)合圖形觀察此命令使用44/229圖形窗口創(chuàng)建和分割

subplot(m,n,k)命令。在圖形區(qū)域中顯示多個圖形窗口。m為上下分割數(shù)，n為左右分割數(shù)，k為第k子圖編號。例：將一個圖形分為9個子圖，在第k個子圖畫sin(kx)圖象.程序：

clf，b=2*pi;x=linspace(0,b,50);fork=1:9y=sin(k*

x);subplot(3,3,k),plot(x,y),axis([0,2*pi,-1,1])end45/229若干有用指令clf：去除圖形窗口已經(jīng)有內(nèi)容.shg：顯示圖形窗口。clear、clearx：去除工作空間已經(jīng)有變量。figure(n):打開第n個圖形窗口help:…:續(xù)行號46/229繪制二元函數(shù)基本步驟：1.生成二維網(wǎng)格點2.計算函數(shù)在網(wǎng)格點上值3.繪制函數(shù)圖形47/229三維繪圖（meshgrid指令：生成網(wǎng)格點）觀察meshgrid指令效果。程序：a=-0.98;b=0.98;c=-1;d=1;n=10;x=linspace(a,b,n);y=linspace(c,d,n);[X,Y]=meshgrid(x,y);plot(X,Y,'+')觀察結(jié)果48/229三維繪圖（計算函數(shù)值，定義域淘汰）程序：fori=1:nforj=1:nif(1-X(i,j))<eps1|X(i,j)-Y(i,j)<eps1z(i,j)=NaN;elsez(i,j)=1000*sqrt(1-X(i,j))^-1.*log(X(i,j)-Y(i,j));endendend49/229三維繪圖（繪圖指令）mesh(X,Y,z)：在三維空間中繪出由(X,Y,z)表示曲面；meshz(X,Y,z)：除了含有mesh功效外，還畫出上下高度線，meshc(X,Y,z)：除了含有mesh功效外，還在曲面下方畫出函數(shù)z=f(x,y)等值線圖，surf(X,Y,z)：也是三維繪圖指令，與mesh區(qū)分在于mesh繪出彩色線，surf繪出彩色面，運行exp2_1,觀察效果50/229三維繪圖（等值線指令）表現(xiàn)二維函數(shù)圖形另一個方式是繪制等值線圖。contour(X,Y,z,n)：n條等高線，n可缺省；contourf(X,Y,z,n)：等值線間用不一樣顏色填滿，有更加好視覺效果；contour3(X,Y,z,n)：在三維空間畫出等值線圖colorbar：將顏色與函數(shù)值對應(yīng)起來顯示在圖中。51/229三維繪圖（等值線指令，繼續(xù)exp2_2顯示效果）clf,contour(X,Y,z,40),colorbarcontourf(X,Y,z,40),colorbarcontour3(X,Y,z,40),colormap([0,0,0])為等值線標上函數(shù)值:可套用下面程序格式.[cs,h]=contour(X,Y,z,15);clabel(cs,h,'labelspacing',244)labelspace是數(shù)值標識之間相隔寬度，默認值為144,這里取了244，52/229空間曲線和運動方向表現(xiàn)一條空間曲線能夠用矢量函數(shù)表示為它速度矢量表現(xiàn)為曲線切矢量:53/229空間曲線和運動方向表現(xiàn)很顯然飛行曲線方程為：

54/229繪制空間曲線（指令）

plot3(x,y,z)：繪制三維空間曲線，使用方法和plot類似。quiver(X,Y,u,v)：繪制二維矢量，在坐標矩陣點[X,Y]處繪制矢量[u,v],其中u為矢量x坐標，v為矢量y坐標，其維數(shù)大于2。quiver3（X,Y,Z,u,v,w）：繪制三維矢量，使用方法與quiver類似。Gradient：[Fx,Fy,Fz]=gradient(F)為函數(shù)F數(shù)值梯度55/229繪制空間曲線（程序講解exp2_3）exp2_3.mclf,t=linspace(0,1.5,20);x=t.^2;y=(2/3)*t.^3;z=(6/4)*t.^4-(1/3)*t.^3;plot3(x,y,z,'r.-,'linewidth',1,'markersize',10),holdonVx=gradient(x);Vy=gradient(y);Vz=gradient(z);h=quiver3(x,y,z,Vx,Vy,Vz),set(h,'linewidth',1),gridonaxis([01.501.5040])xlabel('x'),ylabel('y'),zlabel('z')，boxon運行程序56/229應(yīng)用、思索和練習(xí)用平行截面法討論由曲面z=x^2-y^2組成馬鞍面形狀。對于二重積分，積分指示線方法是很有用，當然你首先得了解一下什么是積分指示線法，請查閱高等數(shù)學(xué)相關(guān)內(nèi)容，然后設(shè)計一個數(shù)學(xué)試驗，然后用Matlab繪圖工具表現(xiàn)這一方法。

57/229應(yīng)用、思索和練習(xí)

繪制微分方程y/dx=xy,y(0)=0.4斜率場，并將解曲線畫在圖中，觀察斜率場和解曲線關(guān)系。

58/229應(yīng)用、思索和練習(xí)地球表面氣溫差異很大，而且隨時間改變，要繪制氣溫分布圖絕不是一件輕易事情。不過赤道溫度最高，而在兩級最冷，在中間地帶則是過渡帶。所以可粗略將這種氣溫分布情況用圖形表現(xiàn)出來，試繪制地球表面氣溫分布示意圖。

59/229應(yīng)用、思索和練習(xí)60/229應(yīng)用、思索和練習(xí)(若干特殊圖形)x=[1:10];y=[56348110356];subplot(2,3,1),bar(x,y),axis([110111])subplot(2,3,2),hist(y,x),axis([11014])subplot(2,3,3),stem(x,y,'k'),axis([110111])subplot(2,3,4),stairs(x,y,'k'),axis([110111])subplot(2,3,5),x=[130.55];explode=[0001];pie(x,explode)subplot(2,3,6),z=0:0.1:100;x=sin(z);y=cos(z).*10;comet3(x,y,z)

61/229應(yīng)用、思索和練習(xí)(周期函數(shù)推廣）62/229應(yīng)用、思索和練習(xí)(線性P周期函數(shù)）63/229應(yīng)用、思索和練習(xí)(線性P周期函數(shù))64/229應(yīng)用、思索和練習(xí)(線性P周期函數(shù)）p周期函數(shù)是一個特殊線性p周期函數(shù)，為M=0情形。線性p周期函數(shù)是不是一個線性函數(shù)同一個p周期函數(shù)之和呢？怎樣證實這一點？65/229應(yīng)用、思索和練習(xí)（循環(huán)比賽名次）有若干支球隊參加循環(huán)比賽，他們兩兩交鋒，每場比賽只計勝敗，不允許平局，循環(huán)賽結(jié)束后要依據(jù)他們比賽成績排列名次。一個方法是計算得分，得分是每支球隊獲勝場次，依據(jù)各隊得分排知名次，決定冠軍隊。

66/229應(yīng)用、思索和練習(xí)（循環(huán)比賽名次即使計算各隊得分很輕易，但有時按得分排名方法并不一定合理。假定有4支球隊,記做v1~v4。在一次循環(huán)賽中，v1得分為2，v2得分為2，v3得分為1，v4得分為1，能夠把得分寫成4維向量形式：s=[2211]T，在這種情況下應(yīng)該怎樣決定[v1,v2,v3,v4]名次呢？67/229函數(shù)式－直接確定型模型

試驗三68/229從系統(tǒng)分析觀點了解函數(shù)y=f（x）x：自變量，y：因變量，f：映射規(guī)則函數(shù)不是枯燥數(shù)學(xué)符號函數(shù)是系統(tǒng)函數(shù)是數(shù)學(xué)模型是描述自然現(xiàn)象有力工具69/229黑箱模型和經(jīng)驗函數(shù)白箱：映射規(guī)則f已知；灰箱：映射規(guī)則f部分已知；黑箱：映射規(guī)則f未知。對于黑箱模型，只知道輸入輸出數(shù)據(jù)，需依據(jù)這些數(shù)據(jù)近似決定映射規(guī)則f70/229經(jīng)驗函數(shù)(機床加工問題）用程控銑床加工機翼斷面下輪廓線時每一刀只能沿x方向和y方向走非常小一步。表3-1給出了下輪廓線上部分數(shù)據(jù)但工藝要求銑床沿x方向每次只能移動0.1單位.這時需求出當x坐標每改變0.1單位時y坐標。試完成加工所需數(shù)據(jù)，畫出曲線.71/229交通事故調(diào)查（司機有責任嗎？）一輛汽車在拐彎時急剎車，結(jié)果沖進路邊溝里．警察聞訊趕到現(xiàn)場，對汽車留在路上剎車痕跡進行了測量，確定剎車痕跡近似為一圓弧。訊問司機時，他說，當車進入彎道時剎車失靈，且進入彎道時汽車時速為40英里／小時。此速度系該路速度上限。經(jīng)過驗車證實該車制動器在事故發(fā)生時確實失靈，但司機所說車速是否真實呢？72/229交通事故調(diào)查（司機有責任嗎？）通常，作一條基準線來測量剎車痕跡.（水平）距離x沿基準線測得，（垂直）距離y與x垂直.表3-6給出了剎車痕跡測量相關(guān)值.警察測量了路坡度，發(fā)覺這段路是平．請給出一個使警察能夠查對速度計算方法．73/229

航行區(qū)域警示線某海域上頻繁地有各種噸位船只經(jīng)過。為確保船只航行安全，相關(guān)機構(gòu)在低潮時對水深進行了測量,表3-8是他們提供測量數(shù)據(jù)：表3-8.水道水深測量數(shù)據(jù)x 129.0140.0103.588.0185.5195.0105.5 y 7.5141.523.0147.022.5137.585.5 z 4868688 x 157.5107.577.081.0162.0162.0117.5 y -6.5-81.03.056.5-66.584.0-33.5 z 9988949 74/229航行區(qū)域警示線其中（x,y）為測量點，z為（x,y）處水深(英尺)。船噸位能夠用其吃水深度來反應(yīng)，分為4英尺、4.5英尺、5英尺和5.5英尺4檔。航運部門要在矩形海域（75，200）×（－50，150）上為不一樣噸位航船設(shè)置警示標識。請依據(jù)測量數(shù)據(jù)描述該海域地貌，并繪制不一樣噸位警示線，供航運部門使用。水深z是區(qū)域坐標（x,y）函數(shù)z=z(x,y），測量數(shù)據(jù)只是它部分取值?？衫L制函數(shù)圖象和等值線圖，將不一樣吃水線標識圖上75/229插值與擬合（基本原理和區(qū)分）已知有n+1個節(jié)點（xj,yj），j=0,1,…,n其中xj互不相同節(jié)點(xj,yj)可看成由某個函數(shù)y=f（x）產(chǎn)生f解析表示式可能十分復(fù)雜或不存在封閉形式,也能夠是未知76/229插值與擬合（基本原理和區(qū)分）插值：結(jié)構(gòu)一個相對簡單函數(shù)y=g(x)使g經(jīng)過全部節(jié)點即使g(xj)=yj，j=0,1,…,n用g(x)作為函數(shù)f(x)近似。插值指令yi=interp1(x1,y1,xi,'method')對應(yīng)于插值函數(shù)yi=g(xi)其中x1,y1為節(jié)點向量method＝四個選項：‘nearest’為近鄰插值；‘linear’為線性插值；‘spline’為樣條插值；'cubic'為立方函數(shù)插值。77/229插值與擬合（基本原理和區(qū)分）多項式擬合對給定數(shù)據(jù)（xj,yj），j=0,1,…,n選取適當階數(shù)多項式(也可采取其它形式函數(shù))比如二次多項式g(x)=ax^2+bx+c使g(x)盡可能迫近（擬合）這些數(shù)據(jù)擬合指令polyfit、polyval用p=polyfit(x1,y1,m)做m次多項式擬合擬合數(shù)據(jù)向量為x1,y1多項式系數(shù)為p=[p(1),…,p(m),p(m+1)]即g(x)=p(1)x^m+…p(m)x+p(m+1)用y=polyval(p,x)計算在x處多項式值y78/229觀察插值、擬合效果運行觀察程序exp3_1.m選取一個已知函數(shù)作為參考，并將這一函數(shù)圖象用虛線顯示在圖中。觀察程序允許用鼠標選取節(jié)點按鼠標左鍵選點，按右鍵選最終一個點觀察不一樣選點方式對各種插值和擬合效果影響79/229程序注解（inline指令）定義內(nèi)聯(lián)函數(shù)：inline指令g=inline('x^2-x^4');程序80/229程序注解（ginput）[x,y,button]=ginput(n)用鼠標在屏幕選n個點，返回這n個點，存于x,y中。button統(tǒng)計了選點時使用鼠標鍵方式：1為左鍵、2為中間鍵、3為右鍵。81/229程序注解（插值擬合）xx=linspace(a,b,n);%定義自變量xxynearest=interp1(x1,y1,xx,'nearest');ylinear=interp1(x1,y1,xx,'linear');yspline=interp1(x1,y1,xx,'spline');[p,c]=polyfit(x1,y1,4);ypolyfit=polyval(p,xx);

82/229程序注解（插值擬合）subplot(2,2,1),h=plot(xx,ynearest,'r-');set(h,'linewidth',2)subplot(2,2,2),h=plot(xx,ylinear,'r-');set(h,'linewidth',2);subplot(2,2,3),h=plot(xx,yspline,'r-');set(h,'linewidth',2)subplot(2,2,4),h=plot(xx,ypolyfit,'r-');set(h,'linewidth',2)83/229插值擬合效果觀察沿曲線選?。硞€節(jié)點，保持等間隔。當節(jié)點較少時，插值效果怎樣？加密節(jié)點，共８個等距節(jié)點，觀察插值效果，假如去掉中間一個節(jié)點，插值效果又會怎樣？有意偏離原來曲線，假如誤差較大，將會怎樣呢？84/229

微分、積分和微分方程試驗四85/229定積分--連續(xù)求和86/229定積分--連續(xù)求和87/229三種方法計算數(shù)值積分（1）定義法，取近似和極限。高等數(shù)學(xué)中不是重點內(nèi)容但數(shù)值積分各種算法卻是基于定義建立

（2）用不定積分計算定積分。不定積分是求導(dǎo)逆運算，而定積分是連續(xù)變量求和（曲邊梯形面積）表面上看是兩個完全不一樣概念，經(jīng)過牛頓－萊布尼茲公式聯(lián)絡(luò)在一起，（3）解微分方程計算定積分88/229微積分學(xué)基本定理尤其，F(xiàn)(b)-F(a)就是所需定積分.在高等數(shù)學(xué)中總是期望求出不定積分封閉解.但數(shù)值積分是更有用工具。牛頓－萊布尼茲公式不愧為微積分“基本定理”。89/229基本定理推廣（解微分方程計算定積分）90/229基本定理推廣（解微分方程計算定積分）91/229解微分方程Eular折線法92/229解微分方程Eular折線法將區(qū)間n=4等分（共有5個分點)；計算分點和對應(yīng)函數(shù)值（x(1),x(2),x(3)x(4)x(5))（f(1),f(2),f(3),f(4),f(5))在第一個子區(qū)間[x(1),x(2)]上，畫出折線段y(2)=y(1)+f(1)*(x-x(1))代替解曲線段y(x)，這里y(1)=y0=0折線段起點為[x(1),y(1)]，終點為[x(2),y(2)].運行exp4_1.m，觀察第二、三、四子區(qū)間情況。93/229符號微積分用Matlab符號工具箱（SymbolicToolbox）能夠進行符號演算94/229符號微積分(創(chuàng)建符號變量）

symvar創(chuàng)建單個符號變量；symsvar1var2…

創(chuàng)建多個符號變量；f=sym(‘符號表示式’)創(chuàng)建符號表示式，賦予f；

equ=sym('equation')創(chuàng)建符號方程

。95/229符號微積分(極限）limit(‘表示式’,var,a)：求當var→a，表示式極限例：求極限：symsxaI1=limit(‘(sin(x)-sin(3*x))/sin(x)’,x,0)運行結(jié)果96/229符號微積分(求導(dǎo)）diff(f,‘var’,n)求f對變量varn階導(dǎo)數(shù)缺省n時為求一階導(dǎo)數(shù)缺省變量'var'時，默認變量為x可用來求單變量函數(shù)導(dǎo)數(shù)多變量函數(shù)偏導(dǎo)數(shù)還能夠求抽象函數(shù)導(dǎo)數(shù)97/229符號微積分(求導(dǎo)）例：求symsxyf=sym('exp(-2*x)*

cos(3*

x^(1/2))')diff(f,x)運行98/229符號微積分(求導(dǎo)）symsxyg=sym('g(x,y)')f=sym('f(x,y,g(x,y))')diff(f,x)diff(f,x,2)運行例：求99/229符號微積分(積分）int(f,var)：求函數(shù)f不定積分；int(f,var,積分下限，積分上限)：求函數(shù)f定積分或廣義積分例：求不定積分symsxyzI1=int(sin(x*y+z),z)100/229符號微積分(積分）symsxyzI2=int(1/(3+2*x+x^2),x,0,1)I3=int(1/(3+2*x+x^2),x,-inf,inf)101/229符號微積分(化簡、提取和代入)符號運算結(jié)果比較繁瑣，使用化簡指令可對其進行化簡。不過不能指望機器能夠完成一切，人推理往往必須。慣用化簡指令以下展開指令：expand(表示式)；因式分解：factor(表示式)降冪排列：collect（表示式，var）；普通化簡：simplify(A)；102/229符號微積分(化簡、提取和代入)觀察：將展開(a+x)^6-(a-x)^6，然后作因式分解。t_expand=expand(t)t_factor=factor(t_expand)t_simplify=simplify(t)觀察結(jié)果103/229數(shù)值微積分（梯形公式和辛普森公式）trapz(x,y)，按梯形公式計算近似積分；其中步長x=[x0x1…xn]和函數(shù)值y=[f0f1…fn]為同維向量，q=quad('fun',a,b,tol,trace,P1,P2,...)（低階方法，辛普森自適應(yīng)遞歸法求積）q=quad8('fun',a,b,tol,trace,P1,P2,...)（高階方法，自適應(yīng)法Cotes求積）在一樣精度下高階方法quad8要求節(jié)點較少。104/229［x,y］=ode23('fun',tspan,y0,option)（低階龍格－庫塔函數(shù)）［x,y］=ode45('fun',tspan,y0,option)（高階龍格－庫塔函數(shù)）105/229應(yīng)用、思索和練習(xí)(追擊問題)我緝私雷達發(fā)覺，距離d處有一走私船正以勻速a沿直線行駛，緝私艦馬上以最大速度（勻速v）追趕。若用雷達進行跟蹤，保持船瞬時速度方向一直指向走私船，緝私艦運動軌跡是怎樣？是否能夠追上走私船？假如能追上，需要用多長時間？106/229應(yīng)用、思索和練習(xí)(追擊問題)107/229應(yīng)用、思索和練習(xí)(追擊問題)r=dsolve(‘eq1,eq2,…’,‘cond1,cond2,…’,‘v’)方程符號解symsydrxs1=dsolve('D2x=-r*sqrt(1+Dx^2)/y','x(20)=0','Dx(20)=0','y')

xs=simplify(xs1)運行結(jié)果，畫彗星圖108/229應(yīng)用、思索和練習(xí)(追擊問題)r=dsolve(‘eq1,eq2,…’,‘cond1,cond2,…’,‘v’)方程符號解symsydrxs1=dsolve('D2x=-r*sqrt(1+Dx^2)/y','x(20)=0','Dx(20)=0','y')

xs=simplify(xs1)運行結(jié)果，畫彗星圖109/229應(yīng)用、思索和練習(xí)(追擊問題，假如雷達失效)當緝私艦雷達發(fā)覺d處有一走私船后，雷達突然損壞若假定走私船作勻速直線運動（但不知方向），且緝私艦艇速度v大于走私船速度a，則緝私艦應(yīng)采取什么樣航行路線，不論走私船從哪個方向逃跑，都能追捕上它？110/229實時動畫制作（見試驗10）觀察：模擬彈簧振動討論最簡單情形，一彈簧系統(tǒng)作橫向運動，其位移由u=2+cos(t)所決定，仿真彈簧振動111/229實時動畫制作(初始化、見試驗10）程序講解animinit('onecart1Animation')axis([-26-1010]);holdon;u=2;

xy=[0000uuu+1u+1uu;-1.201.2001.21.2-1.2-1.20];x=xy(1,:);y=xy(2,:);plot([-1020],[-1.4-1.4],'k-','LineWidth',2);hndl=plot(x,y,'k-','EraseMode','XOR','LineWidth',2)

112/229實時動畫制作(初始化、見試驗10）zxy10-2set(gca,'UserData',hndl);fort=1:0.025:1000;u=2+exp(-0.00*t)*cos(t);x=[0000uuu+1u+1uu];

hndl=get(gca,'UserData');set(hndl,'XData',x,'YData',y);drawnowend113/229電影動畫制作（zxy7_3)moviein、getframe、movie指令x=-8:0.5:8;[XX,YY]=meshgrid(x);r=sqrt(XX.^2+YY.^2)+eps;Z=sin(r)./r;surf(Z);%畫出禎theAxes=axis;%保留坐標值，使得全部幀都在同一坐標系中114/229電影動畫制作fmat=moviein(20);%創(chuàng)建動畫矩陣，保留20禎forj=1:20;%循環(huán)創(chuàng)建動畫數(shù)據(jù)surf(sin(2*pi*j/20)*Z,Z)%畫出每一步曲面axis(theAxes)%使用相同坐標系fmat(:,j)=getframe;%拷貝禎到矩陣fmat中endmovie(fmat,10)%演示動畫10次115/229應(yīng)用、思索和練習(xí)(槍支設(shè)計)槍支發(fā)怒后，氣體壓強隨子彈在膛內(nèi)運動而改變。槍管長度x單位為m。壓強p是距離x函數(shù)，經(jīng)過實測得到了一批數(shù)據(jù)，子彈射出槍管時出口速度是設(shè)計者關(guān)心問題，假如一只槍管長0.6096m，其膛孔面積4.56×10-5m2，子彈重量0.956N，試決定這種型號槍支出口速度；更普通，確定出口速度和槍管長度關(guān)系曲線，繪制這一曲線，并作出適當標識。這么問題和你在高等數(shù)學(xué)中處理積分有什么區(qū)分嗎？116/229應(yīng)用、思索和練習(xí)(天然氣井開采量)東方天然氣企業(yè)鉆了一口新氣井，他們希望研究一下將這口井于供氣管路聯(lián)接經(jīng)濟性計算此井壓強隨時間改變曲線，由此得到流量Q與時間t關(guān)系，以此預(yù)計此井總開采量。117/229

最優(yōu)化試驗試驗六118/229最正確水槽斷面問題（矩形斷面）用寬l=24cm長方鐵板折成一個斷面為矩形水槽，問怎樣折法可使水槽斷面面積到達最大119/229最正確水槽斷面問題（梯形斷面）將問題１推廣等腰梯形水槽，問怎樣折法可使水槽斷面面積到達最大?120/229最正確水槽斷面問題（對稱五邊形斷面）將鐵板折成如圖所表示對稱五邊形，問怎樣折法可使水槽斷面面積到達最大？121/229最正確水槽斷面問題（五邊）122/229最正確水槽斷面問題（五邊）運行zxy6_6[s,fval]=fmincon('zxy6_6S',x0,A,b,[],[],lb,ub)求解123/229最正確水槽斷面問題（多邊和無限邊）優(yōu)化變量數(shù)與最大斷面面積關(guān)系斷面形狀 優(yōu)化變量數(shù) 最大斷面積cm2

矩形斷面 172 梯形斷面283.14 對稱五邊形488.637

將鐵板折成對稱7邊形、9邊形，普通為對稱2n+1邊形能夠期望最大斷面面積得到深入增加隨之而來是計算代價也隨之增加。124/229最正確水槽斷面問題（無限邊）125/229最正確水槽斷面問題（無限邊）126/229最正確水槽斷面問題（無限邊）127/229微分法求最大和最?。ǜ叩葦?shù)學(xué)）128/229微分法求最大和最小（高等數(shù)學(xué)）129/229微分法求最大和最?。ǜ叩葦?shù)學(xué)）運行zxy6_1.msymsx1x2%定義符號變量。f=x1^3-x2^3+3*x1^2+3*x2^2-9*x1;%函數(shù)z。v=[x1x2];df=jacobian(f,v)%計算雅可比。[X,Y]=solve(df(1),df(2))%用指令solve求駐點。運行zxy6_2畫等值線圖并將點標注在圖上130/229微分法求最大和最?。ǜ叩葦?shù)學(xué)）jacobian(f,v)：計算函數(shù)符號梯度、雅可比矩陣如：若f(v),v=[v1v2]則df=[df/dv1df/dv2]如：若f(v)=[f1(v)f2(v)],v=[v1v2]則df=[df1/dv1df1/dv1

df2/dv1df2/dv2]131/229微分法求最大和最?。ǜ叩葦?shù)學(xué)）solve指令：solve('eqn1','eqn2',…,'eqnn')求n個方程eqn1，eqn2，…，eqnn所組成方程組根（符號解）132/229盲人下山與迭代尋優(yōu)一個盲人處于山上某一點x0要下到谷底，他應(yīng)怎樣做？133/229盲人下山與迭代尋優(yōu)134/229Matlab優(yōu)化工具箱介紹多元函數(shù)無約束優(yōu)化指令fminunc、fminsearch剖析135/229Matlab優(yōu)化工具箱介紹觀察：在命令窗口鍵入bandemo選擇不一樣方法觀察對香蕉函數(shù)優(yōu)化結(jié)果和過程。136/229Matlab優(yōu)化工具箱介紹[x,fval,exitflag,output,grad,hessian]=fminunc(fun,x0,options,P1,…)其中有些項能夠缺省，如exitflag,output,grad,hessian,options,P1,P2,…等等。x0是初始點；fun是目標函數(shù)，能夠用inline指令或建立M文件方法生成目標函數(shù)；137/229Matlab優(yōu)化工具箱介紹參數(shù)傳遞：使用變量P1,P2,…可在目標函數(shù)和主程序之間需要傳遞一些參數(shù)也可使用全局變量Gobal說明來進行傳遞。輸出：x,fval,exitflag,output,grad,hessian為輸出信息，最優(yōu)點、最優(yōu)函數(shù)值、算法結(jié)束狀態(tài)（exitflag>0算法收斂；=0到達最大步驟而停頓；<0算法不收斂）、算法結(jié)束后一些信息（如迭代次數(shù)、所使用優(yōu)化方法等，可在命令窗口查看output詳細內(nèi)容）、梯度值和海森矩陣，除x之外均可缺省。138/229Matlab優(yōu)化工具箱介紹OPTIONS（控制參數(shù)）OPTIONS是一個數(shù)組，有多項內(nèi)容使用optimset對它進行修改和設(shè)置。optimset(‘屬性’，‘屬性值’,…)修改OPTIONS默認值，如默認屬性‘LargeScale’＝‘on’，(使用“信賴域算法”)。假如要使用其它方法，就要修改此項設(shè)置。139/229Matlab優(yōu)化工具箱介紹內(nèi)聯(lián)函數(shù)inlineinline('函數(shù)表示式','變量1','變量2',…)，不混同情形下變量項能夠缺省。f=inline('100*(x(2)-x(1)^2)^2+(1-x(1))^2'),x=[2,2],y=f(x),備注140/229Matlab優(yōu)化工具箱介紹

用M文件生成目標函數(shù)（套用以下格式）myfun.mfunction[f,g,H]=myfun(x)％關(guān)鍵詞function不可省，函數(shù)myfun和M文件同名。f=…%計算目標函數(shù)x。ifnargout>1%如有兩個輸出量（[目標函數(shù)，梯度]）。g=…%計算g為函數(shù)x點解析梯度（可?。?。ifnargout>2%如有三個輸出量（[目標函數(shù)，梯度，海森陣]）。H=…%H為函數(shù)在x點海森陣，（可省）。end141/229Matlab優(yōu)化工具箱介紹（zxy6_4講解運行）bandemo.m簡化和剖析程序zxy6_4.m是對bandemo.m簡化基本結(jié)構(gòu)為：（1）繪制香蕉函數(shù)等值線圖，并將StartPoint和Solution標在圖形上。（2）用Switch語句結(jié)構(gòu)，允許程序選取BFGS、DFP、最速下降法和單純形法等四種優(yōu)化方法。142/229Matlab優(yōu)化工具箱介紹多變量約束優(yōu)化指令fmincon[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options,P1,P2,…)上面命令等價于143/229Matlab優(yōu)化工具箱介紹線性規(guī)劃linprog指令[x,fval,exitflag,output,lambda]=linprog(f,A,b,Aeq,beq,lb,ub,x0,options)

算法選擇：options=optimset(‘largescale’,‘off’)，單純形方法；options=optimset('largescale','on')，內(nèi)點法(默認)。144/229Matlab優(yōu)化工具箱介紹一元函數(shù)尋優(yōu)fminbnd指令[x,fval,exitflag,output]=fminbnd(fun,x1,x2,options,P1,P2,...)此時x,x1,x2是標量，f(x)為標量函數(shù)。145/229Matlab優(yōu)化工具箱介紹Quadprog：解二次規(guī)劃lsqnonlin：解非線性最小二乘lsqcurvefit：非線性數(shù)據(jù)擬合lsqnonneg：非負系數(shù)最小二乘法。lsqlin：約束最小二乘146/229應(yīng)用思索與練習(xí)（計算最正確水槽斷面面積）試推導(dǎo)對稱2n+1邊形面積普通公式選擇一系列n值，仿照zxy6_6.m計算它們最大斷面面積，觀察計算結(jié)果規(guī)律性。在工程實踐中并不能確保每一次計算都能夠成功，不過本問題即使不成功，你是否也能洞察結(jié)果？147/229對盲人下山問題，引入一個有界約束區(qū)域，試用圖形表現(xiàn)函數(shù)在區(qū)域邊界上圖象。能夠用等值線或用曲頂柱體曲面顯示函數(shù)在區(qū)域改變情況。不過提議單獨用二維繪圖指令plot畫出它們曲線圖，觀察函數(shù)在邊界極值情況。應(yīng)用思索與練習(xí)（盲人約束下山）148/229應(yīng)用思索與練習(xí)（盲人約束下山）結(jié)合高等數(shù)學(xué)知識，假如可能，用Matlab符號演算指令求出函數(shù)在不一樣約束下極值點和最值點（比如可用Lagrange函數(shù)方法處理這些問題）。你也能夠在盲人下山模擬中對有約束情況進行討論，這時盲人應(yīng)該怎樣前進呢？149/229應(yīng)用思索與練習(xí)（啤酒配方問題）某啤酒廠希望用原料摻水方法生產(chǎn)一個復(fù)合標準低成本啤酒。其標準要求為：酒精含量為3.1％；發(fā)酵前平均比重在1.034～1.040之間；顏色在8～10EBC單位之間；每升混合物中，蛇麻子脂含量在20～25mg之間。請依據(jù)相關(guān)數(shù)據(jù)算出最優(yōu)配方。150/229應(yīng)用思索與練習(xí)（

儲能飛輪設(shè)計）下面表示式用于設(shè)計儲能用飛輪。準則是儲備能量最大。用約束條件限定了重量、直徑、轉(zhuǎn)速和厚度，試計算最優(yōu)解。你能確定算出解是最優(yōu)嗎？151/229應(yīng)用思索與練習(xí)（齒輪減速器設(shè)計）抽去各變量物理意義，齒輪減速器最優(yōu)設(shè)計模型以下：這是一個含有7個變量、23個約束優(yōu)化問題。試對其進行計算。你可能會碰到很大困難，你能想方法處理這些困難嗎？152/229應(yīng)用思索與練習(xí)（齒輪減速器設(shè)計）153/229方程求根、不動點和迭代試驗七154/229隱函數(shù)存在定理可視化155/229隱函數(shù)存在定理可視化選擇特殊例子運行zxy7_1.m，畫出曲面z=F(x,y)、x-y平面圖像和它們交線。畫出曲線z=F(x0,y)；（備注）156/229隱函數(shù)存在定理可視化157/229隱函數(shù)存在定理可視化確定隱函數(shù)，方程求根[x,fval,exitflag,output]=fzero(fun,x0,options)zxy7-2.m158/229蛛網(wǎng)圖與不動點迭代稱滿足方程f(x)=x點x為函數(shù)f不動點求函數(shù)f不動點。能夠從一個初始點x0出發(fā)，以格式xn+1=f(xn)進行迭代；得到x0,x1,x2,……,xn,…..假如該數(shù)列是收斂，則159/229蛛網(wǎng)圖與不動點迭代160/229蛛網(wǎng)圖與不動點迭代運行觀察程序zxy7_3，了解蛛網(wǎng)圖原理161/229簡單和復(fù)雜：二次函數(shù)迭代和混沌觀察對二次函數(shù)f(x)=rx(1-x)進行迭代，其中0<r<1是一個可變參數(shù)。1）固定若干個不一樣值，觀察迭代序列極限；迭代N次，略去前n個迭代值，并將后N–n個迭代值畫在r-x坐標系中（zxy7_4)2）用蛛網(wǎng)圖觀察三種不一樣類型迭代。(zxy7_5)3）加密r取值，得到加密Feigenbaum圖。（zxy7_4改變參數(shù))162/229線性代數(shù)試驗試驗八163/229向量組線性關(guān)系164/229向量組線性關(guān)系（產(chǎn)生相關(guān)向量，運行zxy8_1）產(chǎn)生向量：產(chǎn)生m個n維向量，且各向量分量均在[a,b]之間。clearn=3;m=2;a=-10;b=10;rand('seed',32),A=unifrnd(a,b,[n,m])組合向量：產(chǎn)生m=2個組合系數(shù),將組合得到新向量并入矩陣A中：x=unifrnd(-1,1,[1,m]),A(:,3)=x(1)*A(:,1)+x(2)*A(:,2)165/229向量組線性關(guān)系（產(chǎn)生相關(guān)向量，運行zxy8_1)運行zxy8_1A=-1.7383-9.17070.02568.70646.6219-9.08424.70967.5495-0.0246-8.9245-3.53318.3272-5.74703.9105-0.0038-0.2352-6.61972.1934166/229向量組線性關(guān)系（產(chǎn)生相關(guān)性判別）167/229Gauss消元法（運行rref(A),rrefmovie(A)）

rref(A)將矩陣A做行初等變換階梯化。B=rref(A)B=1.00000-0.0011-0.53600.58510.259001.0000-0.0026-0.8478-0.83300.9415000000rrefmovie(A)：觀察到行初等變換過程

168/229Gauss消元法(同解方程）Rref(A)169/229Gauss消元法(解）向量形式170/229Gauss消元法(基礎(chǔ)解系）Ax=0基礎(chǔ)解系可由下式計算：X=[-B(1:r,r+1:m+s);eye(l)]其中r=rank(A)，l=m+s-r為基礎(chǔ)解系個數(shù)。r=2;m=2;s=4;B=[-B(1:2,r+1:m+s);eye(m+s-r)]問題：怎樣用Matlab解普通非線性齊次方程組，如A(:,1:4)X=A(:,7)？171/229應(yīng)用練習(xí)與思索(平面四桿機構(gòu)設(shè)計)172/229應(yīng)用練習(xí)與思索(平面四桿機構(gòu)設(shè)計)某操縱裝置采取四桿鉸鏈機構(gòu)。已知兩連架桿（L1,L3）輸入角和輸出角滿足下表數(shù)據(jù)所表示對應(yīng)關(guān)系，機架長度L4=50mm，試確定其余三桿長度。173/229應(yīng)用練習(xí)與思索(平面四桿機構(gòu)設(shè)計)174/229應(yīng)用練習(xí)與思索(平面四桿機構(gòu)設(shè)計)確定四桿長度，并用Matlab繪圖指令用圖形表示你設(shè)計結(jié)果。你需要設(shè)計一個表現(xiàn)方案，使人能夠很明白看出你設(shè)計結(jié)果是正確。假如只利用表前兩組對應(yīng)角度值，設(shè)計方案還是唯一嗎？計算一下結(jié)果。一樣給出直觀表示。體會到自由變量含義？假如表中值為4組對應(yīng)角度值，你就碰到超定方程了。它沒有準確解，只有近似解。你愿意用Matlab去解它嗎？試一試。175/229應(yīng)用練習(xí)與思索(用Matlab做線性代數(shù)題)176/229應(yīng)用練習(xí)與思索(用Matlab做線性代數(shù)題)symsaa1=[1;4;0;2];a2=[2;7;1;3];a3=[0;1;-1;1];a4=[3;10;a;4];A=[a1,a2,a3,a4]fori=2:4%行初等變換A(i,:)=A(i,:)-A(1,:)*A(i,1);endA(2,:)=A(2,:)/A(2,2);fori=3:4A(i,:)=A(i,:)-A(2,:)*A(i,2);end177/229矩陣相同化簡178/229矩陣相同化簡179/229矩陣相同化簡選擇方陣A，如二階方陣

A=[1/5,99/100;1,0];選擇一個初始點(二維列向量)，按下面公式進行迭代:xk+1=Axk觀察這些迭代點位置和趨向180/229矩陣相同化簡(程序zxy8_2.m迭代部分)Clear,clf,n=40;a=-20*100;b=-a;c=a;d=b;p=0.1;A=[1/5,99/100;1,0];axis([abcd]),grid,holdonbutton=1whilebutton==1[xi,yi,button]=ginput(1);plot(xi,yi,'ko'),holdon,X0=[xi;yi];X=X0;fori=1:nX=[A*X,X0];h=plot(X(1,1),X(2,1),'k.',X(1,1:2),X(2,1:2),'k:');set(h,'MarkerSize',6),grid,holdonquiver([X(1,2),1]',[X(2,2),1]',[X(1,1)-X(1,2),0]',[X(2,1)-X(2,2),0]',p)endend181/229矩陣相同化簡(程序zxy8_2.m畫直線)p=60;x=linspace(a,b,30);[pc,lamda]=eig(A),pc=-pc;z1=pc(2,1)/pc(1,1)*x;z2=pc(2,2)/pc(1,2)*x;plot(x,z1,'linewidth',2),h=quiver([500,501]',[-1000,-1001]',[pc(1,1),0]',[pc(2,1),0]',p)set(h,'linewidth',2,'color','red'),182/229矩陣相同化簡(高維線性離散動力系統(tǒng))動力系統(tǒng)理論基本目標是了解迭代過程最終或漸進性態(tài)，假如迭代過程是一時間為自變量微分方程，則該理論試圖預(yù)言微分方程解在遙遠未來或遙遠過去最終性態(tài)。假如過程是像函數(shù)迭代離散過程，則這種理論希望了解伴隨n變大，迭代點最終性態(tài)。就是說，動力系統(tǒng)提出了一個聽起來非數(shù)學(xué)問題：這些點跑到哪里去？當它們抵達那里又在干些什么？183/229矩陣相同化簡(高維線性離散動力系統(tǒng))記線性動力系統(tǒng)L(x)，它是L(x)=Lx其中L是適當維數(shù)方陣，我們只考慮2～3維情形。184/229矩陣相同化簡(高維線性離散動力系統(tǒng))185/229矩陣相同化簡(高維線性離散動力系統(tǒng))考查前面在2維情況例子中標準形A=[2,0;0,1/2]它含有混合特征值運行zxy8_3.m進行觀察迭代情形。觀察B＝PAP^-1迭代情況，比較異同。186/229主成份分析和線性變換氣象分析預(yù)報中需要分析很多變量指標。何抓住主要特點，用較少指標代替原來較多指標，又能綜合反應(yīng)原來較多指標信息，就是實際工作者所關(guān)心問題。（降維）主成份分析方法為此提供了一個有效伎倆。187/229主成份分析和線性變換設(shè)有兩個變量指標：X1：代表某地10月副高強度指數(shù)；X2：代表該地10月副高面積指數(shù)。188/229主成份分析和線性變換運行zxy8_4.m畫出散點圖。怎樣找到適當坐標軸，使信息損失盡可能??？189/229主成份分析和線性變換190/229主成份分析和線性變換求協(xié)方差矩陣求正交矩陣P，滿足運行觀察191/229主成份分析和線性變換(運行zxy8_4)192/229線性變換運行zxy8_5193/229數(shù)學(xué)試驗與Matlab/~zhou/index.html制作人周曉陽華中科技大學(xué)數(shù)學(xué)系zxyhust@194/229

Galton釘板和二項分布

分布列意義195/229Galton釘板模擬英國生物統(tǒng)計學(xué)家Galton設(shè)計了Galton板右邊是一個5層Galton釘板196/229Galton釘板模擬(原理)在一板上釘有n排釘子自頂端扔進一小球任其自由下落在下落過程中小球碰到釘子，左右落下機會相等最終小球落入底板中某一個格子圖中用012345表示這6個格子197/229Galton釘板模擬(博彩問題)在每一格子中放上適當價值獎品如依次為1010.20.218(元）扔一次小球你要付1元給莊主假如小球落入某個格子你將取得對應(yīng)價值獎品你合算嗎？莊主會盈利嗎？198/229Galton釘板模擬(扔1萬個小球)小球落入哪一個格子是不確定所以要計算落入每一個格子可能性試想向Galton板中扔10000個小球這些小球?qū)⒍逊e起來小球堆積形狀告訴了我們什么呢？199/229Galton釘板模擬(程序zxy9_1.m)(1)確定釘子位置：將釘子橫、縱坐標存放在兩個矩陣X和Y之中。(2)選取0<p<1,將[0,1]區(qū)間分成兩段:[0,p)和[p，1]。(3)產(chǎn)生隨機數(shù)r=rand(1，1),假如r<p，讓小球向右落下；若r>p，讓小球向左落下。(見備注）(4)將這一過程重復(fù)n次，并用直線連接小球落下時所經(jīng)過點，這么就模擬了小球從頂端隨機地落入某一格子過程。

200/229Galton釘板模擬(程序zxy9_1.m)(5)模擬小球堆積形狀。輸入扔球次數(shù)m（比如m=100)計算落在第i個格子小球數(shù)在總球數(shù)m中所占百分比f(i)當模擬結(jié)束時，就得到了頻率:f(0),f(1),f(2),f(3),f(4)畫出它們圖形。就是小球堆積形狀201/229Galton釘板模擬(程序zxy9_1.m)（6）動畫指令結(jié)構(gòu)moviein(n)：創(chuàng)建動畫矩陣；

制作動畫矩陣數(shù)據(jù)；getframe：拷貝動畫矩陣；movie(Mat,m)：播放動畫矩陣m次，（zxy7_6演示、講解，備注）202/229Galton釘板模擬(程序zxy9_1.m)運行zxy9_1.m203/229一個模擬結(jié)果扔100個小球向右概率p=0.5要改變參數(shù)觀察一下不一樣模擬結(jié)果嗎？這很輕易.自己動手試試吧204/229隨機變量及其分布當你扔小球時，你和莊家關(guān)心什么？？？？？？？？？？？？？對，是小球落入格子編號數(shù)X（有些繞口，但很主要）在投球前，你不能說你小球會落在第0個格子。但你能夠說小球?qū)⒙湓诘赬個格子X是一個隨機數(shù)是概率論中主要討論對象-----隨機變量?。?！205/229隨機變量及其分布實際上,更應(yīng)該關(guān)心是X分布列分布列是小球落在各格子里概率：P(X=0),P(X=1),P(X=2),P(X=3),P(X=4)，P(X=5)想一想，它是不是表現(xiàn)了大量投球后小球堆積極限形狀呢備注（比較頻率和概率）206/229Bernoulli試驗和二項分布不要把Galton釘板簡單地看成消遣它是一個有用概率模型當你學(xué)習(xí)了概率論，你將知道Bernoulli試驗?zāi)Ｐ蚽重Bernoulli試驗成功次數(shù)X服從二項分布B(n,p).上面模擬對應(yīng)于n=5,p=0.5情形207/229二項分布列

隨機變量X～B(n,p)，則它分布列為208/229統(tǒng)計工具箱用指令f=binopdf(x,n,p)可計算二項分布分布列用F=binocdf(x,n,p)可計算二項分布分