clear all;
A=[21  7  -1
   5   7   7
   4  -4  20];
k=1; x=A(:,k); 
[n nn]=size(x);
w=zeros(n,1);
g=norm(x(k:n));
p=sign(x(k));
s=sqrt(2*g*(g+p*x(k)));
w(k)=(x(k)+p*g)/s;
w(k+1:n)=x(k+1:n)/s;
H=eye(n)-2*w*w',
B=H*A,

k=2; x=B(:,k); 
[n nn]=size(x);
w=zeros(n,1);
g=norm(x(k:n));
p=sign(x(k));
s=sqrt(2*g*(g+p*x(k)));
w(k)=(x(k)+p*g)/s;
w(k+1:n)=x(k+1:n)/s;
H=eye(n)-2*w*w',
C=H*B,


