Home > matpower7.1 > mp-opt-model > lib > nlps_fmincon.m

nlps_fmincon

PURPOSE ^

NLPS_FMINCON Nonlinear programming (NLP) Solver based on FMINCON.

SYNOPSIS ^

function [x, f, eflag, output, lambda] = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt)

DESCRIPTION ^

NLPS_FMINCON  Nonlinear programming (NLP) Solver based on FMINCON.
   [X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
       NLPS_FMINCON(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT)
   [X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_FMINCON(PROBLEM)
   A wrapper function providing a standardized interface for using
   FMINCON from the Optimization Toolbox 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
           fmincon_opt - options struct for FMINCON, value in verbose
                   overrides these options
               opts - struct of other values to be passed directly to
                   OPTIMSET or OPTIMOPTIONS (overridden by fields below)
               alg - (4) algorithm codes for FMINCON
                   1 - active-set (not suitable for large problems)
                   2 - interior-point, w/default 'bfgs' Hessian approx
                   3 - interior-point, w/ 'lbfgs' Hessian approx
                   4 - interior-point, w/exact user-supplied Hessian
                   5 - interior-point, w/Hessian via finite differences
                   6 - sqp (not suitable for large problems)
               tol_x - (1e-4) termination tol on x
               tol_f - (1e-4) termination tol on f
               max_it - 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 : FMINCON exit flag
           (see FMINCON documentation for details)
       OUTPUT : FMINCON output struct
           (see FMINCON documentation for details)
       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_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);

       x = nlps_fmincon(f_fcn, x0);
       x = nlps_fmincon(f_fcn, x0, A, l);
       x = nlps_fmincon(f_fcn, x0, A, l, u);
       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin);
       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax);
       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn);
       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn);
       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
       x = nlps_fmincon(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_fmincon(...);
       [x, f] = nlps_fmincon(...);
       [x, f, exitflag] = nlps_fmincon(...);
       [x, f, exitflag, output] = nlps_fmincon(...);
       [x, f, exitflag, output, lambda] = nlps_fmincon(...);

   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_fmincon(problem);

   See also NLPS_MASTER, FMINCON.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [x, f, eflag, output, lambda] = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt)
0002 %NLPS_FMINCON  Nonlinear programming (NLP) Solver based on FMINCON.
0003 %   [X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
0004 %       NLPS_FMINCON(F_FCN, X0, A, L, U, XMIN, XMAX, GH_FCN, HESS_FCN, OPT)
0005 %   [X, F, EXITFLAG, OUTPUT, LAMBDA] = NLPS_FMINCON(PROBLEM)
0006 %   A wrapper function providing a standardized interface for using
0007 %   FMINCON from the Optimization Toolbox to solve the following NLP
0008 %   (nonlinear programming) problem:
0009 %
0010 %   Minimize a function F(X) beginning from a starting point X0, subject to
0011 %   optional linear and nonlinear constraints and variable bounds.
0012 %
0013 %       min F(X)
0014 %        X
0015 %
0016 %   subject to
0017 %
0018 %       G(X) = 0            (nonlinear equalities)
0019 %       H(X) <= 0           (nonlinear inequalities)
0020 %       L <= A*X <= U       (linear constraints)
0021 %       XMIN <= X <= XMAX   (variable bounds)
0022 %
0023 %   Inputs (all optional except F_FCN and X0):
0024 %       F_FCN : handle to function that evaluates the objective function,
0025 %           its gradients and Hessian for a given value of X. If there
0026 %           are nonlinear constraints, the Hessian information is
0027 %           provided by the HESS_FCN function passed in the 9th argument
0028 %           and is not required here. Calling syntax for this function:
0029 %               [F, DF, D2F] = F_FCN(X)
0030 %       X0 : starting value of optimization vector X
0031 %       A, L, U : define the optional linear constraints. Default
0032 %           values for the elements of L and U are -Inf and Inf,
0033 %           respectively.
0034 %       XMIN, XMAX : optional lower and upper bounds on the
0035 %           X variables, defaults are -Inf and Inf, respectively.
0036 %       GH_FCN : handle to function that evaluates the optional
0037 %           nonlinear constraints and their gradients for a given
0038 %           value of X. Calling syntax for this function is:
0039 %               [H, G, DH, DG] = GH_FCN(X)
0040 %           where the columns of DH and DG are the gradients of the
0041 %           corresponding elements of H and G, i.e. DH and DG are
0042 %           transposes of the Jacobians of H and G, respectively.
0043 %       HESS_FCN : handle to function that computes the Hessian of the
0044 %           Lagrangian for given values of X, lambda and mu, where
0045 %           lambda and mu are the multipliers on the equality and
0046 %           inequality constraints, g and h, respectively. The calling
0047 %           syntax for this function is:
0048 %               LXX = HESS_FCN(X, LAM)
0049 %           where lambda = LAM.eqnonlin and mu = LAM.ineqnonlin.
0050 %       OPT : optional options structure with the following fields,
0051 %           all of which are also optional (default values shown in
0052 %           parentheses)
0053 %           verbose (0) - controls level of progress output displayed
0054 %               0 = no progress output
0055 %               1 = some progress output
0056 %               2 = verbose progress output
0057 %           fmincon_opt - options struct for FMINCON, value in verbose
0058 %                   overrides these options
0059 %               opts - struct of other values to be passed directly to
0060 %                   OPTIMSET or OPTIMOPTIONS (overridden by fields below)
0061 %               alg - (4) algorithm codes for FMINCON
0062 %                   1 - active-set (not suitable for large problems)
0063 %                   2 - interior-point, w/default 'bfgs' Hessian approx
0064 %                   3 - interior-point, w/ 'lbfgs' Hessian approx
0065 %                   4 - interior-point, w/exact user-supplied Hessian
0066 %                   5 - interior-point, w/Hessian via finite differences
0067 %                   6 - sqp (not suitable for large problems)
0068 %               tol_x - (1e-4) termination tol on x
0069 %               tol_f - (1e-4) termination tol on f
0070 %               max_it - maximum number of iterations
0071 %       PROBLEM : The inputs can alternatively be supplied in a single
0072 %           PROBLEM struct with fields corresponding to the input arguments
0073 %           described above: f_fcn, x0, A, l, u, xmin, xmax,
0074 %                            gh_fcn, hess_fcn, opt
0075 %
0076 %   Outputs:
0077 %       X : solution vector
0078 %       F : final objective function value
0079 %       EXITFLAG : FMINCON exit flag
0080 %           (see FMINCON documentation for details)
0081 %       OUTPUT : FMINCON output struct
0082 %           (see FMINCON documentation for details)
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_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
0101 %
0102 %       x = nlps_fmincon(f_fcn, x0);
0103 %       x = nlps_fmincon(f_fcn, x0, A, l);
0104 %       x = nlps_fmincon(f_fcn, x0, A, l, u);
0105 %       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin);
0106 %       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax);
0107 %       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn);
0108 %       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn);
0109 %       x = nlps_fmincon(f_fcn, x0, A, l, u, xmin, xmax, gh_fcn, hess_fcn, opt);
0110 %       x = nlps_fmincon(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_fmincon(...);
0115 %       [x, f] = nlps_fmincon(...);
0116 %       [x, f, exitflag] = nlps_fmincon(...);
0117 %       [x, f, exitflag, output] = nlps_fmincon(...);
0118 %       [x, f, exitflag, output, lambda] = nlps_fmincon(...);
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_fmincon(problem);
0149 %
0150 %   See also NLPS_MASTER, FMINCON.
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 %%----- input argument handling  -----
0161 %% gather inputs
0162 if nargin == 1 && isstruct(f_fcn)       %% problem struct
0163     p = f_fcn;
0164     f_fcn = p.f_fcn;
0165     x0 = p.x0;
0166     nx = size(x0, 1);       %% number of optimization variables
0167     if isfield(p, 'opt'),       opt = p.opt;            else,   opt = [];       end
0168     if isfield(p, 'hess_fcn'),  hess_fcn = p.hess_fcn;  else,   hess_fcn = '';  end
0169     if isfield(p, 'gh_fcn'),    gh_fcn = p.gh_fcn;      else,   gh_fcn = '';    end
0170     if isfield(p, 'xmax'),      xmax = p.xmax;          else,   xmax = [];      end
0171     if isfield(p, 'xmin'),      xmin = p.xmin;          else,   xmin = [];      end
0172     if isfield(p, 'u'),         u = p.u;                else,   u = [];         end
0173     if isfield(p, 'l'),         l = p.l;                else,   l = [];         end
0174     if isfield(p, 'A'),         A = p.A;                else,   A=sparse(0,nx); end
0175 else                                    %% individual args
0176     nx = size(x0, 1);       %% number of optimization variables
0177     if nargin < 10
0178         opt = [];
0179         if nargin < 9
0180             hess_fcn = '';
0181             if nargin < 8
0182                 gh_fcn = '';
0183                 if nargin < 7
0184                     xmax = [];
0185                     if nargin < 6
0186                         xmin = [];
0187                         if nargin < 5
0188                             u = [];
0189                             if nargin < 4
0190                                 l = [];
0191                                 A = sparse(0,nx);
0192                             end
0193                         end
0194                     end
0195                 end
0196             end
0197         end
0198     end
0199 end
0200 
0201 %% set default argument values if missing
0202 if isempty(A) || (~isempty(A) && (isempty(l) || all(l == -Inf)) && ...
0203                                  (isempty(u) || all(u == Inf)))
0204     A = sparse(0,nx);           %% no limits => no linear constraints
0205 end
0206 nA = size(A, 1);                %% number of original linear constraints
0207 if isempty(u)                   %% By default, linear inequalities are ...
0208     u = Inf(nA, 1);             %% ... unbounded above and ...
0209 end
0210 if isempty(l)
0211     l = -Inf(nA, 1);            %% ... unbounded below.
0212 end
0213 if isempty(xmin)                %% By default, optimization variables are ...
0214     xmin = -Inf(nx, 1);         %% ... unbounded below and ...
0215 end
0216 if isempty(xmax)
0217     xmax = Inf(nx, 1);          %% ... unbounded above.
0218 end
0219 if isempty(gh_fcn)
0220     nonlinear = false;          %% no nonlinear constraints present
0221 else
0222     nonlinear = true;           %% we have some nonlinear constraints
0223 end
0224 
0225 %% default options
0226 if ~isempty(opt) && isfield(opt, 'verbose') && ~isempty(opt.verbose)
0227     verbose = opt.verbose;
0228 else
0229     verbose = 0;
0230 end
0231 
0232 %% split l <= A*x <= u into less than, equal to, greater than, and
0233 %% doubly-bounded sets
0234 ieq = find( abs(u-l) <= eps );          %% equality
0235 igt = find( u >=  1e10 & l > -1e10 );   %% greater than, unbounded above
0236 ilt = find( l <= -1e10 & u <  1e10 );   %% less than, unbounded below
0237 ibx = find( (abs(u-l) > eps) & (u < 1e10) & (l > -1e10) );
0238 Af  = [ A(ilt, :); -A(igt, :); A(ibx, :); -A(ibx, :) ];
0239 bf  = [ u(ilt);   -l(igt);     u(ibx);    -l(ibx)];
0240 Afeq = A(ieq, :);
0241 bfeq = u(ieq);
0242 
0243 %% basic optimset options needed for fmincon
0244 if isfield(opt, 'fmincon_opt')
0245     fmincon_opt = opt.fmincon_opt;
0246 else
0247     fmincon_opt = struct();
0248 end
0249 if isfield(fmincon_opt, 'opts')
0250     fopts = fmincon_opt.opts;
0251 else
0252     fopts = struct();
0253 end
0254 if isfield(fmincon_opt, 'tol_x')
0255     fopts.TolX = fmincon_opt.tol_x;
0256 end
0257 if isfield(fmincon_opt, 'tol_f')
0258     fopts.TolFun = fmincon_opt.tol_f;
0259 end
0260 if isfield(fmincon_opt, 'max_it') && fmincon_opt.max_it ~= 0
0261     fopts.MaxIter = fmincon_opt.max_it;
0262     fopts.MaxFunEvals = 4 * fmincon_opt.max_it;
0263 end
0264 fopts.GradObj = 'on';
0265 fopts.GradConstr = 'on';
0266 if verbose < 1
0267     fopts.Display = 'off';
0268 elseif verbose < 2
0269     fopts.Display = 'final';
0270 elseif verbose < 3
0271     fopts.Display = 'iter';
0272 else
0273     fopts.Display = 'testing';
0274 end
0275 
0276 %% select algorithm
0277 if have_feature('fmincon_ipm')
0278     %% set default algorithm, if not specified
0279     if ~isfield(fmincon_opt, 'alg') || fmincon_opt.alg == 0
0280         fmincon_opt.alg = 4;  %% interior-point, w/ exact user-supplied Hessian
0281     end
0282     if fmincon_opt.alg == 4 && isempty(hess_fcn)
0283         fmincon_opt.alg = 2;  %% interior-point, w/ default 'bfgs' Hessian approx
0284     end
0285     switch fmincon_opt.alg
0286         case 1      %% active-set (does not use sparse matrices, not suitable for large problems)
0287             fopts.Algorithm   = 'active-set';
0288             Af = full(Af);
0289             Afeq = full(Afeq);
0290         case 2      %% interior-point, w/ default 'bfgs' Hessian approx
0291             fopts.Algorithm   = 'interior-point';
0292         case 3      %% interior-point, w/ 'lbfgs' Hessian approx
0293             fopts.Algorithm   = 'interior-point';
0294             fopts.Hessian     = 'lbfgs';
0295         case 4      %% interior-point, w/ exact user-supplied Hessian
0296             fopts.Algorithm   = 'interior-point';
0297             fopts.Hessian     = 'user-supplied';
0298             fmc_hessian = @(x, lambda)hess_fcn(x, lambda, 1);
0299             fopts.HessFcn     = fmc_hessian;
0300         case 5      %% interior-point, w/ finite-diff Hessian
0301             fopts.Algorithm = 'interior-point';
0302             fopts.Hessian = 'fin-diff-grads';
0303             fopts.SubproblemAlgorithm = 'cg';
0304 %             fopts.SubProblem = 'cg';
0305 %             fopts.HessianApproximation = 'finite-difference';
0306         case 6      %% sqp (does not use sparse matrices, not suitable for large problems)
0307             fopts.Algorithm = 'sqp';
0308             Af = full(Af);
0309             Afeq = full(Afeq);
0310         otherwise
0311             error('nlps_fmincon: unknown algorithm specified in ''fmincon.alg'' option');
0312     end
0313 else
0314     fopts.LargeScale = 'off';
0315     Af = full(Af);
0316     Afeq = full(Afeq);
0317 end
0318 if have_feature('optimoptions')
0319     fmoptions = optimoptions('fmincon');
0320     fmoptions = nested_struct_copy(fmoptions, fopts);
0321 else
0322     fmoptions = optimset(fopts);
0323 end
0324 % fmoptions = optimset(fmoptions, 'DerivativeCheck', 'on', 'FinDiffType', 'central', 'FunValCheck', 'on');
0325 % fmoptions = optimset(fmoptions, 'Diagnostics', 'on');
0326 
0327 %%-----  run solver  -----
0328 [x, f, eflag, output, Lambda] = ...
0329   fmincon(f_fcn, x0, Af, bf, Afeq, bfeq, xmin, xmax, gh_fcn, fmoptions);
0330 success = (eflag > 0);
0331 
0332 %% fix Lambdas
0333 %% (shadow prices on equality variable bounds can come back on wrong limit)
0334 kl = find(Lambda.lower < 0 & Lambda.upper == 0);
0335 Lambda.upper(kl) = -Lambda.lower(kl);
0336 Lambda.lower(kl) = 0;
0337 ku = find(Lambda.upper < 0 & Lambda.lower == 0);
0338 Lambda.lower(ku) = -Lambda.upper(ku);
0339 Lambda.upper(ku) = 0;
0340 
0341 %% package up results
0342 nlt = length(ilt);
0343 ngt = length(igt);
0344 nbx = length(ibx);
0345 
0346 %% extract multipliers for linear constraints
0347 kl = find(Lambda.eqlin < 0);
0348 ku = find(Lambda.eqlin > 0);
0349 
0350 mu_l = zeros(size(u));
0351 mu_l(ieq(kl)) = -Lambda.eqlin(kl);
0352 mu_l(igt) = Lambda.ineqlin(nlt+(1:ngt));
0353 mu_l(ibx) = Lambda.ineqlin(nlt+ngt+nbx+(1:nbx));
0354 
0355 mu_u = zeros(size(u));
0356 mu_u(ieq(ku)) = Lambda.eqlin(ku);
0357 mu_u(ilt) = Lambda.ineqlin(1:nlt);
0358 mu_u(ibx) = Lambda.ineqlin(nlt+ngt+(1:nbx));
0359 
0360 lambda = struct( ...
0361     'lower', Lambda.lower, ...
0362     'upper', Lambda.upper, ...
0363     'eqnonlin', Lambda.eqnonlin, ...
0364     'ineqnonlin', Lambda.ineqnonlin, ...
0365     'mu_l', mu_l, ...
0366     'mu_u', mu_u );

Generated on Fri 09-Oct-2020 11:21:31 by m2html © 2005