COMPUTE_COST Computes a user-defined cost.
F_U = COMPUTE_COST(OM, X)
F_U = COMPUTE_COST(OM, X, NAME)
Computes the value of a user defined cost, either for all user
defined costs or for a named set of costs. Requires calling
BUILD_COST_PARAMS first to build the full set of parameters.
Let X be the full set of optimization variables and F_U(X, CP) be the
user-defined cost at X, corresponding to the set of cost parameters in
the CP struct returned by GET_COST_PARAMS, where CP is a struct with the
following fields:
N - nw x nx sparse matrix
Cw - nw x 1 vector
H - nw x nw sparse matrix (optional, all zeros by default)
dd, mm - nw x 1 vectors (optional, all ones by default)
rh, kk - nw x 1 vectors (optional, all zeros by default)
These parameters are used as follows to compute F_U(X, CP)
R = N*x - rh
/ kk(i), R(i) < -kk(i)
K(i) = < 0, -kk(i) <= R(i) <= kk(i)
\ -kk(i), R(i) > kk(i)
RR = R + K
U(i) = / 0, -kk(i) <= R(i) <= kk(i)
\ 1, otherwise
DDL(i) = / 1, dd(i) = 1
\ 0, otherwise
DDQ(i) = / 1, dd(i) = 2
\ 0, otherwise
Dl = diag(mm) * diag(U) * diag(DDL)
Dq = diag(mm) * diag(U) * diag(DDQ)
w = (Dl + Dq * diag(RR)) * RR
F_U(X, CP) = 1/2 * w'*H*w + Cw'*w
See also OPF_MODEL, ADD_COSTS, BUILD_COST_PARAMS, GET_COST_PARAMS.0001 function f = compute_cost(om, x, name) 0002 %COMPUTE_COST Computes a user-defined cost. 0003 % F_U = COMPUTE_COST(OM, X) 0004 % F_U = COMPUTE_COST(OM, X, NAME) 0005 % 0006 % Computes the value of a user defined cost, either for all user 0007 % defined costs or for a named set of costs. Requires calling 0008 % BUILD_COST_PARAMS first to build the full set of parameters. 0009 % 0010 % Let X be the full set of optimization variables and F_U(X, CP) be the 0011 % user-defined cost at X, corresponding to the set of cost parameters in 0012 % the CP struct returned by GET_COST_PARAMS, where CP is a struct with the 0013 % following fields: 0014 % N - nw x nx sparse matrix 0015 % Cw - nw x 1 vector 0016 % H - nw x nw sparse matrix (optional, all zeros by default) 0017 % dd, mm - nw x 1 vectors (optional, all ones by default) 0018 % rh, kk - nw x 1 vectors (optional, all zeros by default) 0019 % 0020 % These parameters are used as follows to compute F_U(X, CP) 0021 % 0022 % R = N*x - rh 0023 % 0024 % / kk(i), R(i) < -kk(i) 0025 % K(i) = < 0, -kk(i) <= R(i) <= kk(i) 0026 % \ -kk(i), R(i) > kk(i) 0027 % 0028 % RR = R + K 0029 % 0030 % U(i) = / 0, -kk(i) <= R(i) <= kk(i) 0031 % \ 1, otherwise 0032 % 0033 % DDL(i) = / 1, dd(i) = 1 0034 % \ 0, otherwise 0035 % 0036 % DDQ(i) = / 1, dd(i) = 2 0037 % \ 0, otherwise 0038 % 0039 % Dl = diag(mm) * diag(U) * diag(DDL) 0040 % Dq = diag(mm) * diag(U) * diag(DDQ) 0041 % 0042 % w = (Dl + Dq * diag(RR)) * RR 0043 % 0044 % F_U(X, CP) = 1/2 * w'*H*w + Cw'*w 0045 % 0046 % See also OPF_MODEL, ADD_COSTS, BUILD_COST_PARAMS, GET_COST_PARAMS. 0047 0048 % MATPOWER 0049 % $Id: compute_cost.m,v 1.6 2011/12/01 17:03:59 cvs Exp $ 0050 % by Ray Zimmerman, PSERC Cornell 0051 % Copyright (c) 2008-2010 by Power System Engineering Research Center (PSERC) 0052 % 0053 % This file is part of MATPOWER. 0054 % See http://www.pserc.cornell.edu/matpower/ for more info. 0055 % 0056 % MATPOWER is free software: you can redistribute it and/or modify 0057 % it under the terms of the GNU General Public License as published 0058 % by the Free Software Foundation, either version 3 of the License, 0059 % or (at your option) any later version. 0060 % 0061 % MATPOWER is distributed in the hope that it will be useful, 0062 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0063 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0064 % GNU General Public License for more details. 0065 % 0066 % You should have received a copy of the GNU General Public License 0067 % along with MATPOWER. If not, see <http://www.gnu.org/licenses/>. 0068 % 0069 % Additional permission under GNU GPL version 3 section 7 0070 % 0071 % If you modify MATPOWER, or any covered work, to interface with 0072 % other modules (such as MATLAB code and MEX-files) available in a 0073 % MATLAB(R) or comparable environment containing parts covered 0074 % under other licensing terms, the licensors of MATPOWER grant 0075 % you additional permission to convey the resulting work. 0076 0077 if nargin < 3 0078 cp = get_cost_params(om); 0079 else 0080 cp = get_cost_params(om, name); 0081 end 0082 0083 [N, Cw, H, dd, rh, kk, mm] = deal(cp.N, cp.Cw, cp.H, cp.dd, ... 0084 cp.rh, cp.kk, cp.mm); 0085 nw = size(N, 1); 0086 r = N * x - rh; %% Nx - rhat 0087 iLT = find(r < -kk); %% below dead zone 0088 iEQ = find(r == 0 & kk == 0); %% dead zone doesn't exist 0089 iGT = find(r > kk); %% above dead zone 0090 iND = [iLT; iEQ; iGT]; %% rows that are Not in the Dead region 0091 iL = find(dd == 1); %% rows using linear function 0092 iQ = find(dd == 2); %% rows using quadratic function 0093 LL = sparse(iL, iL, 1, nw, nw); 0094 QQ = sparse(iQ, iQ, 1, nw, nw); 0095 kbar = sparse(iND, iND, [ ones(length(iLT), 1); 0096 zeros(length(iEQ), 1); 0097 -ones(length(iGT), 1)], nw, nw) * kk; 0098 rr = r + kbar; %% apply non-dead zone shift 0099 M = sparse(iND, iND, mm(iND), nw, nw); %% dead zone or scale 0100 diagrr = sparse(1:nw, 1:nw, rr, nw, nw); 0101 0102 %% linear rows multiplied by rr(i), quadratic rows by rr(i)^2 0103 w = M * (LL + QQ * diagrr) * rr; 0104 0105 f = full((w' * H * w) / 2 + Cw' * w);