Skip to content

Commit 1cf19d4

Browse files
rdzmanrdzman
committed
Use opt_model/get_soln() to extract results
rather than doing the indexing manually.
1 parent a43f587 commit 1cf19d4

3 files changed

Lines changed: 54 additions & 49 deletions

File tree

‎CHANGES.md‎

Lines changed: 7 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -5,10 +5,15 @@ Change history for MOST
55
Changes since 1.0.2
66
-------------------
77

8+
*Requires MATPOWER with 7.1 or later (for MP-Opt-Model 2.2 or later).*
9+
10+
#### 9/9/20
11+
- Use `@opt_model/get_soln()` to extract variable and shadow price
12+
results, rather than doing the indexing manually.
13+
814
#### 3/19/20
915
- Convert to using `@opt_model/solve()` method rather than calling
10-
`miqps_matpower()` or `qps_matpower()` directly. Requires MATPOWER
11-
7.1 or later.
16+
`miqps_matpower()` or `qps_matpower()` directly.
1217
- **INCOMPATIBLE CHANGE**: Update objective function value returned in
1318
`mdo.QP.f` to include the previously missing constant term.
1419

‎docs/src/MOST-manual/MOST-manual.tex‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3353,6 +3353,7 @@ \subsubsection*{Other Changes}
33533353
\item Update \code{most\_summary()} to include sections for fixed loads
33543354
and storage expected stored energy.
33553355
\item Convert to using \code{@opt\_model/solve()} method rather than calling \code{miqps\_matpower()} or \code{qps\_matpower()} directly. Requires \matpower{} \hl{7.1} or later.
3356+
\item Use \code{@opt\_model/get\_soln()} to extract variable and shadow price results, rather than doing the indexing manually. Requires \matpower{} \hl{7.1} or later.
33563357
\end{itemize}
33573358

33583359
\subsubsection*{Incompatible Changes}

‎lib/most.m‎

Lines changed: 46 additions & 47 deletions
Original file line numberDiff line numberDiff line change
@@ -1980,7 +1980,6 @@
19801980
mdo = mdi; %% initialize output
19811981
mdo.om = om.copy(); %% make copy of opt_model object, so changes to
19821982
%% output obj (mdo) don't modify input obj (mdi)
1983-
[vv, ll] = mdo.om.get_idx();
19841983
if mpopt.most.solve_model
19851984
%% check consistency of model options (in case mdi was built in previous call)
19861985
if mdi.DCMODEL ~= mo.DCMODEL
@@ -2037,9 +2036,10 @@
20372036
fprintf('- Post-processing results.\n');
20382037
end
20392038
if success
2039+
om = mdo.om;
20402040
for t = 1:nt
20412041
if UC
2042-
mdo.UC.CommitSched(:, t) = mdo.QP.x(vv.i1.u(t):vv.iN.u(t));
2042+
mdo.UC.CommitSched(:, t) = om.get_soln('var', 'u', {t});
20432043
end
20442044
for j = 1:mdi.idx.nj(t)
20452045
for k = 1:mdi.idx.nc(t,j)+1
@@ -2049,7 +2049,7 @@
20492049
mpc.bus(:, VM) = 1;
20502050
end
20512051
% Injections and shadow prices
2052-
mpc.gen(:, PG) = baseMVA * mdo.QP.x(vv.i1.Pg(t,j,k):vv.iN.Pg(t,j,k));
2052+
mpc.gen(:, PG) = baseMVA * om.get_soln('var', 'Pg', {t,j,k});
20532053
%% need to update Qg for loads consistent w/constant power factor
20542054
Pmin = mpc.gen(:, PMIN);
20552055
Qmin = mpc.gen(:, QMIN);
@@ -2059,10 +2059,11 @@
20592059
mpc.gen(ivl, QG) = mpc.gen(ivl, PG) .* Qlim ./ Pmin(ivl);
20602060
if mdo.DCMODEL
20612061
%% bus angles
2062-
mpc.bus(:, VA) = (180/pi) * mdo.QP.x(vv.i1.Va(t,j,k):vv.iN.Va(t,j,k));
2062+
Va = om.get_soln('var', 'Va', {t,j,k});
2063+
mpc.bus(:, VA) = (180/pi) * Va;
20632064

20642065
%% nodal prices
2065-
price = (mdo.QP.lambda.mu_u(ll.i1.Pmis(t,j,k):ll.iN.Pmis(t,j,k))-mdo.QP.lambda.mu_l(ll.i1.Pmis(t,j,k):ll.iN.Pmis(t,j,k))) / baseMVA;
2066+
price = (om.get_soln('lin', 'mu_u', 'Pmis', {t,j,k})-om.get_soln('lin', 'mu_l', 'Pmis', {t,j,k})) / baseMVA;
20662067
mpc.bus(:, LAM_P) = price;
20672068

20682069
%% line flows and line limit shadow prices
@@ -2073,50 +2074,48 @@
20732074
mpc.branch(:, MU_SF) = 0;
20742075
mpc.branch(:, MU_ST) = 0;
20752076
ion = find(mpc.branch(:, BR_STATUS));
2076-
rows = ll.i1.Pf(t,j,k):ll.iN.Pf(t,j,k);
2077-
cols = vv.i1.Va(t,j,k):vv.iN.Va(t,j,k);
2078-
lf = baseMVA * (mdo.QP.A(rows,cols) * mdo.QP.x(cols) + mdo.flow(t,j,k).PLsh);
2077+
AA = om.params_lin_constraint('Pf', {t,j,k});
2078+
lf = baseMVA * (AA * Va + mdo.flow(t,j,k).PLsh);
20792079
mpc.branch(ion, PF) = lf;
20802080
mpc.branch(ion, PT) = -lf;
2081-
mpc.branch(ion, MU_SF) = mdo.QP.lambda.mu_u(rows) / baseMVA;
2082-
mpc.branch(ion, MU_ST) = mdo.QP.lambda.mu_l(rows) / baseMVA;
2081+
mpc.branch(ion, MU_SF) = om.get_soln('lin', 'mu_u', 'Pf', {t,j,k}) / baseMVA;
2082+
mpc.branch(ion, MU_ST) = om.get_soln('lin', 'mu_l', 'Pf', {t,j,k}) / baseMVA;
20832083
else
20842084
%% system price
2085-
price = (mdo.QP.lambda.mu_l(ll.i1.Pmis(t,j,k):ll.iN.Pmis(t,j,k))-mdo.QP.lambda.mu_u(ll.i1.Pmis(t,j,k):ll.iN.Pmis(t,j,k))) / baseMVA;
2085+
price = (om.get_soln('lin', 'mu_l', 'Pmis', {t,j,k})-om.get_soln('lin', 'mu_u', 'Pmis', {t,j,k})) / baseMVA;
20862086
mpc.bus(:, LAM_P) = price;
20872087
end
20882088
if UC
20892089
% igenon does not contain gens ousted because of a contingency or
20902090
% a forced-off UC.CommitKey
20912091
igenon = find(mpc.gen(:, GEN_STATUS));
2092-
u = mdo.QP.x(vv.i1.u(t):vv.iN.u(t));
2093-
mpc.gen(igenon, GEN_STATUS) = u(igenon);
2092+
mpc.gen(igenon, GEN_STATUS) = mdo.UC.CommitSched(igenon, t);
20942093
gs = mpc.gen(igenon, GEN_STATUS) > 0; % gen status
20952094
mpc.gen(:, MU_PMAX) = 0;
20962095
mpc.gen(:, MU_PMIN) = 0;
20972096
mpc.gen(igenon, MU_PMAX) = gs .* ...
2098-
mdo.QP.lambda.mu_u(ll.i1.uPmax(t,j,k):ll.iN.uPmax(t,j,k)) / baseMVA;
2097+
om.get_soln('lin', 'mu_u', 'uPmax', {t,j,k}) / baseMVA;
20992098
mpc.gen(igenon, MU_PMIN) = gs .* ...
2100-
mdo.QP.lambda.mu_u(ll.i1.uPmin(t,j,k):ll.iN.uPmin(t,j,k)) / baseMVA;
2099+
om.get_soln('lin', 'mu_u', 'uPmin', {t,j,k}) / baseMVA;
21012100
if mdo.QCoordination
21022101
mpc.gen(:, MU_QMAX) = 0;
21032102
mpc.gen(:, MU_QMIN) = 0;
21042103
mpc.gen(igenon, MU_QMAX) = gs .* ...
2105-
mdo.QP.lambda.mu_u(ll.i1.uQmax(t,j,k):ll.iN.uQmax(t,j,k)) / baseMVA;
2104+
om.get_soln('lin', 'mu_u', 'uQmax', {t,j,k}) / baseMVA;
21062105
mpc.gen(igenon, MU_QMIN) = gs .* ...
2107-
mdo.QP.lambda.mu_u(ll.i1.uQmin(t,j,k):ll.iN.uQmin(t,j,k)) / baseMVA;
2106+
om.get_soln('lin', 'mu_u', 'uQmin', {t,j,k}) / baseMVA;
21082107
end
21092108
else
21102109
gs = mpc.gen(:, GEN_STATUS) > 0; % gen status
21112110
mpc.gen(:, MU_PMAX) = gs .* ...
2112-
mdo.QP.lambda.upper(vv.i1.Pg(t,j,k):vv.iN.Pg(t,j,k)) / baseMVA;
2111+
om.get_soln('var', 'mu_u', 'Pg', {t,j,k}) / baseMVA;
21132112
mpc.gen(:, MU_PMIN) = gs .* ...
2114-
mdo.QP.lambda.lower(vv.i1.Pg(t,j,k):vv.iN.Pg(t,j,k)) / baseMVA;
2113+
om.get_soln('var', 'mu_l', 'Pg', {t,j,k}) / baseMVA;
21152114
if mdo.QCoordination
21162115
mpc.gen(:, MU_QMAX) = gs .* ...
2117-
mdo.QP.lambda.upper(vv.i1.Qg(t,j,k):vv.iN.Qg(t,j,k)) / baseMVA;
2116+
om.get_soln('var', 'mu_u', 'Qg', {t,j,k}) / baseMVA;
21182117
mpc.gen(:, MU_QMIN) = gs .* ...
2119-
mdo.QP.lambda.lower(vv.i1.Qg(t,j,k):vv.iN.Qg(t,j,k)) / baseMVA;
2118+
om.get_soln('var', 'mu_l', 'Qg', {t,j,k}) / baseMVA;
21202119
end
21212120
end
21222121
if mdi.IncludeFixedReserves
@@ -2125,34 +2124,34 @@
21252124
r.R = z;
21262125
r.prc = z;
21272126
r.mu = struct('l', z, 'u', z, 'Pmax', z);
2128-
r.totalcost = sum(mdo.om.eval_quad_cost(mdo.QP.x, 'Rcost', {t,j,k}));
2129-
r.R(r.igr) = mdo.QP.x(vv.i1.R(t,j,k):vv.iN.R(t,j,k)) * baseMVA;
2127+
r.totalcost = sum(om.get_soln('qdc', 'Rcost', {t,j,k}));
2128+
r.R(r.igr) = om.get_soln('var', 'R', {t,j,k}) * baseMVA;
2129+
R_mu_l = om.get_soln('lin', 'mu_l', 'Rreq', {t,j,k});
21302130
for gg = r.igr
21312131
iz = find(r.zones(:, gg));
2132-
kk = ll.i1.Rreq(t,j,k):ll.iN.Rreq(t,j,k);
2133-
r.prc(gg) = sum(mdo.QP.lambda.mu_l(kk(iz))) / baseMVA;
2132+
r.prc(gg) = sum(R_mu_l(iz)) / baseMVA;
21342133
end
2135-
r.mu.l(r.igr) = mdo.QP.lambda.lower(vv.i1.R(t,j,k):vv.iN.R(t,j,k)) / baseMVA;
2136-
r.mu.u(r.igr) = mdo.QP.lambda.upper(vv.i1.R(t,j,k):vv.iN.R(t,j,k)) / baseMVA;
2137-
r.mu.Pmax(r.igr) = mdo.QP.lambda.mu_u(ll.i1.Pg_plus_R(t,j,k):ll.iN.Pg_plus_R(t,j,k)) / baseMVA;
2134+
r.mu.l(r.igr) = om.get_soln('var', 'mu_l', 'R', {t,j,k}) / baseMVA;
2135+
r.mu.u(r.igr) = om.get_soln('var', 'mu_u', 'R', {t,j,k}) / baseMVA;
2136+
r.mu.Pmax(r.igr) = om.get_soln('lin', 'mu_u', 'Pg_plus_R', {t,j,k}) / baseMVA;
21382137
mpc.reserves = r;
21392138
end
21402139
mdo.flow(t,j,k).mpc = mpc; %% stash modified mpc in output struct
21412140
end
21422141
end
21432142
% Contract, contingency reserves, energy limits
2144-
mdo.results.Pc(:,t) = baseMVA * mdo.QP.x(vv.i1.Pc(t):vv.iN.Pc(t));
2145-
mdo.results.Rpp(:,t) = baseMVA * mdo.QP.x(vv.i1.Rpp(t):vv.iN.Rpp(t));
2146-
mdo.results.Rpm(:,t) = baseMVA * mdo.QP.x(vv.i1.Rpm(t):vv.iN.Rpm(t));
2143+
mdo.results.Pc(:,t) = baseMVA * om.get_soln('var', 'Pc', {t});
2144+
mdo.results.Rpp(:,t) = baseMVA * om.get_soln('var', 'Rpp', {t});
2145+
mdo.results.Rpm(:,t) = baseMVA * om.get_soln('var', 'Rpm', {t});
21472146
if ns
2148-
mdo.results.Sm(:,t) = baseMVA * mdo.QP.x(vv.i1.Sm(t):vv.iN.Sm(t));
2149-
mdo.results.Sp(:,t) = baseMVA * mdo.QP.x(vv.i1.Sp(t):vv.iN.Sp(t));
2147+
mdo.results.Sm(:,t) = baseMVA * om.get_soln('var', 'Sm', {t});
2148+
mdo.results.Sp(:,t) = baseMVA * om.get_soln('var', 'Sp', {t});
21502149
end
21512150
end
21522151
% Ramping reserves
21532152
for t = 1:mdo.idx.ntramp
2154-
mdo.results.Rrp(:,t) = baseMVA * mdo.QP.x(vv.i1.Rrp(t):vv.iN.Rrp(t));
2155-
mdo.results.Rrm(:,t) = baseMVA * mdo.QP.x(vv.i1.Rrm(t):vv.iN.Rrm(t));
2153+
mdo.results.Rrp(:,t) = baseMVA * om.get_soln('var', 'Rrp', {t});
2154+
mdo.results.Rrm(:,t) = baseMVA * om.get_soln('var', 'Rrm', {t});
21562155
end
21572156
% Expected energy prices for generators, per generator and per period,
21582157
% both absolute and conditional on making it to that period
@@ -2172,13 +2171,13 @@
21722171
mdo.results.RppPrices = zeros(ng, nt);
21732172
mdo.results.RpmPrices = zeros(ng, nt);
21742173
for t = 1:nt
2175-
mdo.results.RppPrices(:, t) = mdo.QP.lambda.lower(vv.i1.Rpp(t):vv.iN.Rpp(t)) / baseMVA;
2176-
mdo.results.RpmPrices(:, t) = mdo.QP.lambda.lower(vv.i1.Rpm(t):vv.iN.Rpm(t)) / baseMVA;
2174+
mdo.results.RppPrices(:, t) = om.get_soln('var', 'mu_l', 'Rpp', {t}) / baseMVA;
2175+
mdo.results.RpmPrices(:, t) = om.get_soln('var', 'mu_l', 'Rpm', {t}) / baseMVA;
21772176
for j = 1:mdi.idx.nj(t);
21782177
for k = 1:mdi.idx.nc(t,j)+1
21792178
ii = find(mdi.flow(t,j,k).mpc.gen(:, GEN_STATUS) > 0);
2180-
mdo.results.RppPrices(ii, t) = mdo.results.RppPrices(ii, t) + mdo.QP.lambda.mu_l(ll.i1.dPpRp(t,j,k):ll.iN.dPpRp(t,j,k)) / baseMVA;
2181-
mdo.results.RpmPrices(ii, t) = mdo.results.RpmPrices(ii, t) + mdo.QP.lambda.mu_l(ll.i1.dPmRm(t,j,k):ll.iN.dPmRm(t,j,k)) / baseMVA;
2179+
mdo.results.RppPrices(ii, t) = mdo.results.RppPrices(ii, t) + om.get_soln('lin', 'mu_l', 'dPpRp', {t,j,k}) / baseMVA;
2180+
mdo.results.RpmPrices(ii, t) = mdo.results.RpmPrices(ii, t) + om.get_soln('lin', 'mu_l', 'dPmRm', {t,j,k}) / baseMVA;
21822181
end
21832182
end
21842183
end
@@ -2190,17 +2189,17 @@
21902189
for j1 = 1:mdo.idx.nj(t)
21912190
for j2 = 1:mdo.idx.nj(t+1)
21922191
if mdi.tstep(t+1).TransMask(j2,j1)
2193-
mdo.results.RrpPrices(:, t) = mdo.results.RrpPrices(:, t) + mdo.QP.lambda.mu_l(ll.i1.Rrp(t,j1,j2):ll.iN.Rrp(t,j1,j2)) / baseMVA;
2194-
mdo.results.RrmPrices(:, t) = mdo.results.RrmPrices(:, t) + mdo.QP.lambda.mu_l(ll.i1.Rrm(t,j1,j2):ll.iN.Rrm(t,j1,j2)) / baseMVA;
2192+
mdo.results.RrpPrices(:, t) = mdo.results.RrpPrices(:, t) + om.get_soln('lin', 'mu_l', 'Rrp', {t,j1,j2}) / baseMVA;
2193+
mdo.results.RrmPrices(:, t) = mdo.results.RrmPrices(:, t) + om.get_soln('lin', 'mu_l', 'Rrm', {t,j1,j2}) / baseMVA;
21952194
end
21962195
end
21972196
end
21982197
end
21992198
% then last period only if specified for with terminal state
22002199
if ~mdo.OpenEnded
22012200
for j1 = 1:mdo.idx.nj(nt)
2202-
mdo.results.RrpPrices(:, nt) = mdo.results.RrpPrices(:, nt) + mdo.QP.lambda.mu_l(ll.i1.Rrp(nt,j1,1):ll.iN.Rrp(nt,j1,1)) / baseMVA;
2203-
mdo.results.RrmPrices(:, nt) = mdo.results.RrmPrices(:, nt) + mdo.QP.lambda.mu_l(ll.i1.Rrm(nt,j1,1):ll.iN.Rrm(nt,j1,1)) / baseMVA;
2201+
mdo.results.RrpPrices(:, nt) = mdo.results.RrpPrices(:, nt) + om.get_soln('lin', 'mu_l', 'Rrp', {nt,j1,1}) / baseMVA;
2202+
mdo.results.RrmPrices(:, nt) = mdo.results.RrmPrices(:, nt) + om.get_soln('lin', 'mu_l', 'Rrm', {nt,j1,1}) / baseMVA;
22042203
end
22052204
end
22062205
% Expected wear and tear costs per gen and period
@@ -2240,7 +2239,7 @@
22402239
end
22412240
% If Cyclic storage, pull InitialStorage value out of x
22422241
if ns && mdo.Storage.ForceCyclicStorage
2243-
mdo.Storage.InitialStorage = baseMVA * mdo.QP.x(vv.i1.S0:vv.iN.S0);
2242+
mdo.Storage.InitialStorage = baseMVA * om.get_soln('var', 'S0');
22442243
end
22452244
% Compute expected storage state trajectory
22462245
mdo.Storage.ExpectedStorageState = zeros(ns,nt);
@@ -2267,14 +2266,14 @@
22672266
if nzds
22682267
mdo.results.Z = zeros(nzds, ntds);
22692268
for t = 1:ntds
2270-
mdo.results.Z(:,t) = mdo.QP.x(vv.i1.Z(t):vv.iN.Z(t));
2269+
mdo.results.Z(:,t) = om.get_soln('var', 'Z', {t});
22712270
end
22722271
end
22732272
mdo.results.Y = zeros(nyds, ntds);
22742273
if nyds
22752274
for t = 1:ntds
2276-
mdo.results.Y(:, t) = ...
2277-
mdo.QP.A(ll.i1.DSy(t):ll.iN.DSy(t), :) * mdo.QP.x;
2275+
AA = om.params_lin_constraint('DSy', {t});
2276+
mdo.results.Y(:, t) = AA * mdo.QP.x;
22782277
end
22792278
end
22802279
end

0 commit comments

Comments
 (0)