hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_glopt_bnd_mcs_solve (e05jb)

Purpose

nag_glopt_bnd_mcs_solve (e05jb) is designed to find the global minimum or maximum of an arbitrary function, subject to simple bound-constraints using a multi-level coordinate search method. Derivatives are not required, but convergence is only guaranteed if the objective function is continuous in a neighbourhood of a global optimum. It is not intended for large problems.
The initialization function nag_glopt_bnd_mcs_init (e05ja) must have been called before calling nag_glopt_bnd_mcs_solve (e05jb).

Syntax

[bl, bu, list, numpts, initpt, x, obj, comm, user, ifail] = e05jb(objfun, ibound, bl, bu, list, numpts, initpt, monit, comm, 'n', n, 'iinit', iinit, 'sdlist', sdlist, 'user', user)
[bl, bu, list, numpts, initpt, x, obj, comm, user, ifail] = nag_glopt_bnd_mcs_solve(objfun, ibound, bl, bu, list, numpts, initpt, monit, comm, 'n', n, 'iinit', iinit, 'sdlist', sdlist, 'user', user)
nag_glopt_bnd_mcs_init (e05ja) must be called before calling nag_glopt_bnd_mcs_solve (e05jb), or any of the option-setting or option-getting functions nag_glopt_bnd_mcs_optset_string (e05jd), nag_glopt_bnd_mcs_optset_char (e05je), nag_glopt_bnd_mcs_optset_int (e05jf), nag_glopt_bnd_mcs_optset_real (e05jg), nag_glopt_bnd_mcs_option_check (e05jh), nag_glopt_bnd_mcs_optget_char (e05jj), nag_glopt_bnd_mcs_optget_int (e05jk) or nag_glopt_bnd_mcs_optget_real (e05jl).
You must not alter the number of non-fixed variables in your problem or the contents of the array comm between calls of the functions nag_glopt_bnd_mcs_init (e05ja), nag_glopt_bnd_mcs_solve (e05jb), nag_glopt_bnd_mcs_optset_string (e05jd), nag_glopt_bnd_mcs_optset_char (e05je), nag_glopt_bnd_mcs_optset_int (e05jf), nag_glopt_bnd_mcs_optset_real (e05jg), nag_glopt_bnd_mcs_option_check (e05jh), nag_glopt_bnd_mcs_optget_char (e05jj), nag_glopt_bnd_mcs_optget_int (e05jk) or nag_glopt_bnd_mcs_optget_real (e05jl).
Note: the interface to this routine has changed since earlier releases of the toolbox:
Mark 24: iinit optional
.

Description

nag_glopt_bnd_mcs_solve (e05jb) is designed to solve modestly sized global optimization problems having simple bound-constraints only; it finds the global optimum of a nonlinear function subject to a set of bound constraints on the variables. Without loss of generality, the problem is assumed to be stated in the following form:
minimize F(x)  subject to  xu  and  u,
xRn
minimize xRn F(x)   subject to   x u   and   u ,
where F(x)F(x) (the objective function) is a nonlinear scalar function (assumed to be continuous in a neighbourhood of a global minimum), and the bound vectors are elements of RnR-n, where RR- denotes the extended reals R{,}R{-,}. Relational operators between vectors are interpreted elementwise.
The optional parameter Maximize should be set if you wish to solve maximization, rather than minimization, problems.
If certain bounds are not present, the associated elements of  or uu can be set to special values that will be treated as - or + +. See the description of the optional parameter Infinite Bound Size. Phrases in this document containing terms like ‘unbounded values’ should be understood to be taken relative to this optional parameter.
Fixing variables (that is, setting li = uili=ui for some ii) is allowed in nag_glopt_bnd_mcs_solve (e05jb).
A typical excerpt from a function calling nag_glopt_bnd_mcs_solve (e05jb) is:
[comm, ifail] = e05ja(n_r);
[comm, ifail] = e05jd(optstr, comm);
[..., ifail] = e05jb(objfun, ...);
where nag_glopt_bnd_mcs_optset_string (e05jd) sets the optional parameter and value specified in optstr.
The initialization function nag_glopt_bnd_mcs_init (e05ja) does not need to be called before each invocation of nag_glopt_bnd_mcs_solve (e05jb). You should be aware that a call to the initialization function will reset each optional parameter to its default value, and, if you are using repeatable randomized initialization lists (see the description of the parameter iinit), the random state stored in the array comm will be destroyed.
You must supply a function that evaluates F(x)F(x); derivatives are not required.
The method used by nag_glopt_bnd_mcs_solve (e05jb) is based on MCS, the Multi-level Coordinate Search method described in Huyer and Neumaier (1999), and the algorithm it uses is described in detail in Section [Algorithmic Details].

References

Huyer W and Neumaier A (1999) Global optimization by multi-level coordinate search Journal of Global Optimization 14 331–355

Parameters

Compulsory Input Parameters

1:     objfun – function handle or string containing name of m-file
objfun must evaluate the objective function F(x)F(x) for a specified nn-vector xx.
[f, user, inform] = objfun(n, x, nstate, user)

Input Parameters

1:     n – int64int32nag_int scalar
nn, the number of variables.
2:     x(n) – double array
xx, the vector at which the objective function is to be evaluated.
3:     nstate – int64int32nag_int scalar
If nstate = 1nstate=1 then nag_glopt_bnd_mcs_solve (e05jb) is calling objfun for the first time. This parameter setting allows you to save computation time if certain data must be read or calculated only once.
4:     user – Any MATLAB object
objfun is called from nag_glopt_bnd_mcs_solve (e05jb) with the object supplied to nag_glopt_bnd_mcs_solve (e05jb).

Output Parameters

1:     f – double scalar
Must be set to the value of the objective function at xx, unless you have specified termination of the current problem using inform.
2:     user – Any MATLAB object
3:     inform – int64int32nag_int scalar
Must be set to a value describing the action to be taken by the solver on return from objfun. Specifically, if the value is negative the solution of the current problem will terminate immediately; otherwise, computations will continue.
2:     ibound – int64int32nag_int scalar
Indicates whether the facility for dealing with bounds of special forms is to be used. ibound must be set to one of the following values.
ibound = 0ibound=0
You will supply  and uu individually.
ibound = 1ibound=1
There are no bounds on xx.
ibound = 2ibound=2
There are semi-infinite bounds 0x0x.
ibound = 3ibound=3
There are constant bounds = 1=1 and u = u1u=u1.
Note that it only makes sense to fix any components of xx when ibound = 0ibound=0.
Constraint: ibound = 0ibound=0, 11, 22 or 33.
3:     bl(n) – double array
4:     bu(n) – double array
n, the dimension of the array, must satisfy the constraint n > 0n>0.
blbl is , the array of lower bounds. bubu is uu, the array of upper bounds.
If ibound = 0ibound=0, you must set bl(i)bli to ii and bu(i)bui to uiui, for i = 1,2,,ni=1,2,,n. If a particular xixi is to be unbounded below, the corresponding bl(i)bli should be set to infbnd-infbnd, where infbndinfbnd is the value of the optional parameter Infinite Bound Size. Similarly, if a particular xixi is to be unbounded above, the corresponding bu(i)bui should be set to infbndinfbnd.
If ibound = 1ibound=1 or 22, arrays bl and bu need not be set on input.
If ibound = 3ibound=3, you must set bl(1)bl1 to 11 and bu(1)bu1 to u1u1. The remaining elements of bl and bu will then be populated by these initial values.
Constraints:
  • if ibound = 0ibound=0, bl(i)bu(i)blibui, for i = 1,2,,ni=1,2,,n;
  • if ibound = 3ibound=3, bl(1) < bu(1)bl1<bu1.
5:     list(n,sdlist) – double array
n, the first dimension of the array, must satisfy the constraint n > 0n>0.
This parameter need not be set on entry if you wish to use one of the preset initialization methods (iinit3iinit3).
list is the ‘initialization list’: whenever a sub-box in the algorithm is split for the first time (either during the initialization procedure or later), for each non-fixed coordinate ii the split is done at the values list(i,1 : numpts(i))listi1:numptsi, as well as at some adaptively chosen intermediate points. The array sections list(i,1 : numpts(i)) listi1:numptsi , for i = 1,2,,ni=1,2,,n, must be in ascending order with each entry being distinct. In this context, ‘distinct’ should be taken to mean relative to the safe-range parameter (see nag_machine_real_safe (x02am)).
Constraint: if x(i)xi is not fixed, list(i,1 : numpts(i))listi1:numptsi is in ascending order with each entry being distinct, for i = 1,2,,ni=1,2,,nbl(i)list(i,j)bu(i)blilistijbui, for i = 1,2,,ni=1,2,,n and j = 1,2,,numpts(i)j=1,2,,numptsi.
6:     numpts(n) – int64int32nag_int array
n, the dimension of the array, must satisfy the constraint n > 0n>0.
This parameter need not be set on entry if you wish to use one of the preset initialization methods (iinit3iinit3).
numpts encodes the number of splitting points in each non-fixed dimension.
Constraints:
  • if x(i)xi is not fixed, numpts(i) sdlist numptsi sdlist ;
  • numpts(i) 3 numptsi 3 , for i = 1,2,,ni=1,2,,n.
7:     initpt(n) – int64int32nag_int array
n, the dimension of the array, must satisfy the constraint n > 0n>0.
This parameter need not be set on entry if you wish to use one of the preset initialization methods (iinit3iinit3).
You must designate a point stored in list that you wish nag_glopt_bnd_mcs_solve (e05jb) to consider as an ‘initial point’ for the purposes of the splitting procedure. Call this initial point x*x*. The coordinates of x*x* correspond to a set of indices JiJi, for i = 1,2,,ni=1,2,,n, such that xi * xi* is stored in list(i,Ji) listiJi , for i = 1,2,,ni=1,2,,n. You must set initpt(i) = Ji initpti = Ji , for i = 1,2,,ni=1,2,,n.
Constraint: if x(i)xi is not fixed, 1 initpt(i) sdlist 1 initpti sdlist , for i = 1,2,,ni=1,2,,n.
8:     monit – function handle or string containing name of m-file
monit may be used to monitor the optimization process. It is invoked upon every successful completion of the procedure in which a sub-box is considered for splitting. It will also be called just before nag_glopt_bnd_mcs_solve (e05jb) exits if that splitting procedure was not successful.
If no monitoring is required, monit may be string 'e05jbk'
[user, inform] = monit(n, ncall, xbest, icount, ninit, list, numpts, initpt, nbaskt, xbaskt, boxl, boxu, nstate, user)

Input Parameters

1:     n – int64int32nag_int scalar
nn, the number of variables.
2:     ncall – int64int32nag_int scalar
The cumulative number of calls to objfun.
3:     xbest(n) – double array
The current best point.
4:     icount(66) – int64int32nag_int array
An array of counters.
icount(1)icount1
nboxesnboxes, the current number of sub-boxes.
icount(2)icount2
nclocncloc, the cumulative number of calls to objfun made in local searches.
icount(3)icount3
nlocnloc, the cumulative number of points used as start points for local searches.
icount(4)icount4
nsweepnsweep, the cumulative number of sweeps through levels.
icount(5)icount5
mm, the cumulative number of splits by initialization list.
icount(6)icount6
ss, the current lowest level containing non-split boxes.
5:     ninit – int64int32nag_int scalar
The maximum over ii of the number of points in coordinate ii at which to split according to the initialization list list. See also the description of the parameter numpts.
6:     list(n,ninit) – double array
The initialization list.
7:     numpts(n) – int64int32nag_int array
The number of points in each coordinate at which to split according to the initialization list list.
8:     initpt(n) – int64int32nag_int array
A pointer to the ‘initial point’ in list. Element initpt(i)initpti is the column index in list of the iith coordinate of the initial point.
9:     nbaskt – int64int32nag_int scalar
The number of points in the ‘shopping basket’ xbaskt.
10:   xbaskt(n,nbaskt) – double array
Note: the jjth candidate minimum has its iith coordinate stored in xbaskt(j,i)xbasktji, for i = 1,2,,ni=1,2,,n and j = 1,2,,nbasktj=1,2,,nbaskt.
The ‘shopping basket’ of candidate minima.
11:   boxl(n) – double array
The array of lower bounds of the current search box.
12:   boxu(n) – double array
The array of upper bounds of the current search box.
13:   nstate – int64int32nag_int scalar
Is set by nag_glopt_bnd_mcs_solve (e05jb) to indicate at what stage of the minimization monit was called.
nstate = 1nstate=1
This is the first time that monit has been called.
nstate = -1nstate=-1
This is the last time monit will be called.
nstate = 0nstate=0
This is the first and last time monit will be called.
14:   user – Any MATLAB object
monit is called from nag_glopt_bnd_mcs_solve (e05jb) with the object supplied to nag_glopt_bnd_mcs_solve (e05jb).

Output Parameters

1:     user – Any MATLAB object
2:     inform – int64int32nag_int scalar
Must be set to a value describing the action to be taken by the solver on return from monit. Specifically, if the value is negative the solution of the current problem will terminate immediately; otherwise, computations will continue.
9:     comm(lcomm) – double array

Optional Input Parameters

1:     n – int64int32nag_int scalar
Default: The dimension of the arrays bl, bu, numpts, initpt and the first dimension of the array list. (An error is raised if these dimensions are not equal.)
nn, the number of variables.
Constraint: n > 0n>0.
2:     iinit – int64int32nag_int scalar
Selects which initialization method to use.
iinit = 0iinit=0
Simple initialization (boundary and midpoint), with
numpts(i) = 3 numptsi=3 , initpt(i) = 2 initpti=2  and
list(i,j) = (bl(i), (bl(i) + bu(i)) / 2 ,bu(i)) listij = (bli, ( bli + bui ) / 2 ,bui) ,
for i = 1,2,,n i=1,2,,n  and j = 1,2,3j=1,2,3.
iinit = 1iinit=1
Simple initialization (off-boundary and midpoint), with
numpts(i) = 3 numptsi=3 , initpt(i) = 2 initpti=2  and
list(i,j) = ( (5bl(i) + bu(i)) / 6 , (bl(i) + bu(i)) / 2 , (bl(i) + 5bu(i)) / 6 ) listij = ( ( 5bli+ bui) / 6 , ( bli + bui ) / 2 , ( bli + 5bui ) / 6 ) ,
for i = 1,2,,n i=1,2,,n  and j = 1,2,3j=1,2,3.
iinit = 2iinit=2
Initialization using linesearches.
iinit = 3iinit=3
You are providing your own initialization list.
iinit = 4iinit=4
Generate a random initialization list.
For more information on methods iinit = 2iinit=2, 33 or 44 see Section [Initialization and Sweeps].
If ‘infinite’ values (as determined by the value of the optional parameter Infinite Bound Size) are detected by nag_glopt_bnd_mcs_solve (e05jb) when you are using a simple initialization method (iinit = 0iinit=0 or 11), a safeguarded initialization procedure will be attempted, to avoid overflow.
Default: iinit = 0iinit=0
Constraint: iinit = 0iinit=0, 11, 22, 33 or 44.
3:     sdlist – int64int32nag_int scalar
Default: The second dimension of the array list.
The second dimension of the array list as declared in the (sub)program from which nag_glopt_bnd_mcs_solve (e05jb) is called. sdlist is, at least, the maximum over ii of the number of points in coordinate ii at which to split according to the initialization list list; that is, sdlistmaxi numpts(i)sdlistmaxinumptsi.
Internally, nag_glopt_bnd_mcs_solve (e05jb) uses list to determine sets of points along each coordinate direction to which it fits quadratic interpolants. Since fitting a quadratic requires at least three distinct points, this puts a lower bound on sdlist. Furthermore, in the case of initialization by linesearches (iinit = 2iinit=2) internal storage considerations require that sdlist be at least 192192, but not all of this space may be used.
Constraints:
4:     user – Any MATLAB object
user is not used by nag_glopt_bnd_mcs_solve (e05jb), but is passed to objfun and monit. Note that for large objects it may be more efficient to use a global variable which is accessible from the m-files than to use user.

Input Parameters Omitted from the MATLAB Interface

lcomm iuser ruser

Output Parameters

1:     bl(n) – double array
2:     bu(n) – double array
Unless ifail = 1ifail=1 or 22 on exit, bl and bu are the actual arrays of bounds used by nag_glopt_bnd_mcs_solve (e05jb).
3:     list(n,sdlist) – double array
Unless ifail = 1ifail=1, 22 or 999-999 on exit, the actual initialization data used by nag_glopt_bnd_mcs_solve (e05jb). If you wish to monitor the contents of list you are advised to do so solely through monit, not through the output value here.
4:     numpts(n) – int64int32nag_int array
Unless ifail = 1ifail=1, 22 or 999-999 on exit, the actual initialization data used by nag_glopt_bnd_mcs_solve (e05jb).
5:     initpt(n) – int64int32nag_int array
Unless ifail = 1ifail=1, 22 or 999-999 on exit, the actual initialization data used by nag_glopt_bnd_mcs_solve (e05jb).
6:     x(n) – double array
If ifail = 0ifail=0, contains an estimate of the global optimum (see also Section [Accuracy]).
7:     obj – double scalar
If ifail = 0ifail=0, contains the function value at x.
If you request early termination of nag_glopt_bnd_mcs_solve (e05jb) using inform in objfun or the analogous inform in monit, there is no guarantee that the function value at x equals obj.
8:     comm(lcomm) – double array
lcomm = 100lcomm=100.
9:     user – Any MATLAB object
10:   ifail – int64int32nag_int scalar
ifail = 0ifail=0 unless the function detects an error (see [Error Indicators and Warnings]).
nag_glopt_bnd_mcs_solve (e05jb) returns with ifail = 0ifail=0 if your termination criterion has been met: either a target value has been found to the required relative error (as determined by the values of the optional parameters Target Objective Value, Target Objective Error and Target Objective Safeguard), or the best function value was static for the number of sweeps through levels given by the optional parameter Static Limit. The latter criterion is the default.

Error Indicators and Warnings

Errors or warnings detected by the function:

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

  ifail = 1ifail=1
Constraint: lcomm100lcomm100.
Initialization function nag_glopt_bnd_mcs_init (e05ja) has not been called.
  ifail = 2ifail=2
A value of Splits Limit (smaxsmax) smaller than nr + 3nr+3 was set.
Constraint: ibound = 0ibound=0, 11, 22 or 33.
Constraint: if ibound = 0ibound=0 then bl(i) bu(i) bli bui , for i = 1,2,,ni=1,2,,n; if ibound = 3ibound=3 then bl(1) < bu(1)bl1<bu1.
Constraint: if ibound = 3ibound=3 then bl(1) < bu(1)bl1<bu1.
Constraint: if iinit = 2iinit=2 then sdlist192sdlist192.
Constraint: if iinit2iinit2 then sdlist3sdlist3.
Constraint: if x(i)xi is not fixed then initpt(i) sdlist initpti sdlist , for i = 1,2,,ni=1,2,,n.
Constraint: if x(i)xi is not fixed then initpt(i)1initpti1, for i = 1,2,,ni=1,2,,n.
Constraint: if x(i)xi is not fixed then list(i,j) bu(i) listij bui , for i = 1,2,,ni=1,2,,n and j = 1,2,,numpts(i)j=1,2,,numptsi.
Constraint: if x(i)xi is not fixed then list(i,j) bl(i) listij bli , for i = 1,2,,ni=1,2,,n and j = 1,2,,numpts(i)j=1,2,,numptsi.
Constraint: if x(i)xi is not fixed then numpts(i)sdlistnumptsisdlist, for i = 1,2,,ni=1,2,,n.
Constraint: if x(i)xi is not fixed then numpts(i)3numptsi3, for i = 1,2,,ni=1,2,,n.
Constraint: iinit = 0iinit=0, 11, 22, 33 or 44.
Constraint: n > 0n>0.
On entry, user-supplied section list(i,1 : numpts(i))listi1:numptsi contained ndistndist distinct elements, and ndist < numpts(i)ndist<numptsi.
On entry, user-supplied section list(i,1 : numpts(i))listi1:numptsi was not in ascending order.
The number of non-fixed variables nr = 0nr=0.
Constraint: nr > 0nr>0.
  ifail = 3ifail=3
A finite initialization list could not be computed internally. Consider reformulating the bounds on the problem, try providing your own initialization list, use the randomization option (iinit = 4iinit=4) or vary the value of Infinite Bound Size.
The user-supplied initialization list contained infinite values, as determined by the optional parameter Infinite Bound Size.
  ifail = 4ifail=4
The division procedure completed but your target value could not be reached.
Despite every sub-box being processed Splits Limit times, the target value you provided in Target Objective Value could not be found to the tolerances given in Target Objective Error and Target Objective Safeguard. You could try reducing Splits Limit or the objective tolerances.
  ifail = 5ifail=5
The function evaluations limit was exceeded.
Approximately Function Evaluations Limit function calls have been made without your chosen termination criterion being satisfied.
W ifail = 6ifail=6
User-supplied monitoring function requested termination.
User-supplied objective function requested termination.
  ifail = 7ifail=7
An error occurred during initialization. It is likely that points from the initialization list are very close together. Try relaxing the bounds on the variables or use a different initialization method.
An error occurred during linesearching.
  ifail = 999ifail=-999
Dynamic memory allocation failed.

Accuracy

If ifail = 0ifail=0 on exit, then the vector returned in the array x is an estimate of the solution xx whose function value satisfies your termination criterion: the function value was static for Static Limit sweeps through levels, or
F(x) objval max ( objerr × |objval| ,objsfg) ,
F(x) - objval max( objerr × |objval| ,objsfg) ,
where objvalobjval is the value of the optional parameter Target Objective Value, objerrobjerr is the value of the optional parameter Target Objective Error, and objsfgobjsfg is the value of the optional parameter Target Objective Safeguard.

Further Comments

For each invocation of nag_glopt_bnd_mcs_solve (e05jb), local workspace arrays of fixed length are allocated internally. The total size of these arrays amounts to 13nr + smax113nr+smax-1 integer elements, where smaxsmax is the value of the optional parameter Splits Limit and nrnr is the number of non-fixed variables, and (2 + nr)sdlist + 2n + 21nr + 3nr2 + 1(2+nr)sdlist+2n+21nr+3nr2+1 double elements. In addition, if you are using randomized initialization lists (see the description of the parameter iinit), a further 2121 integer elements are allocated internally.
In order to keep track of the regions of the search space that have been visited while looking for a global optimum, nag_glopt_bnd_mcs_solve (e05jb) internally allocates arrays of increasing sizes depending on the difficulty of the problem. Two of the main factors that govern the amount allocated are the number of sub-boxes (call this quantity nboxesnboxes) and the number of points in the ‘shopping basket’ (the parameter nbaskt on entry to monit). Safe, pessimistic upper bounds on these two quantities are so large as to be impractical. In fact, the worst-case number of sub-boxes for even the most simple initialization list (when ninit = 3ninit=3 on entry to monit) grows like nrnr nr nr. Thus nag_glopt_bnd_mcs_solve (e05jb) does not attempt to estimate in advance the final values of nboxesnboxes or nbaskt for a given problem. There are a total of 55 integer arrays and 4 + nr + ninit4+nr+ninit double arrays whose lengths depend on nboxesnboxes, and there are a total of 22 integer arrays and 3 + n + nr3+n+nr double arrays whose lengths depend on nbaskt. nag_glopt_bnd_mcs_solve (e05jb) makes a fixed initial guess that the maximum number of sub-boxes required will be 1000010000 and that the maximum number of points in the ‘shopping basket’ will be 10001000. If ever a greater amount of sub-boxes or more room in the ‘shopping basket’ is required, nag_glopt_bnd_mcs_solve (e05jb) performs reallocation, usually doubling the size of the inadequately-sized arrays. Clearly this process requires periods where the original array and its extension exist in memory simultaneously, so that the data within can be copied, which compounds the complexity of nag_glopt_bnd_mcs_solve (e05jb)'s memory usage. It is possible (although not likely) that if your problem is particularly difficult to solve, or of a large size (hundreds of variables), you may run out of memory.
One array that could be dynamically resized by nag_glopt_bnd_mcs_solve (e05jb) is the ‘shopping basket’ (xbaskt on entry to monit). If the initial attempt to allocate 1000nr1000nr doubles for this array fails, monit will not be called on exit from nag_glopt_bnd_mcs_solve (e05jb).
nag_glopt_bnd_mcs_solve (e05jb) performs better if your problem is well-scaled. It is worth trying (by guesswork perhaps) to rescale the problem if necessary, as sensible scaling will reduce the difficulty of the optimization problem, so that nag_glopt_bnd_mcs_solve (e05jb) will take less computer time.

Example

function nag_glopt_bnd_mcs_solve_example
% Problem data for peaks function
prob = 'peaks';
xres = 100;
yres = 100;

bl = [-3; -3];
bu = -bl;
fglob = -6.55; % Approx.
xglob = [0.23; -1.63]; % Approx.



% Initialize nag_glopt_bnd_mcs_solve
n = length(bl);
[comm, ifail] = nag_glopt_bnd_mcs_init;

if (ifail == 0)

    % Vanilla call.
    disp(sprintf('\n'));
    disp('Solve with no options or init.-list data');

    ibound = int64(0);             % All bounds will be given;
    list = zeros(n,3);             % Only need to _declare_ the init.-list
    numpts = zeros(n, 1, 'int64'); % data: these will be _set_ internally.
    initpt = zeros(n, 1, 'int64');

    [bl, bu, listOut, numptsOut, initptOut, ...
     xbest, obj, comm, userOut, ifail] = ...
        nag_glopt_bnd_mcs_solve(@objective, ibound, bl, bu, ...
                           list, numpts, initpt, @monitor, comm);

    disp(['nag_glopt_bnd_mcs_solve (no options) exited with ifail = ' num2str(ifail)]);

    if (ifail == 0)
        disp('xbest:');
        disp(xbest);
        disp(['obj = ' num2str(obj)]);
    end

    % Set some options. No need to reinitialize: n hasn't changed, and we
    % didn't set any options above.
    disp(sprintf('\n'));
    disp('Solve with options and init.-list data');

    % Echo the setting of opt. params.
    comm = nag_glopt_bnd_mcs_optset_string('List', comm);

    comm = nag_glopt_bnd_mcs_optset_string('Function Evaluations Limit = 100000', comm);
    comm = nag_glopt_bnd_mcs_optset_int('Static Limit', int64(3*n), comm);

    % Increase infbnd by factor of 10.
    infbnd = nag_glopt_bnd_mcs_optget_real('Infinite Bound Size', comm);
    comm = nag_glopt_bnd_mcs_optset_real('Infinite Bound Size', 10*infbnd, comm);

    comm = nag_glopt_bnd_mcs_optset_char('Local Searches', 'on', comm);

    % Set the initialization-list data.
    iinit = int64(3); % We're providing the data this time:
    list = zeros(n, 3);
    list(:, 1) = bl; list(:, 3) = bu;
    list(:, 2) = [-1; 0];
    numpts = repmat(int64(3), n, 1); % 3 splitting points for each dim;
    initpt = repmat(int64(2), n, 1); % 2nd pt in each row to be the 'init.' pt.

    [bl, bu, listOut, numptsOut, initptOut, ...
     xbest, obj, comm, userOut, ifail] = ...
        nag_glopt_bnd_mcs_solve(@objective, ibound, bl, bu, list, ...
              numpts, initpt, @monitor, comm, 'iinit', iinit);

    disp(['nag_glopt_bnd_mcs_solve (options) exited with ifail = ' num2str(ifail)]);

    if (ifail == 0)
        disp('xbest:');
        disp(xbest);
        disp(['obj = ' num2str(obj)]);
    end

end

function [f,user,inform] = objective(n,x,nstate,user)

if (n==2)
    inform = int64(0);
else
    inform = int64(-1);
end

if (inform >= 0)

    % We're prepared to evaluate the objective at this point

    if (nstate == 1)
        disp(sprintf('\n'));
        disp('(OBJFUN was just called for the first time)');
    end

    f = peaks(x(1), x(2));
end
function [user,inform] = ...
    monitor(n,ncall,xbest,icount,ninit,list,numpts,initpt,nbaskt,xbaskt,...
            boxl,boxu,nstate,user)

inform = int64(0);

if (nstate == 0 || nstate == 1)
    disp(sprintf('\n'))
    disp('*** Begin monitoring information ***');
    disp(sprintf('\n'))
end

if (nstate <= 0)
    disp(['Total sub-boxes = ' num2str(icount(1))]);
    disp(['Total function evaluations = ' num2str(ncall)]);
    disp(['Total function evaluations used in local searches = ' num2str(icount(2))]);
    disp(['Total points used in local search = ' num2str(icount(3))]);
    disp(['Total sweeps through levels = ' num2str(icount(4))]);
    disp(['Total splits by init. list = ' num2str(icount(5))]);
    disp(['Lowest level with nonsplit boxes = ' num2str(icount(6))]);
    disp(['Number of candidate minima in the ''shopping basket'' = ' num2str(nbaskt)]);
    disp('Shopping basket:');
    disp(xbaskt);
    disp(sprintf('\n'))
    disp('*** End monitoring information ***');
    disp(sprintf('\n'))
end
 


Solve with no options or init.-list data


(OBJFUN was just called for the first time)


*** Begin monitoring information ***


Total sub-boxes = 228
Total function evaluations = 197
Total function evaluations used in local searches = 88
Total points used in local search = 13
Total sweeps through levels = 12
Total splits by init. list = 5
Lowest level with nonsplit boxes = 7
Number of candidate minima in the 'shopping basket' = 2
Shopping basket:
   -1.3474    0.2283
    0.2045   -1.6255



*** End monitoring information ***


nag_glopt_bnd_mcs_solve (no options) exited with ifail = 0
xbest:
    0.2283
   -1.6255

obj = -6.5511


Solve with options and init.-list data
      FUNCTION EVALUATIONS LIMIT =           100000
      STATIC LIMIT =                6
      INFINITE BOUND SIZE =     1.1579208923731620E+78
      LOCAL SEARCHES = on


(OBJFUN was just called for the first time)


*** Begin monitoring information ***


Total sub-boxes = 146
Total function evaluations = 170
Total function evaluations used in local searches = 103
Total points used in local search = 7
Total sweeps through levels = 7
Total splits by init. list = 5
Lowest level with nonsplit boxes = 4
Number of candidate minima in the 'shopping basket' = 2
Shopping basket:
    0.2283   -1.3474
   -1.6255    0.2045



*** End monitoring information ***


nag_glopt_bnd_mcs_solve (options) exited with ifail = 0
xbest:
    0.2283
   -1.6255

obj = -6.5511

function e05jb_example
% Problem data for peaks function
prob = 'peaks';
xres = 100;
yres = 100;

bl = [-3; -3];
bu = -bl;
fglob = -6.55; % Approx.
xglob = [0.23; -1.63]; % Approx.



% Initialize e05jb
n = length(bl);
[comm, ifail] = e05ja;

if (ifail == 0)

    % Vanilla call.
    disp(sprintf('\n'));
    disp('Solve with no options or init.-list data');

    ibound = int64(0);             % All bounds will be given;
    list = zeros(n,3);             % Only need to _declare_ the init.-list
    numpts = zeros(n, 1, 'int64'); % data: these will be _set_ internally.
    initpt = zeros(n, 1, 'int64');

    [bl, bu, listOut, numptsOut, initptOut, ...
     xbest, obj, comm, userOut, ifail] = ...
        e05jb(@objective, ibound, bl, bu, ...
                           list, numpts, initpt, @monitor, comm);

    disp(['e05jb (no options) exited with ifail = ' num2str(ifail)]);

    if (ifail == 0)
        disp('xbest:');
        disp(xbest);
        disp(['obj = ' num2str(obj)]);
    end

    % Set some options. No need to reinitialize: n hasn't changed, and we
    % didn't set any options above.
    disp(sprintf('\n'));
    disp('Solve with options and init.-list data');

    % Echo the setting of opt. params.
    comm = e05jd('List', comm);

    comm = e05jd('Function Evaluations Limit = 100000', comm);
    comm = e05jf('Static Limit', int64(3*n), comm);

    % Increase infbnd by factor of 10.
    infbnd = e05jl('Infinite Bound Size', comm);
    comm = e05jg('Infinite Bound Size', 10*infbnd, comm);

    comm = e05je('Local Searches', 'on', comm);

    % Set the initialization-list data.
    iinit = int64(3); % We're providing the data this time:
    list = zeros(n, 3);
    list(:, 1) = bl; list(:, 3) = bu;
    list(:, 2) = [-1; 0];
    numpts = repmat(int64(3), n, 1); % 3 splitting points for each dim;
    initpt = repmat(int64(2), n, 1); % 2nd pt in each row to be the 'init.' pt.

    [bl, bu, listOut, numptsOut, initptOut, ...
     xbest, obj, comm, userOut, ifail] = ...
        e05jb(@objective, ibound, bl, bu, list, ...
              numpts, initpt, @monitor, comm, 'iinit', iinit);

    disp(['e05jb (options) exited with ifail = ' num2str(ifail)]);

    if (ifail == 0)
        disp('xbest:');
        disp(xbest);
        disp(['obj = ' num2str(obj)]);
    end

end

function [f,user,inform] = objective(n,x,nstate,user)

if (n==2)
    inform = int64(0);
else
    inform = int64(-1);
end

if (inform >= 0)

    % We're prepared to evaluate the objective at this point

    if (nstate == 1)
        disp(sprintf('\n'));
        disp('(OBJFUN was just called for the first time)');
    end

    f = peaks(x(1), x(2));
end
function [user,inform] = ...
    monitor(n,ncall,xbest,icount,ninit,list,numpts,initpt,nbaskt,xbaskt,...
            boxl,boxu,nstate,user)

inform = int64(0);

if (nstate == 0 || nstate == 1)
    disp(sprintf('\n'))
    disp('*** Begin monitoring information ***');
    disp(sprintf('\n'))
end

if (nstate <= 0)
    disp(['Total sub-boxes = ' num2str(icount(1))]);
    disp(['Total function evaluations = ' num2str(ncall)]);
    disp(['Total function evaluations used in local searches = ' num2str(icount(2))]);
    disp(['Total points used in local search = ' num2str(icount(3))]);
    disp(['Total sweeps through levels = ' num2str(icount(4))]);
    disp(['Total splits by init. list = ' num2str(icount(5))]);
    disp(['Lowest level with nonsplit boxes = ' num2str(icount(6))]);
    disp(['Number of candidate minima in the ''shopping basket'' = ' num2str(nbaskt)]);
    disp('Shopping basket:');
    disp(xbaskt);
    disp(sprintf('\n'))
    disp('*** End monitoring information ***');
    disp(sprintf('\n'))
end
 


Solve with no options or init.-list data


(OBJFUN was just called for the first time)


*** Begin monitoring information ***


Total sub-boxes = 228
Total function evaluations = 197
Total function evaluations used in local searches = 88
Total points used in local search = 13
Total sweeps through levels = 12
Total splits by init. list = 5
Lowest level with nonsplit boxes = 7
Number of candidate minima in the 'shopping basket' = 2
Shopping basket:
   -1.3474    0.2283
    0.2045   -1.6255



*** End monitoring information ***


e05jb (no options) exited with ifail = 0
xbest:
    0.2283
   -1.6255

obj = -6.5511


Solve with options and init.-list data
      FUNCTION EVALUATIONS LIMIT =           100000
      STATIC LIMIT =                6
      INFINITE BOUND SIZE =     1.1579208923731620E+78
      LOCAL SEARCHES = on


(OBJFUN was just called for the first time)


*** Begin monitoring information ***


Total sub-boxes = 146
Total function evaluations = 170
Total function evaluations used in local searches = 103
Total points used in local search = 7
Total sweeps through levels = 7
Total splits by init. list = 5
Lowest level with nonsplit boxes = 4
Number of candidate minima in the 'shopping basket' = 2
Shopping basket:
    0.2283   -1.3474
   -1.6255    0.2045



*** End monitoring information ***


e05jb (options) exited with ifail = 0
xbest:
    0.2283
   -1.6255

obj = -6.5511

the remainder of this document is intended for more advanced users. Section [Algorithmic Details] contains a detailed description of the algorithm. This information may be needed in order to understand Section [Optional Parameters], which describes the optional parameters that can be set by calls to nag_glopt_bnd_mcs_optset_string (e05jd), nag_glopt_bnd_mcs_optset_char (e05je), nag_glopt_bnd_mcs_optset_int (e05jf) and/or nag_glopt_bnd_mcs_optset_real (e05jg).

Algorithmic Details

Here we summarise the main features of the MCS algorithm used in nag_glopt_bnd_mcs_solve (e05jb), and we introduce some terminology used in the description of the function and its arguments. We assume throughout that we will only do any work in coordinates ii in which xixi is free to vary. The MCS algorithm is fully described in Huyer and Neumaier (1999).

Initialization and Sweeps

Each sub-box is determined by a basepoint xx and an opposite point yy. We denote such a sub-box by B[x,y]B[x,y]. The basepoint is allowed to belong to more than one sub-box, is usually a boundary point, and is often a vertex.
An initialization procedure produces an initial set of sub-boxes. Whenever a sub-box is split along a coordinate ii for the first time (in the initialization procedure or later), the splitting is done at three or more user-defined values {xij}j{xij}j at which the objective function is sampled, and at some adaptively chosen intermediate points. At least four children are generated. More precisely, we assume that we are given
i xi1 < xi2 < < xiLi ui ,  Li 3 ,  for ​ i = 1,2,,n
i xi1 < xi2 < < xiLi ui ,  Li 3 ,  for ​ i=1,2,,n
and a vector pp that, for each ii, locates within {xij}j{xij}j the iith coordinate of an initial point x0x0; that is, if xi0 = xijxi0=xij for some j = 1,2,,Lij=1,2,,Li, then pi = jpi=j. A good guess for the global optimum can be used as x0x0.
The initialization points and the vectors  and pp are collectively called the initialization list (and sometimes we will refer to just the initialization points as ‘the initialization list’, whenever this causes no confusion). The initialization data may be input by you, or they can be set to sensible default values by nag_glopt_bnd_mcs_solve (e05jb): if you provide them yourself, list(i,j)listij should contain xijxij, numpts(i)numptsi should contain LiLi, and initpt(i)initpti should contain pipi, for i = 1,2,,ni=1,2,,n and j = 1,2,,Lij=1,2,,Li; if you wish nag_glopt_bnd_mcs_solve (e05jb) to use one of its preset initialization methods, you could choose one of two simple, three-point methods (see Figure 1). If the list generated by one of these methods contains infinite values, attempts are made to generate a safeguarded list using the function subint(x,y)subint(x,y) (which is also used during the splitting procedure, and is described in Section [Splitting]). If infinite values persist, nag_glopt_bnd_mcs_solve (e05jb) exits with ifail = 3ifail=3. There is also the option to generate an initialization list with the aid of linesearches (by setting iinit = 2iinit=2). Starting with the absolutely smallest point in the root box, linesearches are made along each coordinate. For each coordinate, the local minimizers found by the linesearches are put into the initialization list. If there were fewer than three minimizers, they are augmented by close-by values. The final preset initialization option (iinit = 4iinit=4) generates a randomized list, so that independent multiple runs may be made if you suspect a global optimum has not been found. Each call to the initialization function nag_glopt_bnd_mcs_init (e05ja) resets the initial-state vector for the Wichmann–Hill base-generator that is used. Depending on whether you set the optional parameter Repeatability to ‘ON’ or ‘OFF’, the random state is initialized to give a repeatable or non-repeatable sequence. Then, a random integer between 33 and sdlist is selected, which is then used to determine the number of points to be generated in each coordinate; that is, numpts becomes a constant vector, set to this value. The components of list are then generated, from a uniform distribution on the root box if the box is finite, or else in a safeguarded fashion if any bound is infinite. The array initptinitpt is set to point to the best point in list.
Given an initialization list (preset or otherwise), nag_glopt_bnd_mcs_solve (e05jb) evaluates FF at x0x0, and sets the initial estimate of the global minimum, x*x*, to x0x0. Then, for i = 1,2,,ni=1,2,,n, the objective function FF is evaluated at Li1Li-1 points that agree with x*x* in all but the iith coordinate. We obtain pairs (j,fij) ( x^ j , f i j ) , for j = 1,2,,Lij=1,2,,Li, with: x* = j1x*=x^j1, say; with, for jj1jj1,
x̂kj =
{ xk * if ​k ≠ i; xkj otherwise;
x^kj = { xk* if ​ki; xkj otherwise;
and with
fij = F (j) .
fij = F ( x^ j ) .
The point having the smallest function value is renamed x*x* and the procedure is repeated with the next coordinate.
Once nag_glopt_bnd_mcs_solve (e05jb) has a full set of initialization points and function values, it can generate an initial set of sub-boxes. Recall that the root box is B[x,y] = [,u]B[x,y]=[,u], having basepoint x = x0x=x0. The opposite point yy is a corner of [,u][,u] farthest away from xx, in some sense. The point xx need not be a vertex of [,u][,u], and yy is entitled to have infinite coordinates. We loop over each coordinate ii, splitting the current box along coordinate ii into 2Li22Li-2, 2Li12Li-1 or 2Li2Li sub-intervals with exactly one of the ijx^ij as endpoints, depending on whether two, one or none of the ijx^ij are on the boundary. Thus, as well as splitting at ij x^ i j , for j = 1,2,,Lij=1,2,,Li, we split at additional points zij z i j , for j = 2,3,,Lij=2,3,,Li. These additional zijzij are such that
zij = ij1 + qm (ijij1) ,   j = 2,,Li ,
zij = x^ i j-1 + qm ( x^ i j - x^ i j-1 ) ,   j=2,,Li ,
where qq is the golden-section ratio (sqrt(5)1) / 2(5-1)/2, and the exponent mm takes the value 11 or 22, chosen so that the sub-box with the smaller function value gets the larger fraction of the interval. Each child sub-box gets as basepoint the point obtained from x*x* by changing xi * xi* to the xijxij that is a boundary point of the corresponding iith coordinate interval; this new basepoint therefore has function value fijfij. The opposite point is derived from yy by changing yiyi to the other end of that interval.
nag_glopt_bnd_mcs_solve (e05jb) can now rank the coordinates based on an estimated variability of FF. For each ii we compute the union of the ranges of the quadratic interpolant through any three consecutive ijx^ij, taking the difference between the upper and lower bounds obtained as a measure of the variability of FF in coordinate ii. A vector ππ is populated in such a way that coordinate ii has the πiπith highest estimated variability. For tiebreaks, when the x*x* obtained after splitting coordinate ii belongs to two sub-boxes, the one that contains the minimizer of the quadratic models is designated the current sub-box for coordinate i + 1i+1.
Boxes are assigned levels in the following manner. The root box is given level 11. When a sub-box of level ss is split, the child with the smaller fraction of the golden-section split receives level s + 2s+2; all other children receive level s + 1s+1. The box with the better function value is given the larger fraction of the splitting interval and the smaller level because then it is more likely to be split again more quickly. We see that after the initialization procedure the first level is empty and the non-split boxes have levels 2,,nr + 22,,nr+2, so it is meaningful to choose smaxsmax much larger than nrnr. Note that the internal structure of nag_glopt_bnd_mcs_solve (e05jb) demands that smaxsmax be at least nr + 3nr+3.
Examples of initializations in two dimensions are given in Figure 1. In both cases the initial point is x0 = ( + u) / 2x0=(+u)/2; on the left the initialization points are
x1 = ,  x2 = ( + u) / 2 ,  x3 = u ,
x1 = ,  x2 = (+u) / 2 ,  x3 = u ,
while on the right the points are
x1 = (5 + u) / 6 ,  x2 = ( + u) / 2 ,  x3 = ( + 5u) / 6 .
x1 = ( 5 + u ) / 6 ,  x2 = ( + u ) / 2 ,  x3 = ( + 5 u ) / 6 .
In Figure 1, basepoints and levels after initialization are displayed. Note that these initialization lists correspond to iinit = 0iinit=0 and iinit = 1iinit=1, respectively.
Examples of the initialization procedure
Figure 1: Examples of the initialization procedure
After initialization, a series of sweeps through levels is begun. A sweep is defined by three steps:
(i) scan the list of non-split sub-boxes. Fill a record list bb according to bs = 0bs=0 if there is no box at level ss, and with bsbs pointing to a sub-box with the lowest function value among all sub-boxes with level ss otherwise, for 0 < s < smax0<s<smax;
(ii) the sub-box with label bsbs is a candidate for splitting. If the sub-box is not to be split, according to the rules described in Section [Splitting], increase its level by 11 and update bs + 1bs+1 if necessary. If the sub-box is split, mark it so, insert its children into the list of sub-boxes, and update bb if any child with level ss yields a strict improvement of FF over those sub-boxes at level ss;
(iii) increment ss by 11. If s = smaxs=smax then displaying monitoring information and start a new sweep; else if bs = 0bs=0 then repeat this step; else display monitoring information and go to the previous step.
Clearly, each sweep ends after at most smax1smax-1 visits of the third step.

Splitting

Each sub-box is stored by nag_glopt_bnd_mcs_solve (e05jb) as a set of information about the history of the sub-box: the label of its parent, a label identifying which child of the parent it is, etc. Whenever a sub-box B[x,y]B[x,y] of level s < smaxs<smax is a candidate for splitting, as described in Section [Initialization and Sweeps], we recover xx, yy, and the number, njnj, of times coordinate jj has been split in the history of BB. Sub-box BB could be split in one of two ways.
(i) Splitting by rank
If s > 2nr (min nj + 1) s > 2nr ( minnj+1 ) , the box is always split. The splitting index is set to a coordinate ii such that ni = min njni=minnj.
(ii) Splitting by expected gain
If s 2nr (min nj + 1) s 2nr ( minnj+1 ) , the sub-box could be split along a coordinate where a maximal gain in function value is expected. This gain is estimated according to a local separable quadratic model obtained by fitting to 2nr + 12nr+1 function values. If the expected gain is too small the sub-box is not split at all, and its level is increased by 11.
Eventually, a sub-box that is not eligible for splitting by expected gain will reach level 2nr (min nj + 1) + 1 2nr ( minnj+1 ) +1  and then be split by rank, as long as smaxsmax is large enough. As smaxsmax, the rule for splitting by rank ensures that each coordinate is split arbitrarily often.
Before describing the details of each splitting method, we introduce the procedure for correctly handling splitting at adaptive points and for dealing with unbounded intervals. Suppose we want to split the iith coordinate interval {xi,yi}{xi,yi}, where we define {xi,yi} = [min (xi,yi),max (xi,yi)]{xi,yi}=[min(xi,yi),max(xi,yi)], for xiRxiR and yiRyiR-, and where xx is the basepoint of the sub-box being considered. The descendants of the sub-box should shrink sufficiently fast, so we should not split too close to xixi. Moreover, if yiyi is large we want the new splitting value to not be too large, so we force it to belong to some smaller interval {ξ,ξ}{ξ,ξ}, determined by
ξ = subint (xi,yi) ,  ξ = xi + (ξxi) / 10 ,
ξ = subint (xi,yi) ,  ξ = xi + ( ξ - xi ) / 10 ,
where the function subintsubint is defined by
subint (x,y) =
{ sign(y) if ​ 1000|x| < 1 ​ and ​ |y| > 1000 ; 10sign(y)|x| if ​ 1000|x| ≥ 1 ​ and ​ |y| > 1000|x| ; y otherwise.
subint (x,y) = { sign(y) if ​ 1000|x|<1 ​ and ​ |y|>1000 ; 10sign(y)|x| if ​ 1000|x|1 ​ and ​ |y|>1000|x| ; y otherwise.

Splitting by rank

Consider a sub-box BB with level s > 2nr (min nj + 1) s > 2nr ( minnj+1 ) . Although the sub-box has reached a high level, there is at least one coordinate along which it has not been split very often. Among the ii such that ni = min nj ni = minnj  for BB, select the splitting index to be the coordinate with the lowest πiπi (and hence highest variability rank). ‘Splitting by rank’ refers to the ranking of the coordinates by nini and πiπi.
If ni = 0ni=0, so that BB has never been split along coordinate ii, the splitting is done according to the initialization list and the adaptively chosen golden-section split points, as described in Section [Initialization and Sweeps]. Also as covered there, new basepoints and opposite points are generated. The children having the smaller fraction of the golden-section split (that is, those with larger function values) are given level min {s + 2,smax}min{s+2,smax}. All other children are given level s + 1s+1.
Otherwise, BB ranges between xixi and yiyi in the iith coordinate direction. The splitting value is selected to be zi = xi + 2(subint(xi,yi)xi) / 3zi=xi+2(subint(xi,yi)-xi)/3; we are not attempting to split based on a large reduction in function value, merely in order to reduce the size of a large interval, so zizi may not be optimal. Sub-box BB is split at zizi and the golden-section split point, producing three parts and requiring only one additional function evaluation, at the point xx obtained from xx by changing the iith coordinate to zizi. The child with the smaller fraction of the golden-section split is given level min {s + 2,smax}min{s+2,smax}, while the other two parts are given level s + 1s+1. Basepoints are assigned as follows: the basepoint of the first child is taken to be xx, and the basepoint of the second and third children is the point xx. Opposite points are obtained by changing yiyi to the other end of the iith coordinate-interval of the corresponding child.

Splitting by expected gain

When a sub-box BB has level s 2nr (min nj + 1) s 2nr ( minnj+1 ) , we compute the optimal splitting index and splitting value from a local separable quadratic used as a simple local approximation of the objective function. To fit this curve, for each coordinate we need two additional points and their function values. Such data may be recoverable from the history of BB: whenever the iith coordinate was split in the history of BB, we obtained values that can be used for the current quadratic interpolation in coordinate ii.
We loop over ii; for each coordinate we pursue the history of BB back to the root box, and we take the first two points and function values we find, since these are expected to be closest to the current basepoint xx. If the current coordinate has not yet been split we use the initialization list. Then we generate a local separable model e(ξ)e(ξ) for F(ξ)F(ξ) by interpolation at xx and the 2nr2nr additional points just collected:
n
e(ξ) = F(x) + ei(ξi).
i = 1
e(ξ) = F(x) + i=1 n ei (ξi) .
We define the expected gain ie^i in function value when we evaluate at a new point obtained by changing coordinate ii in the basepoint, for each ii, based on two cases:
(i) ni = 0ni=0. We compute the expected gain as
i = min {fij}fipi.
1j Li
e^i = min 1j Li { fij } - f i pi .
Again, we split according to the initialization list, with the new basepoints and opposite points being as before.
(ii) ni > 0ni>0. Now, the iith component of our sub-box ranges from xixi to yiyi. Using the quadratic partial correction function
ei (ξi) = αi (ξixi) + βi (ξixi) 2
ei (ξi) = αi ( ξi - xi ) + βi ( ξi - xi ) 2
we can approximate the maximal gain expected when changing xixi only. We will choose the splitting value from {ξ,ξ}{ξ,ξ}. We compute
i = min ei(ξi)
ξi {ξ,ξ}
e^i = min ξi {ξ,ξ} ei (ξi)
and call zizi the minimizer in {ξ,ξ}{ξ,ξ}.
If the expected best function value fexpfexp satisfies
fexp = F(x) + min i < fbest,
1in
fexp = F(x) + min 1in e^i < fbest ,
(1)
where fbestfbest is the current best function value (including those function values obtained by local optimization), we expect the sub-box to contain a better point and so we split it, using as splitting index the component with minimal ie^i. Equation (1) prevents wasting function calls by avoiding splitting sub-boxes whose basepoints have bad function values. These sub-boxes will eventually be split by rank anyway.
We now have a splitting index and a splitting value zizi. The sub-box is split at zizi as long as ziyiziyi, and at the golden-section split point; two or three children are produced. The larger fraction of the golden-section split receives level s + 1s+1, while the smaller fraction receives level min {s + 2,smax}min{s+2,smax}. If it is the case that ziyiziyi and the third child is larger than the smaller of the two children from the golden-section split, the third child receives level s + 1s+1. Otherwise it is given the level min {s + 2,smax}min{s+2,smax}. The basepoint of the first child is set to xx, and the basepoint of the second (and third if it exists) is obtained by changing the iith coordinate of xx to zizi. The opposite points are again derived by changing yiyi to the other end of the iith coordinate interval of BB.
If equation (1) does not hold, we expect no improvement. We do not split, and we increase the level of BB by 11.

Local Search

The local optimization algorithm used by nag_glopt_bnd_mcs_solve (e05jb) uses linesearches along directions that are determined by minimizing quadratic models, all subject to bound constraints. Triples of vectors are computed using coordinate searches based on linesearches. These triples are used in triple search procedures to build local quadratic models for FF. A trust-region-type approach to minimize these models is then carried out, and more information about the coordinate search and the triple search can be found in Huyer and Neumaier (1999).
The local search starts by looking for better points without being too local, by making a triple search using points found by a coordinate search. This yields a new point and function value, an approximation of the gradient of the objective, and an approximation of the Hessian of the objective. Then the quadratic model for FF is minimized over a small box, with the solution to that minimization problem then being used as a linesearch direction to minimize the objective. A measure rr is computed to quantify the predictive quality of the quadratic model.
The third stage is the checking of termination criteria. The local search will stop if more than loclimloclim visits to this part of the local search have occurred, where loclimloclim is the value of the optional parameter Local Searches Limit. If that is not the case, it will stop if the limit on function calls has been exceeded (see the description of the optional parameter Function Evaluations Limit). The final criterion checks if no improvement can be made to the function value, or whether the approximated gradient gg is small, in the sense that
|g|T max (|x|,|xold|) < loctol (f0f) .
|g|T max(|x|,| x old |) < loctol (f0-f) .
The vector xoldxold is the best point at the start of the current loop in this iterative local-search procedure, the constant loctolloctol is the value of the optional parameter Local Searches Tolerance, ff is the objective value at xx, and f0f0 is the smallest function value found by the initialization procedure.
Next, nag_glopt_bnd_mcs_solve (e05jb) attempts to move away from the boundary, if any components of the current point lie there, using linesearches along the offending coordinates. Local searches are terminated if no improvement could be made.
The fifth stage carries out another triple search, but this time it does not use points from a coordinate search, rather points lying within the trust-region box are taken.
The final stage modifies the trust-region box to be bigger or smaller, depending on the quality of the quadratic model, minimizes the new quadratic model on that box, and does a linesearch in the direction of the minimizer. The value of rr is updated using the new data, and then we go back to the third stage (checking of termination criteria).
The Hessians of the quadratic models generated by the local search may not be positive definite, so nag_glopt_bnd_mcs_solve (e05jb) uses the general nonlinear optimizer nag_opt_nlp2_sparse_solve (e04vh) to minimize the models.

Optional Parameters

Several optional parameters in nag_glopt_bnd_mcs_solve (e05jb) define choices in the problem specification or the algorithm logic. In order to reduce the number of formal parameters of nag_glopt_bnd_mcs_solve (e05jb) these optional parameters have associated default values that are appropriate for most problems. Therefore, you need only specify those optional parameters whose values are to be different from their default values.
The remainder of this section can be skipped if you wish to use the default values for all optional parameters.
The following is a list of the optional parameters available. A full description of each optional parameter is provided in Section [Description of the optional parameters].
Optional parameters may be specified by calling one, or more, of the functions nag_glopt_bnd_mcs_optset_string (e05jd), nag_glopt_bnd_mcs_optset_char (e05je), nag_glopt_bnd_mcs_optset_int (e05jf) and nag_glopt_bnd_mcs_optset_real (e05jg) before a call to nag_glopt_bnd_mcs_solve (e05jb).
All optional parameters not specified by you are set to their default values. Valid values of optional parameters specified by you are unaltered by nag_glopt_bnd_mcs_solve (e05jb) and so remain in effect for subsequent calls to nag_glopt_bnd_mcs_solve (e05jb), unless you explicitly change them.

Description of the Optional Parameters

For each option, we give a summary line, a description of the optional parameter and details of constraints.
The summary line contains:
Option names are case-insensitive and must be provided in full; abbreviations are not recognized.
Defaults  
This special keyword is used to reset all optional parameters to their default values, and any random state stored in the array comm will be destroyed.
Any option value given with this keyword will be ignored. This optional parameter cannot be queried or got.
Function Evaluations Limit  ii
Default = 100nr2=100nr2
This puts an approximate limit on the number of function calls allowed. The total number of calls made is checked at the top of an internal iteration loop, so it is possible that a few calls more than nfnf may be made.
Constraint: nf > 0nf>0.
Infinite Bound Size  rr
Default = rmax(1/4)=rmax14
This defines the ‘infinite’ bound infbndinfbnd in the definition of the problem constraints. Any upper bound greater than or equal to infbndinfbnd will be regarded as  (and similarly any lower bound less than or equal to infbnd-infbnd will be regarded as -).
Constraint: rmax(1/4)infbndrmax(1/2)rmax14infbndrmax12.
Local Searches  aa
Default = 'ON'='ON'
If you want to try to accelerate convergence of nag_glopt_bnd_mcs_solve (e05jb) by starting local searches from candidate minima, you will require lcsrchlcsrch to be ‘ON’.
Constraint: lcsrch = 'ON'​ or ​'OFF'lcsrch='ON'​ or ​'OFF'.
Local Searches Limit  ii
Default = 50=50
This defines the maximal number of iterations to be used in the trust-region loop of the local-search procedure.
Constraint: loclim > 0loclim>0.
Local Searches Tolerance  rr
Default = 2ε=2ε
The value of loctolloctol is the multiplier used during local searches as a stopping criterion for when the approximated gradient is small, in the sense described in Section [Local Search].
Constraint: loctol2εloctol2ε.
Minimize  
Default
Maximize  
These keywords specify the required direction of optimization. Any option value given with these keywords will be ignored.
Nolist  
Default
List  
These options control the echoing of each optional parameter specification as it is supplied. List turns printing on, Nolist turns printing off. The output is sent to the current advisory message unit (as defined by nag_file_set_unit_advisory (x04ab)).
Any option value given with these keywords will be ignored. This optional parameter cannot be queried or got.
Repeatability  aa
Default = 'OFF'='OFF'
For use with random initialization lists (iinit = 4iinit=4). When set to ‘ON’, an internally-initialized random state is stored in the array comm for use in subsequent calls to nag_glopt_bnd_mcs_solve (e05jb).
Constraint: repeat = 'ON'​ or ​'OFF'repeat='ON'​ or ​'OFF'.
Splits Limit  ii
Default = d(nr + 2) / 3=d(nr+2)/3
Along with the initialization list list, this defines a limit on the number of times the root box will be split along any single coordinate direction. If Local Searches is ‘OFF’ you may find the default value to be too small.
Constraint: smax > nr + 2smax>nr+2.
Static Limit  ii
Default = 3nr=3nr
As the default termination criterion, computation stops when the best function value is static for stclimstclim sweeps through levels. This parameter is ignored if you have specified a target value to reach in Target Objective Value.
Constraint: stclim > 0stclim>0.
Target Objective Error  rr
Default = ε(1/4)=ε14
If you have given a target objective value to reach in objvalobjval (the value of the optional parameter Target Objective Value), objerrobjerr sets your desired relative error (from above if Minimize is set, from below if Maximize is set) between obj and objvalobjval, as described in Section [Accuracy]. See also the description of the optional parameter Target Objective Safeguard.
Constraint: objerr2εobjerr2ε.
Target Objective Safeguard  rr
Default = ε(1/2)=ε12
If you have given a target objective value to reach in objvalobjval (the value of the optional parameter Target Objective Value), objsfgobjsfg sets your desired safeguarded termination tolerance, for when objvalobjval is close to zero.
Constraint: objsfg2εobjsfg2ε.
Target Objective Value  rr
This parameter may be set if you wish nag_glopt_bnd_mcs_solve (e05jb) to use a specific value as the target function value to reach during the optimization. Setting objvalobjval overrides the default termination criterion determined by the optional parameter Static Limit.

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013