clear;
clc;
yb=zeros(5,5);
yb(1,1)=(1.37874-6.26166i)/2;yb(1,2)=-0.62402+3.90015i;yb(1,3)=-0.75471+2.64150i;
yb(2,2)=(1.45390-66.98082i)/2;yb(2,3)=-0.82987+3.11203i;yb(2,4)=63.49206i;
yb(3,3)=(1.58459-35.73786i)/2;yb(3,5)=31.74603i;
yb(4,4)=(-66.66667i)/2;
yb(5,5)=(-33.33333i)/2;
yb=yb+conj(yb');
G=real(yb);
B=imag(yb);
s=[-1.6-0.8i;-2-i;-3.7-1.3i;5];
U0=[1;1;1;1.05;1.05];
U1=zeros(4,1);
deta_U=zeros(4,1);
eps1=10^-4;
eps_U=1;
theta=0;
k=0;
while (eps_U>eps1)
q=0;
for j=1:5
q=q+U0(4)*conj(yb(4,j))*conj(U0(j));
end
q=imag(q);
s(4)=5+sqrt(-1)*q;
p=0;
for j=1:5
p=p+yb(4,j)*U0(j);
end
p=p-yb(4,4)*U0(4);
U1(4)=(conj(s(4))/conj(U0(4))-p)/yb(4,4);
theta=angle(U1(4));
U1(4)=1.05*(cos(theta)+sqrt(-1)*sin(theta));
for i=1:3
a=0;
b=0;
for j=1:i-1
a=a+yb(i,j)*U1(j);
end
for j=i+1:4
b=b+yb(i,j)*U0(j);
end
U1(i)=(conj(s(i))/conj(U0(i))-yb(i,5)*U0(5)-a-b)/yb(i,i);
end
for i=1:4
deta_U=U1(i)-U0(i);
U0(i)=U1(i);
end
eps_U=max(abs(deta_U));
k=k+1;
end
disp(abs(U0));
disp(U0);
theta=zeros(5,1);
for i=1:5
theta(i)=angle(U0(i));
end
disp(theta);
y0=[-4i;-8i;-4i;0;0];
s=zeros(5,5);
for i=1:5
for j=1:5
s(i,j)=U0(i)*(conj(U0(i))*conj(y0(i))+(conj(U0(i))-conj(U0(j)))*conj(-yb(i,j)));
end
end
for i=1:5
s(i,i)=0;
end
deta_s=zeros(5,5);
for i=1:5
for j=1:5
deta_s(i,j)=s(i,j)+s(j,i);
end
end
disp(s);
disp(deta_s);