原文来自http://hi.baidu.com/zzz700/blog/item/51f05207998fe1134bfb5126.html

http://www.mathworks.cn/help/optim/ug/fmincon.html

fmincon

Find minimum of constrained nonlinear multivariable function

http://www.mathworks.cn/help/optim/ug/_brhkghv-11.html#br9p_ry

Writing Constraints

...............

Equation

Finds the minimum of a problem specified by

x, b, beq, lb, and ub are vectors, A and Aeq are matrices, c(x) and ceq(x) are functions that return vectors, and f(x) is a function that returns a scalar. f(x), c(x), and ceq(x) can be nonlinear functions.

Syntax

x = fmincon(fun,x0,A,b) x = fmincon(fun,x0,A,b,Aeq,beq) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options) x = fmincon(problem) [x,fval] = fmincon(...) [x,fval,exitflag] = fmincon(...) [x,fval,exitflag,output] = fmincon(...) [x,fval,exitflag,output,lambda] = fmincon(...) [x,fval,exitflag,output,lambda,grad] = fmincon(...) [x,fval,exitflag,output,lambda,grad,hessian] = fmincon(...)

............

Writing Constraints

On this page…
Types of Constraints Iterations Can Violate Constraints Bound Constraints Linear Inequality Constraints Linear Equality Constraints Nonlinear Constraints Example: Using All Types of Constraints

On this page…

Types of Constraints

Iterations Can Violate Constraints

Bound Constraints

Linear Inequality Constraints

Linear Equality Constraints

Nonlinear Constraints

Example: Using All Types of Constraints

Types of Constraints

Optimization Toolbox solvers have special forms for constraints:

Bound Constraints — Lower and upper bounds on individual components: x ≥ l and x ≤ u.
Linear Inequality Constraints — A·x ≤ b. A is an m-by-n matrix, which represents m constraints for an n-dimensional vector x. b is m-dimensional.
Linear Equality Constraints — Aeq·x = beq. Equality constraints have the same form as inequality constraints.
Nonlinear Constraints — c(x) ≤ 0 and ceq(x) = 0. Both c and ceq are scalars or vectors representing several constraints.

Optimization Toolbox functions assume that inequality constraints are of the form c_i(x) ≤ 0 or A x ≤ b. Express greater-than constraints as less-than constraints by multiplying them by –1. For example, a constraint of the form c_i(x) ≥ 0 is equivalent to the constraint –c_i(x) ≤ 0. A constraint of the form A·x ≥ b is equivalent to the constraint –A·x ≤ –b. For more information, see Linear Inequality Constraints and Nonlinear Constraints.

You can sometimes write constraints in several ways. For best results, use the lowest numbered constraints possible:

Bounds
Linear equalities
Linear inequalities
Nonlinear equalities
Nonlinear inequalities

For example, with a constraint 5 x ≤ 20, use a bound x ≤ 4 instead of a linear inequality or nonlinear inequality.

For information on how to pass extra parameters to constraint functions, see Passing Extra Parameters.

Iterations Can Violate Constraints

Be careful when writing your objective and constraint functions. Intermediate iterations can lead to points that are infeasible (do not satisfy constraints). If you write objective or constraint functions that assume feasibility, these functions can error or give unexpected results.

For example, if you take a square root or logarithm of x, and x < 0, the result is not real. You can try to avoid this error by setting 0 as a lower bound on x. Nevertheless, an intermediate iteration can violate this bound.

Algorithms That Satisfy Bound Constraints

Some solver algorithms satisfy bound constraints at every iteration:

fmincon interior-point, sqp, and trust-region-reflective algorithms
lsqcurvefit trust-region-reflective algorithm
lsqnonlin trust-region-reflective algorithm
fminbnd

Solvers and Algorithms That Can Violate Bound Constraints

The following solvers and algorithms can violate bound constraints at intermediate iterations:

fmincon active-set algorithm
fgoalattain solver
fminimax solver
fseminf solver

Bound Constraints

Lower and upper bounds limit the components of the vector x. You need not give gradients for this type of constraint; solvers calculate them automatically. Bounds do not affect Hessians.

If you know bounds on the location of an optimum, you can obtain faster and more reliable solutions by explicitly including these bounds in your problem formulation.

Give bounds as vectors with the same length as x.

If a particular component has no lower bound, use –Inf as the bound; similarly, use Inf if a component has no upper bound.
If you have only bounds of one type (upper or lower), you do not need to write the other type. For example, if you have no upper bounds, you do not need to supply a vector of Infs. Also, if only the first m out of n components have bounds, then you need only supply a vector of length m containing bounds.

For example, suppose your bounds are:

x₃ ≥ 8
x₂ ≤ 3.

Write the constraint vectors as

l = [-Inf; -Inf; 8]
u = [Inf; 3] or u = [Inf; 3; Inf].

For a more complex example of bounds, see Example: Linear Programming.

Linear Inequality Constraints

Linear inequality constraints have the form A·x ≤ b. When A is m-by-n, there are m constraints on a variable x with ncomponents. You supply the m-by-n matrix A and the m-component vector b. You do not need to give gradients for this type of constraint; solvers calculate them automatically. Linear inequalities do not affect Hessians.

For example, suppose that you have the following linear inequalities as constraints:

x₁ + x₃ ≤ 4,
2x₂ – x₃ ≥ –2,
x₁ – x₂ + x₃ – x₄ ≥ 9.

Here m = 3 and n = 4.

Write these using the following matrix A and vector b:

Notice that the "greater than" inequalities were first multiplied by –1 in order to get them into "less than" inequality form. In MATLAB syntax:

A = [1 0 1 0;
    0 -2 1 0;
    -1 1 -1 1];
b = [4;2;-9];

For a more complex example of linear constraints, see Example: Linear Programming.

Linear Equality Constraints

Linear equalities have the form Aeq·x = beq, which represents m equations with n-component vector x. You supply the m-by-n matrix Aeq and the m-component vector beq. You do not need to give a gradient for this type of constraint; solvers calculate them automatically. Linear equalities do not affect Hessians. The form of this type of constraint is the same as for Linear Inequality Constraints.

Nonlinear Constraints

Nonlinear inequality constraints have the form c(x) ≤ 0, where c is a vector of constraints, one component for each constraint. Similarly, nonlinear equality constraints are of the form ceq(x) = 0.

Nonlinear constraint functions must return both inequality and equality constraints, even if they do not both exist. Return empty [] for a nonexistent constraint.

If you provide gradients for c and ceq, your solver can run faster and give more reliable results.

Providing a gradient has another advantage. A solver can reach a point x such that x is feasible, but finite differences around x always lead to an infeasible point. In this case, a solver can fail or halt prematurely. Providing a gradient allows a solver to proceed.

For example, suppose that you have the following inequalities as constraints:

Write these constraints in a function file as follows:

function [c,ceq]=ellipseparabola(x)
c(1) = (x(1)^2)/9 + (x(2)^2)/4 - 1;
c(2) = x(1)^2 - x(2) - 1;
ceq = [];
end

ellipseparabola returns empty [] as the nonlinear equality function. Also, both inequalities were put into ≤ 0 form.

Including Gradients in Constraint Functions

To include gradient information, write a conditionalized function as follows:

function [c,ceq,gradc,gradceq]=ellipseparabola(x)
c(1) = x(1)^2/9 + x(2)^2/4 - 1;
c(2) = x(1)^2 - x(2) - 1;
ceq = [];

if nargout > 2
    gradc = [2*x(1)/9, 2*x(1); ...
             x(2)/2, -1];
    gradceq = [];
end

See Writing Scalar Objective Functions for information on conditionalized functions. The gradient matrix has the form

gradc_{i, j} = [∂c(j)/∂x_i].

The first column of the gradient matrix is associated with c(1), and the second column is associated with c(2). This is the transpose of the form of Jacobians.

To have a solver use gradients of nonlinear constraints, indicate that they exist by using optimset:

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

Make sure to pass the options structure to your solver:

[x,fval] = fmincon(@myobj,x0,A,b,Aeq,beq,lb,ub, ...
           @ellipseparabola,options)

If you have a Symbolic Math Toolbox license, you can calculate gradients and Hessians automatically, as described inExample: Using Symbolic Math Toolbox Functions to Calculate Gradients and Hessians.

Anonymous Nonlinear Constraint Functions

For information on anonymous objective functions, see Anonymous Function Objectives.

Nonlinear constraint functions must return two outputs. The first output corresponds to nonlinear inequalities, and the second corresponds to nonlinear equalities.

Anonymous functions return just one output. So how can you write an anonymous function as a nonlinear constraint?

The deal function distributes multiple outputs. For example, suppose your nonlinear inequalities are

Suppose that your nonlinear equality is

x₂ = tanh(x₁).

Write a nonlinear constraint function as follows:

c = @(x)[x(1)^2/9 + x(2)^2/4 - 1;
        x(1)^2 - x(2) - 1];
ceq = @(x)tanh(x(1)) - x(2);
nonlinfcn = @(x)deal(c(x),ceq(x));

To minimize the function cosh(x₁) + sinh(x₂) subject to the constraints in nonlinfcn, use fmincon:

obj = @(x)cosh(x(1))+sinh(x(2));
opts = optimset('Algorithm','sqp');
z = fmincon(obj,[0;0],[],[],[],[],[],[],nonlinfcn,opts)

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is 
non-decreasing in feasible directions, to within the default 
value of the function tolerance, and constraints are satisfied 
to within the default value of the constraint tolerance.

z =
   -0.6530
   -0.5737

To check how well the resulting point z satisfies the constraints, use nonlinfcn:

[cout,ceqout] = nonlinfcn(z)

cout =
   -0.8704
         0

ceqout =
     1.1102e-016

z indeed satisfies all the constraints to within the default value of the TolCon constraint tolerance, 1e-6.

Example: Using All Types of Constraints

This section contains an example of a nonlinear minimization problem with all possible types of constraints. The objective function is in the subfunction myobj(x). The nonlinear constraints are in the subfunction myconstr(x). This example does not use gradients.

function [x fval exitflag] = fullexample
x0 = [1; 4; 5; 2; 5];
lb = [-Inf; -Inf;  0; -Inf;   1];
ub = [ Inf;  Inf; 20];
Aeq = [1 -0.3 0 0 0];
beq = 0;
A = [0 0  0 -1  0.1
     0 0  0  1 -0.5
     0 0 -1  0  0.9];
b = [0; 0; 0];
opts = optimset('Algorithm','sqp');
     
[x,fval,exitflag]=fmincon(@myobj,x0,A,b,Aeq,beq,lb,ub,...
                                  @myconstr,opts)

%---------------------------------------------------------
function f = myobj(x)

f = 6*x(2)*x(5) + 7*x(1)*x(3) + 3*x(2)^2;

%---------------------------------------------------------
function [c, ceq] = myconstr(x)

c = [x(1) - 0.2*x(2)*x(5) - 71
     0.9*x(3) - x(4)^2 - 67];
ceq = 3*x(2)^2*x(5) + 3*x(1)^2*x(3) - 20.875;

Calling fullexample produces the following display in the Command Window:

[x fval exitflag] = fullexample;

Local minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the default value of the function tolerance,
and constraints are satisfied to within the default value of the constraint tolerance.

x =
    0.6114
    2.0380
    1.3948
    0.1572
    1.5498

fval =
   37.3806

exitflag =
     1

.................

Examples

Find values of x that minimize f(x) = –x₁x₂x₃, starting at the point x = [10;10;10], subject to the constraints:

0 ≤ x₁ + 2x₂ + 2x₃ ≤ 72.

Write a file that returns a scalar value f of the objective function evaluated at x:
```
function f = myfun(x)
f = -x(1) * x(2) * x(3);
```
Rewrite the constraints as both less than or equal to a constant,

–x₁–2x₂–2x₃ ≤ 0
x₁ + 2x₂ + 2x₃≤ 72
Since both constraints are linear, formulate them as the matrix inequality A·x ≤ b, where
```
A = [-1 -2 -2; ...
      1  2  2];
b = [0;72];
```

Supply a starting point and invoke an optimization routine:

x0 = [10;10;10];    % Starting guess at the solution
[x,fval] = fmincon(@myfun,x0,A,b);

After 11 iterations, the solution is
```
x
x =
    24.0000
    12.0000
    12.0000
```
where the function value is
```
fval
fval =
    -3.4560e+03
```
and linear inequality constraints evaluate to be less than or equal to 0:
```
A*x-b
ans =
    -72
      0
```

Notes

Trust-Region-Reflective Optimization

To use the trust-region-reflective algorithm, you must

Supply the gradient of the objective function in fun.
Set GradObj to 'on' in options.
Specify the feasible region using one, but not both, of the following types of constraints:
- Upper and lower bounds constraints
- Linear equality constraints, in which the equality constraint matrix Aeq cannot have more rows than columns

You cannot use inequality constraints with the trust-region-reflective algorithm. If the preceding conditions are not met,fmincon reverts to the active-set algorithm.

fmincon returns a warning if you do not provide a gradient and the Algorithm option is trust-region-reflective.fmincon permits an approximate gradient to be supplied, but this option is not recommended; the numerical behavior of most optimization methods is considerably more robust when the true gradient is used.

The trust-region-reflective method in fmincon is in general most effective when the matrix of second derivatives, i.e., the Hessian matrix H(x), is also computed. However, evaluation of the true Hessian matrix is not required. For example, if you can supply the Hessian sparsity structure (using the HessPattern option in options), fmincon computes a sparse finite-difference approximation to H(x).

If x0 is not strictly feasible, fmincon chooses a new strictly feasible (centered) starting point.

If components of x have no upper (or lower) bounds, fmincon prefers that the corresponding components of ub (orlb) be set to Inf (or -Inf for lb) as opposed to an arbitrary but very large positive (or negative in the case of lower bounds) number.

Take note of these characteristics of linearly constrained minimization:

A dense (or fairly dense) column of matrix Aeq can result in considerable fill and computational cost.
fmincon removes (numerically) linearly dependent rows in Aeq; however, this process involves repeated matrix factorizations and therefore can be costly if there are many dependencies.
Each iteration involves a sparse least-squares solution with matrix

where R^T is the Cholesky factor of the preconditioner. Therefore, there is a potential conflict between choosing an effective preconditioner and minimizing fill in .

Active-Set Optimization

If equality constraints are present and dependent equalities are detected and removed in the quadratic subproblem,'dependent' appears under the Procedures heading (when you ask for output by setting the Display option to'iter'). The dependent equalities are only removed when the equalities are consistent. If the system of equalities is not consistent, the subproblem is infeasible and 'infeasible' appears under the Procedures heading.

Algorithms

Trust-Region-Reflective Optimization

The trust-region-reflective algorithm is a subspace trust-region method and is based on the interior-reflective Newton method described in [3] and [4]. Each iteration involves the approximate solution of a large linear system using the method of preconditioned conjugate gradients (PCG). See the trust-region and preconditioned conjugate gradient method descriptions in fmincon Trust Region Reflective Algorithm.

Active-Set Optimization

fmincon uses a sequential quadratic programming (SQP) method. In this method, the function solves a quadratic programming (QP) subproblem at each iteration. fmincon updates an estimate of the Hessian of the Lagrangian at each iteration using the BFGS formula (see fminunc and references [7] and [8]).

fmincon performs a line search using a merit function similar to that proposed by [6], [7], and [8]. The QP subproblem is solved using an active set strategy similar to that described in [5]. fmincon Active Set Algorithm describes this algorithm in detail.

See also SQP Implementation for more details on the algorithm used.

Interior-Point Optimization

This algorithm is described in fmincon Interior Point Algorithm. There is more extensive description in [1], [41], and [9].

SQP Optimization

The fmincon 'sqp' algorithm is similar to the 'active-set' algorithm described in Active-Set Optimization. fmincon SQP Algorithm describes the main differences. In summary, these differences are:

Limitations

fmincon is a gradient-based method that is designed to work on problems where the objective and constraint functions are both continuous and have continuous first derivatives.

When the problem is infeasible, fmincon attempts to minimize the maximum constraint value.

The trust-region-reflective algorithm does not allow equal upper and lower bounds. For example, if lb(2)==ub(2),fmincon gives this error:

Equal upper and lower bounds not permitted in this
large-scale method.
Use equality constraints and the medium-scale
method instead.

There are two different syntaxes for passing a Hessian, and there are two different syntaxes for passing a HessMultfunction; one for trust-region-reflective, and another for interior-point.

For trust-region-reflective, the Hessian of the Lagrangian is the same as the Hessian of the objective function. You pass that Hessian as the third output of the objective function.

For interior-point, the Hessian of the Lagrangian involves the Lagrange multipliers and the Hessians of the nonlinear constraint functions. You pass the Hessian as a separate function that takes into account both the position x and the Lagrange multiplier structure lambda.

Trust-Region-Reflective Coverage and Requirements

Additional Information Needed	For Large Problems
Must provide gradient for `f(x)` in `fun`.	Provide sparsity structure of the Hessian or compute the Hessian in `fun`. The Hessian should be sparse. Aeq should be sparse.

References

[1] Byrd, R.H., J. C. Gilbert, and J. Nocedal, "A Trust Region Method Based on Interior Point Techniques for Nonlinear Programming," Mathematical Programming, Vol 89, No. 1, pp. 149–185, 2000.

[2] Byrd, R.H., Mary E. Hribar, and Jorge Nocedal, "An Interior Point Algorithm for Large-Scale Nonlinear Programming, SIAM Journal on Optimization," SIAM Journal on Optimization, Vol 9, No. 4, pp. 877–900, 1999.

[3] Coleman, T.F. and Y. Li, "An Interior, Trust Region Approach for Nonlinear Minimization Subject to Bounds," SIAM Journal on Optimization, Vol. 6, pp. 418–445, 1996.

[4] Coleman, T.F. and Y. Li, "On the Convergence of Reflective Newton Methods for Large-Scale Nonlinear Minimization Subject to Bounds," Mathematical Programming, Vol. 67, Number 2, pp. 189–224, 1994.

[5] Gill, P.E., W. Murray, and M.H. Wright, Practical Optimization, London, Academic Press, 1981.

[6] Han, S.P., "A Globally Convergent Method for Nonlinear Programming," Vol. 22, Journal of Optimization Theory and Applications, p. 297, 1977.

[7] Powell, M.J.D., "A Fast Algorithm for Nonlinearly Constrained Optimization Calculations," Numerical Analysis, ed. G.A. Watson, Lecture Notes in Mathematics, Springer Verlag, Vol. 630, 1978.

[8] Powell, M.J.D., "The Convergence of Variable Metric Methods For Nonlinearly Constrained Optimization Calculations,"Nonlinear Programming 3 (O.L. Mangasarian, R.R. Meyer, and S.M. Robinson, eds.), Academic Press, 1978.

[9] Waltz, R. A., J. L. Morales, J. Nocedal, and D. Orban, "An interior algorithm for nonlinear optimization that combines line search and trust region steps," Mathematical Programming, Vol 107, No. 3, pp. 391–408, 2006.

posted on 2011-10-22 21:20 emy_yu 阅读(954) 评论(0) 收藏举报