COMPUTE_COST Computes a user-defined cost. F_U = COMPUTE_COST(OM, X) F_U = COMPUTE_COST(OM, X, NAME) F_U = COMPUTE_COST(OM, X, NAME, IDX) 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 OPT_MODEL, ADD_COSTS, BUILD_COST_PARAMS, GET_COST_PARAMS.
0001 function f = compute_cost(om, x, name, idx) 0002 %COMPUTE_COST Computes a user-defined cost. 0003 % F_U = COMPUTE_COST(OM, X) 0004 % F_U = COMPUTE_COST(OM, X, NAME) 0005 % F_U = COMPUTE_COST(OM, X, NAME, IDX) 0006 % 0007 % Computes the value of a user defined cost, either for all user 0008 % defined costs or for a named set of costs. Requires calling 0009 % BUILD_COST_PARAMS first to build the full set of parameters. 0010 % 0011 % Let X be the full set of optimization variables and F_U(X, CP) be the 0012 % user-defined cost at X, corresponding to the set of cost parameters in 0013 % the CP struct returned by GET_COST_PARAMS, where CP is a struct with the 0014 % following fields: 0015 % N - nw x nx sparse matrix 0016 % Cw - nw x 1 vector 0017 % H - nw x nw sparse matrix (optional, all zeros by default) 0018 % dd, mm - nw x 1 vectors (optional, all ones by default) 0019 % rh, kk - nw x 1 vectors (optional, all zeros by default) 0020 % 0021 % These parameters are used as follows to compute F_U(X, CP) 0022 % 0023 % R = N*x - rh 0024 % 0025 % / kk(i), R(i) < -kk(i) 0026 % K(i) = < 0, -kk(i) <= R(i) <= kk(i) 0027 % \ -kk(i), R(i) > kk(i) 0028 % 0029 % RR = R + K 0030 % 0031 % U(i) = / 0, -kk(i) <= R(i) <= kk(i) 0032 % \ 1, otherwise 0033 % 0034 % DDL(i) = / 1, dd(i) = 1 0035 % \ 0, otherwise 0036 % 0037 % DDQ(i) = / 1, dd(i) = 2 0038 % \ 0, otherwise 0039 % 0040 % Dl = diag(mm) * diag(U) * diag(DDL) 0041 % Dq = diag(mm) * diag(U) * diag(DDQ) 0042 % 0043 % w = (Dl + Dq * diag(RR)) * RR 0044 % 0045 % F_U(X, CP) = 1/2 * w'*H*w + Cw'*w 0046 % 0047 % See also OPT_MODEL, ADD_COSTS, BUILD_COST_PARAMS, GET_COST_PARAMS. 0048 0049 % MATPOWER 0050 % Copyright (c) 2008-2015 by Power System Engineering Research Center (PSERC) 0051 % by Ray Zimmerman, PSERC Cornell 0052 % 0053 % $Id: compute_cost.m 2644 2015-03-11 19:34:22Z ray $ 0054 % 0055 % This file is part of MATPOWER. 0056 % Covered by the 3-clause BSD License (see LICENSE file for details). 0057 % See http://www.pserc.cornell.edu/matpower/ for more info. 0058 0059 if nargin < 3 0060 cp = get_cost_params(om); 0061 elseif nargin < 4 0062 cp = get_cost_params(om, name); 0063 else 0064 cp = get_cost_params(om, name, idx); 0065 end 0066 0067 [N, Cw, H, dd, rh, kk, mm] = deal(cp.N, cp.Cw, cp.H, cp.dd, ... 0068 cp.rh, cp.kk, cp.mm); 0069 nw = size(N, 1); 0070 r = N * x - rh; %% Nx - rhat 0071 iLT = find(r < -kk); %% below dead zone 0072 iEQ = find(r == 0 & kk == 0); %% dead zone doesn't exist 0073 iGT = find(r > kk); %% above dead zone 0074 iND = [iLT; iEQ; iGT]; %% rows that are Not in the Dead region 0075 iL = find(dd == 1); %% rows using linear function 0076 iQ = find(dd == 2); %% rows using quadratic function 0077 LL = sparse(iL, iL, 1, nw, nw); 0078 QQ = sparse(iQ, iQ, 1, nw, nw); 0079 kbar = sparse(iND, iND, [ ones(length(iLT), 1); 0080 zeros(length(iEQ), 1); 0081 -ones(length(iGT), 1)], nw, nw) * kk; 0082 rr = r + kbar; %% apply non-dead zone shift 0083 M = sparse(iND, iND, mm(iND), nw, nw); %% dead zone or scale 0084 diagrr = sparse(1:nw, 1:nw, rr, nw, nw); 0085 0086 %% linear rows multiplied by rr(i), quadratic rows by rr(i)^2 0087 w = M * (LL + QQ * diagrr) * rr; 0088 0089 f = full((w' * H * w) / 2 + Cw' * w);