For the four bus power system in power world shown below, which can be downloade
ID: 2315154 • Letter: F
Question
For the four bus power system in power world shown below, which can be downloaded from Blackboard: Reset it to Flat Start, then using Newton -Raphson power flow method to solve the power flow (determine the voltage magnitude and angle at each bus); Using Fast Decoupled Newton -Raphson method to solve the same problem; Using DC power flow to solve the same problem. Please use table to summarize the result for each question (such as voltage magnitude, voltage angle, real power, reactive power). Please add the screenshot from Power world for each question. Using DC power flow method and the data from the Power world model to solve the same problem manually. Compare the results with Power world solution.Explanation / Answer
newton's
clc
clear;
n=3;
v=[1.04 1.0 1.04];
y=[5.88228-j*23.50514 -2.9427+j*11.7676 -2.9427+j*11.7676
-2.9427+j*11.7676 5.88228-j*23.50514 -2.9427+j*11.7676
-2.9427+j*11.7676 -2.9427+j*11.7676 5.88228-j*23.50514];
bus=ones(n,1);
qlmx=zeros(n,1);
qlmn=zeros(n,1);
vmg=zeros(n,1);
qlmx(3)=1.5;
qlmn(3)=0;
vmg(2)=1.04;
diff=10;
nf=1;
ps=[inf 0.5 -1.5];
qs=[inf 1 0];
vp=v;
while (diff>0.0001 || nf==1),
eq=1;
for i=2:n,
bus(3)=2;
sc(i)=0;
sv=0;
for k=1:n,
sv=sv+y(i,k)*v(k);
end
sc(i)=v(i)*conj(sv);
p(i)=real(sc(i));
q(i)=imag(sc(i));
if bus(i)==2,
if (q(i)>qlmx(i) || q(i)<qlmn(i)),
if(q(i)<qlmn(i)),
q(i)=qlmn(i);
else
q(i)=qlmx(i);
end
bus(i)=1;
else
bus(i)=2;
end
end
if bus(i)==1,
er(eq)='p';
es(eq)=i;
mm(eq)=ps(i)-p(i);
er(eq+1)='q';
es(eq+1)=i;
mm(eq+1)=qs(i)-q(i);
cr(eq)='d';
cs(eq)=i;
cr(eq+1)='v';
cs(eq+1)=i;
eq=eq+2;
else
er(eq)='p';
es(eq)=i;
cr(eq)='d';
cs(eq)=i;
mm(eq)=ps(i)-p(i);
eq=eq+1;
end
end
mm
eq=eq-1;
nfq=eq;
up=zeros(eq,1);
abs(v)
abs(vp)
vp=v;
pause
for ceq=1:eq,
for cc=1:eq,
am=real(y(es(ceq),cs(cc))*v(cs(cc)));
bm=imag(y(es(ceq),cs(cc))*v(cs(cc)));
ei=real(v(es(ceq)));
fi=imag(v(es(ceq)));
if er(ceq)=='p' && cr(cc)=='d',
if es(ceq)~=cs(cc),
H=am*fi-bm*ei;
else
H=-q(es(ceq))-imag(y(es(ceq),cs(ceq)))*abs(v(es(ceq))^2);
end
jacob(ceq,cc)=H;
end
if er(ceq)=='p' && cr(cc)=='v',
if es(ceq)~=cs(cc),
N=am*ei+bm*fi;
else
N=p(es(ceq))+real(y(es(ceq),cs(ceq)))*abs(v(es(ceq))^2);
end
jacob(ceq,cc)=N;
end
if er(ceq)=='q' && cr(cc)=='d',
if es(ceq)~=cs(cc),
J=-(am*ei+bm*fi);
else
J=p(es(ceq))-real(y(es(ceq),cs(ceq)))*abs(v(es(ceq))^2);
end
jacob(ceq,cc)=J;
end
if er(ceq)=='q' && cr(cc)=='v',
if es(ceq)~=cs(cc),
L=am*fi-bm*ei;
else
L=q(es(ceq))-imag(y(es(ceq),cs(ceq)))*abs(v(es(ceq))^2);
end
jacob(ceq,cc)=L;
end
end
end
jacob
pause
up=inv(jacob)*mm';
nfq=1;
for i=2:n,
if bus(i)==1,
deldif=up(nfq);
newdelta=angle(v(i))+deldif;
vdif=up(nfq+1)
newvtg=abs(v(i))+vdif;
v(i)=newvtg*(cos(newdelta)+j*sin(newdelta));
nfq=nfq+2;
else
deldif=up(nfq);
newdelta=angle(v(i))+deldif;
v(i)=abs(v(i))*(cos(newdelta)+j*sin(newdelta));
nfq=nfq+1;
end
end
diff=max(abs(abs(v(2:n))-abs(vp(2:n))));
nf=nf+1;
diff
newdelta
q
nf
end
fast decoupled
clc
clear;
n=3;
v=[1.04 1.0 1.04];
y=[5.88228-j*23.50514 -2.9427+j*11.7676 -2.9427+j*11.7676
-2.9427+j*11.7676 5.88228-j*23.50514 -2.9427+j*11.7676
-2.9427+j*11.7676 -2.9427+j*11.7676 5.88228-j*23.50514];
bus=ones(n,1);
qlmx=zeros(n,1);
qlmn=zeros(n,1);
vmg=zeros(n,1);
qlmx(3)=1.5;
qlmn(3)=0;
vmg(2)=1.04;
diff=10;nf=1;
ps=[inf 0.5 -1.5];
qs=[inf 1 0];
vp=v;
while (diff>0.0001 || nf==1),
eq=1;
abs(v)
abs(vp)
vp=v;
for i=2:n,
bus(3)=2;
sc(i)=0;sv=0;
for k=1:n,
sv=sv+y(i,k)*v(k);
end
sc(i)=v(i)*conj(sv);
p(i)=real(sc(i));
q(i)=imag(sc(i));
if bus(i)==2,
if (q(i)>qlmx(i) || q(i)<qlmn(i)),
if(q(i)<qlmn(i)),
q(i)=qlmn(i);
else
q(i)=qlmx(i);
end
bus(i)=1;
else
bus(i)=2;
end
end
if bus(i)==1,
er(eq)='p';
es(eq)=i;
mm(eq)=(ps(i)-p(i))/abs(v(i));
er(eq+1)='q';
es(eq+1)=i;
mm(eq+1)=(qs(i)-q(i))/abs(v(i));
cr(eq)='d';
cs(eq)=i;
cr(eq+1)='v';
cs(eq+1)=i;
eq=eq+2;
else
er(eq)='p';
es(eq)=i;
cr(eq)='d';
cs(eq)=i;
mm(eq)=(ps(i)-p(i))/abs(v(i));
eq=eq+1;
end
end
mm
eq=eq-1;
nfq=eq;
up=zeros(eq,1);
for ceq=1:eq,
for cc=1:eq,
bm=imag(y(es(ceq),cs(cc)));
if er(ceq)=='p' && cr(cc)=='d',
if es(ceq)~=cs(cc),
H=-bm;
else
H=-imag(y(es(ceq),cs(ceq)));
end
jacob(ceq,cc)=H;
end
if er(ceq)=='q' && cr(cc)=='v',
if es(ceq)~=cs(cc),
L=-bm;
else
L=-imag(y(es(ceq),cs(ceq)));
end
jacob(ceq,cc)=L;
end
end
end
jacob
pause
up=inv(jacob)*mm';
nfq=1;
for i=2:n,
if bus(i)==1,
deldif=up(nfq);
newdelta=angle(v(i))+deldif;
vdif=up(nfq+1)
newvtg=abs(v(i))+vdif;
v(i)=newvtg*(cos(newdelta)+j*sin(newdelta));
nfq=nfq+2;
else
deldif=up(nfq);
newdelta=angle(v(i))+deldif;
v(i)=abs(v(i))*(cos(newdelta)+j*sin(newdelta));
nfq=nfq+1;
end
end
diff=max(abs(abs(v(2:n))-abs(vp(2:n))));
nf=nf+1;
newdelta
q
nf
end
Related Questions
drjack9650@gmail.com
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.