MATLAB数值域绘图演算法.docVIP

第四章 数值域的绘图演算法 设A为 n阶复数矩阵, 其数值域(numerical range)为复数子集合 W(A)={x*Ax: x?Cn, x*x=1}. 为人所知, 此集合为包含A的固有值的凸集合. 当n=2时, W(A)为一椭圆盘, 其焦距为两个固有值?, ?, 短轴长为 |A|22-|?|2- |?|2 之根号. 本章将讨论W(A)的绘图演算法, 重点是方法二与方法三的MATLAB演算法. 方法一 : Theorem 1 ([7]) Let A be an n-square complex matrix. Then W(A) = ? W(Axy), x,y are o.n. where Axy = [x*Ax x*Ay y*Ax y*Ay] ? M2 应用 Theorem 1, 只需用BASIC程式即可绘W(A), 参阅[7]的BASIC程式: step1 随机取两个orthonormal向量 x, y step2 绘椭圆 W(Axy) step3 goto step1重复直到图形成型. 方法二 : Theorem 2 ([6]) Let A be an n-square complex matrix. Then the right vertical supporting line of W(ei?A) is x = ?1(?), the maximal eigenvalue of the matrix H?=(ei?A+e-i?A*)/2 应用 Theorem 2绘W(A)边界切线的MATLAB演算法: step1 input A step2 initialize ? from 0 to 2? step3 form H? step4 compute the maximal eigenvalue of H? step5 draw the line ?1(?)e-i? which is the supporting line of e-i? W(ei?A)=W(A) step6 goto step step2 下面我们以矩阵 a= [0 2 0 ; 0 0 1 ; 0 0 3 ] 为例 , 画50条直线,可以看出一曲线与50条直线相切,这个曲线称为50条直线的包络(envelope), 亦即W(A)的边界. clear; a=[0 2 0;0 0 1;0 0 3] %a=input(matrix A=) i=sqrt(-1); H=(a+a)/2;K=(a-a)/(2*i); for k=1:50 t=(k/50)*pi; Ht=cos(t)*H-sin(t)*K;Kt=cos(t)*K+sin(t)*H; a1=max(real(eig(Ht))); b1=max(real(eig(Kt))); a2=min(real(eig(Ht))); b2=min(real(eig(Kt))); A1=[a1 a1 a1];B1=[b1 b2 b2]; X(:,2*k-1)=cos(t)*A1+sin(t)*B1; Y(:,2*k-1)=-sin(t)*A1+cos(t)*B1; A2=[a1 a2 a2];B2=[b1 b1 b2]; X(:,2*k)=cos(t)*A2+sin(t)*B2; Y(:,2*k)=-sin(t)*A2+cos(t)*B2; end %for plot(X,Y) a = 0 2 0 0 0 1 0 0 3 方法三 : Theorem 3 ([3]) Let A be an n-square complex matrix. Then the point