Experiment No.: – 10
Aim of the Experiment: –
Study on AC load flow using Gauss-seidel method.
Apparatus Required: –
Computer System with MATLAB installed.
Diagram: –
Theory: –
MATLAB Program: –
Line data:
% Line Data for Y-Bus Formation.
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;];
y bus:
% Program to form Admittance And Impedance Bus Formation….
function ybus = ybusppg(); % Returns ybus
linedata = linedata6(); % Calling “linedata6.m” for Line Data…
fb = linedata(:,1); % From bus number…
tb = linedata(:,2); % To bus number…
r = linedata(:,3); % Resistance, R…
x = linedata(:,4); % Reactance, X…
b = linedata(:,5); % Ground Admittance, B/2…
z = r + i*x; % Z matrix…
y = 1./z; % To get inverse of each element…
b = i*b; % Make B imaginary…
nbus = max(max(fb),max(tb)); % no. of buses…
nbranch = length(fb); % no. of branches…
ybus = zeros(nbus,nbus); % Initialise YBus…
% Formation of the Off Diagonal Elements…
for k=1:nbranch
ybus(fb(k),tb(k)) = -y(k);
ybus(tb(k),fb(k)) = ybus(fb(k),tb(k));
end
% Formation of Diagonal Elements….
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; % Bus Admittance Matrix
zbus = inv(ybus); % Bus Impedance Matrix
disp(y)
Polar to rectangle:
% Polar to Rectangular Conversion
function rect = pol2rect(r,o) % r = magnitude, o = angle in radians.
rect = rcos(o) + jrsin(o); % rect = real + jimag
Gauss-seidel: –
% Program for Gauss – Seidel Load Flow Analysis
% Assumption, Bus 1 is considered as Slack bus.
clc;
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.00001) % 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)-jQ(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.