MATLAB 第四章 MATLAB的數值計算功能

第四章 MATLAB 的數值計算功能Chapter 4: Numerical computation of MATLAB一、多項式(Polynomial)1多項式的表達與創(chuàng)建(Expression and Creating of polynomial)(1) 多項式的表達（expression of polynomial）_Matlab用行矢量表達多項式系數(Coefficient)，各元素按變量的降冪順序排列，如多項式為：P(x)=a0xn+a1xn-1+a2xn-2an-1x+an則其系數矢量(Vector of coefficient)為：P=a0 a1 an-1 an如將根矢量(Vector of root)表示為：ar= ar1 ar2 arn則根矢量與系數矢量之間關系為：（x-ar1）(x- ar2) (x- arn)= a0xn+a1xn-1+a2xn-2an-1x+an(2)多項式的創(chuàng)建（polynomial creating）a)系數矢量的直接輸入法利用poly2sym函數直接輸入多項式的系數矢量，就可方便的建立符號形式的多項式。例：創(chuàng)建多項式x3-4x2+3x+2poly2sym(1 -4 3 2)ans =x3-4*x2+3*x+2POLY Convert roots to polynomial. POLY(A), when A is an N by N matrix, is a row vector with N+1 elements which are the coefficients of the characteristic polynomial, DET(lambda*EYE(SIZE(A) - A) . POLY(V), when V is a vector, is a vector whose elements are the coefficients of the polynomial whose roots are the elements of V . For vectors, ROOTS and POLY are inverse functions of each other, up to ordering, scaling, and roundoff error.b) 由根矢量創(chuàng)建多項式通過調用函數 p=poly(ar)產生多項式的系數矢量, 再利用poly2sym函數就可方便的建立符號形式的多項式。注：（1）根矢量元素為n ，則多項式系數矢量元素為n+1；（2）函數poly2sym(pa) 把多項式系數矢量表達成符號形式的多項式，缺省情況下自變量符號為x，可以指定自變量。（3）使用簡單繪圖函數ezplot可以直接繪制符號形式多項式的曲線。例 1：由根矢量創(chuàng)建多項式。將多項式(x-6)(x-3)(x-8)表示為系數形式 a=6 3 8 %根矢量pa=poly(a) %求系數矢量ppa=poly2sym(pa) %以符號形式表示原多項式ezplot(ppa,-50,50)pa = 1 -17 90 -144ppa =x3-17*x2+90*x-144注：含復數根的根矢量所創(chuàng)建的多項式要注意：（1）要形成實系數多項式，根矢量中的復數根必須共軛成對； （2）含復數根的根矢量所創(chuàng)建的多項式系數矢量中，可能帶有很小的虛部，此時可采用取實部的命令(real)把虛部濾掉。進行多項式的求根運算時，有兩種方法，一是直接調用求根函數roots，poly和 roots 互為逆函數。另一種是先把多項式轉化為伴隨矩陣，然后再求其特征值，該特征值即是多項式的根。例 3： 由給定復數根矢量求多項式系數矢量。r=-0.5 -0.3+0.4i -0.3-0.4i;p=poly(r)pr=real(p)ppr=poly2sym(pr)p = 1.0000 1.1000 0.5500 0.1250pr = 1.0000 1.1000 0.5500 0.1250ppr =x3+11/10*x2+11/20*x+1/8c) 特征多項式輸入法用poly函數可實現由矩陣的特征多項式系數創(chuàng)建多項式。條件：特征多項式系數矢量的第一個元素必須為一。例 2： 求三階方陣A的特征多項式系數，并轉換為多項式形式。a=6 3 8;7 5 6; 1 3 5Pa=poly(a) %求矩陣的特征多項式系數矢量Ppa=poly2sym(pa)Pa = 1.0000 -16.0000 38.0000 -83.0000Ppa =x3-17*x2+90*x-144注：n 階方陣的特征多項式系數矢量一定是n +1階的。注：（1）要形成實系數多項式，根矢量中的復數根必須共軛成對； （2）含復數根的根矢量所創(chuàng)建的多項式系數矢量中，可能帶有很小的虛部，此時可采用取實部的命令(real)把虛部濾掉。進行多項式的求根運算時，有兩種方法，一是直接調用求根函數roots，poly和 roots 互為逆函數。另一種是先把多項式轉化為伴隨矩陣，然后再求其特征值，該特征值即是多項式的根。例 4： 將多項式的系數表示形式轉換為根表現形式。求 x3-6x2-72x-27的根a=1 -6 -72 -27r=roots(a)r = 12.1229 -5.7345 -0.3884MATLAB約定，多項式系數矢量用行矢量表示，根矢量用列矢量表示。1. 多項式的乘除運算(Multiplication and division of polynomial)多項式乘法用函數conv(a,b)實現， 除法用函數deconv(a,b)實現。例1：a(s)=s2+2s+3, b(s)=4s2+5s+6,計算 a(s)與 b(s)的乘積。a=1 2 3; b=4 5 6;c=conv(a,b)cs=poly2sym(c,s)c = 4 13 28 27 18cs =4*s4+13*s3+28*s2+27*s+18例2： 展開(s2+2s+2)(s+4)(s+1) （多個多項式相乘）c=conv(1,2,2,conv(1,4,1,1)cs=poly2sym(c,s) %（指定變量為s）c = 1 7 16 18 8cs =s4+7*s3+16*s2+18*s+8例2：求多項式s4+7*s3+16*s2+18*s+8分別被(s+4),(s+3)除后的結果。c=1 7 16 18 8;q1,r1=deconv(c,1,4) %q商矢量， r余數矢量q2,r2=deconv(c,1,3)cc=conv(q2,1,3) %對除(s+3)結果檢驗test=(c-r2)=cc)q1 = 1 3 4 2r1 = 0 0 0 0 0q2 = 1 4 4 6r2 = 0 0 0 0 -10cc = 1 7 16 18 18test = 1 1 1 1 11. 其他常用的多項式運算命令(Other computation command of polynomial)pa=polyval(p,s) 按數組運算規(guī)則計算給定s時多項式p的值。pm=polyvalm(p,s) 按矩陣運算規(guī)則計算給定s時多項式p的值。r,p,k=residue(b,a) 部分分式展開，b,a分別是分子分母多項式系數矢量，r,p,k分別是留數、極點和直項矢量p=polyfit(x,y,n) 用n階多項式擬合x，y矢量給定的數據。polyder(p) 多項式微分。注： 對于多項式b(s)與不重根的n階多項式a(s)之比，其部分分式展開為： 式中：p1,p2,pn稱為極點(poles)，r1,r2,rn 稱為留數(residues)，k(s)稱為直項(direct terms)，假如a(s)含有m重根pj，則相應部分應寫成：RESIDUE Partial-fraction expansion (residues). R,P,K = RESIDUE(B,A) finds the residues, poles and direct term of a partial fraction expansion of the ratio of two polynomials B(s)/A(s). If there are no multiple roots,B(s) R(1) R(2) R(n) - = - + - + . + - + K(s) A(s) s - P(1) s - P(2) s - P(n)Vectors B and A specify the coefficients of the numerator and denominator polynomials in descending powers of s. The residuesare returned in the column vector R, the pole locations in column vector P, and the direct terms in row vector K. The number of poles is n = length(A)-1 = length(R) = length(P). The direct term coefficient vector is empty if length(B) n 可求出最小二乘解*欠定系統(tǒng)：(Underdetermind system) m n, then only the first n columns of Q are computed.4. 特征值與特征矢量(Eigenvalues and eigenvectors).MATLAB中使用函數eig計算特征值和 特征矢量，有兩種調用方法：*e=eig(a), 其中e是包含特征值的矢量；*v,d=eig(a), 其中v是一個與a相同的nn階矩陣，它的每一列是矩陣a的一個特征值所對應的特征矢量，d為對角陣，其對角元素即為矩陣a的特征值。例：計算特征值和特征矢量。a=34 25 15; 18 35 9; 41 21 9e=eig(a)v,d=eig(a)a = 34 25 15 18 35 9 41 21 9e = 68.5066 15.5122 -6.0187v = -0.6227 -0.4409 -0.3105 -0.4969 0.6786 -0.0717 -0.6044 -0.5875 0.9479d = 68.5066 0 0 0 15.5122 0 0 0 -6.0187EIG Eigenvalues and eigenvectors.E = EIG(X) is a vector containing the eigenvalues of a square matrix X.V,D = EIG(X) produces a diagonal matrix D of eigenvalues and a full matrix V whose columns are the corresponding eigenvectors so that X*V = V*D.V,D = EIG(X,nobalance) performs the computation with balancing disabled, which sometimes gives more accurate results for certain problems with unusual scaling. If X is symmetric, EIG(X,nobalance) is ignored since X is already balanced.5. 奇異值分解.( Singular value decomposition).如存在兩個矢量u,v及一常數c,使得矩陣A滿足：Av=cu, Au=cv稱c為奇異值，稱u,v為奇異矢量。 將奇異值寫成對角方陣，而相對應的奇異矢量作為列矢量則可寫成兩個正交矩陣U，V， 使得： AV=U， AU=V 因為U，V正交，所以可得奇異值表達式： A=UV。一個m行n列的矩陣A經奇異值分解，可求得m行m列的U, m行n列的矩陣和n行n列的矩陣V.。奇異值分解用svd函數實現，調用格式為；u,s,v=svd(a) SVD Singular value decomposition.U,S,V = SVD(X) produces a diagonal matrix S, of the same dimension as X and with nonnegative diagonal elements in decreasing order, and unitary matrices U and V so that X = U*S*V.S = SVD(X) returns a vector containing the singular values.U,S,V = SVD(X,0) produces the economy size decomposition. If X is m-by-n with m n, then only the first n columns of U are computed and S is n-by-n.例： 奇異值分解。a=8 5; 7 3;4 6;u,s,v=svd(a) % s為奇異值對角方陣u = -0.6841 -0.1826 -0.7061 -0.5407 -0.5228 0.6591 -0.4895 0.8327 0.2589s = 13.7649 0 0 3.0865 0 0v = -0.8148 -0.5797 -0.5797 0.8148 五 數據分析(Data Analyaia)MATLAB對數據分析有兩條約定：（1） 若輸入量X是矢量，則不論是行矢量還是列矢量，運算是對整個矢量進行的； （2）若輸入量X是數組，（或稱矩陣），則命令運算是按列進行的。即默認每個列是有一個變量的不同“觀察“所得的數據組成。 1. 基本統(tǒng)計命令 （表4-1）例： 做各種基本統(tǒng)計運算。A=5 -10 -6 0;2 6 3 -3;-9 5 -10 11;-22 17 0 -19;-1 6 -4 4Amax=max(A) %找A各列的最大元素Amin=min(A) %找A各列的最小元素Amed=median(A) %找A各列的中位元素Amean=mean(A) %找A各列的平均值Astd=std(A) %求A各列的標準差Aprod=prod(A) %求A各列元素的積Asum=sum(A) %求A各列元素的和S=cumsum(A) %求A各列元素的累積和P=cumprod(A) %求A各列元素的累積j積I=sort(A) %使A的各列元素按遞增排列A = 5 -10 -6 0 2 6 3 -3 -9 5 -10 11 -22 17 0 -19 -1 6 -4 4Amax = 5 17 3 11Amin = -22 -10 -10 -19Amed = -1 6 -4 0Amean = -5.0000 4.8000 -3.4000 -1.4000Astd = 10.8397 9.6281 5.0794 11.1490Aprod = -1980 -30600 0 0Asum = -25 24 -17 -7S = 5 -10 -6 0 7 -4 -3 -3 -2 1 -13 8 -24 18 -13 -11 -25 24 -17 -7P = 5 -10 -6 0 10 -60 -18 0 -90 -300 180 0 1980 -5100 0 0 -1980 -30600 0 0I = -22 -10 -10 -19 -9 5 -6 -3 -1 6 -4 0 2 6 0 4 5 17 3 11求矩陣元素的最大值、最小值可用： Amax=max(maxA) 或 Amax=max(A(:)， Amin=min(min(A) 或 Amin=min(A(:) 2協方差陣和相關陣(Covariance matrix and Correlation coefficients).（表 42）例： 計算協方差和相關陣。x=rand(10,3); y=rand(10,3);cx=cov(x) %求協方差陣cy=cov(y)cxy=cov(x,y) %求兩隨機變量的協方差px=corrcoef(x) %求相關陣pxy=corrcoef(x,y) %求兩隨機變量的（22）相關系數cx = 0.0483 -0.0066 0.0146 -0.0066 0.0283 0.0154 0.0146 0.0154 0.0978cy = 0.1177 0.0073 -0.0127 0.0073 0.0239 -0.0230 -0.0127 -0.0230 0.0772cxy = 0.0550 0.0023 0.0023 0.0697px = 1.0000 -0.1783 0.2118 -0.1783 1.0000 0.2934 0.2118 0.2934 1.0000pxy = 1.0000 0.0372 0.0372 1.0000COV Covariance matrix.COV(X), if X is a vector, returns the variance. For matrices, where each row is an observation, and each column a variable, COV(X) is the covariance matrix. DIAG(COV(X) is a vector of variances for each column, and SQRT(DIAG(COV(X) is a vector of standard deviations. COV(X,Y), where X and Y are vectors of equal length, is equivalent to COV(X(:) Y(:). COV(X) or COV(X,Y) normalizes by (N-1) where N is the number of observations. This makes COV(X) the best unbiased estimate of the covariance matrix if the observations are from a normal distribution.CORRCOEF Correlation coefficie