※ 引述《snow3804 (snow3804)》之铭言：
: matlab有hess(A)可以将矩阵A转换成Hessenberg矩阵
: 但我想知道这中间是怎么计算的
: 搜寻电脑里有hess.m档只是打开却只看到注解而已
: http://tinyurl.com/m6v7p96
: 有没有办法可以知道是怎么实作的
推 iHakka: matlab是用lapack去实作的，可以追回去看看 08/13 19:27
推 infernodimon: lapack里有 印象中是用householder 08/14 19:50
我讲一下原由,我在这里查到Hessenberg的matlab程式码
http://www.auburn.edu/~tamtiny/lecture%2010.pdf
function [H,Q]=houshess(A)
n=max(size(A));
Q=eye(n);
H=A;
for k=1:(n-2)
[v,beta]=vhouse(H(k+1:n,k));
I=eye(k);
N=zeros(k,n-k);
m=length(v);
R=eye(m)-beta*v*v';
HH=H(k+1:n,k:n);
H(k+1:n,k:n)=R*H(k+1:n,k:n);
HH=H(1:n,k+1:n);
H(1:n,k+1:n)=H(1:n,k+1:n)*R;
P=[I, N; N', R];
Q=Q*P;
end
return
function [v,beta]=vhouse(x)
n=length(x);
x=x/norm(x);
s=x(2:n)'*x(2:n);
v=[1; x(2:n)];
if (s==0), beta=0;
else
mu=sqrt(x(1)^2+s);
if (x(1) <= 0)
v(1)=x(1)-mu;
else
v(1)=-s/(x(1)+mu);
end
beta=2*v(1)^2/(s+v(1)^2);
v=v/v(1);
end
return
但执行的结果和matla内建的hess()结果不一样
http://www.mathworks.com/help/matlab/ref/hess.html
A=[-149 -50 -154;537 180 546;-27 -9 -25];
H=houshess(A)
H =
-149.0000 -42.2037 156.3165
537.6783 152.5511 -554.9272
0 0.0728 2.4489
H=hess(A)
H =
-149.0000 42.2037 -156.3165
-537.6783 152.5511 -554.9272
0 0.0728 2.4489
看得出来其实只是有些地方差个负号而已
当然这两个都是Hessenberg矩阵
但我还是想知道程式码是差在哪里造成正负号不同
lapack我努力去找到程式码,只是lapack是用fortran写的,但我看不懂程式码
http://www.netlib.org/lapack/explore-3.1.1-html/sgehrd.f.html