Home > matpower7.1 > lib > dcopf_solver.m

dcopf_solver

PURPOSE ^

DCOPF_SOLVER Solves a DC optimal power flow.

SYNOPSIS ^

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

DESCRIPTION ^

DCOPF_SOLVER  Solves a DC optimal power flow.

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

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

   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
           .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, OPT_MODEL/SOLVE.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [results, success, raw] = dcopf_solver(om, mpopt)
0002 %DCOPF_SOLVER  Solves a DC optimal power flow.
0003 %
0004 %   [RESULTS, SUCCESS, RAW] = DCOPF_SOLVER(OM, MPOPT)
0005 %
0006 %   Inputs are an OPF model object and a MATPOWER options struct.
0007 %
0008 %   Outputs are a RESULTS struct, SUCCESS flag and RAW output struct.
0009 %
0010 %   RESULTS is a MATPOWER case struct (mpc) with the usual baseMVA, bus
0011 %   branch, gen, gencost fields, along with the following additional
0012 %   fields:
0013 %       .order      see 'help ext2int' for details of this field
0014 %       .x          final value of optimization variables (internal order)
0015 %       .f          final objective function value
0016 %       .mu         shadow prices on ...
0017 %           .var
0018 %               .l  lower bounds on variables
0019 %               .u  upper bounds on variables
0020 %           .lin
0021 %               .l  lower bounds on linear constraints
0022 %               .u  upper bounds on linear constraints
0023 %
0024 %   SUCCESS     1 if solver converged successfully, 0 otherwise
0025 %
0026 %   RAW         raw output in form returned by MINOS
0027 %       .xr     final value of optimization variables
0028 %       .pimul  constraint multipliers
0029 %       .info   solver specific termination code
0030 %       .output solver specific output information
0031 %
0032 %   See also OPF, OPT_MODEL/SOLVE.
0033 
0034 %   MATPOWER
0035 %   Copyright (c) 2000-2020, Power Systems Engineering Research Center (PSERC)
0036 %   by Ray Zimmerman, PSERC Cornell
0037 %   and Carlos E. Murillo-Sanchez, PSERC Cornell & Universidad Nacional de Colombia
0038 %
0039 %   This file is part of MATPOWER.
0040 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0041 %   See https://matpower.org for more info.
0042 
0043 %%----- initialization -----
0044 %% define named indices into data matrices
0045 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0046     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0047 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0048     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0049     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0050 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0051     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0052     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0053 [PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;
0054 
0055 %% unpack data
0056 mpc = om.get_mpc();
0057 [baseMVA, bus, gen, branch, gencost] = ...
0058     deal(mpc.baseMVA, mpc.bus, mpc.gen, mpc.branch, mpc.gencost);
0059 cp = om.get_cost_params();
0060 Bf = om.get_userdata('Bf');
0061 Pfinj = om.get_userdata('Pfinj');
0062 [vv, ll] = om.get_idx();
0063 
0064 %% problem dimensions
0065 nb = size(bus, 1);          %% number of buses
0066 nl = size(branch, 1);       %% number of branches
0067 ny = om.getN('var', 'y');   %% number of piece-wise linear costs
0068 
0069 %% options
0070 model = om.problem_type();
0071 opt = mpopt2qpopt(mpopt, model);
0072 if strcmp(opt.alg, 'OSQP')
0073     opt.x0 = [];    %% disable provided starting point for OSQP
0074 end
0075 
0076 %% try to select an interior initial point, unless requested not to
0077 if mpopt.opf.start < 2 && ...
0078         (strcmp(opt.alg, 'MIPS') || strcmp(opt.alg, 'IPOPT'))
0079     [x0, xmin, xmax] = om.params_var();     %% init var & bounds
0080     s = 1;                      %% set init point inside bounds by s
0081     lb = xmin; ub = xmax;
0082     lb(xmin == -Inf) = -1e10;   %% replace Inf with numerical proxies
0083     ub(xmax ==  Inf) =  1e10;
0084     x0 = (lb + ub) / 2;         %% set x0 mid-way between bounds
0085     k = find(xmin == -Inf & xmax < Inf);    %% if only bounded above
0086     x0(k) = xmax(k) - s;                    %% set just below upper bound
0087     k = find(xmin > -Inf & xmax == Inf);    %% if only bounded below
0088     x0(k) = xmin(k) + s;                    %% set just above lower bound
0089     Varefs = bus(bus(:, BUS_TYPE) == REF, VA) * (pi/180);
0090     x0(vv.i1.Va:vv.iN.Va) = Varefs(1);  %% angles set to first reference angle
0091     if ny > 0
0092         ipwl = find(gencost(:, MODEL) == PW_LINEAR);
0093         c = gencost(sub2ind(size(gencost), ipwl, NCOST+2*gencost(ipwl, NCOST)));    %% largest y-value in CCV data
0094         x0(vv.i1.y:vv.iN.y) = max(c) + 0.1 * abs(max(c));
0095     end
0096     opt.x0 = x0;
0097 end
0098 
0099 %%-----  run opf  -----
0100 [x, f, eflag, output, lambda] = om.solve(opt);
0101 success = (eflag == 1);
0102 
0103 %%-----  calculate return values  -----
0104 if ~any(isnan(x))
0105     %% update solution data
0106     Va = x(vv.i1.Va:vv.iN.Va);
0107     Pg = x(vv.i1.Pg:vv.iN.Pg);
0108 
0109     %% update voltages & generator outputs
0110     bus(:, VM) = ones(nb, 1);
0111     bus(:, VA) = Va * 180/pi;
0112     gen(:, PG) = Pg * baseMVA;
0113 
0114     %% compute branch flows
0115     branch(:, [QF, QT]) = zeros(nl, 2);
0116     branch(:, PF) = (Bf * Va + Pfinj) * baseMVA;
0117     branch(:, PT) = -branch(:, PF);
0118 end
0119 
0120 %% package up results
0121 mu_l = lambda.mu_l;
0122 mu_u = lambda.mu_u;
0123 muLB = lambda.lower;
0124 muUB = lambda.upper;
0125 
0126 %% update Lagrange multipliers
0127 il = find(branch(:, RATE_A) ~= 0 & branch(:, RATE_A) < 1e10);
0128 bus(:, [LAM_P, LAM_Q, MU_VMIN, MU_VMAX]) = zeros(nb, 4);
0129 gen(:, [MU_PMIN, MU_PMAX, MU_QMIN, MU_QMAX]) = zeros(size(gen, 1), 4);
0130 branch(:, [MU_SF, MU_ST]) = zeros(nl, 2);
0131 bus(:, LAM_P)       = (mu_u(ll.i1.Pmis:ll.iN.Pmis) - mu_l(ll.i1.Pmis:ll.iN.Pmis)) / baseMVA;
0132 if ~isempty(il)
0133     branch(il, MU_SF)   = mu_u(ll.i1.Pf:ll.iN.Pf) / baseMVA;
0134     branch(il, MU_ST)   = mu_l(ll.i1.Pf:ll.iN.Pf) / baseMVA;
0135 end
0136 gen(:, MU_PMIN)     = muLB(vv.i1.Pg:vv.iN.Pg) / baseMVA;
0137 gen(:, MU_PMAX)     = muUB(vv.i1.Pg:vv.iN.Pg) / baseMVA;
0138 pimul = [
0139   mu_l - mu_u;
0140  -ones(ny>0, 1);    %% dummy entry corresponding to linear cost row in A (in MINOS)
0141   muLB - muUB
0142 ];
0143 
0144 mu = struct( ...
0145   'var', struct('l', muLB, 'u', muUB), ...
0146   'lin', struct('l', mu_l, 'u', mu_u) );
0147 
0148 results = mpc;
0149 [results.bus, results.branch, results.gen, ...
0150     results.om, results.x, results.mu, results.f] = ...
0151         deal(bus, branch, gen, om, x, mu, f);
0152 
0153 raw = struct('xr', x, 'pimul', pimul, 'info', eflag, 'output', output);
Generated on Fri 09-Oct-2020 11:21:31 by m2html © 2005