-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathGridData118bus.m
81 lines (74 loc) · 2 KB
/
GridData118bus.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
74
75
76
77
78
79
80
81
%% mpc data
clear
clc
close all
mpc = loadcase('case118.m');
[Ngen,~] = size(mpc.gen); % Number of generators
Cgen_k = (mpc.gencost(:,5).*mpc.gen(:,9)+mpc.gencost(:,6))*mpc.baseMVA/1e4; % 1e4 $/p.u.
Ingen = mpc.gen(:,1);
Nbus = 118; % Number of buses
Iref = 1; % Ref bus
PD = mpc.bus(:, 3) / mpc.baseMVA; % P demand
Igen = zeros(Nbus, Ngen); % bus of generator
for i = 1: Ngen
Igen(mpc.gen(i, 1), i) = 1;
end
PMAXgen = mpc.gen(1: Ngen, 9) / mpc.baseMVA;
PMINgen = mpc.gen(1: Ngen, 10) / mpc.baseMVA;
Ibranch = mpc.branch(:, 1: 2); % branch: from bus, to bus
[Nbranch, ~] = size(Ibranch);
BR_R = mpc.branch(:,3);
BR_X = mpc.branch(:,4);
BR_G = real(1./(BR_R+1i*BR_X));
BR_B = imag(1./(BR_R+1i*BR_X));
Sbranch = 2*ones(Nbranch,1); % branch capacity
IFrom = zeros(Nbranch, Nbus);
ITo = zeros(Nbranch, Nbus);
for i = 1: Nbranch
IFrom(i, Ibranch(i, 1)) = 1;
ITo(i, Ibranch(i, 2)) = 1;
end
PTDF = zeros(Nbranch+Nbus-1,Nbus);
Arel = IFrom(:,2:Nbus)-ITo(:,2:Nbus);
Mrel = [Arel',zeros(Nbus-1,Nbus-1);diag(BR_X),Arel];
for i = 2:Nbus
erel = zeros(Nbranch+Nbus-1,1);
erel(i-1) = 1;
PTDF(:,i) = Mrel\erel;
end
PTDF = PTDF(1:Nbranch,:);
PTDF_l = [PTDF(1:12,:);PTDF(14:33,:);PTDF(35:41,:)];
Sbranch_l = [Sbranch(1:12);Sbranch(14:33);Sbranch(35:41)];
Num_T = 672;
Delta_t = 1;
eta_c = 0.95;
eta_d = 0.95;
w_l = 0.1;
w_h = 0.9;
Nes = 2; % number of energy storage
Ines = [15;18];
Ies = zeros(Nbus, Nes);
for i = 1:Nes
Ies(Ines(i),i) = 1;
end
ESe = [0.2;0.2];
ESp = 0.2*ESe;
load('Data118.mat','PDdata','PVdata','Winddata','PDu');
PD0 = PDdata;
NR = 2;
InR = [6;15];
IR = zeros(Nbus,NR);
for i = 1:NR
IR(InR(i),i) = 1;
end
PRMAX = [1;1];
PR = zeros(NR,Num_T);
PR(1,:) = PRMAX(1)*PVdata(1:Num_T);
PR(2,:) = PRMAX(2)*Winddata(1:Num_T);
losskappa = 0.005;
losskappa = -losskappa*ones(Nbus,1);
losskappa(Ingen) = -losskappa(Ingen);
losskappa(InR) = -losskappa(InR);
save Grid.mat
load('CCarbon118.mat','CCarbon');
save('CCarbon.mat','CCarbon');