mips

mips(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt)

mips() - MATPOWER Interior Point Solver (MIPS).

[x, f, exitflag, output, lambda] = ...
    mips(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);

x = mips(f_fcn, x0)
x = mips(f_fcn, x0, A, l)
x = mips(f_fcn, x0, A, l, u)
x = mips(f_fcn, x0, A, l, u, xmin)
x = mips(f_fcn, x0, A, l, u, xmin, xmax)
x = mips(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn)
x = mips(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt)
x = mips(problem)
        where problem is a struct with fields:
            f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt
            all fields except f_fcn and x0 are optional
x = mips(...)
[x, f] = mips(...)
[x, f, exitflag] = mips(...)
[x, f, exitflag, output] = mips(...)
[x, f, exitflag, output, lambda] = mips(...)

Primal-dual interior point method for NLP (nonlinear programming). Minimize a function \(f(\x)\) beginning from a starting point x0, subject to optional linear and nonlinear constraints and variable bounds.

(1)\[\min_\x f(\x)\]

subject to

(2)\[\g(\x) = 0\]

(3)\[\h(\x) \le 0\]

(4)\[\param{\rvec{l}} \le \param{\rmat{A}} \x \le \param{\rvec{u}}\]

(5)\[\param{\x}_\mathrm{min} \le \x \le \param{\x}_\mathrm{max}\]

where (1)–(5) are, respectively, the objective function, nonlinear equality constraints, nonlinear inequality constraints, linear constraints, and variable bounds.

Inputs:

(all optional except f_fcn and x0 )
f_fcn (function handle) – handle to function that evaluates the objective function \(f(\x)\), its gradients and Hessian for a given value of \(x\). If there are nonlinear constraints, the Hessian information is provided by the hess_fcn function passed in the 9th argument and is not required here. Calling syntax for this function:
```
[f, df, d2f] = f_fcn(x)
```
x0 (double) – starting value of optimization vector \(\x\)
A, l, u (double) – (optional, respective defaults \([empty, -\infty, +\infty]\) ) matrix \(\param{\cmat{A}}\) and vectors \(\param{\rvec{l}}\) and \(\param{\rvec{u}}\) to define the linear constraints in (4)
xmin, xmax (double) – (optional, respective defaults \([-\infty, +\infty]\) ) lower and upper bounds, \(\param{\x}_\mathrm{min}\) and \(\param{\x}_\mathrm{max}\), on the optimization variables \(\x\)
gh_fcn (function handle) – (optional) handle to function that evaluates the nonlinear constraints \(\g(\x)\) and \(\h(\x)\) and their gradients for a given value of \(x\). Calling syntax for this function is:
```
[h, g, dh, dg] = gh_fcn(x)
```
where the columns of dh and dg are the gradients of the corresponding elements of h and g, i.e. dh and dg are transposes of the Jacobians of h and g, respectively.
hess_fcn (function handle) – (optional) handle to function that computes the Hessian of the Lagrangian for given values of \(\x\), \(\lam\) and \(\muv\), where \(\lam\) and \(\muv\) are the multipliers on the equality and inequality constraints, \(\g\) and \(\h\), respectively. The calling syntax for this function is:
```
Lxx = hess_fcn(x, lam)
```
where \(\lam\) = lam.eqnonlin and \(\muv\) = lam.ineqnonlin.
opt (struct) – (optional) options structure with the following fields, all of which are also optional (default values shown in parentheses)
- verbose (0) - controls level of progress output displayed
  - 0 = no progress output
  - 1 = some progress output
  - 2 = verbose progress output
- linsolver ( '' ) - linear system solver for solving update steps
  - '' = default solver (currently same as '\')
  - '\' = built-in \ operator
  - 'PARDISO' = PARDISO solver (if available)
- feastol (1e-6) - termination tolerance for feasibility condition
- gradtol (1e-6) - termination tolerance for gradient condition
- comptol (1e-6) - termination tolerance for complementarity condition
- costtol (1e-6) - termination tolerance for cost condition
- max_it (150) - maximum number of iterations
- step_control (0) - set to 1 to enable step-size control
- sc.red_it (20) - maximum number of step-size reductions if step-control is on
- cost_mult (1) - cost multiplier used to scale the objective function for improved conditioning. Note: This value is also passed as the 3rd argument to the Hessian evaluation function so that it can appropriately scale the objective function term in the Hessian of the Lagrangian.
- xi (0.99995) - constant \(\xi\) used in \(alpha\) updates [1]
- sigma (0.1) - centering parameter \(\sigma\) [1]
- z0 (1) - used to initialize slack variables [1]
- alpha_min (1e-8) - returns exitflag = -1 if either alpha parameter, \(\alpha_p\) or \(\alpha_d\) becomes smaller than this value [1]
- rho_min (0.95) - lower bound on rho_t [1]
- rho_max (1.05) - upper bound on rho_t [1]
- mu_threshold (1e-5) - KT multipliers smaller than this value for non-binding constraints are forced to zero
- max_stepsize (1e10) - returns exitflag = -1 if the 2-norm of the reduced Newton step exceeds this value [1]
problem (struct) – The inputs can alternatively be supplied in a single problem struct with fields corresponding to the input arguments described above, namely, f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, and opt, where all except f_fcn and x0 are optional

Outputs:

x (double) – solution vector, \(\x\)
f (double) – final objective function value, \(f(\x)\)
exitflag (integer) – exit flag
- 1 = first order optimality conditions satisfied
- 0 = maximum number of iterations reached
- -1 = numerically failed
output (struct) – output struct with fields:
- iterations - number of iterations performed
- hist - struct array with trajectories of the following:
  - feascond, gradcond, compcond, costcond, gamma, stepsize, obj, alphap, alphad
- message - exit message
lambda (struct) – struct containing the Langrange and Kuhn-Tucker multipliers on the constraints, with fields:
- eqnonlin - nonlinear equality constraints
- ineqnonlin - nonlinear inequality constraints
- mu_l - lower (left-hand) limit on linear constraints
- mu_u - upper (right-hand) limit on linear constraints
- lower - lower bound on optimization variables
- upper - upper bound on optimization variables

Note: The calling syntax is almost identical to that of fmincon() from MathWorks’ Optimization Toolbox. The main difference is that the linear constraints are specified with A, l, u instead of A, B, Aeq, Beq. The functions for evaluating the objective function, constraints and Hessian are identical.

Example: (problem from https://en.wikipedia.org/wiki/Nonlinear_programming):

function [f, df, d2f] = f2(x)
f = -x(1)*x(2) - x(2)*x(3);
if nargout > 1           % gradient is required
    df = -[x(2); x(1)+x(3); x(2)];
    if nargout > 2       % Hessian is required
        d2f = -[0 1 0; 1 0 1; 0 1 0];   % actually not used since
    end                                 % 'hess_fcn' is provided
end

function [h, g, dh, dg] = gh2(x)
h = [ 1 -1 1; 1 1 1] * x.^2 + [-2; -10];
dh = 2 * [x(1) x(1); -x(2) x(2); x(3) x(3)];
g = []; dg = [];

function Lxx = hess2(x, lam, cost_mult)
if nargin < 3, cost_mult = 1; end
mu = lam.ineqnonlin;
Lxx = cost_mult * [0 -1 0; -1 0 -1; 0 -1 0] + ...
        [2*[1 1]*mu 0 0; 0 2*[-1 1]*mu 0; 0 0 2*[1 1]*mu];

problem = struct( ...
    'f_fcn',    @(x)f2(x), ...
    'gh_fcn',   @(x)gh2(x), ...
    'hess_fcn', @(x, lam, cost_mult)hess2(x, lam, cost_mult), ...
    'x0',       [1; 1; 0], ...
    'opt',      struct('verbose', 2) ...
);
[x, f, exitflag, output, lambda] = mips(problem);

Ported by Ray Zimmerman from C code written by H. Wang for his PhD dissertation:

H. Wang, On the Computation and Application of Multi-period Security-Constrained Optimal Power Flow for Real-time Electricity Market Operations, Ph.D. thesis, Electrical and Computer Engineering, Cornell University, May 2007.

Please see also:

H. Wang, C. E. Murillo-Sanchez, R. D. Zimmerman, R. J. Thomas, “On Computational Issues of Market-Based Optimal Power Flow”, IEEE Transactions on Power Systems, Vol. 22, No. 3, Aug. 2007, pp. 1185-1193. doi: 10.1109/TPWRS.2007.901301