Home > matpower4.1 > mips6opf_solver.m

mips6opf_solver

PURPOSE ^

------------------------------ deprecated ------------------------------

SYNOPSIS ^

function [results, success, raw] = mips6opf_solver(om, mpopt)

DESCRIPTION ^

------------------------------  deprecated  ------------------------------
   MATLAB 6.x support to be removed in a future version.
--------------------------------------------------------------------------
MIPS6OPF_SOLVER  Solves AC optimal power flow using MIPS (for MATLAB 6.x).

   [RESULTS, SUCCESS, RAW] = MIPS6OPF_SOLVER(OM, MPOPT)

   Inputs are an OPF model object and a MATPOWER options vector.

   Outputs are a RESULTS struct, SUCCESS flag and RAW output struct.

   RESULTS is a MATPOWER case struct (mpc) with the usual baseMVA, bus
   branch, gen, gencost fields, along with the following additional
   fields:
       .order      see 'help ext2int' for details of this field
       .x          final value of optimization variables (internal order)
       .f          final objective function value
       .mu         shadow prices on ...
           .var
               .l  lower bounds on variables
               .u  upper bounds on variables
           .nln
               .l  lower bounds on nonlinear constraints
               .u  upper bounds on nonlinear constraints
           .lin
               .l  lower bounds on linear constraints
               .u  upper bounds on linear constraints

   SUCCESS     1 if solver converged successfully, 0 otherwise

   RAW         raw output in form returned by MINOS
       .xr     final value of optimization variables
       .pimul  constraint multipliers
       .info   solver specific termination code
       .output solver specific output information

   See also OPF, MIPS6.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [results, success, raw] = mips6opf_solver(om, mpopt)
0002 %------------------------------  deprecated  ------------------------------
0003 %   MATLAB 6.x support to be removed in a future version.
0004 %--------------------------------------------------------------------------
0005 %MIPS6OPF_SOLVER  Solves AC optimal power flow using MIPS (for MATLAB 6.x).
0006 %
0007 %   [RESULTS, SUCCESS, RAW] = MIPS6OPF_SOLVER(OM, MPOPT)
0008 %
0009 %   Inputs are an OPF model object and a MATPOWER options vector.
0010 %
0011 %   Outputs are a RESULTS struct, SUCCESS flag and RAW output struct.
0012 %
0013 %   RESULTS is a MATPOWER case struct (mpc) with the usual baseMVA, bus
0014 %   branch, gen, gencost fields, along with the following additional
0015 %   fields:
0016 %       .order      see 'help ext2int' for details of this field
0017 %       .x          final value of optimization variables (internal order)
0018 %       .f          final objective function value
0019 %       .mu         shadow prices on ...
0020 %           .var
0021 %               .l  lower bounds on variables
0022 %               .u  upper bounds on variables
0023 %           .nln
0024 %               .l  lower bounds on nonlinear constraints
0025 %               .u  upper bounds on nonlinear constraints
0026 %           .lin
0027 %               .l  lower bounds on linear constraints
0028 %               .u  upper bounds on linear constraints
0029 %
0030 %   SUCCESS     1 if solver converged successfully, 0 otherwise
0031 %
0032 %   RAW         raw output in form returned by MINOS
0033 %       .xr     final value of optimization variables
0034 %       .pimul  constraint multipliers
0035 %       .info   solver specific termination code
0036 %       .output solver specific output information
0037 %
0038 %   See also OPF, MIPS6.
0039 
0040 %   MATPOWER
0041 %   $Id: mips6opf_solver.m,v 1.16 2010/06/09 14:56:58 ray Exp $
0042 %   by Ray Zimmerman, PSERC Cornell
0043 %   and Carlos E. Murillo-Sanchez, PSERC Cornell & Universidad Autonoma de Manizales
0044 %   Copyright (c) 2000-2010 by Power System Engineering Research Center (PSERC)
0045 %
0046 %   This file is part of MATPOWER.
0047 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0048 %
0049 %   MATPOWER is free software: you can redistribute it and/or modify
0050 %   it under the terms of the GNU General Public License as published
0051 %   by the Free Software Foundation, either version 3 of the License,
0052 %   or (at your option) any later version.
0053 %
0054 %   MATPOWER is distributed in the hope that it will be useful,
0055 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0056 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0057 %   GNU General Public License for more details.
0058 %
0059 %   You should have received a copy of the GNU General Public License
0060 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0061 %
0062 %   Additional permission under GNU GPL version 3 section 7
0063 %
0064 %   If you modify MATPOWER, or any covered work, to interface with
0065 %   other modules (such as MATLAB code and MEX-files) available in a
0066 %   MATLAB(R) or comparable environment containing parts covered
0067 %   under other licensing terms, the licensors of MATPOWER grant
0068 %   you additional permission to convey the resulting work.
0069 
0070 %%----- initialization -----
0071 %% define named indices into data matrices
0072 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0073     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0074 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0075     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0076     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0077 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0078     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0079     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0080 [PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;
0081 
0082 %% options
0083 verbose = mpopt(31);    %% VERBOSE
0084 feastol = mpopt(81);    %% PDIPM_FEASTOL
0085 gradtol = mpopt(82);    %% PDIPM_GRADTOL
0086 comptol = mpopt(83);    %% PDIPM_COMPTOL
0087 costtol = mpopt(84);    %% PDIPM_COSTTOL
0088 max_it  = mpopt(85);    %% PDIPM_MAX_IT
0089 max_red = mpopt(86);    %% SCPDIPM_RED_IT
0090 step_control = (mpopt(11) == 565);  %% OPF_ALG == 565, MIPS-sc
0091 if feastol == 0
0092     feastol = mpopt(16);    %% = OPF_VIOLATION by default
0093 end
0094 opt = struct(   'feastol', feastol, ...
0095                 'gradtol', gradtol, ...
0096                 'comptol', comptol, ...
0097                 'costtol', costtol, ...
0098                 'max_it', max_it, ...
0099                 'max_red', max_red, ...
0100                 'step_control', step_control, ...
0101                 'cost_mult', 1e-4, ...
0102                 'verbose', verbose  );
0103 
0104 %% unpack data
0105 mpc = get_mpc(om);
0106 [baseMVA, bus, gen, branch, gencost] = ...
0107     deal(mpc.baseMVA, mpc.bus, mpc.gen, mpc.branch, mpc.gencost);
0108 [vv, ll, nn] = get_idx(om);
0109 
0110 %% problem dimensions
0111 nb = size(bus, 1);          %% number of buses
0112 nl = size(branch, 1);       %% number of branches
0113 ny = getN(om, 'var', 'y');  %% number of piece-wise linear costs
0114 
0115 %% linear constraints
0116 [A, l, u] = linear_constraints(om);
0117 
0118 %% bounds on optimization vars
0119 [x0, xmin, xmax] = getv(om);
0120 
0121 %% build admittance matrices
0122 [Ybus, Yf, Yt] = makeYbus(baseMVA, bus, branch);
0123 
0124 %% try to select an interior initial point
0125 ll = xmin; uu = xmax;
0126 ll(xmin == -Inf) = -1e10;   %% replace Inf with numerical proxies
0127 uu(xmax ==  Inf) =  1e10;
0128 x0 = (ll + uu) / 2;
0129 Varefs = bus(bus(:, BUS_TYPE) == REF, VA) * (pi/180);
0130 x0(vv.i1.Va:vv.iN.Va) = Varefs(1);  %% angles set to first reference angle
0131 if ny > 0
0132     ipwl = find(gencost(:, MODEL) == PW_LINEAR);
0133 %     PQ = [gen(:, PMAX); gen(:, QMAX)];
0134 %     c = totcost(gencost(ipwl, :), PQ(ipwl));
0135     c = gencost(sub2ind(size(gencost), ipwl, NCOST+2*gencost(ipwl, NCOST)));    %% largest y-value in CCV data
0136     x0(vv.i1.y:vv.iN.y) = max(c) + 0.1 * abs(max(c));
0137 %     x0(vv.i1.y:vv.iN.y) = c + 0.1 * abs(c);
0138 end
0139 
0140 %% find branches with flow limits
0141 il = find(branch(:, RATE_A) ~= 0 & branch(:, RATE_A) < 1e10);
0142 nl2 = length(il);           %% number of constrained lines
0143 
0144 %%-----  run opf  -----
0145 f_fcn = @opf_costfcn;
0146 gh_fcn = @opf_consfcn;
0147 ipm_hessian = @opf_hessfcn;
0148 [x, f, info, Output, Lambda] = ...
0149   mips6(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, ipm_hessian, opt, om, Ybus, Yf(il,:), Yt(il,:), mpopt, il);
0150 success = (info > 0);
0151 
0152 %% update solution data
0153 Va = x(vv.i1.Va:vv.iN.Va);
0154 Vm = x(vv.i1.Vm:vv.iN.Vm);
0155 Pg = x(vv.i1.Pg:vv.iN.Pg);
0156 Qg = x(vv.i1.Qg:vv.iN.Qg);
0157 V = Vm .* exp(1j*Va);
0158 
0159 %%-----  calculate return values  -----
0160 %% update voltages & generator outputs
0161 bus(:, VA) = Va * 180/pi;
0162 bus(:, VM) = Vm;
0163 gen(:, PG) = Pg * baseMVA;
0164 gen(:, QG) = Qg * baseMVA;
0165 gen(:, VG) = Vm(gen(:, GEN_BUS));
0166 
0167 %% compute branch flows
0168 Sf = V(branch(:, F_BUS)) .* conj(Yf * V);  %% cplx pwr at "from" bus, p.u.
0169 St = V(branch(:, T_BUS)) .* conj(Yt * V);  %% cplx pwr at "to" bus, p.u.
0170 branch(:, PF) = real(Sf) * baseMVA;
0171 branch(:, QF) = imag(Sf) * baseMVA;
0172 branch(:, PT) = real(St) * baseMVA;
0173 branch(:, QT) = imag(St) * baseMVA;
0174 
0175 %% line constraint is actually on square of limit
0176 %% so we must fix multipliers
0177 muSf = zeros(nl, 1);
0178 muSt = zeros(nl, 1);
0179 if ~isempty(il)
0180     muSf(il) = 2 * Lambda.ineqnonlin(1:nl2)       .* branch(il, RATE_A) / baseMVA;
0181     muSt(il) = 2 * Lambda.ineqnonlin((1:nl2)+nl2) .* branch(il, RATE_A) / baseMVA;
0182 end
0183 
0184 %% update Lagrange multipliers
0185 bus(:, MU_VMAX)  = Lambda.upper(vv.i1.Vm:vv.iN.Vm);
0186 bus(:, MU_VMIN)  = Lambda.lower(vv.i1.Vm:vv.iN.Vm);
0187 gen(:, MU_PMAX)  = Lambda.upper(vv.i1.Pg:vv.iN.Pg) / baseMVA;
0188 gen(:, MU_PMIN)  = Lambda.lower(vv.i1.Pg:vv.iN.Pg) / baseMVA;
0189 gen(:, MU_QMAX)  = Lambda.upper(vv.i1.Qg:vv.iN.Qg) / baseMVA;
0190 gen(:, MU_QMIN)  = Lambda.lower(vv.i1.Qg:vv.iN.Qg) / baseMVA;
0191 bus(:, LAM_P)    = Lambda.eqnonlin(nn.i1.Pmis:nn.iN.Pmis) / baseMVA;
0192 bus(:, LAM_Q)    = Lambda.eqnonlin(nn.i1.Qmis:nn.iN.Qmis) / baseMVA;
0193 branch(:, MU_SF) = muSf / baseMVA;
0194 branch(:, MU_ST) = muSt / baseMVA;
0195 
0196 %% package up results
0197 nlnN = getN(om, 'nln');
0198 
0199 %% extract multipliers for nonlinear constraints
0200 kl = find(Lambda.eqnonlin < 0);
0201 ku = find(Lambda.eqnonlin > 0);
0202 nl_mu_l = zeros(nlnN, 1);
0203 nl_mu_u = [zeros(2*nb, 1); muSf; muSt];
0204 nl_mu_l(kl) = -Lambda.eqnonlin(kl);
0205 nl_mu_u(ku) =  Lambda.eqnonlin(ku);
0206 
0207 mu = struct( ...
0208   'var', struct('l', Lambda.lower, 'u', Lambda.upper), ...
0209   'nln', struct('l', nl_mu_l, 'u', nl_mu_u), ...
0210   'lin', struct('l', Lambda.mu_l, 'u', Lambda.mu_u) );
0211 
0212 results = mpc;
0213 [results.bus, results.branch, results.gen, ...
0214     results.om, results.x, results.mu, results.f] = ...
0215         deal(bus, branch, gen, om, x, mu, f);
0216 
0217 pimul = [ ...
0218   results.mu.nln.l - results.mu.nln.u;
0219   results.mu.lin.l - results.mu.lin.u;
0220   -ones(ny>0, 1);
0221   results.mu.var.l - results.mu.var.u;
0222 ];
0223 raw = struct('xr', x, 'pimul', pimul, 'info', info, 'output', Output);

Generated on Mon 26-Jan-2015 15:00:13 by m2html © 2005