-
Notifications
You must be signed in to change notification settings - Fork 2
/
loss_allocation.m
73 lines (72 loc) · 2.95 KB
/
loss_allocation.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
function [ds,LA,gen_bus] = loss_allocation(input_file,epsilon,iter_max,save_out)
%LOSS_ALLOCATION Calculates loss allocation for each pair load-generator.
%
% [ds,LA,gen_bus] = loss_allocation(input_file,epsilon,iter_max,save_out)
% [ds,LA,gen_bus] = loss_allocation(input_file)
%
% Inputs:
% input_file : String containing the name of the file with the case data (default is 'case5')
% epsilon : Termination tolerance on per unit voltage mismatch (default is 1e-8)
% iter_max : Maximum number of iterations
% save_out : Logical variable which indicates whether output file should be written
%
% Outputs:
% ds : Structure with the results with the followind fields:
% Uslack, Ubase, Sbase, branch, gen : data from the input case
% U : Bus voltages
% Sg, gen_bus : Generator bus indicies and power outputs
% Sf, St : Branch flows at "from" and "to" ends
% iter, time : Number of iterations and calculation time
% Sd, Yd : Apparent load demand and shunt admittances at buses (includes capacitors
% and branch shunt admittances)
% Ybranch : Branch shunt admittances
% input_file : Input file name
% LA : Matrix of allocated losses for each pair load-generator.
% LA(i,j) is network loss allocated to load at bus "i" due to power supply from generator "j"
% gen_bus : Vector with generator buses
%
% See also DIST_PF, TRACE_FLOW.
%% default arguments
if nargin < 4
save_out = true;
if nargin < 3
iter_max = 20;
if nargin < 2
epsilon = 1e-8;
if nargin < 1
input_file = 'case5';
end
end
end
end
%% power flow
ds = dist_pf(input_file,epsilon,iter_max);
Ibase = ds.Sbase*1000/ds.Ubase;
[U, Sg, gen_bus, Sf, St, Sd, Yd] = deal(ds.U,ds.Sg,ds.gen_bus,ds.Sf,ds.St,ds.Sd,ds.Yd);
[f, t] = deal(ds.branch(:,1),ds.branch(:,2));
If = conj(Sf./U(f)) - ds.Ybranch/2.*U(f); % substract branch charge currents so that If = It
It = conj(St./U(t)) + ds.Ybranch/2.*U(2:end);
NB = length(U);
NG = length(gen_bus);
%% loss allocation
Id = conj(Sd./U); % load currents
Ic = Yd.*U; % branch charge currents
Id = Id + Ic;
Ig = conj(Sg./U(gen_bus)); % generator currents
% calculate components of load currents by matrix method
Jr = trace_flow(real(Id),f,t,real(If),gen_bus,real(Ig)); % real parts
Ji = trace_flow(-imag(Id),f,t,-imag(If),gen_bus,-imag(Ig)); % imaginary parts (we take negative since the imaginary parts are lagging and are negative at loads, therefore we make them positive)
J = Jr - 1j*Ji;
ds.If = If*Ibase;
ds.It = It*Ibase;
ds.Id = Id*Ibase;
ds.Ig = Ig*Ibase;
ds.J = J*Ibase;
% calculate voltage difference between each bus and each generator bus
DU = ones(NB,1)*U(gen_bus).' - repmat(U,1,NG);
% calculate loss allocation for supplying load by each generator
LA = DU .* conj(J(:,gen_bus)) * ds.Sbase * 1000;
%% solution print
if save_out
print_pf(ds,LA);
end