NLPS_KNITRO Nonlinear programming (NLP) Solver based on Artelys Knitro.
[X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
NLPS_KNITRO(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT)
[X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_KNITRO(PROBLEM)
A wrapper function providing a standardized interface for using
Artelys Knitro to solve the following NLP (nonlinear programming) problem:
Minimize a function F(X) beginning from a starting point X0, subject to
optional linear and nonlinear constraints and variable bounds.
min F(X)
X
subject to
G(X) = 0 (nonlinear equalities)
H(X) <= 0 (nonlinear inequalities)
L <= A*X <= U (linear constraints)
XMIN <= X <= XMAX (variable bounds)
Inputs (all optional except F_FCN and X0):
F_FCN : handle to function that evaluates the objective function,
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 : starting value of optimization vector X
A, L, U : define the optional linear constraints. Default
values for the elements of L and U are -Inf and Inf,
respectively.
XMIN, XMAX : optional lower and upper bounds on the
X variables, defaults are -Inf and Inf, respectively.
GH_FCN : handle to function that evaluates the optional
nonlinear constraints 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 : handle to function that computes the Hessian of the
Lagrangian for given values of X, lambda and mu, where
lambda and mu 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 lambda = LAM.eqnonlin and mu = LAM.ineqnonlin.
OPT : 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
knitro_opt - options struct for Artelys Knitro, value in verbose
overrides these options
opts - struct of other values to be passed directly to
Aretelys Knitro via an options file
tol_x - termination tol on x
tol_f - termination tol on f
maxit - maximum number of iterations
PROBLEM : The inputs can alternatively be supplied in a single
PROBLEM struct with fields corresponding to the input arguments
described above: f_fcn, x0, A, l, u, xmin, xmax,
gh_fcn, hess_fcn, opt
Outputs:
X : solution vector
F : final objective function value
EXITFLAG : exit flag
1 = converged
negative values = Artelys Knitro specific failure codes
(see KNITROMATLAB documentation for details)
https://www.artelys.com/docs/knitro/3_referenceManual/knitromatlabReference.html#return-codes-exit-flags
OUTPUT : KNITROMATLAB output struct
(see KNITROMATLAB documentation for details)
https://www.artelys.com/docs/knitro/3_referenceManual/knitromatlabReference.html#output-structure-fields
iterations - number of iterations performed
funcCount - number of function evaluations
constrviolation - maximum of constraint violations
firstorderopt - measure of first-order optimality
LAMBDA : 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.
Calling syntax options:
[x, f, exitflag, output, lambda] = ...
nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
x = nlps_knitro(f_fcn, x0);
x = nlps_knitro(f_fcn, x0, A, l);
x = nlps_knitro(f_fcn, x0, A, l, u);
x = nlps_knitro(f_fcn, x0, A, l, u, xmin);
x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax);
x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn);
x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn);
x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
x = nlps_knitro(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 = nlps_knitro(...);
[x, f] = nlps_knitro(...);
[x, f, exitflag] = nlps_knitro(...);
[x, f, exitflag, output] = nlps_knitro(...);
[x, f, exitflag, output, lambda] = nlps_knitro(...);
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] = nlps_knitro(problem);
See also NLPS_MASTER, KNITROMATLAB.0001 function [x, f, eflag, output, lambda] = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt) 0002 %NLPS_KNITRO Nonlinear programming (NLP) Solver based on Artelys Knitro. 0003 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... 0004 % NLPS_KNITRO(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT) 0005 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_KNITRO(PROBLEM) 0006 % A wrapper function providing a standardized interface for using 0007 % Artelys Knitro to solve the following NLP (nonlinear programming) problem: 0008 % 0009 % Minimize a function F(X) beginning from a starting point X0, subject to 0010 % optional linear and nonlinear constraints and variable bounds. 0011 % 0012 % min F(X) 0013 % X 0014 % 0015 % subject to 0016 % 0017 % G(X) = 0 (nonlinear equalities) 0018 % H(X) <= 0 (nonlinear inequalities) 0019 % L <= A*X <= U (linear constraints) 0020 % XMIN <= X <= XMAX (variable bounds) 0021 % 0022 % Inputs (all optional except F_FCN and X0): 0023 % F_FCN : handle to function that evaluates the objective function, 0024 % its gradients and Hessian for a given value of X. If there 0025 % are nonlinear constraints, the Hessian information is 0026 % provided by the HESS_FCN function passed in the 9th argument 0027 % and is not required here. Calling syntax for this function: 0028 % [F, DF, D2F] = F_FCN(X) 0029 % X0 : starting value of optimization vector X 0030 % A, L, U : define the optional linear constraints. Default 0031 % values for the elements of L and U are -Inf and Inf, 0032 % respectively. 0033 % XMIN, XMAX : optional lower and upper bounds on the 0034 % X variables, defaults are -Inf and Inf, respectively. 0035 % GH_FCN : handle to function that evaluates the optional 0036 % nonlinear constraints and their gradients for a given 0037 % value of X. Calling syntax for this function is: 0038 % [H, G, DH, DG] = GH_FCN(X) 0039 % where the columns of DH and DG are the gradients of the 0040 % corresponding elements of H and G, i.e. DH and DG are 0041 % transposes of the Jacobians of H and G, respectively. 0042 % HESS_FCN : handle to function that computes the Hessian of the 0043 % Lagrangian for given values of X, lambda and mu, where 0044 % lambda and mu are the multipliers on the equality and 0045 % inequality constraints, g and h, respectively. The calling 0046 % syntax for this function is: 0047 % LXX = HESS_FCN(X, LAM) 0048 % where lambda = LAM.eqnonlin and mu = LAM.ineqnonlin. 0049 % OPT : optional options structure with the following fields, 0050 % all of which are also optional (default values shown in 0051 % parentheses) 0052 % verbose (0) - controls level of progress output displayed 0053 % 0 = no progress output 0054 % 1 = some progress output 0055 % 2 = verbose progress output 0056 % knitro_opt - options struct for Artelys Knitro, value in verbose 0057 % overrides these options 0058 % opts - struct of other values to be passed directly to 0059 % Aretelys Knitro via an options file 0060 % tol_x - termination tol on x 0061 % tol_f - termination tol on f 0062 % maxit - maximum number of iterations 0063 % PROBLEM : The inputs can alternatively be supplied in a single 0064 % PROBLEM struct with fields corresponding to the input arguments 0065 % described above: f_fcn, x0, A, l, u, xmin, xmax, 0066 % gh_fcn, hess_fcn, opt 0067 % 0068 % Outputs: 0069 % X : solution vector 0070 % F : final objective function value 0071 % EXITFLAG : exit flag 0072 % 1 = converged 0073 % negative values = Artelys Knitro specific failure codes 0074 % (see KNITROMATLAB documentation for details) 0075 % https://www.artelys.com/docs/knitro/3_referenceManual/knitromatlabReference.html#return-codes-exit-flags 0076 % OUTPUT : KNITROMATLAB output struct 0077 % (see KNITROMATLAB documentation for details) 0078 % https://www.artelys.com/docs/knitro/3_referenceManual/knitromatlabReference.html#output-structure-fields 0079 % iterations - number of iterations performed 0080 % funcCount - number of function evaluations 0081 % constrviolation - maximum of constraint violations 0082 % firstorderopt - measure of first-order optimality 0083 % LAMBDA : struct containing the Langrange and Kuhn-Tucker 0084 % multipliers on the constraints, with fields: 0085 % eqnonlin - nonlinear equality constraints 0086 % ineqnonlin - nonlinear inequality constraints 0087 % mu_l - lower (left-hand) limit on linear constraints 0088 % mu_u - upper (right-hand) limit on linear constraints 0089 % lower - lower bound on optimization variables 0090 % upper - upper bound on optimization variables 0091 % 0092 % Note the calling syntax is almost identical to that of FMINCON 0093 % from MathWorks' Optimization Toolbox. The main difference is that 0094 % the linear constraints are specified with A, L, U instead of 0095 % A, B, Aeq, Beq. The functions for evaluating the objective 0096 % function, constraints and Hessian are identical. 0097 % 0098 % Calling syntax options: 0099 % [x, f, exitflag, output, lambda] = ... 0100 % nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt); 0101 % 0102 % x = nlps_knitro(f_fcn, x0); 0103 % x = nlps_knitro(f_fcn, x0, A, l); 0104 % x = nlps_knitro(f_fcn, x0, A, l, u); 0105 % x = nlps_knitro(f_fcn, x0, A, l, u, xmin); 0106 % x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax); 0107 % x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn); 0108 % x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn); 0109 % x = nlps_knitro(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt); 0110 % x = nlps_knitro(problem); 0111 % where problem is a struct with fields: 0112 % f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt 0113 % all fields except 'f_fcn' and 'x0' are optional 0114 % x = nlps_knitro(...); 0115 % [x, f] = nlps_knitro(...); 0116 % [x, f, exitflag] = nlps_knitro(...); 0117 % [x, f, exitflag, output] = nlps_knitro(...); 0118 % [x, f, exitflag, output, lambda] = nlps_knitro(...); 0119 % 0120 % Example: (problem from https://en.wikipedia.org/wiki/Nonlinear_programming) 0121 % function [f, df, d2f] = f2(x) 0122 % f = -x(1)*x(2) - x(2)*x(3); 0123 % if nargout > 1 %% gradient is required 0124 % df = -[x(2); x(1)+x(3); x(2)]; 0125 % if nargout > 2 %% Hessian is required 0126 % d2f = -[0 1 0; 1 0 1; 0 1 0]; %% actually not used since 0127 % end %% 'hess_fcn' is provided 0128 % end 0129 % 0130 % function [h, g, dh, dg] = gh2(x) 0131 % h = [ 1 -1 1; 1 1 1] * x.^2 + [-2; -10]; 0132 % dh = 2 * [x(1) x(1); -x(2) x(2); x(3) x(3)]; 0133 % g = []; dg = []; 0134 % 0135 % function Lxx = hess2(x, lam, cost_mult) 0136 % if nargin < 3, cost_mult = 1; end 0137 % mu = lam.ineqnonlin; 0138 % Lxx = cost_mult * [0 -1 0; -1 0 -1; 0 -1 0] + ... 0139 % [2*[1 1]*mu 0 0; 0 2*[-1 1]*mu 0; 0 0 2*[1 1]*mu]; 0140 % 0141 % problem = struct( ... 0142 % 'f_fcn', @(x)f2(x), ... 0143 % 'gh_fcn', @(x)gh2(x), ... 0144 % 'hess_fcn', @(x, lam, cost_mult)hess2(x, lam, cost_mult), ... 0145 % 'x0', [1; 1; 0], ... 0146 % 'opt', struct('verbose', 2) ... 0147 % ); 0148 % [x, f, exitflag, output, lambda] = nlps_knitro(problem); 0149 % 0150 % See also NLPS_MASTER, KNITROMATLAB. 0151 0152 % MP-Opt-Model 0153 % Copyright (c) 2010-2020, Power Systems Engineering Research Center (PSERC) 0154 % by Ray Zimmerman, PSERC Cornell 0155 % 0156 % This file is part of MP-Opt-Model. 0157 % Covered by the 3-clause BSD License (see LICENSE file for details). 0158 % See https://github.com/MATPOWER/mp-opt-model for more info. 0159 0160 %% options 0161 use_ktropts_file = 1; %% use a Knitro options file to pass options 0162 create_ktropts_file = 0; %% generate a Knitro options file on the fly 0163 0164 %%----- input argument handling ----- 0165 %% gather inputs 0166 if nargin == 1 && isstruct(f_fcn) %% problem struct 0167 p = f_fcn; 0168 f_fcn = p.f_fcn; 0169 x0 = p.x0; 0170 nx = size(x0, 1); %% number of optimization variables 0171 if isfield(p, 'opt'), opt = p.opt; else, opt = []; end 0172 if isfield(p, 'hess_fcn'), hess_fcn = p.hess_fcn; else, hess_fcn = ''; end 0173 if isfield(p, 'gh_fcn'), gh_fcn = p.gh_fcn; else, gh_fcn = ''; end 0174 if isfield(p, 'xmax'), xmax = p.xmax; else, xmax = []; end 0175 if isfield(p, 'xmin'), xmin = p.xmin; else, xmin = []; end 0176 if isfield(p, 'u'), u = p.u; else, u = []; end 0177 if isfield(p, 'l'), l = p.l; else, l = []; end 0178 if isfield(p, 'A'), A = p.A; else, A=sparse(0,nx); end 0179 else %% individual args 0180 nx = size(x0, 1); %% number of optimization variables 0181 if nargin < 10 0182 opt = []; 0183 if nargin < 9 0184 hess_fcn = ''; 0185 if nargin < 8 0186 gh_fcn = ''; 0187 if nargin < 7 0188 xmax = []; 0189 if nargin < 6 0190 xmin = []; 0191 if nargin < 5 0192 u = []; 0193 if nargin < 4 0194 l = []; 0195 A = sparse(0,nx); 0196 end 0197 end 0198 end 0199 end 0200 end 0201 end 0202 end 0203 end 0204 0205 %% set default argument values if missing 0206 if isempty(A) || (~isempty(A) && (isempty(l) || all(l == -Inf)) && ... 0207 (isempty(u) || all(u == Inf))) 0208 A = sparse(0,nx); %% no limits => no linear constraints 0209 end 0210 nA = size(A, 1); %% number of original linear constraints 0211 if isempty(u) %% By default, linear inequalities are ... 0212 u = Inf(nA, 1); %% ... unbounded above and ... 0213 end 0214 if isempty(l) 0215 l = -Inf(nA, 1); %% ... unbounded below. 0216 end 0217 if isempty(xmin) %% By default, optimization variables are ... 0218 xmin = -Inf(nx, 1); %% ... unbounded below and ... 0219 end 0220 if isempty(xmax) 0221 xmax = Inf(nx, 1); %% ... unbounded above. 0222 end 0223 if isempty(gh_fcn) 0224 nonlinear = false; %% no nonlinear constraints present 0225 else 0226 nonlinear = true; %% we have some nonlinear constraints 0227 end 0228 0229 %% default options 0230 if ~isempty(opt) && isfield(opt, 'verbose') && ~isempty(opt.verbose) 0231 verbose = opt.verbose; 0232 else 0233 verbose = 0; 0234 end 0235 0236 %% split l <= A*x <= u into less than, equal to, greater than, and 0237 %% doubly-bounded sets 0238 ieq = find( abs(u-l) <= eps ); %% equality 0239 igt = find( u >= 1e10 & l > -1e10 ); %% greater than, unbounded above 0240 ilt = find( l <= -1e10 & u < 1e10 ); %% less than, unbounded below 0241 ibx = find( (abs(u-l) > eps) & (u < 1e10) & (l > -1e10) ); 0242 Af = [ A(ilt, :); -A(igt, :); A(ibx, :); -A(ibx, :) ]; 0243 bf = [ u(ilt); -l(igt); u(ibx); -l(ibx)]; 0244 Afeq = A(ieq, :); 0245 bfeq = u(ieq); 0246 0247 %%----- set up problem ----- 0248 %% build Jacobian and Hessian structure 0249 randx = rand(size(x0)); 0250 nonz = 1e-20; 0251 if nonlinear 0252 [h, g, dhs, dgs] = gh_fcn(randx); 0253 if issparse(dhs) 0254 dhs(dhs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0255 dgs(dgs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0256 Js = [dhs'; dgs']; 0257 else 0258 dhs = []; dgs = []; 0259 Js = []; 0260 end 0261 else 0262 g = []; h = []; 0263 dhs = []; dgs = dhs; Js = dhs; 0264 end 0265 neq = length(g); 0266 niq = length(h); 0267 if isempty(hess_fcn) 0268 [f_, df_, Hs] = f_fcn(randx); %% cost 0269 else 0270 lam = struct('eqnonlin', rand(neq, 1), 'ineqnonlin', rand(niq, 1)); 0271 Hs = hess_fcn(randx, lam, 1); 0272 end 0273 if issparse(Hs) 0274 Hs(Hs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0275 else 0276 Hs = []; 0277 end 0278 0279 %% function definitions 0280 ktr_gh_fcn = @(x)constraints(x, gh_fcn, dgs, dhs); 0281 ktr_hess_fcn = @(x, lambda)hessian(x, lambda, f_fcn, hess_fcn, Hs); 0282 0283 %% basic optimset options needed for ktrlink 0284 if isfield(opt, 'knitro_opt') 0285 knitro_opt = opt.knitro_opt; 0286 else 0287 knitro_opt = struct(); 0288 end 0289 kopts = struct(); 0290 kopts.GradObj = 'on'; 0291 kopts.GradConstr = 'on'; 0292 kopts.HessFcn = ktr_hess_fcn; 0293 kopts.Hessian = 'user-supplied'; 0294 if issparse(Js) 0295 kopts.JacobPattern = Js; 0296 end 0297 if issparse(Hs) 0298 kopts.HessPattern = Hs; 0299 end 0300 0301 if use_ktropts_file 0302 if isfield(knitro_opt, 'opt_fname') && ~isempty(knitro_opt.opt_fname) 0303 opt_fname = knitro_opt.opt_fname; 0304 elseif isfield(knitro_opt, 'opt') && knitro_opt.opt 0305 opt_fname = sprintf('knitro_user_options_%d.txt', knitro_opt.opt); 0306 else 0307 %% create ktropts file 0308 ktropts.algorithm = 1; 0309 ktropts.bar_directinterval = 0; 0310 ktropts.outlev = verbose; 0311 if isfield(knitro_opt, 'opts') %% raw Knitro options for options file 0312 ktropts = nested_struct_copy(ktropts, knitro_opt.opts); 0313 end 0314 if isfield(knitro_opt, 'tol_x') && ~isempty(knitro_opt.tol_x) 0315 ktropts.xtol = knitro_opt.tol_x; 0316 end 0317 if isfield(knitro_opt, 'tol_f') && ~isempty(knitro_opt.tol_f) 0318 ktropts.opttol = knitro_opt.tol_f; 0319 end 0320 if isfield(knitro_opt, 'maxit') && knitro_opt.maxit ~= 0 0321 ktropts.maxit = knitro_opt.maxit; 0322 end 0323 0324 opt_fname = write_ktropts(ktropts); 0325 create_ktropts_file = 1; %% make a note that I created it 0326 end 0327 else 0328 if isfield(knitro_opt, 'opts') %% raw Knitro options for optimset() 0329 kopts = nested_struct_copy(kopts, knitro_opt.opts); 0330 end 0331 kopts.Algorithm = 'interior-point'; 0332 if isfield(knitro_opt, 'tol_x') && ~isempty(knitro_opt.tol_x) 0333 kopts.TolX = knitro_opt.tol_x; 0334 end 0335 if isfield(knitro_opt, 'tol_f') && ~isempty(knitro_opt.tol_f) 0336 kopts.TolFun = knitro_opt.tol_f; 0337 end 0338 if isfield(knitro_opt, 'maxit') && knitro_opt.maxit ~= 0 0339 kopts.MaxIter = knitro_opt.maxit; 0340 kopts.MaxFunEvals = 4 * knitro_opt.maxit; 0341 end 0342 if verbose > 1 0343 kopts.Display = 'iter'; 0344 elseif verbose == 1 0345 kopts.Display = 'final'; 0346 else 0347 kopts.Display = 'off'; 0348 end 0349 opt_fname = []; 0350 end 0351 fmoptions = optimset(kopts); 0352 0353 %%----- run solver ----- 0354 if have_feature('knitromatlab') 0355 [x, f, eflag, output, Lambda] = knitromatlab(f_fcn, x0, Af, bf, Afeq, bfeq, ... 0356 xmin, xmax, ktr_gh_fcn, [], fmoptions, opt_fname); 0357 else 0358 [x, f, eflag, output, Lambda] = ktrlink(f_fcn, x0, Af, bf, Afeq, bfeq, ... 0359 xmin, xmax, ktr_gh_fcn, fmoptions, opt_fname); 0360 end 0361 success = (eflag == 0); 0362 if success 0363 eflag = 1; %% success is 1 (not zero), all other values are Knitro return codes 0364 end 0365 0366 %% delete ktropts file 0367 if create_ktropts_file %% ... but only if I created it 0368 delete(opt_fname); 0369 end 0370 0371 %% fix Lambdas 0372 %% (shadow prices on equality variable bounds can come back on wrong limit) 0373 kl = find(Lambda.lower > 0 & Lambda.upper == 0); 0374 Lambda.upper(kl) = Lambda.lower(kl); 0375 Lambda.lower(kl) = 0; 0376 ku = find(Lambda.upper < 0 & Lambda.lower == 0); 0377 Lambda.lower(ku) = Lambda.upper(ku); 0378 Lambda.upper(ku) = 0; 0379 0380 %% package up results 0381 nlt = length(ilt); 0382 ngt = length(igt); 0383 nbx = length(ibx); 0384 0385 %% extract multipliers for linear constraints 0386 kl = find(Lambda.eqlin < 0); 0387 ku = find(Lambda.eqlin > 0); 0388 0389 mu_l = zeros(size(u)); 0390 mu_l(ieq(kl)) = -Lambda.eqlin(kl); 0391 mu_l(igt) = Lambda.ineqlin(nlt+(1:ngt)); 0392 mu_l(ibx) = Lambda.ineqlin(nlt+ngt+nbx+(1:nbx)); 0393 0394 mu_u = zeros(size(u)); 0395 mu_u(ieq(ku)) = Lambda.eqlin(ku); 0396 mu_u(ilt) = Lambda.ineqlin(1:nlt); 0397 mu_u(ibx) = Lambda.ineqlin(nlt+ngt+(1:nbx)); 0398 0399 lambda = struct( ... 0400 'lower', -Lambda.lower, ... 0401 'upper', Lambda.upper, ... 0402 'eqnonlin', Lambda.eqnonlin, ... 0403 'ineqnonlin', Lambda.ineqnonlin, ... 0404 'mu_l', mu_l, ... 0405 'mu_u', mu_u ); 0406 0407 0408 %----- callback functions ----- 0409 function [h, g, dh, dg] = constraints(x, gh_fcn, dgs, dhs) 0410 if isempty(gh_fcn) 0411 nx = length(x); 0412 h = []; g = h; dh = sparse(0, nx); dg = dh; 0413 else 0414 [h, g, dh, dg] = gh_fcn(x); 0415 %% add sparse structure (with tiny values) to current matrices to 0416 %% ensure that sparsity structure matches that supplied 0417 if issparse(dh) 0418 dg = dg + dgs; 0419 dh = dh + dhs; 0420 end 0421 end 0422 0423 function H = hessian(x, lambda, f_fcn, hess_fcn, Hs) 0424 if isempty(hess_fcn) 0425 [f_, df_, H] = f_fcn(x); %% cost 0426 else 0427 %% add sparse structure (with tiny values) to current matrices to 0428 %% ensure that sparsity structure matches that supplied 0429 H = hess_fcn(x, lambda, 1); 0430 if issparse(H) 0431 H = H + Hs; 0432 end 0433 end 0434 0435 %%----- write_ktropts ----- 0436 function fname = write_ktropts(ktropts) 0437 0438 %% generate file name 0439 fname = sprintf('ktropts_%06d.txt', fix(1e6*rand)); 0440 0441 %% open file 0442 [fd, msg] = fopen(fname, 'wt'); %% write options file 0443 if fd == -1 0444 error('could not create %d : %s', fname, msg); 0445 end 0446 0447 %% write options 0448 fields = fieldnames(ktropts); 0449 for k = 1:length(fields) 0450 fprintf(fd, '%s %g\n', fields{k}, ktropts.(fields{k})); 0451 end 0452 0453 %% close file 0454 if fd ~= 1 0455 fclose(fd); 0456 end