0001 function hh = plot_storage(md, idx, varargin)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0052 MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0053 QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0054
0055
0056 my_xlabel = 'Period';
0057 nt = md.idx.nt;
0058 ns = md.idx.ns;
0059 nj_max = max(md.idx.nj);
0060 baseMVA = md.mpc.baseMVA;
0061
0062
0063 if nargin < 2
0064 idx = [];
0065 end
0066 if isempty(idx)
0067 idx = (1:ns)';
0068 end
0069 nidx = length(idx);
0070 if nidx > 1 && size(idx, 1) == 1
0071 idx = idx';
0072 end
0073 g = md.Storage.UnitIdx(idx);
0074 b = md.mpc.gen(g, GEN_BUS);
0075
0076
0077 opt = struct( ...
0078 'saveit', false, ...
0079 'saveall', false, ...
0080 'savepath', '', ...
0081 'savename', 'stored-energy-%s.pdf', ...
0082 'separation', 0.8, ...
0083 'sort_tol', 1e-6, ...
0084 'size_factor', 1, ...
0085 'show_grid', true, ...
0086 'show_adjusted_dispatches', true, ...
0087 'show_expected_initial', true, ...
0088 'show_expected_final', true, ...
0089 'show_dispatches', false );
0090
0091
0092 if mod(length(varargin), 2)
0093 if ~isstruct(varargin{1})
0094 error('plot_storage: Single OPT argument must be a struct');
0095 end
0096 myopt = varargin{1};
0097 k = 2;
0098 else
0099 myopt = struct;
0100 k = 1;
0101 end
0102 while k < length(varargin)
0103 opt_name = varargin{k};
0104 opt_val = varargin{k+1};
0105 if ~isfield(opt, opt_name)
0106 error('plot_name: ''%s'' is not a valid option name', opt_name);
0107 end
0108 myopt.(opt_name) = opt_val;
0109 k = k + 2;
0110 end
0111 fields = fieldnames(myopt);
0112 for f = 1:length(fields)
0113 opt.(fields{f}) = myopt.(fields{f});
0114 end
0115
0116
0117 if opt.saveall && nidx > 1
0118 if isempty(strfind(opt.savename, '%s'))
0119 error('plot_storage: ''savename'' must include a ''%%s'' placeholder when ''saveall'' option is true.');
0120 end
0121 for i = 1:nidx
0122 plot_storage(md, idx(i), opt, 'saveit', true);
0123 end
0124 end
0125
0126
0127 MinStorageLevel = md.Storage.MinStorageLevel;
0128 MaxStorageLevel = md.Storage.MaxStorageLevel;
0129 if size(MinStorageLevel, 1) == 1 && ns > 1
0130 MinStorageLevel = ones(ns, 1) * MinStorageLevel;
0131 end
0132 if size(MinStorageLevel, 2) == 1 && nt > 1
0133 MinStorageLevel = MinStorageLevel * ones(1, nt);
0134 end
0135 if size(MaxStorageLevel, 1) == 1 && ns > 1
0136 MaxStorageLevel = ones(ns, 1) * MaxStorageLevel;
0137 end
0138 if size(MaxStorageLevel, 2) == 1 && nt > 1
0139 MaxStorageLevel = MaxStorageLevel * ones(1, nt);
0140 end
0141 if isempty(md.Storage.InEff)
0142 InEff = 1;
0143 else
0144 InEff = md.Storage.InEff;
0145 end
0146 if size(InEff, 1) == 1 && ns > 1
0147 InEff = ones(ns, 1) * InEff;
0148 end
0149 if size(InEff, 2) == 1 && nt > 1
0150 InEff = InEff * ones(1, nt);
0151 end
0152 if isempty(md.Storage.OutEff)
0153 OutEff = 1;
0154 else
0155 OutEff = md.Storage.OutEff;
0156 end
0157 if size(OutEff, 1) == 1 && ns > 1
0158 OutEff = ones(ns, 1) * OutEff;
0159 end
0160 if size(OutEff, 2) == 1 && nt > 1
0161 OutEff = OutEff * ones(1, nt);
0162 end
0163 if isempty(md.Storage.LossFactor)
0164 LossFactor = 0;
0165 else
0166 LossFactor = md.Storage.LossFactor;
0167 end
0168 if size(LossFactor, 1) == 1 && ns > 1
0169 LossFactor = ones(ns, 1) * LossFactor;
0170 end
0171 if size(LossFactor, 2) == 1 && nt > 1
0172 LossFactor = LossFactor * ones(1, nt);
0173 end
0174 if isempty(md.Storage.rho)
0175 rho = 1;
0176 else
0177 rho = md.Storage.rho;
0178 end
0179 if size(rho, 1) == 1 && ns > 1
0180 rho = ones(ns, 1) * rho;
0181 end
0182 if size(rho, 2) == 1 && nt > 1
0183 rho = rho * ones(1, nt);
0184 end
0185
0186
0187 offset = opt.separation/2;
0188 p = (1:nt)';
0189 pp = 2*p;
0190 ppp(pp-1) = p - offset;
0191 ppp(pp ) = p + offset;
0192 Sp = md.results.Sp(idx, :);
0193 Sm = md.results.Sm(idx, :);
0194 eS = md.Storage.ExpectedStorageState(idx, :);
0195 Smin = MinStorageLevel(idx, :);
0196 Smax = MaxStorageLevel(idx, :);
0197 SSp = NaN(2*nt, 1);
0198 SSm = NaN(2*nt, 1);
0199 eeS = NaN(2*nt, 1);
0200 SSmin(pp) = sum(Smin, 1);
0201 SSmax(pp) = sum(Smax, 1);
0202 SSmin(pp-1) = SSmin(pp);
0203 SSmax(pp-1) = SSmax(pp);
0204
0205 eSi = NaN(nj_max, nt, nidx);
0206 eSf = NaN(nj_max, nt, nidx);
0207 Si_min = NaN(nj_max, nt, nidx);
0208 Si_max = NaN(nj_max, nt, nidx);
0209 dSi = zeros(nj_max, nt, nidx);
0210 Pg = zeros(nj_max, nt, nidx);
0211 jmin = zeros(1, nt);
0212 jmax = zeros(1, nt);
0213 LossCoeff = md.Delta_T * LossFactor/2;
0214 beta1 = (1-LossCoeff) ./ (1+LossCoeff);
0215 beta2 = 1 ./ (1+LossCoeff);
0216 for t = 1:nt
0217 if t == 1
0218 if md.Storage.ForceCyclicStorage
0219 SSm(t) = sum(md.Storage.InitialStorage(idx));
0220 SSp(t) = sum(md.Storage.InitialStorage(idx));
0221 else
0222 SSm(t) = sum(md.Storage.InitialStorageLowerBound(idx));
0223 SSp(t) = sum(md.Storage.InitialStorageUpperBound(idx));
0224 end
0225 eeS(t) = sum(md.Storage.InitialStorage(idx));
0226 else
0227 SSm(2*t-1) = sum(Sm(:, t-1), 1);
0228 SSp(2*t-1) = sum(Sp(:, t-1), 1);
0229 eeS(2*t-1) = sum(eS(:, t-1), 1);
0230 end
0231
0232 SSm(2*t) = sum(Sm(:, t), 1);
0233 SSp(2*t) = sum(Sp(:, t), 1);
0234 eeS(2*t) = sum(eS(:, t), 1);
0235
0236 vv = get_idx(md.om);
0237 for j = 1:md.idx.nj(t)
0238 dS = -md.Delta_T * ...
0239 (InEff(:,t) .* md.QP.x(vv.i1.Psc(t,j,1)-1+(1:ns)) + ...
0240 1./OutEff(:,t) .* md.QP.x(vv.i1.Psd(t,j,1)-1+(1:ns)) );
0241 dSi(j,t,:) = dS(idx) * baseMVA;
0242 Pg(j,t,:) = md.flow(t,j,1).mpc.gen(g, PG);
0243 Lij = md.tstep(t).Li( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :);
0244 Lfj = md.tstep(t).Lf( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :);
0245 Mj = md.tstep(t).Mg( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :) + ...
0246 md.tstep(t).Mh( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :);
0247 Nj = md.tstep(t).Ng( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :) + ...
0248 md.tstep(t).Nh( j:md.idx.nj(t):(ns-1)*md.idx.nj(t)+j, :);
0249 temp_eSi = Lij * md.Storage.InitialStorage/baseMVA + Mj * md.QP.x;
0250 temp_eSf = Lfj * md.Storage.InitialStorage/baseMVA + Nj * md.QP.x;
0251 eSi(j,t,:) = baseMVA * temp_eSi(idx);
0252 eSf(j,t,:) = baseMVA * temp_eSf(idx);
0253 if t == 1
0254 if md.Storage.ForceCyclicStorage
0255 Si_min(j,t,:) = md.Storage.InitialStorage(idx);
0256 Si_max(j,t,:) = md.Storage.InitialStorage(idx);
0257 else
0258 Si_min(j,t,:) = md.Storage.InitialStorageLowerBound(idx);
0259 Si_max(j,t,:) = md.Storage.InitialStorageUpperBound(idx);
0260 end
0261 else
0262 tmp_eSi = reshape(eSi(j,t,:), nidx, 1);
0263 Si_min(j,t,:) = rho(idx,t) .* Sm(:,t-1) + (1-rho(idx,t)) .* tmp_eSi;
0264 Si_max(j,t,:) = rho(idx,t) .* Sp(:,t-1) + (1-rho(idx,t)) .* tmp_eSi;
0265 end
0266 end
0267
0268
0269
0270
0271
0272
0273
0274
0275 tmpS = reshape(Si_min(:,t,:), nj_max, nidx);
0276 rows = opt.sort_tol * round(opt.sort_tol^-1 * [tmpS * beta1(idx,t) + reshape(dSi(:,t,:), nj_max, nidx) * beta2(idx,t) tmpS]);
0277 [junk, i] = sortrows(rows, [1 2]);
0278 jmin(t) = i(1);
0279 tmpS = reshape(Si_max(:,t,:), nj_max, nidx);
0280 rows = opt.sort_tol * round(opt.sort_tol^-1 * [tmpS * beta1(idx,t) + reshape(dSi(:,t,:), nj_max, nidx) * beta2(idx,t) tmpS]);
0281 [junk, i] = sortrows(rows, [-1 -2]);
0282 jmax(t) = i(1);
0283
0284 end
0285 SSi_min = zeros(nt,1);
0286 SSi_max = zeros(nt,1);
0287 for t = 1:nt
0288 SSi_min(t) = sum(Si_min(jmin(t), t, :), 3);
0289 SSi_max(t) = sum(Si_max(jmax(t), t, :), 3);
0290 end
0291
0292
0293
0294 plot(ppp, SSmax+0.001, 'LineStyle', 'none')
0295 hold on
0296 plot(ppp, SSmin-0.001, 'LineStyle', 'none')
0297
0298
0299 v = axis;
0300 if opt.show_grid
0301 m = v(3) - 100*(v(4) - v(3));
0302 M = v(4) + 100*(v(4) - v(3));
0303 for k = 0:nt
0304 line([k; k], [m; M], 'LineWidth', 0.25, 'Color', 0.9*[1 1 1]);
0305 end
0306 end
0307 axis(v);
0308
0309 plot(ppp, SSmin, ':k', 'LineWidth', 1);
0310 plot(ppp, SSmax, ':k', 'LineWidth', 1);
0311 for t = 1:nt
0312 for j = 1:md.idx.nj(t)
0313 if opt.show_adjusted_dispatches
0314 Si = sum(Si_min(j,t,:), 3);
0315 Sf = reshape(Si_min(j,t,:), 1, nidx) * beta1(idx, t) + ...
0316 reshape(dSi(j,t,:), 1, nidx) * beta2(idx,t);
0317 line(t+offset*[-1;1], [Si; Sf], 'Color', 0.8*[1 1 1]);
0318 Si = sum(Si_max(j,t,:), 3);
0319 Sf = reshape(Si_max(j,t,:), 1, nidx) * beta1(idx, t) + ...
0320 reshape(dSi(j,t,:), 1, nidx) * beta2(idx,t);
0321 line(t+offset*[-1;1], [Si; Sf], 'Color', 0.8*[1 1 1]);
0322 end
0323
0324 if opt.show_dispatches
0325 Si = sum(Si_min(j,t,:), 3);
0326 Sf = Si - sum(Pg(j,t,:), 3);
0327 line(t+offset*[-1;1], [Si; Sf], 'Color', 0.8*[1 0.7 1]);
0328 Si = sum(Si_max(j,t,:), 3);
0329 Sf = Si - sum(Pg(j,t,:), 3);
0330 line(t+offset*[-1;1], [Si; Sf], 'Color', 0.8*[1 0.7 1]);
0331 end
0332 end
0333 end
0334 plot(ppp, SSp, '--', 'LineWidth', 1, 'Color', 0.8*[0 1 0]);
0335 plot(ppp, SSm, '--r', 'LineWidth', 1, 'Color', 0.9*[1,0,0]);
0336 plot(p+offset, SSp(pp), 'v', 'LineWidth', 1, 'Color', 0.8*[0 1 0], 'MarkerSize', 6*opt.size_factor);
0337 plot(p+offset, SSm(pp), '^r', 'LineWidth', 1, 'Color', 0.9*[1,0,0], 'MarkerSize', 6*opt.size_factor);
0338 plot(ppp, eeS, 'Color', [0 0 1], 'LineWidth', 2);
0339 if opt.show_expected_initial
0340 if any(any(rho(idx,:) > 0))
0341 plot(p-offset, sum(Si_min, 3), '+', 'LineWidth', 1, 'Color', 0.9*[1,0,0], 'MarkerSize', 4*opt.size_factor);
0342 plot(p-offset, sum(Si_max, 3), '+', 'LineWidth', 1, 'Color', 0.8*[0,1,0], 'MarkerSize', 4*opt.size_factor);
0343 end
0344 plot(p-offset, SSi_min, '.', 'LineWidth', 1, 'MarkerSize', 13*opt.size_factor, 'Color', 0.9*[1 0 0]);
0345 plot(p-offset, SSi_max, '.', 'LineWidth', 1, 'MarkerSize', 13*opt.size_factor, 'Color', 0.8*[0 1 0]);
0346 plot(p-offset, sum(eSi, 3), 'o', 'LineWidth', 0.5, 'MarkerSize', 5*opt.size_factor);
0347 end
0348 if opt.show_expected_final
0349 plot(p+offset, sum(eSf, 3), 'x', 'LineWidth', 0.5, 'MarkerSize', 5*opt.size_factor);
0350 end
0351
0352 hold off;
0353
0354 if nidx == 1
0355 title(sprintf('Stored Energy for Storage Unit %d (Gen %d) @ Bus %d', idx, g, b), 'FontSize', 18*opt.size_factor);
0356 else
0357 if nidx == ns && all(idx == (1:ns)')
0358 txt = 'All';
0359 else
0360 txt = 'Selected';
0361 end
0362 title(sprintf('Stored Energy for %s Storage Units', txt), 'FontSize', 18*opt.size_factor);
0363 end
0364 ylabel('Energy, MWh', 'FontSize', 16*opt.size_factor);
0365 xlabel(my_xlabel, 'FontSize', 16*opt.size_factor);
0366
0367 set(gca, 'FontSize', 12*opt.size_factor);
0368 h = gcf;
0369 set(h, 'PaperOrientation', 'landscape');
0370
0371 set(h, 'PaperPosition', [0 0 11 8.5]);
0372 if opt.saveit || opt.saveall
0373 if nidx == 1
0374 txt = sprintf('%d', idx);
0375 elseif nidx == ns && all(idx == (1:ns)')
0376 txt = 'all';
0377 else
0378 txt = 'selected';
0379 end
0380 if isempty(strfind(opt.savename, '%s'))
0381 pdf_name = opt.savename;
0382 else
0383 pdf_name = sprintf(opt.savename, txt);
0384 end
0385 print('-dpdf', fullfile(opt.savepath, pdf_name));
0386 end
0387 if nargout
0388 hh = h;
0389 end