NLPS_IPOPT Nonlinear programming (NLP) Solver based on IPOPT.
[X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
NLPS_IPOPT(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT)
[X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_IPOPT(PROBLEM)
A wrapper function providing a standardized interface for using
IPOPT 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
ipopt_opt - options struct for IPOPT, value in verbose
overrides these options
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
0 = failed to converge
OUTPUT : output struct with the following fields:
status - see IPOPT documentation for INFO.status
https://coin-or.github.io/Ipopt/IpReturnCodes__inc_8h_source.html
iterations - number of iterations performed (INFO.iter)
cpu - see IPOPT documentation for INFO.cpu
eval - see IPOPT documentation for INFO.eval
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_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
x = nlps_ipopt(f_fcn, x0);
x = nlps_ipopt(f_fcn, x0, A, l);
x = nlps_ipopt(f_fcn, x0, A, l, u);
x = nlps_ipopt(f_fcn, x0, A, l, u, xmin);
x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax);
x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn);
x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn);
x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
x = nlps_ipopt(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_ipopt(...);
[x, f] = nlps_ipopt(...);
[x, f, exitflag] = nlps_ipopt(...);
[x, f, exitflag, output] = nlps_ipopt(...);
[x, f, exitflag, output, lambda] = nlps_ipopt(...);
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_ipopt(problem);
See also NLPS_MASTER, IPOPT.0001 function [x, f, eflag, output, lambda] = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt) 0002 %NLPS_IPOPT Nonlinear programming (NLP) Solver based on IPOPT. 0003 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... 0004 % NLPS_IPOPT(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT) 0005 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_IPOPT(PROBLEM) 0006 % A wrapper function providing a standardized interface for using 0007 % IPOPT 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 % ipopt_opt - options struct for IPOPT, value in verbose 0057 % overrides these options 0058 % PROBLEM : The inputs can alternatively be supplied in a single 0059 % PROBLEM struct with fields corresponding to the input arguments 0060 % described above: f_fcn, x0, A, l, u, xmin, xmax, 0061 % gh_fcn, hess_fcn, opt 0062 % 0063 % Outputs: 0064 % X : solution vector 0065 % F : final objective function value 0066 % EXITFLAG : exit flag 0067 % 1 = converged 0068 % 0 = failed to converge 0069 % OUTPUT : output struct with the following fields: 0070 % status - see IPOPT documentation for INFO.status 0071 % https://coin-or.github.io/Ipopt/IpReturnCodes__inc_8h_source.html 0072 % iterations - number of iterations performed (INFO.iter) 0073 % cpu - see IPOPT documentation for INFO.cpu 0074 % eval - see IPOPT documentation for INFO.eval 0075 % LAMBDA : struct containing the Langrange and Kuhn-Tucker 0076 % multipliers on the constraints, with fields: 0077 % eqnonlin - nonlinear equality constraints 0078 % ineqnonlin - nonlinear inequality constraints 0079 % mu_l - lower (left-hand) limit on linear constraints 0080 % mu_u - upper (right-hand) limit on linear constraints 0081 % lower - lower bound on optimization variables 0082 % upper - upper bound on optimization variables 0083 % 0084 % Note the calling syntax is almost identical to that of FMINCON 0085 % from MathWorks' Optimization Toolbox. The main difference is that 0086 % the linear constraints are specified with A, L, U instead of 0087 % A, B, Aeq, Beq. The functions for evaluating the objective 0088 % function, constraints and Hessian are identical. 0089 % 0090 % Calling syntax options: 0091 % [x, f, exitflag, output, lambda] = ... 0092 % nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt); 0093 % 0094 % x = nlps_ipopt(f_fcn, x0); 0095 % x = nlps_ipopt(f_fcn, x0, A, l); 0096 % x = nlps_ipopt(f_fcn, x0, A, l, u); 0097 % x = nlps_ipopt(f_fcn, x0, A, l, u, xmin); 0098 % x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax); 0099 % x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn); 0100 % x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn); 0101 % x = nlps_ipopt(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt); 0102 % x = nlps_ipopt(problem); 0103 % where problem is a struct with fields: 0104 % f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt 0105 % all fields except 'f_fcn' and 'x0' are optional 0106 % x = nlps_ipopt(...); 0107 % [x, f] = nlps_ipopt(...); 0108 % [x, f, exitflag] = nlps_ipopt(...); 0109 % [x, f, exitflag, output] = nlps_ipopt(...); 0110 % [x, f, exitflag, output, lambda] = nlps_ipopt(...); 0111 % 0112 % Example: (problem from https://en.wikipedia.org/wiki/Nonlinear_programming) 0113 % function [f, df, d2f] = f2(x) 0114 % f = -x(1)*x(2) - x(2)*x(3); 0115 % if nargout > 1 %% gradient is required 0116 % df = -[x(2); x(1)+x(3); x(2)]; 0117 % if nargout > 2 %% Hessian is required 0118 % d2f = -[0 1 0; 1 0 1; 0 1 0]; %% actually not used since 0119 % end %% 'hess_fcn' is provided 0120 % end 0121 % 0122 % function [h, g, dh, dg] = gh2(x) 0123 % h = [ 1 -1 1; 1 1 1] * x.^2 + [-2; -10]; 0124 % dh = 2 * [x(1) x(1); -x(2) x(2); x(3) x(3)]; 0125 % g = []; dg = []; 0126 % 0127 % function Lxx = hess2(x, lam, cost_mult) 0128 % if nargin < 3, cost_mult = 1; end 0129 % mu = lam.ineqnonlin; 0130 % Lxx = cost_mult * [0 -1 0; -1 0 -1; 0 -1 0] + ... 0131 % [2*[1 1]*mu 0 0; 0 2*[-1 1]*mu 0; 0 0 2*[1 1]*mu]; 0132 % 0133 % problem = struct( ... 0134 % 'f_fcn', @(x)f2(x), ... 0135 % 'gh_fcn', @(x)gh2(x), ... 0136 % 'hess_fcn', @(x, lam, cost_mult)hess2(x, lam, cost_mult), ... 0137 % 'x0', [1; 1; 0], ... 0138 % 'opt', struct('verbose', 2) ... 0139 % ); 0140 % [x, f, exitflag, output, lambda] = nlps_ipopt(problem); 0141 % 0142 % See also NLPS_MASTER, IPOPT. 0143 0144 % MP-Opt-Model 0145 % Copyright (c) 2010-2020, Power Systems Engineering Research Center (PSERC) 0146 % by Ray Zimmerman, PSERC Cornell 0147 % 0148 % This file is part of MP-Opt-Model. 0149 % Covered by the 3-clause BSD License (see LICENSE file for details). 0150 % See https://github.com/MATPOWER/mp-opt-model for more info. 0151 0152 %%----- input argument handling ----- 0153 %% gather inputs 0154 if nargin == 1 && isstruct(f_fcn) %% problem struct 0155 p = f_fcn; 0156 f_fcn = p.f_fcn; 0157 x0 = p.x0; 0158 nx = size(x0, 1); %% number of optimization variables 0159 if isfield(p, 'opt'), opt = p.opt; else, opt = []; end 0160 if isfield(p, 'hess_fcn'), hess_fcn = p.hess_fcn; else, hess_fcn = ''; end 0161 if isfield(p, 'gh_fcn'), gh_fcn = p.gh_fcn; else, gh_fcn = ''; end 0162 if isfield(p, 'xmax'), xmax = p.xmax; else, xmax = []; end 0163 if isfield(p, 'xmin'), xmin = p.xmin; else, xmin = []; end 0164 if isfield(p, 'u'), u = p.u; else, u = []; end 0165 if isfield(p, 'l'), l = p.l; else, l = []; end 0166 if isfield(p, 'A'), A = p.A; else, A=sparse(0,nx); end 0167 else %% individual args 0168 nx = size(x0, 1); %% number of optimization variables 0169 if nargin < 10 0170 opt = []; 0171 if nargin < 9 0172 hess_fcn = ''; 0173 if nargin < 8 0174 gh_fcn = ''; 0175 if nargin < 7 0176 xmax = []; 0177 if nargin < 6 0178 xmin = []; 0179 if nargin < 5 0180 u = []; 0181 if nargin < 4 0182 l = []; 0183 A = sparse(0,nx); 0184 end 0185 end 0186 end 0187 end 0188 end 0189 end 0190 end 0191 end 0192 0193 %% set default argument values if missing 0194 if isempty(A) || (~isempty(A) && (isempty(l) || all(l == -Inf)) && ... 0195 (isempty(u) || all(u == Inf))) 0196 A = sparse(0,nx); %% no limits => no linear constraints 0197 end 0198 nA = size(A, 1); %% number of original linear constraints 0199 if isempty(u) %% By default, linear inequalities are ... 0200 u = Inf(nA, 1); %% ... unbounded above and ... 0201 end 0202 if isempty(l) 0203 l = -Inf(nA, 1); %% ... unbounded below. 0204 end 0205 if isempty(xmin) %% By default, optimization variables are ... 0206 xmin = -Inf(nx, 1); %% ... unbounded below and ... 0207 end 0208 if isempty(xmax) 0209 xmax = Inf(nx, 1); %% ... unbounded above. 0210 end 0211 if isempty(gh_fcn) 0212 nonlinear = false; %% no nonlinear constraints present 0213 else 0214 nonlinear = true; %% we have some nonlinear constraints 0215 end 0216 0217 %% default options 0218 if ~isempty(opt) && isfield(opt, 'verbose') && ~isempty(opt.verbose) 0219 verbose = opt.verbose; 0220 else 0221 verbose = 0; 0222 end 0223 0224 %% make sure A is sparse 0225 if ~issparse(A) 0226 A = sparse(A); 0227 end 0228 0229 %% replace equality variable bounds with an equality constraint 0230 %% (since IPOPT does not return shadow prices on variables that it eliminates) 0231 kk = find(xmin == xmax); 0232 nk = length(kk); 0233 if nk 0234 A = [ A; sparse((1:nk)', kk, 1, nk, nx) ]; 0235 l = [ l; xmin(kk) ]; 0236 u = [ u; xmax(kk) ]; 0237 xmin(kk) = -Inf; 0238 xmax(kk) = Inf; 0239 nA = size(A, 1); %% updated number of linear constraints 0240 end 0241 0242 %%----- set up problem ----- 0243 %% build Jacobian and Hessian structure 0244 randx = rand(size(x0)); 0245 nonz = 1e-20; 0246 if nonlinear 0247 [h, g, dhs, dgs] = gh_fcn(randx); 0248 dhs(dhs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0249 dgs(dgs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0250 Js = [dgs'; dhs'; A]; 0251 else 0252 g = []; h = []; 0253 dhs = sparse(0, nx); dgs = dhs; Js = A; 0254 end 0255 neq = length(g); 0256 niq = length(h); 0257 if isempty(hess_fcn) 0258 [f_, df_, Hs] = f_fcn(randx); %% cost 0259 else 0260 lam = struct('eqnonlin', rand(neq, 1), 'ineqnonlin', rand(niq, 1)); 0261 Hs = hess_fcn(randx, lam, 1); 0262 end 0263 if ~issparse(Hs), Hs = sparse(Hs); end %% convert to sparse if necessary 0264 Hs(Hs ~= 0) = nonz; %% set non-zero entries to tiny value (for adding later) 0265 0266 %% set options struct for IPOPT 0267 if ~isempty(opt) && isfield(opt, 'ipopt_opt') && ~isempty(opt.ipopt_opt) 0268 options.ipopt = ipopt_options(opt.ipopt_opt); 0269 else 0270 options.ipopt = ipopt_options; 0271 end 0272 if verbose 0273 options.ipopt.print_level = min(12, verbose*2+1); 0274 else 0275 options.ipopt.print_level = 0; 0276 end 0277 0278 %% extra data to pass to functions 0279 options.auxdata = struct( ... 0280 'f_fcn', f_fcn, ... 0281 'gh_fcn', gh_fcn, ... 0282 'hess_fcn', hess_fcn, ... 0283 'A', A, ... 0284 'nx', nx, ... 0285 'nA', nA, ... 0286 'neqnln', neq, ... 0287 'niqnln', niq, ... 0288 'dgs', dgs, ... 0289 'dhs', dhs, ... 0290 'Hs', Hs ); 0291 0292 %% define variable and constraint bounds 0293 options.lb = xmin; 0294 options.ub = xmax; 0295 options.cl = [zeros(neq, 1); -Inf(niq, 1); l]; 0296 options.cu = [zeros(neq, 1); zeros(niq, 1); u]; 0297 0298 %% assign function handles 0299 funcs.objective = @objective; 0300 funcs.gradient = @gradient; 0301 funcs.constraints = @constraints; 0302 funcs.jacobian = @jacobian; 0303 funcs.hessian = @hessian; 0304 funcs.jacobianstructure = @(d) Js; 0305 funcs.hessianstructure = @(d) tril(Hs); 0306 0307 %%----- run solver ----- 0308 %% run the optimization 0309 if have_feature('ipopt_auxdata') 0310 [x, info] = ipopt_auxdata(x0,funcs,options); 0311 else 0312 [x, info] = ipopt(x0,funcs,options); 0313 end 0314 0315 if info.status == 0 || info.status == 1 0316 eflag = 1; 0317 else 0318 eflag = 0; 0319 end 0320 output = struct('status', info.status); 0321 if isfield(info, 'iter') 0322 output.iterations = info.iter; 0323 else 0324 output.iterations = []; 0325 end 0326 if isfield(info, 'cpu') 0327 output.cpu = info.cpu; 0328 end 0329 if isfield(info, 'eval') 0330 output.eval = info.eval; 0331 end 0332 f = f_fcn(x); 0333 0334 %% check for empty results (in case optimization failed) 0335 if isempty(info.lambda) 0336 lam = NaN(nA, 1); 0337 else 0338 lam = info.lambda; 0339 end 0340 if isempty(info.zl) && ~isempty(xmin) 0341 zl = NaN(nx, 1); 0342 else 0343 zl = info.zl; 0344 end 0345 if isempty(info.zu) && ~isempty(xmax) 0346 zu = NaN(nx, 1); 0347 else 0348 zu = info.zu; 0349 end 0350 0351 %% extract shadow prices for equality var bounds converted to eq constraints 0352 %% (since IPOPT does not return shadow prices on variables that it eliminates) 0353 if nk 0354 offset = neq + niq + nA - nk; 0355 lam_tmp = lam(offset+(1:nk)); 0356 kl = find(lam_tmp < 0); %% lower bound binding 0357 ku = find(lam_tmp > 0); %% upper bound binding 0358 zl(kk(kl)) = -lam_tmp(kl); 0359 zu(kk(ku)) = lam_tmp(ku); 0360 lam(offset+(1:nk)) = []; %% remove these shadow prices 0361 nA = nA - nk; %% reduce dimension accordingly 0362 end 0363 0364 %% extract multipliers for linear constraints 0365 lam_lin = lam(neq+niq+(1:nA)); %% lambda for linear constraints 0366 kl = find(lam_lin < 0); %% lower bound binding 0367 ku = find(lam_lin > 0); %% upper bound binding 0368 mu_l = zeros(nA, 1); 0369 mu_l(kl) = -lam_lin(kl); 0370 mu_u = zeros(nA, 1); 0371 mu_u(ku) = lam_lin(ku); 0372 0373 lambda = struct( ... 0374 'lower', zl, ... 0375 'upper', zu, ... 0376 'eqnonlin', lam(1:neq), ... 0377 'ineqnonlin', lam(neq+(1:niq)), ... 0378 'mu_l', mu_l, ... 0379 'mu_u', mu_u ); 0380 0381 0382 %----- callback functions ----- 0383 function f = objective(x, d) 0384 f = d.f_fcn(x); 0385 0386 function df = gradient(x, d) 0387 [f, df] = d.f_fcn(x); 0388 0389 function c = constraints(x, d) 0390 if isempty(d.gh_fcn) 0391 c = d.A*x; 0392 else 0393 [h, g] = d.gh_fcn(x); 0394 c = [g; h; d.A*x]; 0395 end 0396 0397 function J = jacobian(x, d) 0398 if isempty(d.gh_fcn) 0399 J = d.A; 0400 else 0401 [h, g, dh, dg] = d.gh_fcn(x); 0402 %% add sparse structure (with tiny values) to current matrices to 0403 %% ensure that sparsity structure matches that supplied 0404 J = [(dg + d.dgs)'; (dh + d.dhs)'; d.A]; 0405 end 0406 0407 function H = hessian(x, sigma, lambda, d) 0408 if isempty(d.hess_fcn) 0409 [f_, df_, d2f] = d.f_fcn(x); %% cost 0410 if ~issparse(d2f), d2f = sparse(d2f); end %% convert to sparse if necessary 0411 H = tril(d2f * sigma); 0412 else 0413 lam.eqnonlin = lambda(1:d.neqnln); 0414 lam.ineqnonlin = lambda(d.neqnln+(1:d.niqnln)); 0415 H = d.hess_fcn(x, lam, sigma); 0416 if ~issparse(H), H = sparse(H); end %% convert to sparse if necessary 0417 %% add sparse structure (with tiny values) to current matrices to 0418 %% ensure that sparsity structure matches that supplied 0419 H = tril(H + d.Hs); 0420 end