盖世豪侠

标题: 水帖1号（别问我是什么这只是一地水） [打印本页]

作者: songce1026 时间: 2016-6-3 12:17
标题: 水帖1号（别问我是什么这只是一地水）
clc;
clear all;
format shortEng

n=45;
H=hilb(n);
K=H(:,1:n-5);
x=ones(n-5,1);
b=K*x;
cond=cond(K)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Normal equations and the row oriented algorithm

M=K'*K;
B=K'*b;
% L=chol(M)';
[L,U]=lu(M);
X_Oriented=(L*U)\B

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%arguemented system

X=[eye(n),K;K',zeros(n-5,n-5)]\[b;zeros(n-5,1)];
X_Arguemented=X(n+1:n+n-5,

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Gram-Schmidt
r=zeros(n-5,n-5);
q=zeros(n,n-5);
for k=1:n-5
r(k,k)=norm(K(:,k),2);
if r(k,k)==0
      break
else q(:,k)=K(:,k)/r(k,k);
for j=k+1:n-5
      r(k,j)=q(:,k)'*K(:,j);
      K(:,j)=K(:,j)-r(k,j)*q(:,k);
end
end
end

Q=q;
R=r;
X_Gram=R\(Q'*b)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Householder QR Factorization
K=H(:,1:n-5);
v=zeros(n,n-5);
a=K;
beta=zeros(n-5,1);
for k=1:n-5
alpha=zeros(n-5,1);
alpha(k,1)=-sign(a(k,k))*sqrt(sumsqr(a(k:n,k)));
e=zeros(n,1);
e(k,1)=1;
v(:,k)=[zeros(k-1,1);a(k:n,k)]-alpha(k,1)*e;
beta=v(:,k)'*v(:,k);
if beta==0
      break
end
for j=k:n-5
      gama=v(:,k)'*a(:,j);
      a(:,j)=a(:,j)-(2*gama/beta)*v(:,k);
end
end
bb=a/K*b;
triu(a);
x_House=triu(a)\bb

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%singular value decomposition
[U,S,V]=svd(K);

X_svd=zeros(n-5,1);
for i=1:n-5
x=U(:,i)'*b*V(:,i)/S(i,i);
X_svd=X_svd+x;
end
X_svd

Accuracy_Oriented=cond*norm(K,2)*norm(X_Oriented,2)/norm(b,2)
Accuracy_Arguement=cond*norm(K,2)*norm(X_Arguemented,2)/norm(b,2)
Accuracy_Gram=cond*norm(K,2)*norm(X_Gram,2)/norm(b,2)
Accuracy_Householder=cond*norm(K,2)*norm(x_House,2)/norm(b,2)
Accuracy_SVD=cond*norm(K,2)*norm(X_svd,2)/norm(b,2)

作者: xinshou201632 时间: 2016-6-3 12:27

作者已被禁止，内容自动屏蔽

欢迎光临盖世豪侠 (http://www.jinyingxinxi.cn/bbs/)