Constrained Optimization using Matlab's fmincon

A. Basic Calls (without any special options)
Example1              Example 2
B. Calls with Gradients Supplied
Matlab's HELP DESCRIPTION

For constrained minimization of an objective function f(x) (for maximization use -f), Matlab provides the command fmincon. The objective function must be coded in a function file in the same manner as for fminunc. In these notes this file will be called objfun and saved as objfun.m  in the working directory.

A: Basic calls   top

Without any extra options, fmincon  is called as follows:

- with linear inequality constraints Ax£b only (as in linprog):
[x,fval]=fmincon('objfun',x0,A,b)

- with linear inequality constraints and linear equality constraints Aeq·x=beq only:
[x,fval]=fmincon('objfun',x0,A,b,Aeq,beq)

- with linear inequality and equality constraints, and in addition a lower bound of the form x³lb only:
[x,fval]=fmincon('objfun',x0,A,b,Aeq,beq,lb)
If only a subset of the variables has a lower bound, the components of lb corresponding to variables without lower bound are -Inf. For example, if the variables are (x,y), and x³1 but y has no lower bound, then lb=[1;-Inf].

- with linear inequality and equality constraints and lower as well as an upper bound of the form x£ub only:
[x,fval]=fmincon('objfun',x0,A,b,Aeq,beq,lb,ub)
If only a subset of the variables has an upper bound, the components of ub corresponding to variables without upper bound are Inf. For example, if the variables are (x,y) and x£1 but y has no lower bound, then lb=[1;Inf].

- with linear inequality and equality constraints,  lower and upper bounds, and nonlinear inequality and equality constraints:
[x,fval]=fmincon('objfun',x0,A,b,Aeq,beq,lb,ub,'constraint')
The last input argument in this call  is the name of a function file (denoted  constraint in these notes and saved as constraint.m  in the working directory), in which the nonlinear constraints are coded.

Constraint function file:
constraint.m is a function file (any name can be chosen) in which both the inequality functions  c(x) and the equality constraints ceq(x) are coded and provided in the form of column vectors. The function call

[c,ceq]=constraint(x)

must retrieve c(x) and ceq(x) for given input vector x. Examples of constraint function files are given in Examples 1 and 2 below. If only inequality constraints are given, define ceq=[]. Likewise, if only equality constraints are given, define c=[].

Interpretation:
The retrieved ceq(x) is interpreted by fmincon as equality constraint ceq(x)=0. The inequalities associated with c(x) are interpreted as c(x)£0. Thus, if a constraint of the form c(x)³0 is given, rewrite this as -c(x)£0 and code -c(x) in the constraint function file.

Placeholders:
As shown above, the constraints have to passed to fmincon in the following order:
1. Linear inequality constraints
2. Linear equality constraints
3. Lower bounds
4. Upper bounds
5. Nonlinear constraints
If a certain constraint is required, all other constraints appearing before it have to be inputted as well, even if they are not required in the problem. If this is the case, their input argument is replaced by the placeholder [] (empty input).

Examples:
- If lb and (A,b) are given, but there are no other constraints, the syntax is:
[x,fval]=fmincon('objfun',x0,A,b,[],[],lb)

- If ub and (Aeq,beq) are the only constraints:
[x,fval]=fmincon('objfun',x0,[],[],Aeq,beq,[],ub)

- If only nonlinear constraints are given:
[x,fval]=fmincon('objfun',x0,[],[],[],[],[],[],'constraint')
and function file constraint.m must be provided.
 

Example 1:    top

Find the minimum of

f(x,y)=x4-x2+y2-2x+y
subject to
 
linear inequalities
linear equalities
lower bounds
upper bounds
nonlinear constraints
(a)
--
--
x³0
y£0
--
(b)
--
x+y=0
--
x£1,  y£10
--
(c)
x+y£0
--
--
--
x2+y2£1
(d)
--
--
--
--
x2+y2=1
(e)
--
--
--
--
 x2+y2=1,  x2-y2³1
(f)
--
--
--
--
 x2+y2£1,  x2-y2³

Solution: The objective function is coded as for unconstrained minimization:

function f=objfun(x)

f=x(1)^4-x(1)^2+x(2)^2-2*x(1)+x(2);

For (a), (b) we don't need a constraint function file. The calls are (assuming x0=[value1;value2] is already defined):
(a):   [x,fval]=fmincon('objfun',x0,[],[],[],[],[0;-Inf],[Inf;0])
(b):   [x,fval]=fmincon('objfun',x0,[],[],[1,1],0,[],[1;10])

For (c)-(f)  we need a constraint function file. In each case the first line of the file constraint.m  is:

function [c,ceq]=constraint(x)

followed by an empty line. The commands below the 2nd line are:
 

(c)
(d)
(e)
(f)
c=x(1)^2+x(2)^2-1;
ceq=[];
c=[];
ceq=x(1)^2+x(2)^2-1;
c=1-x(1)^2+x(2)^2;
ceq=x(1)^2+x(2)^2-1;
c1=x(1)^2+x(2)^2-1;
c2=1-x(1)^2+x(2)^2;
c=[c1;c2];ceq=[];

For example, for (f) the full constraint function file is:

function [c,ceq]=constraint(x)

c1=x(1)^2+x(2)^2-1;
c2=1-x(1)^2+x(2)^2;
c=[c1;c2];ceq=[];

Function calls for (c)-(f):
(c):       [x,fval]=fmincon('objfun',x0,[1,1],0,[],[],[],[],'constraint')
(d)-(f):  [x,fval]=fmincon('objfun',x0,[],[],[],[],[],[],'constraint')

Approximate solutions found by fmincon:
 

x0
x
y
fval
(a)
[1;-1]
1.00000006131380
-0.50000014164875
-2.24999999999996
(b)
[1;-1]
0.90852417219345
-0.90852417219345
-2.04426066047301
(c)
[0;0]
0.70710678118746
-0.70710678118746
-1.87132034356109
(d)
[1;0]
0.92894844437517
-0.37020912075712
-2.20932198927909
(e)
[.5;.1]
1.00000000003278
-0.00000810106872
-2.00000810100310
(f)
[.1;.1]
1.00000000000009
0.00000001792512
-1.99999998207488

Example 2:    top

Minimize and maximize the objective function

f(x,y,z)=x3+y3+z3
subject to
x³0,      z£0,     x2+y2+z2=1,     y2³2z2.
Objective function file:
For Minimization:
function f=objfun(x)

f=x(1)^3+x(2)^3+x(3)^3;

For Maximization:
function f=objfun(x)

f=x(1)^3+x(2)^3+x(3)^3;f=-f;

Constraint function file:
function [c,ceq]=constraint(x)

c=2*x(3)^2-x(2)^2;
ceq=x(1)^2+x(2)^2+x(3)^2-1;

Function calls (in command window) and answers:
Minimization:
>> x0=[0;1;2];
>> [x,fval]=fmincon('objfun',x0,[],[],[],[],[0;-Inf;-Inf],[Inf;Inf;0],'constraint')

Warning: Large-scale (trust region) method does not currently solve this type of problem,
switching to medium-scale (line search).
> In C:\MATLABR12\toolbox\optim\fmincon.m at line 213
Optimization terminated successfully:
 Magnitude of directional derivative in search direction
  less than 2*options.TolFun and maximum constraint violation
  is less than options.TolCon
Active Constraints:
     1
     3
x =
   0.92898366078939
  -0.37012154351899
                  0
fval =
  -2.20932218190572

Answer for Maximization (same call, only objective function file was changed):
Warning: Large-scale (trust region) method does not currently solve this type of problem,
switching to medium-scale (line search).
> In C:\MATLABR12\toolbox\optim\fmincon.m at line 213
Optimization terminated successfully:
 Search direction less than 2*options.TolX and
  maximum constraint violation is less than options.TolCon
Active Constraints:
     1
x =
     0
     1
     0
fval =
    -1
 

B: Call of fmincon with gradient information provided    top

As for fminunc the performance of fmincon can be improved if gradient information is supplied. This information can be provided for the objective function, the nonlinear constraint functions, or both. Let's consider Example 1(f) again. The objective function file is extended as:

function [f,gradf]=objfun(x)

f=x(1)^4-x(1)^2+x(2)^2-2*x(1)+x(2);
gradf=[4*x(1)^3-2*x(1)-2;2*x(2)+1];

For providing the gradients of the nonlinear constraints, the constraint function file is extended as:

function [c,ceq,gradc,gradceq]=constraint(x)

c1=x(1)^2+x(2)^2-1;
c2=1-x(1)^2+x(2)^2;
c=[c1;c2];ceq=[];
gradc=[2*x(1),-2*x(1);2*x(2),2*x(2)];
gradceq=[];

Note that the the first column of gradc is the gradient-vector of the first constraint, and the second column of gradc is the gradient vector of the second constraint.

As in the unconstrained case we have to set the gradient option. We want to supply  the gradient of the objective function as well as the nonlinear constraints. The follwoing command sets this option:

>> options = optimset('GradObj','on','GradConstr','on');

In the function call these options are passed to fmincon as input argument after the name of the constraint file:

>> x0=[.1;.1];[x,fval]=fmincon('objfun',x0,[],[],[],[],[],[],'constraint',options)

Warning: Large-scale (trust region) method does not currently solve this type of problem,
switching to medium-scale (line search).
> In C:\MATLABR12\toolbox\optim\fmincon.m at line 213
Optimization terminated successfully:
 Search direction less than 2*options.TolX and
  maximum constraint violation is less than options.TolCon
Active Constraints:
     1
     2
x =
   1.00000000000000
  -0.00000171875724
fval =
  -2.00000171875428
 
 

Matlab's HELP DESCRIPTION   top

FMINCON Finds the constrained minimum of a function of several variables.
    FMINCON solves problems of the form:
        min F(X)  subject to:  A*X  <= B, Aeq*X  = Beq (linear constraints)
         X                       C(X) <= 0, Ceq(X) = 0   (nonlinear constraints)
                                 LB <= X <= UB

    X=FMINCON(FUN,X0,A,B) starts at X0 and finds a minimum X to the function
    FUN, subject to the linear inequalities A*X <= B. FUN accepts input X and
    returns a scalar function value F evaluated at X. X0 may be a scalar,
    vector, or matrix.

    X=FMINCON(FUN,X0,A,B,Aeq,Beq) minimizes FUN subject to the linear equalities
    Aeq*X = Beq as well as A*X <= B. (Set A=[] and B=[] if no inequalities exist.)

    X=FMINCON(FUN,X0,A,B,Aeq,Beq,LB,UB) defines a set of lower and upper
    bounds on the design variables, X, so that the solution is in
    the range LB <= X <= UB. Use empty matrices for LB and UB
    if no bounds exist. Set LB(i) = -Inf if X(i) is unbounded below;
    set UB(i) = Inf if X(i) is unbounded above.

    X=FMINCON(FUN,X0,A,B,Aeq,Beq,LB,UB,NONLCON) subjects the minimization to the
    constraints defined in NONLCON. The function NONLCON accepts X and returns
    the vectors C and Ceq, representing the nonlinear inequalities and equalities
    respectively. FMINCON minimizes FUN such that C(X)<=0 and Ceq(X)=0.
    (Set LB=[] and/or UB=[] if no bounds exist.)

    X=FMINCON(FUN,X0,A,B,Aeq,Beq,LB,UB,NONLCON,OPTIONS) minimizes with the
    default optimization parameters replaced by values in the structure OPTIONS,
    an argument created with the OPTIMSET function.  See OPTIMSET for details.  Used
    options are Display, TolX, TolFun, TolCon, DerivativeCheck, Diagnostics, GradObj,
    GradConstr, Hessian, MaxFunEvals, MaxIter, DiffMinChange and DiffMaxChange,
    LargeScale, MaxPCGIter, PrecondBandWidth, TolPCG, TypicalX, Hessian, HessMult,
    HessPattern. Use the GradObj option to specify that FUN also returns a second
    output argument G that is the partial derivatives of the function df/dX, at the
    point X. Use the Hessian option to specify that FUN also returns a third output
    argument H that is the 2nd partial derivatives of the function (the Hessian) at the
    point X.  The Hessian is only used by the large-scale method, not the
    line-search method. Use the GradConstr option to specify that NONLCON also
    returns third and fourth output arguments GC and GCeq, where GC is the partial
    derivatives of the constraint vector of inequalities C, and GCeq is the partial
    derivatives of the constraint vector of equalities Ceq. Use OPTIONS = [] as a
    place holder if no options are set.

    X=FMINCON(FUN,X0,A,B,Aeq,Beq,LB,UB,NONLCON,OPTIONS,P1,P2,...) passes the
    problem-dependent parameters P1,P2,... directly to the functions FUN
    and NONLCON: feval(FUN,X,P1,P2,...) and feval(NONLCON,X,P1,P2,...).  Pass
    empty matrices for A, B, Aeq, Beq, OPTIONS, LB, UB, and NONLCON to use the
    default values.

    [X,FVAL]=FMINCON(FUN,X0,...) returns the value of the objective
    function FUN at the solution X.

    [X,FVAL,EXITFLAG]=FMINCON(FUN,X0,...) returns a string EXITFLAG that
    describes the exit condition of FMINCON.
    If EXITFLAG is:
       > 0 then FMINCON converged to a solution X.
       0   then the maximum number of function evaluations was reached.
       < 0 then FMINCON did not converge to a solution.

    [X,FVAL,EXITFLAG,OUTPUT]=FMINCON(FUN,X0,...) returns a structure
    OUTPUT with the number of iterations taken in OUTPUT.iterations, the number
    of function evaluations in OUTPUT.funcCount, the algorithm used in
    OUTPUT.algorithm, the number of CG iterations (if used) in OUTPUT.cgiterations,
    and the first-order optimality (if used) in OUTPUT.firstorderopt.

    [X,FVAL,EXITFLAG,OUTPUT,LAMBDA]=FMINCON(FUN,X0,...) returns the Lagrange multipliers
    at the solution X: LAMBDA.lower for LB, LAMBDA.upper for UB, LAMBDA.ineqlin is
    for the linear inequalities, LAMBDA.eqlin is for the linear equalities,
    LAMBDA.ineqnonlin is for the nonlinear inequalities, and LAMBDA.eqnonlin
    is for the nonlinear equalities.

    [X,FVAL,EXITFLAG,OUTPUT,LAMBDA,GRAD]=FMINCON(FUN,X0,...) returns the value of
    the gradient of FUN at the solution X.

    [X,FVAL,EXITFLAG,OUTPUT,LAMBDA,GRAD,HESSIAN]=FMINCON(FUN,X0,...) returns the
    value of the HESSIAN of FUN at the solution X.

    Examples
      FUN can be specified using @:
         X = fmincon(@humps,...)
      In this case, F = humps(X) returns the scalar function value F of the HUMPS function
      evaluated at X.

      FUN can also be an inline object:
         X = fmincon(inline('3*sin(x(1))+exp(x(2))'),[1;1],[],[],[],[],[0 0])
      returns X = [0;0].

    See also OPTIMSET, FMINUNC, FMINBND, FMINSEARCH, @, INLINE.   top