The trust-region reflective method for `fmincon`

can
handle linear equality constraints if no other constraints exist.
Suppose you want to minimize

$$f(x)={\displaystyle \sum _{i=1}^{n-1}\left({\left({x}_{i}^{2}\right)}^{\left({x}_{i+1}^{2}+1\right)}+{\left({x}_{i+1}^{2}\right)}^{\left({x}_{i}^{2}+1\right)}\right)},$$

subject to some linear equality constraints. The objective function
is coded in the function `brownfgh.m`

. This example
takes `n = 1000`

. Furthermore, the `browneq.mat`

file
contains matrices `Aeq`

and `beq`

that
represent the linear constraints *Aeq*·*x* = *beq*. *Aeq* has
100 rows representing 100 linear constraints (so *Aeq* is
a 100-by-1000 matrix).

The file is lengthy so is not included here. View the code with the command

type brownfgh

Because `brownfgh`

computes the gradient and
Hessian values as well as the objective function, you need to use `optimoptions`

to indicate that this information
is available in `brownfgh`

, using the `SpecifyObjectiveGradient`

and `Hessian`

options.

The sparse matrix `Aeq`

and vector `beq`

are
available in the file `browneq.mat`

:

load browneq

The linear constraint system is 100-by-1000, has unstructured
sparsity (use `spy`

`(Aeq)`

to
view the sparsity structure), and is not too badly ill-conditioned:

condest(Aeq*Aeq') ans = 2.9310e+006

fun = @brownfgh; load browneq % Get Aeq and beq, the linear equalities n = 1000; xstart = -ones(n,1); xstart(2:2:n) = 1; options = optimoptions('fmincon','SpecifyObjectiveGradient',true,'HessianFcn','objective',... 'Algorithm','trust-region-reflective'); [x,fval,exitflag,output] = ... fmincon(fun,xstart,[],[],Aeq,beq,[],[],[],options);

`fmincon`

prints the following exit message:

Local minimum possible. fmincon stopped because the final change in function value relative to its initial value is less than the default value of the function tolerance.

The `exitflag`

value of `3`

also
indicates that the algorithm terminated because the change in the
objective function value was less than the tolerance `FunctionTolerance`

.
The final function value is given by `fval`

. Constraints
are satisfied, as you see in `output.constrviolation`

exitflag,fval,output.constrviolation exitflag = 3 fval = 205.9313 ans = 2.2071e-13

The linear equalities are satisfied at `x`

.

norm(Aeq*x-beq) ans = 1.1858e-12