Home > matpower5.1 > t > t_miqps_matpower.m

t_miqps_matpower

PURPOSE ^

T_MIQPS_MATPOWER Tests of MIQPS_MATPOWER MIQP solvers.

SYNOPSIS ^

function t_miqps_matpower(quiet)

DESCRIPTION ^

T_MIQPS_MATPOWER  Tests of MIQPS_MATPOWER MIQP solvers.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function t_miqps_matpower(quiet)
0002 %T_MIQPS_MATPOWER  Tests of MIQPS_MATPOWER MIQP solvers.
0003 
0004 %   MATPOWER
0005 %   Copyright (c) 2010-2015 by Power System Engineering Research Center (PSERC)
0006 %   by Ray Zimmerman, PSERC Cornell
0007 %
0008 %   $Id: t_miqps_matpower.m 2644 2015-03-11 19:34:22Z ray $
0009 %
0010 %   This file is part of MATPOWER.
0011 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0012 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0013 
0014 if nargin < 1
0015     quiet = 0;
0016 end
0017 
0018 algs = {'CPLEX', 'MOSEK', 'GUROBI', 'GLPK', 'OT'};
0019 names = {'CPLEX', 'MOSEK', 'Gurobi', 'glpk', 'intlin/lin/quadprog'};
0020 check = {'cplex', 'mosek', 'gurobi', 'glpk', 'intlinprog'};
0021 does_qp = [1 1 1 0 0];
0022 
0023 n = 48;
0024 nqp = 28;
0025 nmiqp = 6;
0026 t_begin(n*length(algs), quiet);
0027 
0028 for k = 1:length(algs)
0029     if ~isempty(check{k}) && ~have_fcn(check{k})
0030         t_skip(n, sprintf('%s not installed', names{k}));
0031     else
0032         opt = struct('verbose', 0, 'alg', algs{k});
0033         if strcmp(names{k}, 'CPLEX')
0034 %           alg = 0;        %% default uses barrier method with NaN bug in lower lim multipliers
0035             alg = 2;        %% use dual simplex
0036             mpopt = mpoption('cplex.lpmethod', alg, 'cplex.qpmethod', min([4 alg]));
0037             opt.cplex_opt = cplex_options([], mpopt);
0038         end
0039         if strcmp(names{k}, 'MOSEK')
0040             mpopt = mpoption;
0041 %             sc = mosek_symbcon;
0042 %             alg = sc.MSK_OPTIMIZER_DUAL_SIMPLEX;    %% use dual simplex
0043 %             mpopt = mpoption(mpopt, 'mosek.lp_alg', alg );
0044             mpopt = mpoption(mpopt, 'mosek.gap_tol', 1e-10);
0045             opt.mosek_opt = mosek_options([], mpopt);
0046         end
0047 
0048         t = sprintf('%s - 3-d LP : ', names{k});
0049         %% based on example from 'doc linprog'
0050         c = [-5; -4; -6];
0051         A = [1 -1  1;
0052              -3  -2  -4;
0053              3  2  0];
0054         l = [-Inf; -42; -Inf];
0055         u = [20; Inf; 30];
0056         xmin = [0; 0; 0];
0057         x0 = [];
0058         [x, f, s, out, lam] = miqps_matpower([], c, A, l, u, xmin, [], [], [], opt);
0059         t_is(s, 1, 12, [t 'success']);
0060         t_is(x, [0; 15; 3], 6, [t 'x']);
0061         t_is(f, -78, 6, [t 'f']);
0062         t_is(lam.mu_l, [0;1.5;0], 9, [t 'lam.mu_l']);
0063         t_is(lam.mu_u, [0;0;0.5], 9, [t 'lam.mu_u']);
0064         t_is(lam.lower, [1;0;0], 9, [t 'lam.lower']);
0065         t_is(lam.upper, zeros(size(x)), 9, [t 'lam.upper']);
0066 
0067         if does_qp(k)
0068             t = sprintf('%s - unconstrained 3-d quadratic : ', names{k});
0069             %% from http://www.akiti.ca/QuadProgEx0Constr.html
0070             H = [5 -2 -1; -2 4 3; -1 3 5];
0071             c = [2; -35; -47];
0072             x0 = [0; 0; 0];
0073             [x, f, s, out, lam] = miqps_matpower(H, c, [], [], [], [], [], [], [], opt);
0074             t_is(s, 1, 12, [t 'success']);
0075             t_is(x, [3; 5; 7], 8, [t 'x']);
0076             t_is(f, -249, 13, [t 'f']);
0077             t_ok(isempty(lam.mu_l), [t 'lam.mu_l']);
0078             t_ok(isempty(lam.mu_u), [t 'lam.mu_u']);
0079             t_is(lam.lower, zeros(size(x)), 13, [t 'lam.lower']);
0080             t_is(lam.upper, zeros(size(x)), 13, [t 'lam.upper']);
0081         
0082             t = sprintf('%s - constrained 2-d QP : ', names{k});
0083             %% example from 'doc quadprog'
0084             H = [   1   -1;
0085                     -1  2   ];
0086             c = [-2; -6];
0087             A = [   1   1;
0088                     -1  2;
0089                     2   1   ];
0090             l = [];
0091             u = [2; 2; 3];
0092             xmin = [0; 0];
0093             x0 = [];
0094             [x, f, s, out, lam] = miqps_matpower(H, c, A, l, u, xmin, [], x0, [], opt);
0095             t_is(s, 1, 12, [t 'success']);
0096             t_is(x, [2; 4]/3, 7, [t 'x']);
0097             t_is(f, -74/9, 6, [t 'f']);
0098             t_is(lam.mu_l, [0;0;0], 13, [t 'lam.mu_l']);
0099             t_is(lam.mu_u, [28;4;0]/9, 7, [t 'lam.mu_u']);
0100             t_is(lam.lower, zeros(size(x)), 7, [t 'lam.lower']);
0101             t_is(lam.upper, zeros(size(x)), 13, [t 'lam.upper']);
0102 
0103             t = sprintf('%s - constrained 4-d QP : ', names{k});
0104             %% from http://www.jmu.edu/docs/sasdoc/sashtml/iml/chap8/sect12.htm
0105             H = [   1003.1  4.3     6.3     5.9;
0106                     4.3     2.2     2.1     3.9;
0107                     6.3     2.1     3.5     4.8;
0108                     5.9     3.9     4.8     10  ];
0109             c = zeros(4,1);
0110             A = [   1       1       1       1;
0111                     0.17    0.11    0.10    0.18    ];
0112             l = [1; 0.10];
0113             u = [1; Inf];
0114             xmin = zeros(4,1);
0115             x0 = [1; 0; 0; 1];
0116             [x, f, s, out, lam] = miqps_matpower(H, c, A, l, u, xmin, [], x0, [], opt);
0117             t_is(s, 1, 12, [t 'success']);
0118             t_is(x, [0; 2.8; 0.2; 0]/3, 5, [t 'x']);
0119             t_is(f, 3.29/3, 6, [t 'f']);
0120             t_is(lam.mu_l, [6.58;0]/3, 6, [t 'lam.mu_l']);
0121             t_is(lam.mu_u, [0;0], 13, [t 'lam.mu_u']);
0122             t_is(lam.lower, [2.24;0;0;1.7667], 4, [t 'lam.lower']);
0123             t_is(lam.upper, zeros(size(x)), 13, [t 'lam.upper']);
0124 
0125             t = sprintf('%s - (struct) constrained 4-d QP : ', names{k});
0126             p = struct('H', H, 'A', A, 'l', l, 'u', u, 'xmin', xmin, 'x0', x0, 'opt', opt);
0127             [x, f, s, out, lam] = miqps_matpower(p);
0128             t_is(s, 1, 12, [t 'success']);
0129             t_is(x, [0; 2.8; 0.2; 0]/3, 5, [t 'x']);
0130             t_is(f, 3.29/3, 6, [t 'f']);
0131             t_is(lam.mu_l, [6.58;0]/3, 6, [t 'lam.mu_l']);
0132             t_is(lam.mu_u, [0;0], 13, [t 'lam.mu_u']);
0133             t_is(lam.lower, [2.24;0;0;1.7667], 4, [t 'lam.lower']);
0134             t_is(lam.upper, zeros(size(x)), 13, [t 'lam.upper']);
0135         else
0136             t_skip(nqp, sprintf('%s does not handle QP problems', names{k}));
0137         end
0138 
0139         t = sprintf('%s - infeasible LP : ', names{k});
0140         p = struct('A', sparse([1 1]), 'c', [1;1], 'u', -1, 'xmin', [0;0], 'opt', opt);
0141         [x, f, s, out, lam] = miqps_matpower(p);
0142         t_ok(s <= 0, [t 'no success']);
0143 
0144 % opt.verbose = 3;
0145         t = sprintf('%s - 2-d MILP : ', names{k});
0146         %% from MOSEK 6.0 Guided Tour, section  7.13.1, http://docs.mosek.com/6.0/toolbox/node009.html#283040944
0147         c = [-2; -3];
0148         A = sparse([195 273; 4 40]);
0149         u = [1365; 140];
0150         xmax = [4; Inf];
0151         vtype = 'I';
0152         p = struct('c', c, 'A', A, 'u', u, 'xmax', xmax, 'vtype', vtype, 'opt', opt);
0153         [x, f, s, out, lam] = miqps_matpower(p);
0154         t_is(s, 1, 12, [t 'success']);
0155         t_is(x, [4; 2], 12, [t 'x']);
0156         t_is(lam.mu_l, [0; 0], 12, [t 'lam.mu_l']);
0157         t_is(lam.mu_u, [0; 0], 12, [t 'lam.mu_u']);
0158         t_is(lam.lower, [0; 0], 12, [t 'lam.lower']);
0159         t_is(lam.upper, [2; 3], 12, [t 'lam.upper']);
0160 
0161         if does_qp(k)
0162             t = sprintf('%s - 4-d MIQP : ', names{k});
0163             %% from cplexmiqpex.m, CPLEX_Studio_Academic124/cplex/examples/src/matlab/cplexmiqpex.m
0164             H = sparse([ 33   6    0    0;
0165                           6  22   11.5  0;
0166                           0  11.5 11    0;
0167                           0   0    0    0]);
0168             c = [-1 -2 -3 -1]';
0169             Aineq = [-1  1  1 10;
0170                1 -3  1  0];
0171             bineq = [20  30]';
0172             Aeq   = [0  1  0 -3.5];
0173             beq   =  0;
0174             xmin    = [ 0;   0;   0; 2];
0175             xmax    = [40; Inf; Inf; 3];
0176             A = sparse([Aeq; Aineq]);
0177             l = [beq; -Inf; -Inf];
0178             u = [beq; bineq];
0179             vtype = 'CCCI';
0180             p = struct('H', H, 'c', c, 'A', A, 'l', l, 'u', u, ...
0181                 'xmin', xmin, 'xmax', xmax, 'vtype', vtype, 'opt', opt);
0182             [x, f, s, out, lam] = miqps_matpower(p);
0183             t_is(s, 1, 12, [t 'success']);
0184             t_is(x, [7; 7; 0; 2], 7, [t 'x']);
0185             t_is(lam.mu_l, [466; 0; 0], 6, [t 'lam.mu_l']);
0186             t_is(lam.mu_u, [0; 272; 0], 6, [t 'lam.mu_u']);
0187             t_is(lam.lower, [0; 0; 349.5; 4350], 5, [t 'lam.lower']);
0188             t_is(lam.upper, [0; 0; 0; 0], 7, [t 'lam.upper']);
0189         else
0190             t_skip(nmiqp, sprintf('%s does not handle MIQP problems', names{k}));
0191         end
0192 % opt.verbose = 0;
0193     end
0194 end
0195 
0196 t_end;

Generated on Fri 20-Mar-2015 18:23:34 by m2html © 2005