第五章數(shù)值計(jì)算

1、 具體知識(shí)點(diǎn)見(jiàn)導(dǎo)圖具體知識(shí)點(diǎn)見(jiàn)導(dǎo)圖5。5.1 矩陣代數(shù)的計(jì)算矩陣代數(shù)的計(jì)算 矩陣代數(shù)計(jì)算包括求矩陣的各種三角分解、線性方程組的求解、矩陣代數(shù)計(jì)算包括求矩陣的各種三角分解、線性方程組的求解、矩陣的逆和廣義逆、矩陣的模及條件數(shù)、矩陣的譜分解和奇異值分矩陣的逆和廣義逆、矩陣的模及條件數(shù)、矩陣的譜分解和奇異值分解。矩陣代數(shù)的計(jì)算是解決數(shù)學(xué)問(wèn)題的基礎(chǔ)性方法，任何的數(shù)學(xué)模解。矩陣代數(shù)的計(jì)算是解決數(shù)學(xué)問(wèn)題的基礎(chǔ)性方法，任何的數(shù)學(xué)模型不論是統(tǒng)計(jì)模型、微分方程模型還是其他的復(fù)雜數(shù)學(xué)模型最后大型不論是統(tǒng)計(jì)模型、微分方程模型還是其他的復(fù)雜數(shù)學(xué)模型最后大部分將進(jìn)行離散化處理，然后歸為矩陣代數(shù)的計(jì)算。部分將進(jìn)行離散化處

2、理，然后歸為矩陣代數(shù)的計(jì)算。5.1.1 矩陣代數(shù)知識(shí)點(diǎn)搜索矩陣代數(shù)知識(shí)點(diǎn)搜索 在參考文獻(xiàn)在參考文獻(xiàn)Reference中進(jìn)入中進(jìn)入MATLAB Function Reference，即進(jìn)入，即進(jìn)入MATLAB的函數(shù)參考列表，然后雙擊的函數(shù)參考列表，然后雙擊Mathematics，我們就進(jìn)入了，我們就進(jìn)入了MATLAB的核心數(shù)學(xué)方法的函數(shù)的核心數(shù)學(xué)方法的函數(shù)集（或稱(chēng)數(shù)學(xué)方法命令集）。集（或稱(chēng)數(shù)學(xué)方法命令集）。 代數(shù)計(jì)算函數(shù)代數(shù)計(jì)算函數(shù)【例【例5.1】MATLAB在數(shù)值計(jì)算方面的強(qiáng)大功能。在數(shù)值計(jì)算方面的強(qiáng)大功能。傳統(tǒng)求逆法和傳統(tǒng)求逆法和MATLAB左除法求線性方程的性能對(duì)比。左除法求線性方程的性能

3、對(duì)比。rand(state,12);%選定隨機(jī)種子。選定隨機(jī)種子。A=rand(100,100)+1.e8; %生成生成100100均勻分布隨機(jī)矩陣。均勻分布隨機(jī)矩陣。x=ones(100,1); %令解向量令解向量 x 為全為全1的的100元列向量。元列向量。b=A*x; %為使為使 Ax=b 方程一致，用方程一致，用A和和 x 生成生成 b 向量。向量。cond(A) %求求 A 陣的條件數(shù)。陣的條件數(shù)。 ans = 1.4426e+012 從矩陣的條件數(shù)來(lái)看，相當(dāng)大。說(shuō)明矩陣從矩陣的條件數(shù)來(lái)看，相當(dāng)大。說(shuō)明矩陣A是嚴(yán)重病態(tài)矩陣，求以是嚴(yán)重病態(tài)矩陣，求以它為系數(shù)矩陣的線性方程的根將是不精確

4、的，或者說(shuō)是不穩(wěn)定的。它為系數(shù)矩陣的線性方程的根將是不精確的，或者說(shuō)是不穩(wěn)定的。% “求逆求逆”法解恰定方程的誤差、殘差、運(yùn)算次數(shù)和所用時(shí)間法解恰定方程的誤差、殘差、運(yùn)算次數(shù)和所用時(shí)間flops(0);tic %計(jì)數(shù)器置計(jì)數(shù)器置0，啟動(dòng)計(jì)時(shí)器，啟動(dòng)計(jì)時(shí)器Stopwatch Timerxi=inv(A)*b; % xi 是用是用“求逆求逆”法解恰定方程所得的解。法解恰定方程所得的解。ti=toc %關(guān)閉計(jì)時(shí)器，并顯示解方程所用的時(shí)間。關(guān)閉計(jì)時(shí)器，并顯示解方程所用的時(shí)間。ci=flops%“求逆求逆”法解方程所用的運(yùn)算次數(shù)法解方程所用的運(yùn)算次數(shù)eri=norm(x-xi) %解向量解向量 xi 與

5、真解向量與真解向量 x 的范的范-2誤差。誤差。rei=norm(A*xi-b)/norm(b) %方程的范方程的范-2相對(duì)殘差相對(duì)殘差 ti = 0.9300 ； 用時(shí)用時(shí)0.93秒秒ci = 2070322 ； 乘除次數(shù)乘除次數(shù)eri = 3.0708e-004 絕對(duì)誤差絕對(duì)誤差rei = 6.6280e-007 相對(duì)誤差相對(duì)誤差% “左除左除”法解恰定方程的誤差、殘差、運(yùn)算次數(shù)和所用時(shí)間法解恰定方程的誤差、殘差、運(yùn)算次數(shù)和所用時(shí)間flops(0);tic;xd=Ab; % 以上是用以上是用“左除左除”法解方程所得的解。法解方程所得的解。td=toc,cd=flops,erd=norm(x

6、-xd),red=norm(A*xd-b)/norm(b) td = 0.2200 計(jì)算時(shí)間計(jì)算時(shí)間cd = 741872 乘除法的次數(shù)乘除法的次數(shù)erd = 3.2243e-004 絕對(duì)誤差絕對(duì)誤差red = 2.0095e-016 相對(duì)誤差相對(duì)誤差對(duì)比結(jié)果：用左除法是逆陣法運(yùn)算次數(shù)的對(duì)比結(jié)果：用左除法是逆陣法運(yùn)算次數(shù)的4/10，時(shí)間是逆陣除法，時(shí)間是逆陣除法的的1/4。左除法的相對(duì)殘差幾乎為。左除法的相對(duì)殘差幾乎為 “機(jī)器零機(jī)器零” ，比逆陣法相對(duì)精，比逆陣法相對(duì)精度要高的多。度要高的多。5.1.2 矩陣分析矩陣分析Matrix Analysis矩陣分析主要是考察一個(gè)系數(shù)矩陣矩陣分析主要是

7、考察一個(gè)系數(shù)矩陣A對(duì)線性方程組解的靈敏度，即對(duì)線性方程組解的靈敏度，即解的穩(wěn)定性問(wèn)題。例如對(duì)于最小二乘估計(jì)的解解的穩(wěn)定性問(wèn)題。例如對(duì)于最小二乘估計(jì)的解 如果協(xié)方差矩陣具有列相關(guān)性時(shí)，或者說(shuō)該矩陣接近不滿(mǎn)秩如果協(xié)方差矩陣具有列相關(guān)性時(shí)，或者說(shuō)該矩陣接近不滿(mǎn)秩時(shí)，我們稱(chēng)該矩陣為病態(tài)矩陣，它將嚴(yán)重影響最小二乘估計(jì)的穩(wěn)定時(shí)，我們稱(chēng)該矩陣為病態(tài)矩陣，它將嚴(yán)重影響最小二乘估計(jì)的穩(wěn)定性。度量其病態(tài)性的指標(biāo)為條件數(shù)，矩陣的條件數(shù)為該矩陣的最大性。度量其病態(tài)性的指標(biāo)為條件數(shù)，矩陣的條件數(shù)為該矩陣的最大特征值比最小特征值，用來(lái)衡量矩陣計(jì)算的精確度，當(dāng)條件數(shù)為特征值比最小特征值，用來(lái)衡量矩陣計(jì)算的精確度，當(dāng)條件數(shù)為

8、1時(shí)是很好的條件，而條件數(shù)非常大時(shí)則表明矩陣不滿(mǎn)秩，計(jì)算時(shí)將時(shí)是很好的條件，而條件數(shù)非常大時(shí)則表明矩陣不滿(mǎn)秩，計(jì)算時(shí)將產(chǎn)生很大的誤差。產(chǎn)生很大的誤差。矩陣分析包括矩陣分析包括:條件數(shù)的計(jì)算、矩陣模的計(jì)算、矩陣的行列式的計(jì)算和矩陣的標(biāo)準(zhǔn)條件數(shù)的計(jì)算、矩陣模的計(jì)算、矩陣的行列式的計(jì)算和矩陣的標(biāo)準(zhǔn)正交基的計(jì)算等等。正交基的計(jì)算等等。 yXXXb1 【例【例5.1.1】求矩陣】求矩陣A的條件數(shù)、行列式、模等的計(jì)算的條件數(shù)、行列式、模等的計(jì)算 7518244518AA=8,-1,-5；-4,4,-2；18,-5,7 % 定義一個(gè)矩陣定義一個(gè)矩陣c=cond(A) % 計(jì)算矩陣計(jì)算矩陣A的條件數(shù)的條件數(shù)c

9、 = 8.2546n = norm(A) % 計(jì)算矩陣的模計(jì)算矩陣的模n = 21.5249d = det(A) % 計(jì)算矩陣的行列式計(jì)算矩陣的行列式d = 412B = orth(A) % 求矩陣的一個(gè)標(biāo)準(zhǔn)正交基求矩陣的一個(gè)標(biāo)準(zhǔn)正交基 BB =B = -0.2987 -0.9498 -0.0933 -0.2987 -0.9498 -0.0933 0.2469 -0.1713 0.9538 0.2469 -0.1713 0.9538 -0.9219 0.2619 0.2856 -0.9219 0.2619 0.2856EYE=B*B % 驗(yàn)證B是標(biāo)準(zhǔn)正交基 EYE =EYE = 1.0000

10、-0.0000 0.0000 1.0000 -0.0000 0.0000 -0.0000 1.0000 0 -0.0000 1.0000 0 0.0000 0 1.0000 0.0000 0 1.00005.1.3 矩陣的分解和三角分解矩陣的分解和三角分解Factorization 從計(jì)算的角度來(lái)看，處理矩陣問(wèn)題的一個(gè)最有效的方法是，將從計(jì)算的角度來(lái)看，處理矩陣問(wèn)題的一個(gè)最有效的方法是，將一個(gè)一般的矩陣分解為幾個(gè)簡(jiǎn)單矩陣的乘積形式。上三角矩陣為例，一個(gè)一般的矩陣分解為幾個(gè)簡(jiǎn)單矩陣的乘積形式。上三角矩陣為例，由其性質(zhì)，兩個(gè)上三角陣的和、積仍為上三角陣，上三角陣的特征由其性質(zhì)，兩個(gè)上三角陣的和、積

11、仍為上三角陣，上三角陣的特征值就是其對(duì)角線元素。系數(shù)矩陣為上三角陣的線性方程組只需用回值就是其對(duì)角線元素。系數(shù)矩陣為上三角陣的線性方程組只需用回代的方法求解，上三角陣的逆陣仍然是上三角陣。代的方法求解，上三角陣的逆陣仍然是上三角陣。LU分解：分解： 將矩陣將矩陣A分解為分解為A = LU，其中，其中L為單位下三角陣，為單位下三角陣，U為上三角陣。為上三角陣。Cholesky分解：將正定對(duì)稱(chēng)陣分解成分解：將正定對(duì)稱(chēng)陣分解成 A=TT的形式，其中的形式，其中T為為正交變換且為上三角陣。正交變換且為上三角陣。QR分解：分解： 將矩陣分解為將矩陣分解為A=QR的形式，其中的形式，其中Q為一個(gè)正交陣，為

12、一個(gè)正交陣，R為一個(gè)上三角陣。為一個(gè)上三角陣?！纠纠?.1.2】對(duì)矩陣】對(duì)矩陣A進(jìn)行進(jìn)行LU分解分解L,U = lu(A) % 對(duì)矩陣對(duì)矩陣A進(jìn)行進(jìn)行LU分解分解A = 8 -1 -5 -4 4 -2 18 -5 -7L = 0.4444 0.4231 1.0000 -0.2222 1.0000 0 1.0000 0 0U = 18.0000 -5.0000 -7.0000 0 2.8889 -3.5556 0 0 -0.3846【例【例5.1.3】對(duì)對(duì)稱(chēng)正定矩陣進(jìn)行】對(duì)對(duì)稱(chēng)正定矩陣進(jìn)行Cholesky分解分解B = 1 -3 2 -3 10 -5 2 -5 6 T = chol(B) %

13、對(duì)對(duì)稱(chēng)矩陣進(jìn)行對(duì)對(duì)稱(chēng)矩陣進(jìn)行Cholesky分解分解C = T*T % 對(duì)對(duì)T進(jìn)行驗(yàn)證進(jìn)行驗(yàn)證B = 1 -3 2B = 1 -3 2 -3 10 -5 -3 10 -5 2 -5 6 2 -5 6T = 1 -3 2T = 1 -3 2 0 1 1 0 1 1 0 0 1 0 0 1C = 1 -3 2 C = 1 -3 2 -3 10 -5 -3 10 -5 2 -5 6 2 -5 6 5.1.4 求解線性方程組求解線性方程組Linear Equations【例【例5.1.4】求解以下線性方程組】求解以下線性方程組A=1,-3,2;-3,10,-5;2,-5,6 % 定義系數(shù)矩陣定義系數(shù)矩

14、陣Ab=27,-78,64 % 定義常數(shù)定義常數(shù)bx=Ab % 求解線性方程組求解線性方程組b1=A*x最后的解為：最后的解為：x=(1 -4 7) 6478276525103231321xxx5.1.5 矩陣的特征值或奇異值及其與特征向量矩陣的特征值或奇異值及其與特征向量Eigenvalues and Singular Values（1）對(duì)稱(chēng)矩陣的譜分解）對(duì)稱(chēng)矩陣的譜分解,即求矩陣的特征值和特征向量即求矩陣的特征值和特征向量設(shè)矩陣設(shè)矩陣A是是NN維方陣，則譜分解定義如下：維方陣，則譜分解定義如下：UUUDUAn1其中其中D為對(duì)角矩陣，為對(duì)角矩陣，U為正交矩陣。對(duì)角矩陣的元素稱(chēng)為矩陣為正交矩陣

15、。對(duì)角矩陣的元素稱(chēng)為矩陣A的的特征值，矩陣特征值，矩陣U的列向量為對(duì)應(yīng)于矩陣的列向量為對(duì)應(yīng)于矩陣A的特征值的特征向量。的特征值的特征向量。（2）對(duì)一般矩陣的奇異值分解）對(duì)一般矩陣的奇異值分解設(shè)矩陣設(shè)矩陣A為為nm階矩陣，則矩陣的奇異值分解定義為：階矩陣，則矩陣的奇異值分解定義為：mmmnnnmnVDUA其中其中 U 和和 V 分別為分別為 n 和和 m 階正交方陣，階正交方陣，D 為為 nm 階對(duì)角陣，階對(duì)角陣，其非對(duì)角元素均為其非對(duì)角元素均為 0，D的對(duì)角線元素稱(chēng)為矩陣的對(duì)角線元素稱(chēng)為矩陣 A 的奇異值。奇的奇異值。奇異值的分解與矩陣的譜分解方法類(lèi)似。異值的分解與矩陣的譜分解方法類(lèi)似。求解矩

16、陣求解矩陣A的特征值和特征向量的語(yǔ)法為：的特征值和特征向量的語(yǔ)法為：V, D = eig(A)其中其中V為矩陣其列向量為矩陣為矩陣其列向量為矩陣A的特征向量，的特征向量，D為對(duì)角矩陣對(duì)角線為對(duì)角矩陣對(duì)角線元素為矩陣元素為矩陣A的特征值。的特征值。 【例【例5.1.5】求矩陣的特征值和奇異值?！壳缶仃嚨奶卣髦岛推娈愔怠＃?）求上面矩陣）求上面矩陣A的特征值和特征向量。的特征值和特征向量。V,D = eig(A) % 求矩陣求矩陣A的特征值和特征向量的特征值和特征向量V = 0.9654 -0.0091 0.2606 0.2211 0.5581 -0.7998 -0.1381 0.8297 0.5

17、408D = 0.0266 0 0 0 2.6151 0 0 0 14.3583（2）求下列矩陣）求下列矩陣B的奇異值及其對(duì)應(yīng)的向量。的奇異值及其對(duì)應(yīng)的向量。B =9 4；6 8； 2 7 % 矩陣矩陣B為長(zhǎng)方形矩陣為長(zhǎng)方形矩陣U,S,V = svd(B) % 求矩陣求矩陣B的奇異值和對(duì)應(yīng)的兩個(gè)向量的奇異值和對(duì)應(yīng)的兩個(gè)向量B = 9 4 6 8 2 7U = -0.6105 0.7174 0.3355 -0.6646 -0.2336 -0.7098 -0.4308 -0.6563 0.6194S = 14.9359 0 0 5.1883 0 0V = -0.6925 0.7214 -0.721

18、4 -0.69255.2 多項(xiàng)式與插值多項(xiàng)式與插值 在數(shù)學(xué)建模競(jìng)賽中或在科學(xué)研究中，使用頻率最高的數(shù)據(jù)加工在數(shù)學(xué)建模競(jìng)賽中或在科學(xué)研究中，使用頻率最高的數(shù)據(jù)加工技術(shù)是對(duì)實(shí)驗(yàn)數(shù)據(jù)的擬合，擬合有很多方法如多項(xiàng)式方法、最技術(shù)是對(duì)實(shí)驗(yàn)數(shù)據(jù)的擬合，擬合有很多方法如多項(xiàng)式方法、最小二乘法和樣條插法。當(dāng)我們對(duì)給定的數(shù)據(jù)利用上面的方法選小二乘法和樣條插法。當(dāng)我們對(duì)給定的數(shù)據(jù)利用上面的方法選擇了擬合方法后，對(duì)任意給定的自變量點(diǎn)（插值點(diǎn)）我們可以擇了擬合方法后，對(duì)任意給定的自變量點(diǎn)（插值點(diǎn)）我們可以計(jì)算該點(diǎn)的函數(shù)值了。計(jì)算該點(diǎn)的函數(shù)值了。MATLAB約定多項(xiàng)式如下：約定多項(xiàng)式如下： 1121)(nnnnaxaxa

19、xaxP 我們可以用行向量來(lái)表示一個(gè)多項(xiàng)式。如，我們可以用行向量來(lái)表示一個(gè)多項(xiàng)式。如，。我們可以進(jìn)行多項(xiàng)式的乘、除、因子的提取、多項(xiàng)式的簡(jiǎn)化、求。我們可以進(jìn)行多項(xiàng)式的乘、除、因子的提取、多項(xiàng)式的簡(jiǎn)化、求多項(xiàng)式的根等等常規(guī)的多項(xiàng)式的計(jì)算和多項(xiàng)式的擬合。多項(xiàng)式的根等等常規(guī)的多項(xiàng)式的計(jì)算和多項(xiàng)式的擬合。注意：如果某系數(shù)為注意：如果某系數(shù)為0，必須寫(xiě)上以占一個(gè)位置。，必須寫(xiě)上以占一個(gè)位置。,121 naaaP5.2.1 多項(xiàng)式的創(chuàng)建方法多項(xiàng)式的創(chuàng)建方法（1）直接用行向量創(chuàng)建，如：）直接用行向量創(chuàng)建，如：P=1，2，3，5，6，7 %這創(chuàng)建了一個(gè)這創(chuàng)建了一個(gè)5階多項(xiàng)式。階多項(xiàng)式。 （2）利用命令：）利用

20、命令：P = ploy（A） % 產(chǎn)生多項(xiàng)式系數(shù)向量。產(chǎn)生多項(xiàng)式系數(shù)向量。這里這里A為方陣，即多項(xiàng)式為方陣，即多項(xiàng)式 P 為為 A 的特征多項(xiàng)式。的特征多項(xiàng)式?！纠纠?.2.1】求三階魔方方陣】求三階魔方方陣M 的特征多項(xiàng)式。的特征多項(xiàng)式。M=magic(3)PM=poly(M)； % A的特征多項(xiàng)式的特征多項(xiàng)式PPM=poly2str(PM,x) % 將多項(xiàng)式轉(zhuǎn)成習(xí)慣的方式將多項(xiàng)式轉(zhuǎn)成習(xí)慣的方式RM=roots(PM) % 求求A或者說(shuō)特征多項(xiàng)式的根或者說(shuō)特征多項(xiàng)式的根PPA = x3 - 15 x2 - 24 x + 360RM = 15.0000 -4.8990 4.89907632)

21、(345 xxxxxP5.2.2 計(jì)算多項(xiàng)式的操作函數(shù)與例題計(jì)算多項(xiàng)式的操作函數(shù)與例題 多項(xiàng)式的運(yùn)算包括，多項(xiàng)式相乘、多項(xiàng)式相除、求多項(xiàng)式的公多項(xiàng)式的運(yùn)算包括，多項(xiàng)式相乘、多項(xiàng)式相除、求多項(xiàng)式的公因式等等運(yùn)算，表因式等等運(yùn)算，表5.2.1 列出一些基本多項(xiàng)式運(yùn)算命令列出一些基本多項(xiàng)式運(yùn)算命令【例【例5.2.2】由給定的根向量求多項(xiàng)式的系數(shù)向量?！坑山o定的根向量求多項(xiàng)式的系數(shù)向量。R=-0.5,-0.3+0.4*i,-0.3-0.4*i; % 根向量根向量P=poly(R) % R的特征多項(xiàng)式的特征多項(xiàng)式PR=real(P) % 求求PR的實(shí)部的實(shí)部PPR=poly2str(PR,x) 計(jì)算結(jié)果

22、如下：計(jì)算結(jié)果如下：P = 1.0000 1.1000 0.5500 0.1250PR = 1.0000 1.1000 0.5500 0.1250PPR = x3 + 1.1 x2 + 0.55 x + 0.125 注意：注意：（1）要形成實(shí)系數(shù)多項(xiàng)式，根向量中的復(fù)數(shù)必須共軛成對(duì)。）要形成實(shí)系數(shù)多項(xiàng)式，根向量中的復(fù)數(shù)必須共軛成對(duì)。（2）含復(fù)數(shù)的根向量產(chǎn)生的多項(xiàng)式系數(shù)向量）含復(fù)數(shù)的根向量產(chǎn)生的多項(xiàng)式系數(shù)向量P 有可能存在截有可能存在截 斷誤差級(jí)的虛部。我們可用斷誤差級(jí)的虛部。我們可用real命令把虛數(shù)過(guò)濾掉。命令把虛數(shù)過(guò)濾掉。【例【例5.2.3】 計(jì)算有理式計(jì)算有理式 的的“商商”及及“余余”多

23、項(xiàng)式多項(xiàng)式p1=conv(1,0,2,conv(1,4,1,1); %計(jì)算分子多項(xiàng)式計(jì)算分子多項(xiàng)式p2=1 0 1 1;%注意缺項(xiàng)補(bǔ)零注意缺項(xiàng)補(bǔ)零q,r=deconv(p1,p2);cq=商多項(xiàng)式為商多項(xiàng)式為 ; cr=余多項(xiàng)式為余多項(xiàng)式為 ;disp(cq,poly2str(q,s)disp(cr,poly2str(r,s) 計(jì)算結(jié)果為：計(jì)算結(jié)果為：商多項(xiàng)式為商多項(xiàng)式為 s + 5余多項(xiàng)式為余多項(xiàng)式為 5 s2 + 4 s + 3 114232 sssss【例【例5.2.4】我們對(duì)中國(guó)】我們對(duì)中國(guó)1952年到年到1998年的國(guó)民總收入年的國(guó)民總收入GDP的的數(shù)據(jù)，用數(shù)據(jù)，用1、2、3和和4階

24、多項(xiàng)式來(lái)擬合，并將結(jié)果用圖形表示出來(lái)。階多項(xiàng)式來(lái)擬合，并將結(jié)果用圖形表示出來(lái)。GDP的數(shù)據(jù)文件為的數(shù)據(jù)文件為gdp.txtclear; AA=load(e:datagdp.txt);x=AA(:,1); y=AA(:,2);p1 = polyfit(x,y,1) % 用一階多項(xiàng)式擬合用一階多項(xiàng)式擬合p2 = polyfit(x,y,2) % 用二階多項(xiàng)式擬合用二階多項(xiàng)式擬合p3 = polyfit(x,y,4) % 用四階多項(xiàng)式擬合用四階多項(xiàng)式擬合p4 = polyfit(x,y,6) % 用六階多項(xiàng)式擬合用六階多項(xiàng)式擬合f1 = polyval(p1,x); % 計(jì)算一階多項(xiàng)式的插值點(diǎn)計(jì)算一

25、階多項(xiàng)式的插值點(diǎn)f2 = polyval(p2,x); % 計(jì)算二階多項(xiàng)式的插值點(diǎn)計(jì)算二階多項(xiàng)式的插值點(diǎn)f3 = polyval(p4,x); % 計(jì)算三階多項(xiàng)式的插值點(diǎn)計(jì)算三階多項(xiàng)式的插值點(diǎn)f4 = polyval(p4,x); % 計(jì)算四階多項(xiàng)式的插值點(diǎn)計(jì)算四階多項(xiàng)式的插值點(diǎn)subplot(2,2,1), plot(x,y,.,x,f1,-),title(一階圖一階圖)subplot(2,2,2), plot(x,y,.,x,f2,-),title(二階圖二階圖)subplot(2,2,3), plot(x,y,.,x,f3,-),title(三階圖三階圖)subplot(2,2,4),

26、plot(x,y,.,x,f4,-),title(四階圖四階圖)從結(jié)果圖上看，四階多項(xiàng)式和六階多項(xiàng)式擬合的好。我們可以取四從結(jié)果圖上看，四階多項(xiàng)式和六階多項(xiàng)式擬合的好。我們可以取四階多項(xiàng)式，要注意的是，并不是階數(shù)越多越好。階多項(xiàng)式，要注意的是，并不是階數(shù)越多越好。5.2.4 插值計(jì)算插值計(jì)算 插值與多項(xiàng)式擬合的不同之處在于，對(duì)給定的點(diǎn)多項(xiàng)式擬合不插值與多項(xiàng)式擬合的不同之處在于，對(duì)給定的點(diǎn)多項(xiàng)式擬合不必點(diǎn)點(diǎn)通過(guò)，而插值建立的函數(shù)是點(diǎn)點(diǎn)通過(guò)給定的數(shù)據(jù)點(diǎn)。必點(diǎn)點(diǎn)通過(guò)，而插值建立的函數(shù)是點(diǎn)點(diǎn)通過(guò)給定的數(shù)據(jù)點(diǎn)。MATLAB提供了豐富的插值方法，并配有樣條插值工具箱。提供了豐富的插值方法，并配有樣條插值工

27、具箱。 5.2.5 一維插值多項(xiàng)式一維插值多項(xiàng)式一維插值多項(xiàng)式的基本語(yǔ)法為：一維插值多項(xiàng)式的基本語(yǔ)法為：yi = interp1(x,y,xi,yi,method)這里這里x，y 為給定的數(shù)據(jù)，利用該數(shù)據(jù)建立一維插值多項(xiàng)式，為給定的數(shù)據(jù)，利用該數(shù)據(jù)建立一維插值多項(xiàng)式，xi為為待計(jì)算的插值點(diǎn)，待計(jì)算的插值點(diǎn)，yi為插值點(diǎn)的函數(shù)值。為插值點(diǎn)的函數(shù)值。method：為建立插值多項(xiàng)式的方法，可?。簽榻⒉逯刀囗?xiàng)式的方法，可取method = nearest：利用最近臨方法建立插值多項(xiàng)式：利用最近臨方法建立插值多項(xiàng)式method = linear： 利用線性方法建立插值多項(xiàng)式利用線性方法建立插值多項(xiàng)式m

28、ethod = spline： 利用樣條函數(shù)方法建立插值多項(xiàng)式利用樣條函數(shù)方法建立插值多項(xiàng)式method = pchip： 三階厄米特方法建立逐段樣條插值多項(xiàng)式三階厄米特方法建立逐段樣條插值多項(xiàng)式method = cubic： 利用三階樣條函數(shù)建立插值多項(xiàng)式利用三階樣條函數(shù)建立插值多項(xiàng)式【例【例5.2.5】對(duì)國(guó)民經(jīng)濟(jì)數(shù)據(jù)從】對(duì)國(guó)民經(jīng)濟(jì)數(shù)據(jù)從1955年以年以5年為間隔建立插值函年為間隔建立插值函數(shù)，并在求數(shù)，并在求1955年、年、1986年和預(yù)測(cè)年和預(yù)測(cè)1996年的數(shù)據(jù)。年的數(shù)據(jù)。B=ones(3,2)A=load(e:datazggdp.txt)x=A(:,1)y=A(:,2)ly=log(y

29、)x1=x(4,9,14,19,24,29,34,39,44)y1=y(4,9,14,19,24,29,34,39,44)xx=1955,1986,1996yy=interp1(x,y,xx,spline)B(:,1)=xx;B(:,2)=yy;B結(jié)果為：結(jié)果為：1955 910 1966 10201 1996 66851【例【例5.2.6】 在極坐標(biāo)下插值的例子，給定橢圓上的四個(gè)點(diǎn)，用在極坐標(biāo)下插值的例子，給定橢圓上的四個(gè)點(diǎn)，用樣條插值并作圖。樣條插值并作圖。theta=0:0.5:2*pi; % 產(chǎn)生自變量四個(gè)樣點(diǎn)產(chǎn)生自變量四個(gè)樣點(diǎn)y=-0.5 1 -0.5 -1 0.5 1 -0.5 0

30、.5 1 0.5 -1 -0.5 1 0.5; % 給出相應(yīng)的因變量的點(diǎn)和初值給出相應(yīng)的因變量的點(diǎn)和初值theta2=linspace(theta(1),theta(end),100); % 參量稠密化參量稠密化yy=spline(theta,y,theta2); %求稠密點(diǎn)上的插值求稠密點(diǎn)上的插值plot(yy(1,:),yy(2,:),b,y(1,:),y(2,:),or)5.2.6 高維插值高維插值我們以二維插值為例，其語(yǔ)法為：我們以二維插值為例，其語(yǔ)法為：ZI = interp2 (X, Y, Z, XI,YI, method)其中的參數(shù)和一維的類(lèi)似。其中的參數(shù)和一維的類(lèi)似。【例【例5

31、.2.7】以】以peaks函數(shù)為例，首先取較少的數(shù)據(jù)。然后對(duì)其求函數(shù)為例，首先取較少的數(shù)據(jù)。然后對(duì)其求插值函數(shù)，并比較各方法的差異。插值函數(shù)，并比較各方法的差異。x,y = meshgrid(-3:1:3);z = peaks(x,y);surf(x,y,z);title(較少點(diǎn)的表面圖較少點(diǎn)的表面圖) xi,yi = meshgrid(-3:0.25:3);zi1 = interp2(x,y,z,xi,yi,nearest);subplot(2,2,1),surf(xi,yi,zi1),title(nearest方法方法)zi2 = interp2(x,y,z,xi,yi,bilinear)

32、;subplot(2,2,2),surf(xi,yi,zi2),title(bilinear方法方法)zi3 = interp2(x,y,z,xi,yi,cubic);subplot(2,2,3),surf(xi,yi,zi3),title(bicubict方法方法)zi4 = interp2(x,y,z,xi,yi);subplot(2,2,4),surf(xi,yi,zi4),title(系統(tǒng)內(nèi)定方法系統(tǒng)內(nèi)定方法)5.3 微積分微積分 在大學(xué)生數(shù)學(xué)建模競(jìng)賽中微積分知識(shí)出現(xiàn)的比率是很高的，也在大學(xué)生數(shù)學(xué)建模競(jìng)賽中微積分知識(shí)出現(xiàn)的比率是很高的，也是很基礎(chǔ)知識(shí)。比如是很基礎(chǔ)知識(shí)。比如2002年大

33、學(xué)生數(shù)模競(jìng)賽中，年大學(xué)生數(shù)模競(jìng)賽中，A題就是一個(gè)題就是一個(gè)有關(guān)微積分知識(shí)的應(yīng)用題，即有關(guān)微積分知識(shí)的應(yīng)用題，即車(chē)燈線光源的優(yōu)化設(shè)計(jì)車(chē)燈線光源的優(yōu)化設(shè)計(jì)的問(wèn)的問(wèn)題。我們知道車(chē)燈是一個(gè)拋物線型結(jié)構(gòu)，光源在旋轉(zhuǎn)拋物面的題。我們知道車(chē)燈是一個(gè)拋物線型結(jié)構(gòu)，光源在旋轉(zhuǎn)拋物面的焦點(diǎn)處且形狀為小長(zhǎng)條形。我們要研究光源條的長(zhǎng)短對(duì)燈光照焦點(diǎn)處且形狀為小長(zhǎng)條形。我們要研究光源條的長(zhǎng)短對(duì)燈光照射的影響等等。射的影響等等。5.3.1 函數(shù)的數(shù)值導(dǎo)數(shù)與切平面函數(shù)的數(shù)值導(dǎo)數(shù)與切平面函數(shù)曲面的法線函數(shù)曲面的法線 用用i，j，k 表示直角坐標(biāo)系的單位向量，曲面表示直角坐標(biāo)系的單位向量，曲面在（在（x0, y0）處的法線（）處

34、的法線（Normal）n 是是在在 MATLAB 中計(jì)算與繪制法線的命令為：中計(jì)算與繪制法線的命令為：NX，NY，NZ = surfnorm（X，Y，Z） 122 yfxfkjyfixfn在在 MATLAB 中計(jì)算與繪制法線的命令為：中計(jì)算與繪制法線的命令為：NX，NY，NZ = surfnorm（X，Y，Z） 【例【例5.3.1】曲面法線演示?！壳娣ň€演示。y=-1:0.1:1;x=2*cos(asin(y); % 旋轉(zhuǎn)曲面的旋轉(zhuǎn)曲面的“母線母線”X,Y,Z=cylinder(x,20); % 形成旋轉(zhuǎn)曲面形成旋轉(zhuǎn)曲面surfnorm(X(:,11:21),Y(:,11:21),Z(:,1

35、1:21); %在曲面上畫(huà)法線在曲面上畫(huà)法線view(120,18) %控制觀察角控制觀察角 5.3.2 偏導(dǎo)數(shù)和梯度偏導(dǎo)數(shù)和梯度 函數(shù)函數(shù) 的偏導(dǎo)數(shù)定義如下：的偏導(dǎo)數(shù)定義如下：),(yxfxyxfyxxfyxfxx ),(),(lim),(0全微分可借助梯度（全微分可借助梯度（Gradient）定義為：）定義為： dydxfdyyfdxxfyxdf),(數(shù)值差分命令的基本語(yǔ)法：數(shù)值差分命令的基本語(yǔ)法：DX=diff（X） 求求 X 相鄰元素的一階差分相鄰元素的一階差分DX=diff（X，n） 求求X相鄰元素的相鄰元素的 n 階差分階差分DX=diff（X，n，dim）命令維）命令維 dim

36、上上, 求相鄰元素的求相鄰元素的n階差分階差分差分的含義為：差分的含義為：設(shè)序列設(shè)序列則一階差分的定義為：則一階差分的定義為：n階差分依次類(lèi)推。階差分依次類(lèi)推。 nxxxX,21 12312, nnxxxxxxX【例【例5.3.2】差分方法應(yīng)用十分廣泛，例如在經(jīng)濟(jì)數(shù)據(jù)中一階差分表示發(fā)展速度，】差分方法應(yīng)用十分廣泛，例如在經(jīng)濟(jì)數(shù)據(jù)中一階差分表示發(fā)展速度，二階差分表示經(jīng)濟(jì)發(fā)展的加速度。對(duì)我國(guó)二階差分表示經(jīng)濟(jì)發(fā)展的加速度。對(duì)我國(guó)5298年國(guó)民經(jīng)濟(jì)數(shù)據(jù)年國(guó)民經(jīng)濟(jì)數(shù)據(jù)GDP發(fā)展趨勢(shì)發(fā)展趨勢(shì)進(jìn)行分析，求其一、二階差分。并作圖。進(jìn)行分析，求其一、二階差分。并作圖。A=load(e:datagdp.txt)B

37、=A(:,2)B1=diff(B);B2=diff(B,2)subplot(1,3,1),plot(B),title(發(fā)展趨勢(shì)發(fā)展趨勢(shì))subplot(1,3,2),plot(B1),title(發(fā)展速度發(fā)展速度)subplot(1,3,3),plot(B2),title(發(fā)展加速度發(fā)展加速度)梯度梯度:基本語(yǔ)法：基本語(yǔ)法：FX，F(xiàn)Y=gradient（F，h） 求二元函數(shù)的梯度求二元函數(shù)的梯度FX，F(xiàn)Y，F(xiàn)Z， = gradient（F，h1, h2, ） 【例【例5.3.3】 對(duì)函數(shù)對(duì)函數(shù) 作梯度圖作梯度圖 v = -2:0.2:2; x,y = meshgrid(v);z = x .*

38、exp(-x.2 - y.2);subplot(1,2,1);surf(x,y,z),subplot(1,2,2); px,py = gradient(z,.2,.2);contour(x,y,z), hold on, quiver(x,y,px,py) 22yxxe 方向法線的繪圖命令：基本語(yǔ)法為：方向法線的繪圖命令：基本語(yǔ)法為：quiver3(Z,U,V,W)quiver3(X,Y,Z,U,V,W)quiver3(.,scale)quiver3(.,LineSpec)這里：這里：X，Y，Z： 為曲線上的點(diǎn)。為曲線上的點(diǎn)。 U，V，W： 為每一點(diǎn)上切面的法線方向。為每一點(diǎn)上切面的法線方向。

39、scale： 為實(shí)數(shù)值，可正可負(fù)。表示方向箭頭的指向。為實(shí)數(shù)值，可正可負(fù)。表示方向箭頭的指向?！纠纠?.3.4】求函數(shù)的曲面圖，且作出每點(diǎn)的切面法向量圖形。】求函數(shù)的曲面圖，且作出每點(diǎn)的切面法向量圖形。X,Y = meshgrid(-2:0.25:2,-1:0.2:1); % 產(chǎn)生網(wǎng)格點(diǎn)產(chǎn)生網(wǎng)格點(diǎn)Z = X.* exp(-X.2 - Y.2); % 計(jì)算每個(gè)網(wǎng)格點(diǎn)的函數(shù)值計(jì)算每個(gè)網(wǎng)格點(diǎn)的函數(shù)值U,V,W = surfnorm(X,Y,Z); % 計(jì)算每個(gè)點(diǎn)的法向量分量計(jì)算每個(gè)點(diǎn)的法向量分量quiver3(X,Y,Z,U,V,W,0.5); % 作法向量的圖形作法向量的圖形hold onsur

40、f(X,Y,Z);colormap hsvview(-35,45)axis (-2 2 -1 1 -.6 .6)hold off【例【例5.3.5】綜合應(yīng)用。我們來(lái)作一個(gè)旋轉(zhuǎn)拋物面上某點(diǎn)的切面和】綜合應(yīng)用。我們來(lái)作一個(gè)旋轉(zhuǎn)拋物面上某點(diǎn)的切面和該切面的法線。這里用到了作方向法線的命令。該切面的法線。這里用到了作方向法線的命令。x=0:0.1:1;y=sqrt(2*x); % 旋轉(zhuǎn)曲面的旋轉(zhuǎn)曲面的“母線母線”X,Y,Z=cylinder(y,20); % 形成旋轉(zhuǎn)曲面形成旋轉(zhuǎn)曲面U,V,W = surfnorm(X,Y,Z);Z1=Z(4)-(U(4)*(X-X(4)+V(4)*(Y-Y(4)/W

41、(4) % 求點(diǎn)求點(diǎn)X(4),Y(4),Z(4)的切平面的切平面surf(X,Y,Z); % 在曲面上畫(huà)法線在曲面上畫(huà)法線hold on;surf(X,Y,Z1)；hold on 在上面點(diǎn)處畫(huà)切平面在上面點(diǎn)處畫(huà)切平面plot3(X(4),Y(4),Z(4),r.,MarkerSize,30) % 將該點(diǎn)突出放大并用紅色表示將該點(diǎn)突出放大并用紅色表示U(1:3)=NaN;U(5:end)=NaNhold on;quiver3(X,Y,Z,U,V,W,-4) % 作該點(diǎn)處的法線，長(zhǎng)度為作該點(diǎn)處的法線，長(zhǎng)度為4,反方向反方向axis equal % 該命令是必須的，使得視覺(jué)正確該命令是必須的，使得視

42、覺(jué)正確5.3.3 數(shù)值積分?jǐn)?shù)值積分 對(duì)于形如的積分對(duì)于形如的積分積分的一般命令為：積分的一般命令為：q = quad(fun,a,b)q = quad(fun,a,b,tol)q = quad(fun,a,b,tol,trace)q = quad(fun,a,b,tol,trace,p1,p2,.)q,fcnt = quadl(fun,a,b,.)其中其中fun表示被積函數(shù)，系數(shù)表示被積函數(shù)，系數(shù)a，b為積分的上下限，為積分的上下限，tol為用戶(hù)給為用戶(hù)給定的誤差限，如不給出則系統(tǒng)內(nèi)定的值為定的誤差限，如不給出則系統(tǒng)內(nèi)定的值為0.0000001?！纠纠?.3.6】以下是一個(gè)比較著名的例子，計(jì)

43、算積分】以下是一個(gè)比較著名的例子，計(jì)算積分Q = quad(1./(1+25*x.2),-1,1) ans = 0.5494 1122511dxx也可以使用內(nèi)聯(lián)函數(shù)定義被積函數(shù)也可以使用內(nèi)聯(lián)函數(shù)定義被積函數(shù)F = inline(1./(1+25*x.2); % 定義一個(gè)內(nèi)聯(lián)函數(shù)定義一個(gè)內(nèi)聯(lián)函數(shù)Q = quad(F,-1，1); ans = 0.5494【例【例5.3.7】對(duì)空間曲線求積分，設(shè)空間曲線的參數(shù)方程為：】對(duì)空間曲線求積分，設(shè)空間曲線的參數(shù)方程為：x = t.*sin(t), y = t.*cos(t), z=t t = 0:pi/50:10*pi;plot3(t.*sin(t),t.

44、*cos(t),t,LineWidth,2)f = inline(sqrt(4*(t.*sin(t).2 + (t.*cos(t).2 + 1);len = quad(f,0,pi)計(jì)算結(jié)果為：計(jì)算結(jié)果為：len = 8.51915.3.4 二重積分二重積分二重積分的原理與一重積分一樣，使用也很方便。其語(yǔ)法為：二重積分的原理與一重積分一樣，使用也很方便。其語(yǔ)法為：q = dblquad(fun,xmin,xmax,ymin,ymax)q = dblquad(fun,xmin,xmax,ymin,ymax,tol)q = dblquad(fun,xmin,xmax,ymin,ymax,tol,m

45、ethod)【例【例3.3.8】求積分】求積分q=dblquad(exp(-x.2.-y2),-1,1,-1,1) q = 2.2310 1111)(22dxeyx5.4 求函數(shù)的根求函數(shù)的根 函數(shù)求根或稱(chēng)求函數(shù)的根，是初等數(shù)學(xué)就學(xué)習(xí)過(guò)的。然而對(duì)函數(shù)求根或稱(chēng)求函數(shù)的根，是初等數(shù)學(xué)就學(xué)習(xí)過(guò)的。然而對(duì)于復(fù)雜的函數(shù)、方程組求根的運(yùn)算是用搜索的方法根據(jù)函數(shù)的于復(fù)雜的函數(shù)、方程組求根的運(yùn)算是用搜索的方法根據(jù)函數(shù)的特性逐步接近某個(gè)根，在使用特性逐步接近某個(gè)根，在使用MATLAB命令時(shí)應(yīng)該對(duì)函數(shù)有命令時(shí)應(yīng)該對(duì)函數(shù)有所了解，給出某個(gè)根附近的初始點(diǎn)或根所在的某一個(gè)區(qū)間。所了解，給出某個(gè)根附近的初始點(diǎn)或根所在的某

46、一個(gè)區(qū)間。 主要求根命令為：主要求根命令為：roots：求多項(xiàng)式的根：求多項(xiàng)式的根fzero：求一般一元函數(shù)的根：求一般一元函數(shù)的根fsolve：求非線性方程組的根：求非線性方程組的根5.4.1 求多項(xiàng)式和一般一元函數(shù)的根求多項(xiàng)式和一般一元函數(shù)的根【例【例5.4.1】求多項(xiàng)式】求多項(xiàng)式 的根的根clear;p = 1,-6,-72,-27r = roots(p) p = 1 -6 -72 -27r = 12.1229 -5.7345 -0.3884 2772623 xxx求一元函數(shù)的求一元函數(shù)的0點(diǎn)命令點(diǎn)命令 fzero其基本語(yǔ)法為：其基本語(yǔ)法為：x = fzero(fun,x0)x = fz

47、ero(fun,x0,options)x = fzero(fun,x0,options,P1,P2,.)x,fval = fzero(.)x,fval,exitflag = fzero(.)x,fval,exitflag,output = fzero(.)這里這里fun： 被求解的函數(shù)，可用被求解的函數(shù)，可用inline或字符串定義或字符串定義x0： 為初初始值，可以是一個(gè)數(shù)，也可以是一個(gè)區(qū)間為初初始值，可以是一個(gè)數(shù)，也可以是一個(gè)區(qū)間fval： 目標(biāo)函數(shù)值目標(biāo)函數(shù)值options：輸出顯示項(xiàng)，可選擇：輸出顯示項(xiàng)，可選擇display： off : 沒(méi)有輸出沒(méi)有輸出iter： 將每一步的迭代結(jié)果

48、都輸出來(lái)將每一步的迭代結(jié)果都輸出來(lái)notify： 只有當(dāng)?shù)皇諗繒r(shí)才輸出警告信息（系統(tǒng)內(nèi)定值）只有當(dāng)?shù)皇諗繒r(shí)才輸出警告信息（系統(tǒng)內(nèi)定值）TolX： 迭代終止精度。迭代終止精度。【例【例5.4.2】在】在1，2范圍范圍sinx2的根并顯示詳細(xì)迭代步驟。的根并顯示詳細(xì)迭代步驟。options=optimset(Display,iter)x = fzero(sin(x*x),1 2,options)結(jié)果為：結(jié)果為： Func-count x f(x) Procedure 1 1 0.841471 initial 2 2 -0.756802 initial 3 1.52649 0.725271

49、interpolation 4 1.75821 0.0502804 interpolation 5 1.77439 -0.00687537 interpolation 6 1.77245 3.0226e-005 interpolation 7 1.77245 1.62826e-008 interpolation 8 1.77245 -1.2098e-015 interpolation 9 1.77245 1.89882e-015 interpolationZero found in the interval: 1, 2.x = 1.7725對(duì)于任意對(duì)于任意f(x) 可能有可能有0點(diǎn)，可能沒(méi)有，

50、可能有多個(gè)點(diǎn)，可能沒(méi)有，可能有多個(gè)0點(diǎn)，甚至可能點(diǎn)，甚至可能有無(wú)窮的有無(wú)窮的0點(diǎn)。因此找不到一個(gè)通解，我們可以通過(guò)人機(jī)對(duì)話(huà)，將點(diǎn)。因此找不到一個(gè)通解，我們可以通過(guò)人機(jī)對(duì)話(huà)，將函數(shù)顯示在屏幕上，在可能的點(diǎn)上用鼠標(biāo)采樣得到初始點(diǎn)，然后讓函數(shù)顯示在屏幕上，在可能的點(diǎn)上用鼠標(biāo)采樣得到初始點(diǎn)，然后讓計(jì)算機(jī)根據(jù)這些點(diǎn)自動(dòng)分別地求出每一個(gè)根。計(jì)算機(jī)根據(jù)這些點(diǎn)自動(dòng)分別地求出每一個(gè)根?！纠纠?.4.3】編制程序，對(duì)任何函數(shù)求根，在運(yùn)行該程序時(shí)將顯示】編制程序，對(duì)任何函數(shù)求根，在運(yùn)行該程序時(shí)將顯示函數(shù)與零點(diǎn)的圖形，用戶(hù)在可能的零點(diǎn)附近用鼠標(biāo)點(diǎn)擊采點(diǎn)，程序函數(shù)與零點(diǎn)的圖形，用戶(hù)在可能的零點(diǎn)附近用鼠標(biāo)點(diǎn)擊采點(diǎn)，程序

51、通過(guò)循環(huán)計(jì)算每個(gè)采集點(diǎn)的根，并輸出結(jié)果。通過(guò)循環(huán)計(jì)算每個(gè)采集點(diǎn)的根，并輸出結(jié)果。以函數(shù)以函數(shù) 求零點(diǎn)為例，綜合敘述相關(guān)命令的用法。求零點(diǎn)為例，綜合敘述相關(guān)命令的用法。-10-50510-5-4-3-2-101ty(t) tbettfat )(sin)(2y=inline(sin(t)2*exp(-a*t)-b*abs(t),t,a,b); % 建立內(nèi)聯(lián)函數(shù)建立內(nèi)聯(lián)函數(shù)a=0.1;b=0.5;t=-10:0.01:10; % 對(duì)自變量采樣，采樣步長(zhǎng)不宜太大。對(duì)自變量采樣，采樣步長(zhǎng)不宜太大。y_char=vectorize(y); Y=feval(y_char,t,a,b);clf% 作人機(jī)對(duì)話(huà)界

52、面作人機(jī)對(duì)話(huà)界面plot(t,Y,r);hold on,plot(t,zeros(size(t),k); %畫(huà)坐標(biāo)橫軸畫(huà)坐標(biāo)橫軸xlabel(t);ylabel(y(t),hold off title(用鼠標(biāo)在函數(shù)零點(diǎn)附近連續(xù)點(diǎn)五點(diǎn)用鼠標(biāo)在函數(shù)零點(diǎn)附近連續(xù)點(diǎn)五點(diǎn))% 循環(huán)使用求解命令循環(huán)使用求解命令fzero，求鼠標(biāo)采集點(diǎn)附近的函數(shù)根。，求鼠標(biāo)采集點(diǎn)附近的函數(shù)根。for i=1:5t1,y1=ginput(1);ttt(i),yyy(i),exitflag=fzero(y,t1,a,b) % 在在ttt（i）附近搜索）附近搜索0點(diǎn)點(diǎn)endA(:,1)=ttt;A(:,2)=yyy;A 則計(jì)算結(jié)

53、果為：則計(jì)算結(jié)果為：A = x f（x） -2.0074 0.0000 -2.0074 0.0000 -0.5198 0.0000 1.6738 0.0000 1.6738 -0.00005.4.2 求多元非線性方程組的根求多元非線性方程組的根求多元非線性方程組的根，算法較為復(fù)雜，我們可以利用求多元非線性方程組的根，算法較為復(fù)雜，我們可以利用fsolve命令來(lái)求解，其命令語(yǔ)法為：命令來(lái)求解，其命令語(yǔ)法為：x = fsolve(fun,x0)x = fsolve(fun,x0,options)x = fsolve(fun,x0,options,P1,P2, . )x,fval = fsolve(

54、.)x,fval,exitflag = fsolve(.)x,fval,exitflag,output = fsolve(.)x,fval,exitflag,output,jacobian = fsolve(.)這里這里fun：多元非線性方程組：多元非線性方程組x0： 一向量，初始值一向量，初始值fval： 目標(biāo)函數(shù)值目標(biāo)函數(shù)值output: 輸出內(nèi)容控制輸出內(nèi)容控制iterations：迭代次數(shù)輸出：迭代次數(shù)輸出funcCount： 計(jì)算函數(shù)值計(jì)算函數(shù)值algorithm： 算法使用算法使用cgiterations： PCG迭代次數(shù)迭代次數(shù)（僅對(duì)大型問(wèn)題）（僅對(duì)大型問(wèn)題）stepsize：

55、最終步長(zhǎng)的選擇最終步長(zhǎng)的選擇（僅對(duì)中型問(wèn)題）（僅對(duì)中型問(wèn)題）firstorderopt options： 輸出顯示項(xiàng)，可選擇輸出顯示項(xiàng)，可選擇Diagnostics： 輸出診斷信息輸出診斷信息display： off : 沒(méi)有輸出沒(méi)有輸出iter： 將每一步的迭代結(jié)果都輸出來(lái)將每一步的迭代結(jié)果都輸出來(lái)notify： 只有當(dāng)?shù)皇諗繒r(shí)才輸出警告信息（系統(tǒng)內(nèi)定值）只有當(dāng)?shù)皇諗繒r(shí)才輸出警告信息（系統(tǒng)內(nèi)定值）Jacobian: 雅克比方法雅克比方法on : 精確雅克比方法精確雅克比方法off： 近似雅克比方法近似雅克比方法MaxFunEvals：函數(shù)的最大上限：函數(shù)的最大上限MaxIter： 最

56、大迭代次數(shù)最大迭代次數(shù)TolFun： 終止容忍函數(shù)值終止容忍函數(shù)值TolX： 終止迭代容忍度終止迭代容忍度exitflag： 退出條件退出條件 0 表示函數(shù)收斂到某解，正常退出表示函數(shù)收斂到某解，正常退出 達(dá)到最大給定迭代次數(shù)，退出達(dá)到最大給定迭代次數(shù)，退出 0 求解過(guò)程不收斂求解過(guò)程不收斂 【例【例5.4.4】求解二元函數(shù)方程組】求解二元函數(shù)方程組 在在-5，-5為初始值的零點(diǎn)。為初始值的零點(diǎn)。F = 2*x(1) - x(2) - exp(-x(1), -x(1) + 2*x(2) - exp(-x(2);x0 = -5; -5; % 初始值初始值options=optimset(Disp

57、lay,iter); % 為為Option 設(shè)置值設(shè)置值x,fval = fsolve(F,x0,options) % 求零點(diǎn)求零點(diǎn) 02,02),(21212211xxexxyxfexxyxfx = 0.5671 0.5671fval = 1.0e-008 * -0.5319 -0.53195.5求函數(shù)極值及最優(yōu)化問(wèn)題求函數(shù)極值及最優(yōu)化問(wèn)題 MATLAB提供一些求函數(shù)極值的函數(shù)，如單變量函數(shù)求極小提供一些求函數(shù)極值的函數(shù)，如單變量函數(shù)求極小fmins和多變量函數(shù)求極小和多變量函數(shù)求極小fminsearch等。這兩個(gè)函數(shù)的語(yǔ)法等。這兩個(gè)函數(shù)的語(yǔ)法和參數(shù)與求根函數(shù)類(lèi)似不再綴述。我們也可以打開(kāi)最優(yōu)

58、化工具箱使和參數(shù)與求根函數(shù)類(lèi)似不再綴述。我們也可以打開(kāi)最優(yōu)化工具箱使用里面的更專(zhuān)業(yè)的方法?，F(xiàn)將主要的命令列表如下：用里面的更專(zhuān)業(yè)的方法。現(xiàn)將主要的命令列表如下：5.5.1 無(wú)約束最優(yōu)化問(wèn)題無(wú)約束最優(yōu)化問(wèn)題【例【例5.5.1】求的】求的 最小值，給定初始值最小值，給定初始值-1x,fval,exitflag=fminbnd(x2,-1,x,optimset(TolX,1e-12,Display,off)結(jié)果為：結(jié)果為： x = 0fval = 0exitflag = 12x【例【例5.5.2】一個(gè)典型的函數(shù)】一個(gè)典型的函數(shù)是是Rosenbrock banana 函數(shù)，其典型的最小值坐標(biāo)我們已經(jīng)函

59、數(shù)，其典型的最小值坐標(biāo)我們已經(jīng)知道為（知道為（1，1）。首先將它的圖形作出。）。首先將它的圖形作出。（1）作）作Rosenbrock banana 函數(shù)的三維等高圖函數(shù)的三維等高圖x=-3:0.1:3;y=-2:0.1:4; % -3，3-2，4區(qū)域取點(diǎn)區(qū)域取點(diǎn) X,Y=meshgrid(x,y); % 根據(jù)點(diǎn)在區(qū)域打上網(wǎng)格根據(jù)點(diǎn)在區(qū)域打上網(wǎng)格F=100*(Y-X.2).2+(1-X).2; % 每一點(diǎn)處計(jì)算函數(shù)值每一點(diǎn)處計(jì)算函數(shù)值contour3(X,Y,F,300), % 作三維等高線圖形作三維等高線圖形xlabel(x),ylabel(y),axis(-3,3,-2,4,0,inf),v

60、iew(161,45)hold on,plot3(1,1,0,.r,MarkerSize,20),hold off我們可以看到有一條香蕉狀的谷底，事實(shí)上該函數(shù)的最小值不是對(duì)我們可以看到有一條香蕉狀的谷底，事實(shí)上該函數(shù)的最小值不是對(duì)所有的方法都能找到的，因?yàn)橄憬稜畹墓鹊纵^平坦，因此所有的方法都能找到的，因?yàn)橄憬稜畹墓鹊纵^平坦，因此Rosenbrock banana 函數(shù)成了著名的搜索方法的測(cè)試函數(shù)。函數(shù)成了著名的搜索方法的測(cè)試函數(shù)。（2）建立用于搜索的）建立用于搜索的Rosenbrock banana 函數(shù)的符號(hào)形式函數(shù)的符號(hào)形式ff=inline(100*(x(2)-x(1)2)2+(1-x(