Load flow analysis of a given network using MATLAB coding
We will try to code for the given network.
We are having two tasks:
1) Formation of Y bus and Z Bus matrices
for given network
2) Programming of power flow using Gauss Seidel Method
P1.
Program to form Admittance and Impedance Bus formation
Matlab Coding:
clc
clear all
disp('-----------Y BUS
Formation------------');
x=input('Enter the number of nodes: ');
for i=1:1:x
for
j=1:1:x
if(i==j)
a(i,j)=input(strcat('Enter the valur of admittance Y',int2str(i),int2str(0),':
'));
else
a(i,j)=input(strcat('Enter the valur of admittance Y',int2str(i),int2str(j),':
'));
end
end
end
b=a;y=0;
for i=1:1:x
for
j=1:1:x
if
(i==j)
for
k=1:1:x
y=y+b(i,k);
end
a(i,j)=y;
y=0;
else
a(i,j)= -b(i,j);
end
end
end
b
YBUS=a
zbus=inv(YBUS)
Output:
P2.
Programming of power flow using Gauss Seidel Method
Matlab Coding:
function linedata = linedata6() % Returns linedata.
%
| From | To
| R |
X | B/2 |
%
| Bus | Bus
| | |
|
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;];
function busdata = busdata6() % Returns busdata.
%
|Bus | Type | Vsp | theta | PGi | QGi | PLi | QLi | Qmin | Qmax |
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 ];
ybus = ybusppg(); % Calling
program "ybusppg.m" to get Y-Bus.
busdata = busdata6(); % Calling
"busdata6.m" for bus data.
bus = busdata(:,1); % Bus
number.
type = busdata(:,2); % Type of
Bus 1-Slack, 2-PV, 3-PQ.
V = busdata(:,3); % Initial
Bus Voltages.
th = busdata(:,4); % Initial
Bus Voltage Angles.
GenMW = busdata(:,5); % PGi,
Real Power injected into the buses.
GenMVAR = busdata(:,6); % QGi, Reactive Power
injected into the buses.
LoadMW = busdata(:,7); % PLi, Real Power Drawn
from the buses.
LoadMVAR = busdata(:,8); % QLi, Reactive Power Drawn
from the buses.
Qmin = busdata(:,9); % Minimum
Reactive Power Limit
Qmax = busdata(:,10); % Maximum
Reactive Power Limit
nbus = max(bus); % To get
no. of buses
P = GenMW - LoadMW; % Pi = PGi
- PLi, Real Power at the buses.
Q = GenMVAR - LoadMVAR; % Qi = QGi - QLi, Reactive
Power at the buses.
Vprev = V;
toler = 1; %
Tolerence.
iteration = 1; %
iteration starting
while (toler > 0.000001) % Start of while loop
for i =
2:nbus
sumyv
= 0;
for k =
1:nbus
if i ~=
k
sumyv = sumyv + ybus(i,k)* V(k); %
Vk * Yik
end
end
if
type(i) == 2 %
Computing Qi for PV bus
Q(i) = -imag(conj(V(i))*(sumyv + ybus(i,i)*V(i)));
if
(Q(i) > Qmax(i)) || (Q(i) < Qmin(i))
% Checking for Qi Violation.
if Q(i) < Qmin(i) %
Whether violated the lower limit.
Q(i) = Qmin(i);
else %
No, violated the upper limit.
Q(i) = Qmax(i);
end
type(i) = 3; %
If Violated, change PV bus to PQ bus.
end
end
V(i) =
(1/ybus(i,i))*((P(i)-j*Q(i))/conj(V(i)) - sumyv); % Compute
Bus Voltages.
if
type(i) == 2 % For PV Buses, Voltage Magnitude remains same, but Angle
changes.
V(i) = pol2rect(abs(Vprev(i)), angle(V(i)));
end
end
iteration
= iteration + 1; %
Increment iteration count.
toler =
max(abs(abs(V) - abs(Vprev))); %
Calculate tolerance.
Vprev = V;
% Vprev is required for next iteration, V(i) = pol2rect(abs(Vprev(i)), angle(V(i)));
end
% End of while loop / Iteration
iteration
% Total iterations.
V
% Bus Voltages in Complex form.
Vmag = abs(V)
% Final Bus Voltages.
Ang = 180/pi*angle(V) % Final Bus Voltage Angles
in Degree.
gauss
ans =
1.0000 2.0000 0.1000
0.2000 0.0200
1.0000 4.0000 0.0500
0.2000 0.0200
1.0000 5.0000 0.0800
0.3000 0.0300
2.0000 3.0000 0.0500
0.2500 0.0300
2.0000 4.0000 0.0500
0.1000 0.0100
2.0000
5.0000 0.1000 0.3000
0.0200
2.0000 6.0000 0.0700
0.2000 0.0250
3.0000 5.0000 0.1200
0.2600 0.0250
3.0000 6.0000 0.0200
0.1000 0.0100
4.0000 5.0000 0.2000
0.4000 0.0400
5.0000 6.0000 0.1000
0.3000 0.0300
Comments
Post a Comment