第一節(jié) MATLAB 中的矩陣的輸入 - 湖北民族學院(共31頁)

1、第一節(jié) MATLAB 中的矩陣(j zhn)的輸入1 直接(zhji)輸入一、直接(zhji)在工作窗中輸入：A=2, 4, 6, 8;1 3 5 7; 0 0 0 0;1,0,1,0 其意義是定義了矩陣 二、如果矩陣中的元素是等步長的，可以用下面的方法 A=1:0.2:2;1:6;2:2:12A=1:5 “”號在這里表示為轉(zhuǎn)置，而 1:5 中間少了一個循環(huán)步長，此時將步長自動取為 1 。2 增刪改設(shè)已經(jīng)定義 A=1 2 3 4 5；10 8 6 4 2; B=0 1;1 0; C=1 2;2 4，即已定義A= B= C= 1 2 3 4 5 0 1 1 28 6 4 2 1 0 2 4則命令

2、：A=A(:,1:4);C,B,0 2 0 4 將 A 定義成：A= 而 A(:,3)=； 將刪除 A 的第三列 ，得 2 3 4 0 A= 1 2 4 010 8 6 4 2 10 8 4 22 0 1 0 1 2 1 04 1 0 4 2 4 0 43 命令生成使用(shyng) MATLAB 命令生成矩陣(j zhn)一般使用下面的命令 命令(mng lng) linspace，它有兩個格式：a1=linspace(1,100) %生成一個從1到100的有100 個元素的向量a2=linspace(0,1)%仍然是有 100 個元素但是是從 0 到 1 的向量 a3=linspace(0

3、,-1) %請與上一個向量進行比較上面是第一種格式 linspace(a,b)，它是將 a 到 b 等分成 100份形成的向量。第二種格式 linspace(a,b,n) 中的 n 為一個正整數(shù)，表示是從 a 到 b 等分成 n 份后形成的向量。例如a4=linspace(1,100,11) %從 1 到 100 但只形成11 個元素的向量 a5=linspace(1,100,10) %自己體會這個命令作用 a6=linspace(0,1,11) %加上了“”表示轉(zhuǎn)置a7=linspace(0,-1,10) %自己體會這個命令作用2 命令 ones,zeros 分別形成元素全為 1或全為零的矩

4、陣它也有兩種格式。請觀察它們的作用： ones(6,3) %生成 63 階元素全為 1 的矩陣 ones(5) %生成 5 階元素全為 1 的方陣 zeros(3,6) %生成 36 階元素全為零的矩陣 zeros(4) %生成四階元素全為零的方陣 3 命令 diag 生成對角陣及從矩陣的主對角線生成向量，例如：diag(1 3 5 7) %生成了以 1 3 5 7 為主對角線的方陣： ans= 1 0 0 0 0 3 0 0 0 0 5 00 0 7相反如果先定義了一個三階方陣：A=1 2 3;4 5 6;7 8 9顯示(xinsh)： A= 1 2 3 4 5 6 7 8 9則命令(mng

5、 lng) a8=diag(A) 將用 A 的主對角線生成(shn chn)新的列向量：a8= 1 5 9命令 eye(n) 生成 n 階單位方陣，即主對角線上元素為 1，其余元素為零的方陣。例如鍵入：A=eye(5) 將得到：A= 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 00 0 0 0 1 第二節(jié) MATLAB 文件處理 1 文件編輯如果要在 MATLAB 的工作窗定義矩陣，則用鼠標點擊屏幕左上方的 File 選擇項，再從中選擇 New 中的M-file 項并且用鼠標點擊它，就打開了 MATLAB 文件編輯窗并且可以在此窗中定義 MATLAB 矩陣了（注

6、意對于已有的文件，可以選擇 open 來打開它，然后對其進行修改）。在 MATLAB 文件編輯窗中定義的矩陣與工作窗中定義的方法是完全一樣。并且可以在MATLAB 文件編輯窗的菜單中使用菜單命令直接運行。 可以(ky)在MATLAB中使用菜單中的“File”中的“Set path”將當前工作文件夾定義在你正在(zhngzi)工作的文件夾。2 MATLAB工作窗中變量值的保存(bocn)與調(diào)用MATLAB 工作窗中的變量在退出 MATLAB 工作狀態(tài)后值不能保存，如果需要保存，可以使用命令 save 將其存儲到磁盤上，命令格式有兩種：第一種是用二進制格式來存儲。例如先定義三個矩陣：A1=0:3;

7、2*ones(4);4:1:1 ; A2=1 3 2 4;A3=zeros(3,1);生成下列矩陣與向量：鍵入：save file1 A1 A2 A3 %用二進制格式以文件名 file1.mat 存儲 A1,A2,A3save file2.m A1 A3 ascii%用 ascii 碼以文件名 file2.m 存儲 A1,A3我們還要注意：用二進制格式存儲的文件連變量名一起存儲并可再重新調(diào)入時恢復變量的值，而用 ASCII碼存儲的文件只存儲了變量的值，而變量名是沒有的。用二進制格式存儲的變量，可用命令 load 調(diào)用，調(diào)用格式為： load 例如，前面用 save file1 存儲了所有變量

8、A1,A2,A3，調(diào)用時只要鍵入 load file1 即可。 第三節(jié) MATLAB 中的矩陣運算 1 矩陣運算命令與通常線性代數(shù)(xin xn di sh)命令運算的異同一、MATLAB 在運行時是以矩陣為單位進行運算的。它通常有兩種運算，第一種是矩陣運算，運算時滿足線性代數(shù)中矩陣運算所規(guī)定(gudng)一切運算法則，如加、減、乘，乘方即冪運算(yn sun)(當然運算要符合規(guī)定的條件，例如矩陣 A 與矩陣 B 相乘，必須 A 的列數(shù)等于 B的行數(shù))，運算符號：A+B， A-B， A B (注意“*”不能少) An二、不同之處：1、與通常矩陣運算不同之處在：在線性代數(shù)中矩陣不能與數(shù)相加減，而

9、在MATLAB的矩陣運算中允許矩陣與數(shù)相加減。2、函數(shù)命令可以直接作用到矩陣的每一個元素。3、線性代數(shù)中矩陣沒有除法，而MATLAB中有矩陣除法，例如：輸入 A=1:3;4:6;7:9;b=14,32,50;c=Ab顯示： c=204 4、函數(shù)作用到矩陣的每一個元素，例如如果令 A=1 1/2 1/3; 1/2 1/4 1/8*pi，即定義 則 三、MATLAB中除法運算的規(guī)定與意義：1、運算定義：設(shè)已經(jīng)定義好矩陣 A 與矩陣 B，如果矩陣 A 與矩陣B 的行的維數(shù)相同，則 MATLAB 中可以用矩陣 A 左除矩陣 B，即可以令： X=AB如果(rgu)矩陣 A 與矩陣(j zhn)B 的列的

10、維數(shù)相同(xin tn)，則 MATLAB 中可以用矩陣 A 右除矩陣 B，即可以令： X=B/A 2、矩陣除法的意義 給出線性方程組 AX=B，則 X=A 給出線性方程組的一個解。而 X=B/A 給出了線性方程組 XA=B 的一個解。 目前我們用 MATLAB 求線性方程組的解只有三個方法:當 A 是可逆方陣時，X=inv(A)*B 給出線性方程組AX=B 的唯一解。但是 A 不可逆時方法失效?？捎妹?rref 化線性方程組 AX=B 的增廣矩陣為行最簡階梯矩陣方法來求解，但線性方程組可能出現(xiàn)無解(稱為“超定”)、唯一解(稱為“恰定”)及無窮多解(稱為“欠定”)的情形。無論 A 是否可逆都

11、可用MATLAB 除法求解，并且無論何種情形都是唯一解。當方程存在唯一解時，三種方法求出的解是一樣的。但是用除法作的解一般精度更高。方程為“超定”或者“欠定”時解意義則不同。線性方程組 AX=b 為欠定時有無窮個解，X=Ab得到其中解分量中零元素為最多的一個特解。線性方程組 AX=b 為超定時是無解的, 用 X=Ab得到的是使范數(shù) |AX-b| 為最小的解。我們不詳細說明這個范數(shù)的意義，可理解為使 AX-b 最接近于零的解。例如方程組 通解為：輸入(shr) A=1:3;4:6;7:9;b=14,32,50;c=Ab顯示(xinsh)： c=204是其中(qzhng)有一個零分量的特解。輸入

12、d=0 32 50;g=Ad顯示 g= 1.0e+017 * 0.6305 -1.26100.6305再輸入 h=A*g-d顯示 h=329614因此 g 不滿足 A*g=d，只是使 A*g-d 盡可能接近于零。2 常用的數(shù)學運算命令 以下是一些在MATLAB中最常用的數(shù)學運算命令(均用小寫字母，命令的作用可在MATLAB中用 help 查詢其作用與格式)：一、基本函數(shù)運算命令 1、三角函數(shù)： sin cos tan cot sec csc 2、反三角函數(shù)：asin acos atan acot asec acsc atan2(四象限反正切)3、雙曲函數(shù)(hnsh)： sinh cosh ta

13、nh coth sech csch4、反雙曲函數(shù)(hnsh)：asinh acosh atanh acoth asech acsch5、復數(shù)(fsh)：real imag conj angle6、小數(shù)取整：fix(朝零方向)，ceil(朝正無窮方向)，floor(朝負無窮方向)，round(四舍五入)7、對數(shù)與指數(shù)：log10 log exp8、其它：sqrt abs sign sum prod solve二、線性代數(shù)運算命令： 1、det inv rank rref eig cond行列式求值命令 det， 矩陣求秩命令 rank，求方陣的逆方陣命令，inv 求矩陣行最簡階梯陣命令 rref

14、，求特征值與特征向量命令 eig，求矩陣條件數(shù)命令 cond三、多項式運算命令 MATLAB中用向量表示多項式，如 a= 1 2 3 0 4 表示多項式：常用的多項式運算命令有：1、多項式加減法：在次低的向量前面加零后使其元素個數(shù)相同，再按向量加減法運算命令進行。 2、多項式 a 與多項式 b 相乘：c=conv(a, b)； 3、多項式 a除以多項式 b： q, r=deconv(a, b) q 是商， r 余項) 4、多項式求值： p1=polyval(a, x) (多項式 a 在點 x 處的值) 5、方陣多項式求值： p2=polyvalm(a, A) (矩陣多項式 a 在方陣 A 處的

15、值)6、多項式求根： p3=roots(a) (求多項式 a 的根)7、多項式微分： p4=polyder(a)8、多項式擬合： p5=polyfit(x,y,n)例：x=1 1.2 3 4.2 5; y=x.3-2*x.2+3*x-1; p5=polyfit(x,y,3) 將得出： p5=1.0000 -2.0000 3.0000 -1.0000四、高等數(shù)學中的數(shù)值運算(yn sun)命令在MATLAB數(shù)值計算(j sun)中很少有高等數(shù)學中的運算命令，下面只介紹兩個：1、求數(shù)值積分： trapz2、差分(ch fn)： diffMATLAB中也有很多符號運算命令。但不在本課程中。在那些符號

16、運算命令中能實現(xiàn)許多高等數(shù)學中的運算。五、顯示格式 有 format rat, format long, format short, format long g, format short g, format long e , format short e , format hex等等。 第四節(jié) MATLAB 的數(shù)組操作前面我們談到了 MATLAB 中有數(shù)與矩陣的加減法，這是線性代數(shù)中所沒有的。MATLAB 中還有許多與高等數(shù)學中內(nèi)容不同的運算。這些運算對于求解實際問題很有作用，它使編程比其它語言要簡單方便得多。MATLAB 中的數(shù)組運算就是其它計算機語言中所沒有的使編程變得十分簡單的運算。M

17、ATLAB 中的數(shù)組運算只在同階矩陣中進行。設(shè) A,B 是兩個已經(jīng)定義好的同階矩陣： 則數(shù)組運算 A.*B,A./B,A.B 的值分別是：除了(ch le)前面講過的數(shù)與矩陣的加減法(不算(b sun)數(shù)組運算) 外，數(shù)與矩陣之間還有下面(xi mian)的數(shù)組運算： 從上面敘述可見矩陣與同階矩陣之間(及數(shù)與矩陣之間)的數(shù)組運算的定義特點是： 1 定義在兩個同階矩陣或數(shù)與矩陣之間進行；2 定義方法是前面一項與后面一項用小數(shù)點加運算符連結(jié)；3 定義的實質(zhì)是：當定義的運算在矩陣之間時，是相同位置上的元素進行由運算符規(guī)定的運算。當定義的運算是矩陣與數(shù)之間，則是矩陣中的每一個元素與數(shù)進行運算符規(guī)定的運

18、算。對于其他計算機語言(yyn)，這些運算常常要通過雙重循環(huán)才能完成，而在 MATLAB 中只要(zhyo)簡單地用“.”加運算符即可表示(biosh)，這是 MATLAB 特別方便之處。數(shù)組運算及函數(shù)運算定義容易記憶與使用，它在編程中十分方便，我們看后面的例。如果要畫函數(shù) e-xsinx2 +1 在區(qū)間 0,3 上的圖形，操作過程是：x=0:0.01:3; %定義自變量 x 在區(qū)間 0,3 上的值y=exp(x).*sin(x.2)+1 %用數(shù)組運算計算函數(shù)值 yplot(x,y) %作平面上的曲線圖。只用簡單三行命令完成了在其他計算機語言中要用很長程序才能完成的程序，這是 MATLAB 語

19、言非常簡便好用的原因。 第五節(jié) MATLAB 作圖下面(xi mian)介紹 MATLAB 中的作圖命令(mng lng)。1 二維作圖命令(mng lng) plot一、MATLAB 二維作圖命令最常用的是 plot，除作圖外，還有如下可控制的操作：1、一張圖上畫多條曲線(可以使用hold命令)。2、一個屏幕上開多個窗口作圖(使用subplot命令)。3、曲線選擇線型、顏色(在每條曲線后加,.b 或者 ,color,1,0,0.5 等)。4、標注文字，使用命令：title(字符) 圖形標題xlabel(字符) x 軸標注，ylabel(字符) y 軸標注，zlabel(字符) z 軸標注，l

20、egend(字符, 字符,) 可移動標注text(x,y,(z),字符) 在指定坐標標注gtext(字符) 用鼠標在圖形中選擇地方標注5、坐標網(wǎng)格線顯示與關(guān)閉命令 grid二、坐標軸控制(也適用于三維作圖) axis 1、坐標軸刻度范圍控制axis(xmin, xmax, ymin, ymax (, zmin, zmax)。 2、關(guān)閉坐標軸顯示開關(guān) axis(off) 或者 axis off，axis(on) 或者 axis on 3、坐標軸等刻度命令 axis(equal) 或者 axis equal 4、坐標軸等長命令 axis(square) 或者 axis squal 5、使坐標等刻度

21、與等長命令失效命令 axis(normal) 或者 axis normal三、下面介紹一個在實際應(yīng)用中非常有用的命令：ginput該命令在三維空間中仍然可用，但效果不如二維空間中那么好。四、背景色控制命令figure(color, r, g, b) 三維空間中此命令仍然(rngrn)適用。五、程序(chngx)標注在程序(chngx)中“%”后的部分用作標注。六、還有許多二維作圖命令，請大家參看有關(guān)參考書。我們僅就使用二維等值線 contour命令用于求解線性方程組的例子。 MATLAB 二維作圖命令 contour可以用于求二維超越方程組的數(shù)值解。請看 exmaple1.m:x=0:0.01

22、:2.4995;y=0:0.01:2.4995;A=(0*y+1)*sin(x)+y.2*(0*x+1)+log(5-(0*y+1)*x-y*(0*x+1)-7;contour(x,y,A,0,0)hold onB=3*(0*y+1)*x+2.y*(0*x+1)-(5-(0*y+1)*x-y*(0*x+1).3+1;contour(x,y,B,0,0)hold off此程序用于求方程組：的數(shù)值解。此方法的優(yōu)點是：可以看到圖形全貌，可以求出多解情形，可以獲極高精度。練習 1：求曲線與曲線在區(qū)間0, 上的所有交點。練習 2：求方程組 的數(shù)值(shz)解。要求將解代入方程后每個方程左方減去右方后的絕

23、對值小于 0.00001。2 三維作圖命令(mng lng) 一、MATLAB 三維作圖曲線作圖命令(mng lng)最常用的是 plot3，contour3。除作圖外，還有如下可控制的操作：三維曲線作圖是對參數(shù)式為作圖。例如螺旋線： 作圖命令是 t=0:0.01:10*pi;x=cos(t);y=sin(t);z=t/2; plot3(x,y,z)contour3是三維等值線作圖。二、三維曲面作圖 常用的三維曲面作圖命令有：mesh , surf, meshc, meshz , waterfall, 我們重點介紹 mesh 與 surf。這兩個命令都是作曲面圖形，但是 mesh 是用網(wǎng)格作圖

24、，而 surf 是用曲面片作圖。效果有少許不同，但作圖編程是相同的(這個作圖過程比較復雜)。事實上它們都是用曲面參數(shù)式 作的圖，因此要畫曲面(qmin)圖形，首先得會寫上面的參數(shù)式。然后會用上面的參數(shù)式編程。注意在編程時 z 都是二維矩陣(j zhn)，行數(shù)與向量 u 的維數(shù)相同(xin tn)，列數(shù)與 v 的維數(shù)相同，x, y 可以是與 z 維數(shù)相同的二維矩陣也可以是向量。請看下面馬鞍面 的作圖過程：y= -1 :0.1: 1; x=y; z=(0*y+1)*x.2- y.2* (0*x+1); mesh(x,y,z);這個命令等價于y= -1 :0.1: 1; x=y; z= ones(s

25、ize(y)*x.2 - y.2*ones(size(x); surf(x,y,z);上面作的圖是以 z 的高度作為圖形色彩數(shù)據(jù)的(這是缺省時的數(shù)據(jù))。如果加上顏色、光照、紋理、視角等設(shè)置，可以得出十分漂亮的圖形。請看 example2.m。有時作圖的參數(shù)式要使用比較復雜的數(shù)學運算，例如 2001 年全國數(shù)學建模競賽時要求從血管的 100 張切片中讀出血管的中心與半徑。競賽中要求將血管中心線向三個坐標平面投影。但是我們根據(jù)血管中心線可畫出血管的三維圖象。顯然這更好看更直觀。曲面參數(shù)式的推導見廣西大學學報2003年第3期，而程序請看 example3.m。練習3：畫一個求生圈的圖。提示：求生圈的

26、曲面參數(shù)方程為： 練習4：畫出在二維區(qū)域 x2 + y21 上的馬鞍面 z= x2- y2。練習5：作曲面 z=x2+y2 與平面 4x+2y+5z=20 相交的圖(要求平面只畫出被曲面截出的部分，參考右圖)。練習6：作圓柱面 x2+z2=1 與圓柱面 y2+z2=1 相交的圖(要求平面只畫出在每卦限的部分)。3 MATLAB 動畫 MATLAB 中的動畫功能不是太強，并且需要(xyo)使用 MATLAB 的底層作圖命令(mng lng)，即圖形的句柄操作。MATLAB 中的作圖命令如果使用了賦值語句，即將作圖賦給了某個變量，則該變量即為此張圖形的句柄?？梢詫Υ司浔O(shè)置一些圖形屬性，對于句柄圖

27、形我們在前面已經(jīng)用到了，例如(lr)馬鞍面的圖形中的賦值語句aa=surf(X,Y,Z);就是將曲面的句柄賦給 aa 了，然后用set(aa, Facecolor, 0,0.7,0, bb, Facecolor, 0.7,0,0.7)設(shè)置曲面顏色。MATLAB 動畫實際上就是在句柄中設(shè)置了圖形將是被復蓋還是改寫等屬性，然后在作圖時按不同的位置畫圖(同時擦除前面的圖形或者復蓋前面的圖形以實現(xiàn)動畫)。由于作動畫一般要涉及編程，我們將此部分內(nèi)容放到第四節(jié)中。在圖形的屬性控制中要注意區(qū)分各種圖形對象：figure, axes, line, patch, surface, image, text 等等。

28、4 MATLAB 圖象處理 MATLAB 可對圖象作處理，常用的簡單命令有：imread, image, colormapimread 可以讀出圖象數(shù)據(jù)，而image可以將圖象數(shù)據(jù)顯示為圖象，而colormap可以設(shè)置圖象色彩。練習7：本文件夾中有一張圖“figure0”，是2002年全國數(shù)學建模大賽中的一張圖，要求從圖中讀出圖形的內(nèi)切圓圓心與半徑。 在MATLAB圖象處理中還有光照、色彩、視角等許多處理。由于水平與時間關(guān)系，只能簡單介紹。第六節(jié) MATLAB 編程一、MATLAB 編程 MATLAB 可以單步(dn b)運行，也可以編程運行。 與所有軟件編程一樣，MATLAB中也有條件(ti

29、ojin)、分支、循環(huán)等語句。條件(tiojin)語句： if end 或者 if else end循環(huán)語句： for end 或者 while end需要指出，在其他許多語言中需要循環(huán)的地方在MATLAB中是不需要循環(huán)的。下面我們看一個動畫程序：畫單擺的運動程序：pend.m。figure(color,1 1 1)plot(-0.2;0.2,0;0,color,y,linestyle,-,linewidth,10);%畫懸線的橫梁g=9.8; %重力加速度l=1; %線長theta0=pi/6; %初始角度x0=l*sin(theta0); %初始 x 值y0=-l*cos(theta0);

30、 %初始 y 值axis(-0.75,0.75,-1.25,0); %畫坐標范圍axis(off); %關(guān)閉坐標顯示head=line(x0,y0,color,r,linestyle,.,erasemode,xor,markersize,40); %定義小球格式body=line(0;x0,0;y0,color,b,linestyle,-,erasemode,xor); %定義懸線格式t=0; %時間初值dt=0.0005; %時間增量while 1 %死循環(huán)(關(guān)閉窗口后中止循環(huán)) t=t+dt; theta=theta0*cos(sqrt(g/l)*t); %計算時刻 t 時的幅角 x=l*

31、sin(theta); %計算時刻 t 時的 x 坐標 y=-l*cos(theta); %計算時刻 t 時的 y 坐標 set(head,xdata,x,ydata,y); %計算時刻 t 時小球的位置 set(body,xdata,0;x,ydata,0;y); %計算時刻 t 時懸線的位置 drawnow; %重畫小球(xio qi)與懸線end二、程序(chngx)與函數(shù)將MATLAB運行(ynxng)中的每一句寫在擴展名為“m”的文本文件中，就得到MATLAB程序文件。程序中的變量值都在工作窗中。如果是寫成function u , v, =fun(x , y, )則是MATLAB函數(shù)

32、，函數(shù)參數(shù) x , y 可以是變量、向量或者矩陣，返回值 u, v 也可以是變量、向量或者矩陣，存儲時函數(shù)名必須與存儲文件名相同。 例如，運行 redotree(n, B) (設(shè)的初值為 B=0 0;0 1，取 n=5 或者 n=6)，可以畫出一棵樹的圖形。而運行l(wèi)eafage1(n)可以畫出一張樹葉。練習8：作曲面 z=x2+y2 與平面 4x+2y+5z=c，c 從 40 變化到 20 的動畫圖 第七節(jié) 數(shù)據(jù)擬合MATLAB有數(shù)據(jù)插值與數(shù)據(jù)擬合功能，使用的命令是：interpn(注意這里的“n”必須是一個具體正整數(shù)，一般為 n=1, 2, 3, ，表示)。例如：運行exmaple51.m，

33、對南半球氣旋作可視化顯示圖。%exmaple51.m-南半球氣旋數(shù)據(jù)，程序?qū)崿F(xiàn)氣旋數(shù)據(jù)可視化y=5:10:85;x=1:12;z= 2.4, 1.6, 2.4, 3.2, 1.0, 0.5, 0.4, 0.2, 0.5, 0.8, 2.4, 3.6; 18.7, 21.4, 16.2, 9.2, 2.8, 1.7, 1.4, 2.4, 5.8, 9.2, 10.3, 16; 20.8, 18.5, 18.2, 16.6, 12.9, 10.1, 8.3, 11.2, 12.5, 21.1, 23.9 ,25.5;22.1, 20.1, 20.5, 25.1, 29.2, 32.6, 33.0,

34、 31.0, 28.6, 32.0, 28.1, 25.6; 37.3, 28.8, 27.8, 37.2, 40.3, 41.7, 46.2, 39.9, 35.9, 40.3, 38.2, 43.4; 48.2, 36.6, 35.5, 40, 37.6, 35.4, 35, 34.7, 35.7, 39.5, 40, 41.9; 25.6, 24.2, 25.5, 24.6, 21.1, 22.2, 20.2, 21.2, 22.6, 28.5, 25.3, 24.3; 5.3, 5.3, 5.4, 4.9, 4.9, 7.1, 5.3, 7.3, 7, 8.6, 6.3, 6.6; 0

35、.3, 0, 0, 0.3, 0, 0, 0.1, 0.2, 0.3, 0, 0.1, 0.3; xi,yi=meshgrid(1:12,5:85); zi=interp2(x,y,z,xi,yi,cubic); surf(xi,yi,zi)練習9：一天(y tin)24小時每兩小時(xiosh)的溫度記錄如下：9 9 10 18 24 28 27 25 20 18 15 13請用插值方法(fngf)用出一天24小時的溫度變化圖。用命令 curvefit(fun,x0 , x , y) 可以作曲線擬合，其擬合功能比 MATHMATICA 更強。如某數(shù)據(jù)近似符合規(guī)律y(t)=a+be-kt，已知

36、從100到1000時分別對應(yīng)于cdata=0.00001*454 499 535 565 590 610 626 639 650 659;求 a, b, k 使 盡可能地小。此問題的解可由運行exmaple52.m得到。附錄： 內(nèi)插及曲線擬合我們常常會有從實驗或是對物理現(xiàn)象觀察所收集到的數(shù)據(jù)，分析這些數(shù)據(jù)需要有輔助的工具。通常(tngchng)這些原始 數(shù)據(jù)如果會制成圖，是成點狀分布。我們需要將這些點連成曲線，再進一步將曲線轉(zhuǎn)換成有意義的數(shù)學函數(shù) ，才能對數(shù)據(jù)做比較和分析。上述的過程涉及到內(nèi)插(interpolation)及曲線擬合(n h)(curve-fitting)，我們(w men)會

37、在以下介紹。1 內(nèi)插假設(shè)一組已知的數(shù)據(jù)其型態(tài)為 f(xk), =1, 2, , n, x1=a, x xn=c, 假設(shè)某些點 xi并不屬于上述的 xk但是 axic，如果我們要估計這些點的函數(shù)值 f(xi) 就須要做內(nèi)插 (interpolation)。我們可視原數(shù)據(jù)所描述的函數(shù)復雜程度來選擇不同的數(shù)值內(nèi)插方法。 1 內(nèi)插 1.1 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_1_1.htm 一維內(nèi)插 1.2 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_1_2.htm 二維內(nèi)插 1.3 HYPERLINK

38、 /kejian/matlabcomplex/chapter7/ch7_1_3.htm Spline 內(nèi)插 1.1 一維內(nèi)插線性內(nèi)插是假設(shè)在二個已知數(shù)據(jù)中的變化為線性關(guān)系，因此可由已知二點的座標(a, b)去計算通過這二點的 斜線，公式如下： 其中(qzhng) abc 在上式的 b 點即是代表(dibio)要內(nèi)插的點，f(b) 則是要計算的內(nèi)插函數(shù)值。下圖即是一個(y )以二種內(nèi)插法的比較 pcxfile12cm,5cmfig9_1.pcx caption線性式與 spline 函數(shù)的曲線擬合線性內(nèi)插是最簡單的內(nèi)插方法，但其適用范圍很?。蝗绻瓉頂?shù)據(jù)的函數(shù)f有極大的變化，假設(shè)其數(shù)據(jù)點之 間為

39、線性變化并不合理。所以我們可以用二次、三次方程組或是另一種稱為spline函數(shù)來近似原來數(shù)據(jù)的函 數(shù)。MATLAB的一維內(nèi)插函數(shù)是interp1，其語法為interp1(x,y,xi),interp1(x,y,xi,method)；其中的x,y是原已知的 數(shù)據(jù)的x,y值，而xi則是要內(nèi)插的數(shù)據(jù)點，另外method可以設(shè)定內(nèi)插方法有 linear,cubic,spline，分別是一次、三 次方程組和spline函數(shù)，其中預設(shè)方法是linear。如果數(shù)據(jù)的變化較大，以 spline函數(shù)內(nèi)插所形成的曲線最平滑 ，所以效果最好。而三次方程組所得到的內(nèi)插曲線平滑度，則介于線性與spline函數(shù)之間。我們

40、以下面的例子說明。假設(shè)有一個汽車發(fā)動機在定轉(zhuǎn)速下，溫度與時間（單位為sec）的三次量測值如下timetemp1 temp2temp3 00 00 120 110176 260 180220 368 240349 477 310450 5110 405503 其中(qzhng)溫度的數(shù)據(jù)從 20oC變化(binhu)到 503oC，如果要估計在t=2.6, 4.9 sec 的溫度，可以下列指令(zhlng)計算 x=0 1 2 3 4 5; % 鍵入時間 y=0 20 60 68 77 110; % 鍵入第一組時間 y1=interp1(x,y,2.6) % 要內(nèi)插的數(shù)據(jù)點為 2.6 y1 =

41、% 對應(yīng) 2.6 的函數(shù)值為 64.8 64.8 y1=interp1(x,y,2.6 4.9) %內(nèi)插數(shù)據(jù)點為 2.6, 4.9，注意用 將多個內(nèi)插點放在其中 y1 = 64.8 106.7 y1=interp1(x,y,2.6,cubic) % 以三次方程組對數(shù)據(jù)點 2.6 作內(nèi)插 y1 = % 對應(yīng) 2.6 的函數(shù)值為 66.264 66.264 y1=interp1(x,y,2.6,spline) % 以spline函數(shù)對數(shù)據(jù)點 2.6 作內(nèi)插 y1 = % 對應(yīng) 2.6 的函數(shù)值為 66.368 66.368以下的例子還配合繪圖功能，用以比較不同(b tn)內(nèi)插方法的差異。 h=1:

42、12; temp=5 8 9 15 25 29 31 30 22 25 27 24; % 這組溫度數(shù)據(jù)(shj)變化較大 plot(h,temp,-,h,temp,+) % 將線性內(nèi)插結(jié)果(ji gu)繪圖 h_3=1:0.1:12 % 要每0.1小時估計一次溫度值 t_3=interp1(h,temp,h_3,cubic) % 以三次方程組做內(nèi)插 t_s=interp1(h,temp,h_3,spline) % 以spline函數(shù)做內(nèi)插 hold on subplot(1,2,1) plot(h,temp,-,h,temp,+,h_3,t_3) % 將線性及三次方程組內(nèi)插繪圖 subplot

43、(1,2,2) plot(h,temp,-,h,temp,+,h_3,t_s) % 將線性方程組及spline內(nèi)插繪圖 hold off7.1.2 二維內(nèi)插二維內(nèi)插與一維內(nèi)插的區(qū)別(qbi)是二維內(nèi)插數(shù)據(jù)為二維,語法結(jié)構(gòu)為interp2(X,Y,Z,XI,YI),其中X,Y,Z為已知數(shù)據(jù),Z=Z(X,Y),而XI,YI 為要插值的數(shù)據(jù)(shj)點;如果語法結(jié)構(gòu)為interp2(X,Y,Z,XI,YI,method),其中method可以為linear,cubic表示線形或三次方插值,我們以下例說明:假設(shè)(jish)一汽車的轉(zhuǎn)速(單位為:rpm)、溫度（單位為：oC）、時間（單位為：sec）如下

44、表：time speed 0 2000 rpm 3000 rpm 4000 rpm 1 20 110 176 2 60 180 220 3 68 240 349 4 77 310 450 5 110 405 503 其中溫度的數(shù)據(jù)為20oC到503oC，如果要估計t=2.6, sec, rpm=2500的溫度，可以利用下面的語句：d2(:,1)=0 1 2 3 4 5; % 將時間輸入 d2(:,2)=2000 20 60 68 77 110; % 將 rpm=2000的溫度輸入 d2(:,3)=3000 110 180 240 310 405; % 將 rpm=3000 的溫度輸入 d2(:

45、,4)=4000 176 220 349 450 503; % 將 rpm=4000 的溫度輸入 t=d2(2:6,1); %選擇做內(nèi)插的時間rpm=d2(1,2:4); % 選擇做內(nèi)插的 rpm temp=d2(2:6,2:4); % 選擇做內(nèi)插的溫度temp_i=interp2(rpm,t,temp,2500,2.6) % 以線形內(nèi)插決定 rpm=2500,t=2.6 的溫度 temp_i = 140.4000 1.3 Spline 內(nèi)插關(guān)于spline內(nèi)插我們(w men)在1.1節(jié)已介紹過，它可以用interp1指定內(nèi)插方式為spline來做。另一種方式也可以用 spline(x,y,

46、xi)來做，其中的x,y,xi的用法與interp1中的語法相同。事實上這二種方法采用相同的spline 函數(shù)做 運算，也就是當我們執(zhí)行interp1(x,y,xi,spline)時，MATLAB即呼叫spline(x,y,xi)做運算，再將計算結(jié)果傳回 interp1。 我們還是以在 1.1 節(jié)相同的數(shù)據(jù)(shj)做 spline 內(nèi)插示范， x=0 1 2 3 4 5; y=0 20 60 68 77 110; y1=spline(x,y,2.6) y1 = 67.3 y1=spline(x,y,2.6,4.9) y1 = 67.3 105.23 曲線(qxin)擬合曲線擬合 (curve

47、-fitting) 故名思義就是要將一組離散的數(shù)據(jù)以一個近似的曲線方程組來代表。有了這個解析形 態(tài)的方程組，我們就可以很方便去運用它。曲線擬合與前述的內(nèi)插有許多相似之處，但是這二者最大的區(qū) 別在于曲線擬合要找出一個曲線方程組而內(nèi)插僅是要求出內(nèi)插數(shù)值即可。 要找出近似一組數(shù)據(jù)（也就是最能代表(dibio)這些數(shù)據(jù)）的曲線方程組，有許多選擇(xunz)，從簡單的一階線性方程組 到高階的多項式都有。我們以下就介紹(jisho)二種方法：線性回歸(linear regression)，多項式回歸(polynomial regression)。3 曲線擬合 3.1 HYPERLINK /kejian/m

48、atlabcomplex/chapter7/ch7_3_1.htm 線性回歸 3.2 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_3_2.htm 多項式回歸 3.3 HYPERLINK /kejian/matlabcomplex/chapter7/ch7_3_3.htm 多項式擬合及函數(shù)計算 3.1 線性回歸我們以一簡單數(shù)據(jù)組來說明什么是線性回歸。假設(shè)有一組數(shù)據(jù)型態(tài)為 y=y(x)，其中 x=0, 1, 2, 3, 4, 5, y=0, 20, 60, 68, 77, 110 如果我們要以一個最簡單的方程組來近似這組數(shù)據(jù)，則非一階的線性方程組莫屬。

49、先將這組數(shù)據(jù)繪圖如下 plot(x, y)圖中的斜線是我們隨意假設(shè)一階線性方程組 y=20 x，用以代表這些數(shù)據(jù)的一個方程組。以下將上述繪圖的 MATLAB 指令列出，并計算這個線性方程組的 y 值與原數(shù)據(jù) y 值間誤差平方的總合。 x=0 1 2 3 4 5; y=0 20 60 68 77 110; y1=20*x; % 一階線性方程組的 y1 值 sum_sq = sum(y-y1).2); % 誤差平方總合為 573 axis(-1,6,-20,120) plot(x,y1,x,y,o), title(Linear estimate), grid如此(rc)任意的假設(shè)一個線性方程組并無

50、根據(jù)(gnj)，如果換成其它人來設(shè)定就可能采用不同的線性方程組；所以我們 須要有比較精確(jngqu)方式?jīng)Q定理想的線性方程組。我們可以要求誤差平方的總合為最小，做為決定理想的線性方 程序的準則，這樣的方法就稱為最小平方誤差(least squares error)或是線性回歸。MATLAB的polyfit函數(shù)提供了 從一階到高階多項式的回歸法，其語法為polyfit(x,y,n)，其中x,y為輸入數(shù)據(jù)組n為多項式的階數(shù)，n=1就是一階 的線性回歸法。polyfit函數(shù)所建立的多項式可以寫成 從polyfit函數(shù)得到的輸出值就是上述的各項系數(shù)，以一階線性回歸為例n=1，所以只有 二個輸出值。如

51、果指令為coef=polyfit(x,y,n)，則coef(1)= , coef(2)= ,.,coef(n+1)= 。注意上式對n 階的多 項式會有 n+1 項的系數(shù)。我們來看以下的線性回歸的示范： x=0 1 2 3 4 5; y=0 20 60 68 77 110; coef=polyfit(x,y,1); % coef 代表線性回歸的二個輸出值 a0=coef(1); a1=coef(2); ybest=a1*x+a0; % 由線性回歸產(chǎn)生的一階方程組 sum_sq=sum(y-ybest).2); % 誤差(wch)平方總合為 356.82 axis(-1,6,-20,120) plot(x,ybest,x,y,o), ti