?? programengine.m
字號:
Hqij(i,j)=0;
if j==temp(1)
Hqij(i,j)=v_bus(temp(1))*v_bus(temp(2))*(G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))+B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
end
if j==temp(2)
Hqij(i,j)=v_bus(temp(1))*v_bus(temp(2))*(-G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
end
if j==temp(1)+nbus
Hqij(i,j)=2*v_bus(temp(1))*B(temp(1),temp(2))+v_bus(temp(2))*(G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
end
if j==temp(2)+nbus
Hqij(i,j)=v_bus(temp(1))*(G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
end
end
end
%Hqij
%formation of Hpji matrix corresponding to real power flow between bus j
%and i mismatches
%for i=1:length(mmpji)
% for j=1:nstate
% temp=line_tracker(:,i);
% trash=temp(1);
% temp(1)=temp(2);
% temp(2)=trash;
% Hpji(i,j)=0;
% if j==temp(1)
% Hpji(i,j)=v_bus(temp(1))*v_bus(temp(2))*(-G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))+B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
%end
% if j==temp(2)
% Hpji(i,j)=v_bus(temp(1))*v_bus(temp(2))*(G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
%end
%if j==temp(1)+nbus
% Hpji(i,j)=-2*v_bus(temp(1))*G(temp(1),temp(2))+v_bus(temp(2))*(G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))+B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
%end
%if j==temp(2)+nbus
% Hpji(i,j)=v_bus(temp(1))*(G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))+B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
%end
%end
%end
%Hpji
%formation of Hqji matrix corresponding to reactive power flow between bus
%j and i mismatches
%for i=1:length(mmqji)
% for j=1:nstate
% temp=line_tracker(:,i);
% trash=temp(1);
% temp(1)=temp(2);
% temp(2)=trash;
% Hqji(i,j)=0;
% if j==temp(1)
% Hqji(i,j)=v_bus(temp(1))*v_bus(temp(2))*(G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))+B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
% end
% if j==temp(2)
% Hqji(i,j)=v_bus(temp(1))*v_bus(temp(2))*(-G(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2))));
% end
% if j==temp(1)+nbus
% Hqji(i,j)=2*v_bus(temp(1))*B(temp(1),temp(2))+v_bus(temp(2))*(G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
% end
% if j==temp(2)+nbus
% Hqji(i,j)=v_bus(temp(1))*(G(temp(1),temp(2))*sin(ang_bus(temp(1))-ang_bus(temp(2)))-B(temp(1),temp(2))*cos(ang_bus(temp(1))-ang_bus(temp(2))));
% end
% end
%end
%Hqji
%birth of Jacobian
%H=[Hv; Hp; Hq; Hpij; Hqij; Hpji; Hqji]
H=[Hv; Hp; Hq; Hpij; Hqij];
hsize=size(H);
if count==0
fid=fopen('output.txt','w');
fprintf(fid,'Iteration number %d\n\n',count);
fprintf(fid,'The Jacobian\n\n');
for i=1:hsize(1)
fprintf(fid,'%10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f\n',H(i,:));
end
fclose(fid);
end
if count>0
fid=fopen('output.txt','a');
fprintf(fid,'Iteration number %d\n\n',count);
fprintf(fid,'The Jacobian\n\n');
for i=1:hsize(1)
fprintf(fid,'%10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f\n',H(i,:));
end
fclose(fid);
end
%determination of weight matrix
%v_wts
%p_wts
%q_wts
%pij_wts
%qij_wts
%pji_wts
%qji_wts
%W=eye(nmeas,nmeas);
W=zeros(nmeas,nmeas);
%W
for i=1:nbus
%i
for j=1:nmeas
% j
W(i,j)=0;
if j==i
W(i,j)=v_wts(i);
end
end
%W
end
i_temp=i;
for i=i_temp+1:i_temp+nbus
%i
for j=1:nmeas
W(i,j)=0;
if j==i
W(i,j)=p_wts(i-i_temp);
end
end
%W
end
i_temp=i;
for i=i_temp+1:i_temp+nbus
%i
for j=1:nmeas
W(i,j)=0;
if j==i
W(i,j)=q_wts(i-i_temp);
end
end
%W
end
i_temp=i;
for i=i_temp+1:i_temp+nlines
for j=1:nmeas
W(i,j)=0;
if j==i
W(i,j)=pij_wts(i-i_temp);
end
end
%W
end
i_temp=i;
for i=i_temp+1:i_temp+nlines
for j=1:nmeas
W(i,j)=0;
if j==i
W(i,j)=qij_wts(i-i_temp);
end
end
%W
end
%i_temp=i;
%for i=i_temp+1:i_temp+nlines
% for j=1:nmeas
% W(i,j)=0;
% if j==i
% W(i,j)=pji_wts(i-i_temp);
% end
%end
% W
%end
%i_temp=i;
%for i=i_temp+1:i_temp+nlines
% for j=1:nmeas
% W(i,j)=0;
% if j==i
% W(i,j)=qji_wts(i-i_temp);
% end
%end
%W
%end
% size(H)
% size(H')
% size(W)
%determination of state mismatch vector
%W;
Gain=H'*W*H;
%Gain
Gain(1,1)=10000000;
gsize=size(Gain);
fid=fopen('output.txt','a');
fprintf(fid,'\n\nThe following is the Gain matrix.\n\n');
for i=1:gsize(1);
fprintf(fid,'%10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f %10.2f\n',Gain(i,:));
end
fclose(fid);
niag=inv(Gain);
smm=inv(Gain)*H'*W*mm;
%test=smm*1000
count=count+1;
%v_bus
%ang_bus
smm(1)=0; %to ensure that bus angle 1 remains as reference. no change in it
%smm
smma=smm(1:nbus,:);
smmv=smm(nbus+1:2*nbus,:);
%determination of tolerance
plz_giveta=[];
for i=1:length(smm)
if smm(i)>tolerance
plz_giveta(i)=1;
end
end
if length(plz_giveta)>0
flag=1;
v_bus=v_bus+smmv;
ang_bus=ang_bus+smma;
end
%v_bus
%ang_bus
if length(plz_giveta)==0
flag=0;
% clc;
fid=fopen('output.txt','a');
fprintf(fid,'\n\nConvergence has occured in %d iterations.\n',count);
fprintf(fid,'Estimated bus voltages.\n');
fprintf(fid,'Bus number Voltage\n');
for i=1:nbus
fprintf(fid,'%f %f\n',bus_number(i),v_bus(i));
end
fprintf(fid,'Estimated bus angles\n.');
fprintf(fid,'Bus number Angle\n');
for i=1:nbus
fprintf(fid,'%f %f\n',bus_number(i),ang_bus(i));
end
%Real power injection calculation
p_calc=[];
for i=1:nbus
p_calc(i)=0;
for j=1:nbus
p_calc(i)=p_calc(i)+v_bus(i)*v_bus(j)*(G(i,j)*cos(ang_bus(i)-ang_bus(j))+B(i,j)*sin(ang_bus(i)-ang_bus(j)));
end
end
fprintf(fid,'Estimated real power bus injections.\n');
fprintf(fid,'Bus number P\n');
for i=1:nbus
fprintf(fid,'%f %f\n',bus_number(i),p_calc(i));
end
%reactive power injection calculations
q_calc=[];
for i=1:nbus
q_calc(i)=0;
for j=1:nbus
q_calc(i)=q_calc(i)+v_bus(i)*v_bus(j)*(G(i,j)*sin(ang_bus(i)-ang_bus(j))-B(i,j)*cos(ang_bus(i)-ang_bus(j)));
end
end
fprintf(fid,'Estimated reactive power bus injections.\n');
fprintf(fid,'Bus number Q\n');
for i=1:nbus
fprintf(fid,'%f %f\n',bus_number(i),q_calc(i));
end
start_bus=x;
end_bus=y;
%i to j real flows calculation
pij_calc=[];
length(x);
for i=1:length(x)
q=start_bus(i);
r=end_bus(i);
pij_calc(i)=-((v_bus(q))^2)*G(q,r)+v_bus(q)*v_bus(r)*(cos(ang_bus(q)-ang_bus(r))*G(q,r)+sin(ang_bus(q)-ang_bus(r))*B(q,r));
end
fprintf(fid,'Real line flows\n.');
fprintf(fid,'From_bus To_bus Real power flow\n');
for i=1:length(start_bus)
fprintf(fid,'%f %f %f\n',start_bus(i),end_bus(i),pij_calc(i));
end
%i to j reactive flows calculation
qij_calc=[];
for i=1:length(x)
q=start_bus(i);
r=end_bus(i);
qij_calc(i)=((v_bus(q))^2)*B(q,r)+v_bus(q)*v_bus(r)*(sin(ang_bus(q)-ang_bus(r))*G(q,r)-cos(ang_bus(q)-ang_bus(r))*B(q,r));
end
fprintf(fid,'Reactive line flows\n.');
fprintf(fid,'From_bus To_bus Reactive power flow\n');
for i=1:length(start_bus)
fprintf(fid,'%f %f %f\n',start_bus(i),end_bus(i),qij_calc(i));
end
fclose(fid);
break;
end
%if count==1
% fid=fopen('output.txt','w');
% fprintf(fid,'%6.2f %6.ff %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f\n',smm');
%fclose(fid);
%end
%if count>1
% fid=fopen('output.txt','a');
%fprintf(fid,'%6.2f %6.ff %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f\n',smm');
%fclose(fid);
%end
%count
end
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -