幂法-反幂法求解矩阵最大最小特征值及其对应的特征向量17页.doc

会员专享推荐 · 你可能也用得上更多

开通会员，畅读全站文档

本篇文档已为你保留

本篇文档价值 ~~¥35.00~~会员免费解锁

全文不限次阅读，剩余内容即刻解锁
高速下载通道，免排队、无等待
专属客服 1 对 1 答疑，优先响应
每日赠送下载券，本篇文档免费送

有 86 篇商业计划书正在上传

下载原文档

数值计算解矩阵的按模最大最小特征值及对应的特征向量一.幂法1.幂法简介：当矩阵A满足一定条件时，在工程中可用幂法计算其主特征值(按模最大)及其特征向量。矩阵A需要满足的条件为：(1)(2)存在n个线性无关的特征向量，设为1.1计算过程：不全为0，则有可见，当越小时，收敛越快；且当k充分大时，有，对应的特征向量即是。2算法实现3 matlab程序代码function[t,y]=lpowerA,x0,eps,N)%t为所求特征值，y是对应特征向量k=1;z=0;%z相当于y=x0./max(abs(x0));%规范化初始向量x=A*y;%迭代格式b=max(x);%b相当于if abs(z-b)<eps%判断第一次迭代后是否满足要求t=max(x);return;end while abs(z-b)>eps&&k<N k=k+1;z=b;y=x./max(abs(x));x=A*y;b=max(x);end[m,index]=max(abs(x));%这两步保证取出来的按模最大特征值t=x(index);%是原值，而非其绝对值。end 4举例验证选取一个矩阵A，代入程序，得到结果，并与eig(A)的得到结果比较，再计算A*y-t*y，验证y是否是对应的特征向量。结果如下：结果正确，表明算法和代码正确，然后利用此程序计算15阶Hilb矩阵，与eig(A)的得到结果比较，再计算A*y-t*y，验证y是否是对应的特征向量。设置初始向量为x0=ones(15,1)，结果显示如下可见，结果正确。得到了15阶Hilb矩阵的按模最大特征值和对应的特征向量。二.反幂法1.反幂法简介及其理论在工程计算中，可以利用反幂法计算矩阵按模最小特征值及其对应特征向量。其基本理论如下，与幂法基本相同：，可知，A和A-1的特征值互为倒数，求A按模最小特征值即求A-1的按模最大特征值,取倒数即为A的按模最小特征值所以算法基本相同，区别就是在计算2.算法实现3 matlab程序代码function[s,y]=invpower(A,x0,eps,n)%s为按模最小特征值，y是对应特征向量k=1;r=0;%r相当于y=x0./max(abs(x0));%规范化初始向量[L,U]=lu(A);z=L\y;x=U\z;u=max(x);s=1/u;%按模最小为A-1按模最大的倒数.if abs(u-r)<eps%判断第一次迭代后是否满足终止条件return end while abs(u-r)>eps&&k<n%终止条件.k=k+1;r=u;y=x./max(abs(x));z=L\y;x=U\z;u=max(x);end[m,index]=max(abs(x));%这两步保证取出来的按模最大特征值s=1/x(index);%是原值，而非其绝对值。end 4举例验证同幂法一样，选取一个矩阵A，代入程序，得到结果，并与eig(A)的得到结果比较，再计算A*y-t*y，验证y是否是对应的特征向量。可见结果正确，然后利用此程序计算15阶Hilb矩阵，eig(A)的得到结果比较，再计算A*y-s*y，验证y是否是对应的特征向量。设置初始向量为x0=ones(15,1)，结果显示如下可见，结果真确。得到了15阶Hilb矩阵的按模最大特征值和对应的特征向量。三.计算条件数矩阵A的条件数等于A的范数与A的逆的范数的乘积，即cond(A)=‖A‖·‖A^(-1)‖，对应矩阵的3种范数，可以定义3种条件数。函数cond(A,1)、cond(A)或cond(A inf)是判断矩阵病态与否的一种度量，条件数越大表明矩阵的病态程度越大.,而如果A为对称矩阵，如Hilb矩阵，的最大最小特征值，分别为A的最大最小特征值的平方。所以cond(A)为A的最大最小特征值得比值。对于本例中的15阶Hilb矩阵来说，利用上面计算结果得其条件数(选择第二种条件数）为：3.0934e+017；这与直接利用cond(A)得到的结果：2.5083e+017在同一数量级，再次表明了上述算得得最大最小特征值的正确性，同时又表明Hilb矩阵是病态矩阵。四.Aitken商加速法1.简介与原理同幂法和反幂法计算最大和最小特征值类似，如果计算最大特征值，则迭代格式为；计算最小特征值时，迭代格式为。2.算法实现计算按模最大特征值算法如下：类似幂法和反幂法可以写出按模最小特征值算法，此处不再赘述。3.matlab程序代码function[r,y]=aitken(A,x0,eps,n)%r按模最大特征值,y为对应特征向量k=1;a0=0;%a相当于a1=1;%a1相当于r0=1;%相当于2中的y=x0./max(abs(x0));%规范化初始向量x=A*y;a2=max(abs(x));%a2相当于r=a0-(a1-a0)^2/(a2-2*a1+a0);%相当于if(a2-2*a1+a0)==0%若上式中分母为0，则迭代失败，返回disp"初始向量迭代失败"return;end if abs(r-r0)<eps%判断第一次迭代后是否满足要求，如满足，则返回结果return end while abs(r-r0)>eps&&k<n%终止条件k=k+1;a0=a1;a1=a2;r0=r;y=x./max(abs(x));x=A*y;%迭代格式a2=max(abs(x));if(a2-2*a1+a0)==0%若分母为0，则迭代失败，返回return;end r=a0-(a1-a0)^2/(a2-2*a1+a0);[m,index]=max(abs(eig(A)));%以下代码保证取出来的按模最大特征值aa=eig(A);%是原值，而非其绝对值。if aa(index)>0||aa(index)==0 r=r;else r=-r;end end end类似可得按模最小特征值和特征向量的代码如下：与上面类似，所不同的只是迭代格式不同.function[r,y]=invaitken(A,x0,eps,n)k=1;a0=0;a1=1;r0=1;y=x0./max(abs(x0));[L,U]=lu(A);%迭代格式的不同z=L\y;x=U\z;a2=max(abs(x));r=a0-(a1-a0)^2/(a2-2*a1+a0);if(a2-2*a1+a0)==0 disp"初始向量迭代失败"return;end if abs(r-r0)<eps%判断第一次迭代后是否满足要求，如满足，则返回结果return end while abs(r-r0)>eps&&k<n k=k+1;a=b;b=c;r0=r;y=x./max(abs(x));z=L\y;x=U\z;a2=max(abs(x));if(a2-2*a1+a0)==0 return;end r=a0-(a1-a0)^2/(a2-2*a1+a0);end[m,index]=min(abs(eig(A)));%以下代码保证取出来的按模最大特征值aa=eig(A);%是原值，而非其绝对值。if aa(index)>0||aa(index)==0 r=1/r;else r=-1/r;end end 4.计算Hilb矩阵特征值此处不再举例，而是直接应用于15阶Hilb矩阵，初始向量选为ones(15,1)，结果如下，并将结果与幂法和反幂法得到结果比较这与幂法得到的特征值和特征向量一致，表明算法和代码正确；同理，最小特征值结果如下：这与反幂法得到的结果一致，表明结果正确。五，对称矩阵的Rayleigh商加速法1.简介与原理原理如下：2.算法实现3.Matlab程序代码function[r,y]=rayleigh(A,x0,eps,n)%r是特征值，y是特征向量k=1;r0=0;y=x0./max(abs(x0));x=A*y;%迭代格式计算新的x r=dot(y,x)/dot(y,y);%Reyleigh商if abs(r-r0)<eps return end while abs(r-r0)>eps&&k<n k=k+1;r0=r;y=x./max(abs(x));x=A*y;r=dot(y,x)/dot(y,y);end end类似得计算按模最小特征值的Rayleigh商加速法，如下：function[r,y]=invrayleigh(A,x0,eps,n)k=1;r0=0;y=x0./max(abs(x0));[L,U]=lu(A);%迭代格式不同z=L\y;x=U\z;r=dot(y,x)/dot(y,y);if abs(r-r0)<eps return end while abs(r-r0)>eps&&k<n k=k+1;r0=r;y=x./max(abs(x));z=L\y;x=U\z;r=dot(y,x)/dot(y,y);end r=1/r;end 4.计算Hilb矩阵特征值此处不再举例，而是直接应用于15阶Hilb矩阵，初始向量选为ones(15,1)，结果如下，并将结果与幂法和反幂法得到结果比较这与幂法得到结果一致，表明算法和代码正确。同理，最小特征值如下：与反幂法得到结果一致，表明代码和算法正确。