matlab 对称矩阵特征值为负数,【原创】古典Jaccobi方法计算实对称矩阵A的特征值和向...

帖由Matlab技术论坛原创,原帖参见http://www.matlabsky.cn/thread-5243-1-1.html

转帖请注明出去,更多精彩Matlab资源参见http://www.matlabsky.com

愿您健康美丽更动人!

2010-3-1 18:52:51

function [V,D,iter]=eigcj(A,tol,itermax)

% 古典Jaccobi方法计算实对称矩阵A的特征值和向量

% (eigen value and vector of classical Jaccobi

algorithm)

%

% 参数说明

% A:实对称矩阵

% tol:精度,非对角元素平方和

% itermax:最大迭代次数

% V:特征值

% D:特征向量

% iter:实际迭代次数

%

% 实例说明

% A=[-12 3 3;3 1 -2;3 -2 7];

% [V,D,iter]=eigcj(A,1e-10,20)

% [V2,D2]=eig(A) % 使用MATLAB自带eig函数进行验证

%

% 基础知识

% 1、对于实对称矩阵A,古典雅克比算法必定收敛矩阵的全部特征根

% 2、实对称矩阵的特征值均为实数,且存在标准正交的特征向量系

% 3、对任意实对称矩阵A,总存在正交矩阵Q,使得QAQ'=diag(λ1,λ2,...)

%  λ为A的特征根,Q各列对应其特征向量

% 4、设A为实对称矩阵,Q为正交矩阵,则||A||=||QA||=||AQ||=||QAQ'||=特征根平方和

%  其中||X||表示X的F范数

%

% 算法说明

% 0、赋初值H=I,A0=A

% 1、选非对角线主元素|a(p,q)|=max(|a(i,j)|)

% 2、构造Givens矩阵R

% 3、对主元素a(p,q)进行Givens变换

% 3.1 更新H=R*H

% 3.2 更新A=H*A0*H'

% 4、计算SA,矩阵A的非对角元素平方和

% 5、SA是否到达指定精度,是否达到最大迭代次数

% 6、否,继续循环;是,退出

%

% 参考资料

% 封建湖《数值分析原理》.2005年7月.科学出版社.Page238-243

%

[m,n]=size(A);

H=eye(m,n);

A0=A;

SA=sum(A(:).^2)-sum(diag(A).^2);

iter=0;

while SA>tol &&

iter

iter=iter+1;

%

选主元素a(p,q)

T=triu(A,1);

[tmp,idx]=max(abs(T(:)));

%

计算Givens矩阵中的cos,sin值

[p,q]=ind2sub([m,n],idx);

%

计算非对角元素平方和SA

SA=SA-2*A(p,q).^2;

if

A(p,p)==A(q,q)

cos=sqrt(2)/2;

sin=sign(A(p,q))*cos;

else

d=(A(p,p)-A(q,q))/(2*A(p,q));

tan=sign(d)/(abs(d)+sqrt(d^2+1));

cos=1/sqrt(1+tan^2);

sin=tan*cos;

end

%

计算Givens矩阵

R=eye(m,n);

R(p,p)=cos;

R(p,q)=sin;

R(q,p)=-sin;

R(q,q)=cos;

%

更新H矩阵

H=R*H;

%

更新A矩阵

A=H*A0*H';

end

V=diag(A); % 特征根

D=H'; % 特征向量

  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值