Experiment No. 03
Aim of the Experiment: –
Programming for load flow analysis using Newton-Raphson method.
Apparatus Required: –
- MATLAB installed in a Computer System.
Diagram: –
Theory: –
Procedures: –
MATLAB Program:
Bus Data: –
% Bus data for Load Flow Analysis.
function busdata = busdata6()
busdata = [ 1 1 1.05 0 0.0 0 0 0 0 0;
2 2 1.05 0 0.5 0 0 0 -0.5 1.0;
3 2 1.07 0 0.6 0 0 0 -0.5 1.5;
4 3 1.0 0 0.0 0 0.7 0.7 0 0;
5 3 1.0 0 0.0 0 0.7 0.7 0 0;
6 3 1.0 0 0.0 0 0.7 0.7 0 0 ];
Line Data: –
% Line Data for Y-Bus Formation.
function linedata = linedata6()
linedata = [ 1 2 0.10 0.20 0.02;
1 4 0.05 0.20 0.02;
1 5 0.08 0.30 0.03;
2 3 0.05 0.25 0.03;
2 4 0.05 0.10 0.01;
2 5 0.10 0.30 0.02;
2 6 0.07 0.20 0.025;
3 5 0.12 0.26 0.025;
3 6 0.02 0.10 0.01;
4 5 0.20 0.40 0.04;
5 6 0.10 0.30 0.03;];
Ybus program:
%Y BUS MATRIX FORMATION
function ybus = ybusppg();
linedata = linedata6();
fb = linedata(:,1);
tb = linedata(:,2);
r = linedata(:,3);
x = linedata(:,4);
b = linedata(:,5);
z = r + i*x;
y = 1./z;
b = i*b;
nbus = max(max(fb),max(tb));
nbranch = length(fb);
ybus = zeros(nbus,nbus);
for k=1:nbranch
ybus(fb(k),tb(k)) = -y(k);
ybus(tb(k),fb(k)) = ybus(fb(k),tb(k));
end
for m=1:nbus
for n=1:nbranch
if fb(n) == m | tb(n) == m
ybus(m,m) = ybus(m,m) + y(n) + b(n);
end
end
end
ybus;
zbus = inv(ybus);
disp(y)
Newton Raphson:
%NEWTON RAPHSON LOAD FLOW ANALYSIS
nbus = 6;
ybus = ybusppg();
busd = busdata6();
BMva = 100;
bus = busd(:,1);
type = busd(:,2);
V = busd(:,3);
del = busd(:,4);
Pg = busd(:,5)/BMva;
Qg = busd(:,6)/BMva;
Pl = busd(:,7)/BMva;
Ql = busd(:,8)/BMva;
Qmin = busd(:,9)/BMva;
Qmax = busd(:,10)/BMva;
P = Pg – Pl;
Q = Qg – Ql;
Psp = P;
Qsp = Q;
G = real(ybus);
B = imag(ybus);
pv = find(type == 2 | type == 1);
pq = find(type == 3);
nbus=max(nbus)
npv = length(pv);
npq = length(pq);
Tol = 1;
Iter = 1;
while (Tol > 1e-5)
P = zeros(nbus,1);
Q = zeros(nbus,1);
for i = 1:nbus
for k = 1:nbus
P(i) = P(i) + V(i)* V(k)*(G(i,k)*cos(del(i)-del(k)) + B(i,k)*sin(del(i)-del(k)));
Q(i) = Q(i) + V(i)* V(k)*(G(i,k)*sin(del(i)-del(k)) - B(i,k)*cos(del(i)-del(k)));
end
end
if Iter <= 7 && Iter > 2
for n = 2:nbus
if type(n) == 2
QG = Q(n)+Ql(n);
if QG < Qmin(n)
V(n) = V(n) + 0.01;
elseif QG > Qmax(n)
V(n) = V(n) - 0.01;
end
end
end
end
dPa = Psp-P;
dQa = Qsp-Q;
k = 1;
dQ = zeros(npq,1);
for i = 1:nbus
if type(i) == 3
dQ(k,1) = dQa(i);
k = k+1;
end
end
dP = dPa(2:nbus);
M = [dP; dQ];
J1 = zeros(nbus-1,nbus-1);
for i = 1:(nbus-1)
m = i+1;
for k = 1:(nbus-1)
n = k+1;
if n == m
for n = 1:nbus
J1(i,k) = J1(i,k) + V(m)* V(n)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n)));
end
J1(i,k) = J1(i,k) - V(m)^2*B(m,m);
else
J1(i,k) = V(m)* V(n)*(G(m,n)*sin(del(m)-del(n)) - B(m,n)*cos(del(m)-del(n)));
end
end
end
J2 = zeros(nbus-1,npq);
for i = 1:(nbus-1)
m = i+1;
for k = 1:npq
n = pq(k);
if n == m
for n = 1:nbus
J2(i,k) = J2(i,k) + V(n)*(G(m,n)*cos(del(m)-del(n)) + B(m,n)*sin(del(m)-del(n)));
end
J2(i,k) = J2(i,k) + V(m)*G(m,m);
else
J2(i,k) = V(m)*(G(m,n)*cos(del(m)-del(n)) + B(m,n)*sin(del(m)-del(n)));
end
end
end
J3 = zeros(npq,nbus-1);
for i = 1:npq
m = pq(i);
for k = 1:(nbus-1)
n = k+1;
if n == m
for n = 1:nbus
J3(i,k) = J3(i,k) + V(m)* V(n)*(G(m,n)*cos(del(m)-del(n)) + B(m,n)*sin(del(m)-del(n)));
end
J3(i,k) = J3(i,k) - V(m)^2*G(m,m);
else
J3(i,k) = V(m)* V(n)*(-G(m,n)*cos(del(m)-del(n)) - B(m,n)*sin(del(m)-del(n)));
end
end
end
J4 = zeros(npq,npq);
for i = 1:npq
m = pq(i);
for k = 1:npq
n = pq(k);
if n == m
for n = 1:nbus
J4(i,k) = J4(i,k) + V(n)*(G(m,n)*sin(del(m)-del(n)) - B(m,n)*cos(del(m)-del(n)));
end
J4(i,k) = J4(i,k) - V(m)*B(m,m);
else
J4(i,k) = V(m)*(G(m,n)*sin(del(m)-del(n)) - B(m,n)*cos(del(m)-del(n)));
end
end
end
J = [J1 J2; J3 J4];
X = inv(J)*M;
dTh = X(1:nbus-1);
dV = X(nbus:end);
del(2:nbus) = dTh + del(2:nbus);
k = 1;
for i = 2:nbus
if type(i) == 3
V(i) = dV(k) + V(i);
k = k+1;
end
end
Iter = Iter + 1;
Tol = max(abs(M)); end
Iter
V
del
Conclusion: –
To be written by student.