text:
# NAG FL Interfacee05kbf (handle_solve_mcs)

Note: this routine uses optional parameters to define choices in the problem specification and in the details of the algorithm. If you wish to use default settings for all of the optional parameters, you need only read Sections 1 to 10 of this document. If, however, you wish to reset some or all of the settings please refer to Section 11 for a detailed description of the algorithm, and to Section 12 for a detailed description of the specification of the optional parameters.

## ▸▿ Contents

Settings help

FL Name Style:

FL Specification Language:

## 1Purpose

e05kbf 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.

## 2Specification

Fortran Interface
 Subroutine e05kbf ( nvar, x,
 Integer, Intent (In) :: nvar Integer, Intent (Inout) :: iuser(*), ifail Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Real (Kind=nag_wp), Intent (Out) :: x(nvar), rinfo(100), stats(100) Type (c_ptr), Intent (In) :: handle, cpuser External :: objfun, monit
#include <nag.h>
 void e05kbf_ (void **handle, void (NAG_CALL *objfun)(const Integer *nvar, const double x[], double *fx, Integer *inform, Integer iuser[], double ruser[], void **cpuser),void (NAG_CALL *monit)(const Integer *nvar, const double x[], Integer *inform, const double rinfo[], const double stats[], Integer iuser[], double ruser[], void **cpuser),const Integer *nvar, double x[], double rinfo[], double stats[], Integer iuser[], double ruser[], void **cpuser, Integer *ifail)
The routine may be called by the names e05kbf or nagf_glopt_handle_solve_mcs.

## 3Description

e05kbf is aimed at minimizing a nonlinear objective function subject to bound constraints:
 $minimize x∈ℝn f (x) subject to lx ≤ x ≤ ux .$
Here $f$ is a smooth nonlinear function and ${l}_{x}$ and ${u}_{x}$ are $n$-dimensional vectors defining bounds on the variables.
e05kbf serves as a solver for compatible problems stored as a handle. The handle points to an internal data structure which defines the problem and serves as a means of communication for routines in the NAG optimization modelling suite. To define a compatible problem handle, you must call e04raf to initialize it followed, optionally, by e04rgf, e04ref or e04rff to define the objective function and e04rhf to define bounds on the variables. If e04rhf is not called, all the variables will be considered free by the solver. It should be noted that e05kbf always assumes that the gradient of the objective is dense, therefore, defining a sparse structure for the residuals in the call to e04rgf will have no effect. Additionally, the multi-level coordinate search method used by this solver relies on dividing the feasible space in ‘boxes’ (see Section 11 for a more thorough explanation) so it is advisable to define reasonable bounds for the variables using e04rhf. See Section 3.1 in the E04 Chapter Introduction for more details about the NAG optimization modelling suite.
The algorithm behaviour and solver strategy can be modified by various optional parameters (see Section 12) which can be set by e04zmf and e04zpf at any time between the initialization of the handle by e04raf and a call to the solver. The optional parameters' names specific for this solver start with the prefix MCS (Multi-level Coordinate Search). The default values for these optional parameters are chosen to work well in the general case, but it is recommended that you tune them to your particular problem. In particular, if the objective function is known to be numerically difficult, it could be desirable to define a customized initialization list for the algorithm with the optional parameter MCS Initialization Method. For more details on how to create a custom initialization list, please refer to Section 9.2. Once the solver has finished, options may be modified for the next solve. The solver may be called repeatedly with various initialization lists and/or optional parameters.
The method used by e05kbf 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 11.

## 4References

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

## 5Arguments

1: $\mathbf{handle}$Type (c_ptr) Input
On entry: the handle to the problem. It needs to be initialized (e.g., by e04raf) and to hold a problem formulation compatible with e05kbf. It must not be changed between calls to the NAG optimization modelling suite.
2: $\mathbf{objfun}$Subroutine, supplied by the NAG Library or the user. External Procedure
objfun must calculate the value of the nonlinear objective function $f\left(x\right)$ at a specified point $x$. If there is no nonlinear objective (e.g., e04ref or e04rff was called to define a linear or quadratic objective function), objfun will never be called by e05kbf and may be the dummy routine e04kfv (included in the NAG Library.)
The specification of objfun is:
Fortran Interface
 Subroutine objfun ( nvar, x, fx,
 Integer, Intent (In) :: nvar Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar) Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Real (Kind=nag_wp), Intent (Out) :: fx Type (c_ptr), Intent (In) :: cpuser
 void objfun (const Integer *nvar, const double x[], double *fx, Integer *inform, Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of variable values at which the objective function is to be evaluated.
3: $\mathbf{fx}$Real (Kind=nag_wp) Output
On exit: the value of the objective function at $x$.
4: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
On exit: may be used to indicate that the requested objective value could not be computed. Specifically, it can be set to a negative value:
${\mathbf{inform}}<0$
The solver will cleanly exit with ${\mathbf{ifail}}={\mathbf{20}}$ and return the best available point as well as the solve statistics.
5: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
6: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
7: $\mathbf{cpuser}$Type (c_ptr) User Workspace
objfun is called with the arguments iuser, ruser and cpuser as supplied to e05kbf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to objfun.
objfun must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e05kbf is called. Arguments denoted as Input must not be changed by this procedure.
Note: objfun should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e05kbf. If your code inadvertently does return any NaNs or infinities, e05kbf is likely to produce unexpected results.
3: $\mathbf{monit}$Subroutine, supplied by the NAG Library or the user. External Procedure
monit is provided to enable you to monitor the progress of the optimization and, optionally, to terminate the solver early if necessary. It is invoked at the end of every $i$th step where $i$ is given by the MCS Monitor Frequency (the default is $0$, monit is not called). A step is complete when the procedure in which a sub-box is considered for splitting finishes correctly.
monit may be the dummy subroutine e04kfu (included in the NAG Library).
The specification of monit is:
Fortran Interface
 Subroutine monit ( nvar, x,
 Integer, Intent (In) :: nvar Integer, Intent (Inout) :: inform, iuser(*) Real (Kind=nag_wp), Intent (In) :: x(nvar), rinfo(100), stats(100) Real (Kind=nag_wp), Intent (Inout) :: ruser(*) Type (c_ptr), Intent (In) :: cpuser
 void monit (const Integer *nvar, const double x[], Integer *inform, const double rinfo[], const double stats[], Integer iuser[], double ruser[], void **cpuser)
1: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
2: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Input
On entry: the vector $x$ of decision variables at the current iteration.
3: $\mathbf{inform}$Integer Input/Output
On entry: a non-negative value.
On exit: may be used to request the solver to stop immediately. Specifically, if ${\mathbf{inform}}<0$ the solver will terminate immediately with ${\mathbf{ifail}}={\mathbf{20}}$; otherwise, the solver will proceed normally.
4: $\mathbf{rinfo}\left(100\right)$Real (Kind=nag_wp) array Input
On entry: error measures and various indicators at the end of the current iteration as described in rinfo.
5: $\mathbf{stats}\left(100\right)$Real (Kind=nag_wp) array Input
On entry: solver statistics at the end of the current iteration as described in stats.
6: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
7: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
8: $\mathbf{cpuser}$Type (c_ptr) User Workspace
monit is called with the arguments iuser, ruser and cpuser as supplied to e05kbf. You should use the arrays iuser and ruser, and the data handle cpuser to supply information to monit.
monit must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which e05kbf is called. Arguments denoted as Input must not be changed by this procedure.
Note: monit should not return floating-point NaN (Not a Number) or infinity values, since these are not handled by e05kbf. If your code inadvertently does return any NaNs or infinities, e05kbf is likely to produce unexpected results.
4: $\mathbf{nvar}$Integer Input
On entry: $n$, the current number of decision variables $x$ in the model.
5: $\mathbf{x}\left({\mathbf{nvar}}\right)$Real (Kind=nag_wp) array Output
On exit: the final values of the variables, $x$.
6: $\mathbf{rinfo}\left(100\right)$Real (Kind=nag_wp) array Output
On exit: error measures and various indicators at the end of the final iteration as given in the table below:
 $1$ Objective function value $\mathrm{f}\left(x\right)$. $2$ Number of sweeps. $3$–$100$ Reserved for future use.
7: $\mathbf{stats}\left(100\right)$Real (Kind=nag_wp) array Output
On exit: solver statistics at the end of the final iteration as given in the table below:
 $1$ Number of calls to the objective function. $2$ Total time spent in the solver (including time spent evaluating the objective). $3$ Total time spent evaluating the objective function. $4$ Number of sub-boxes. $5$ Number of splits. $6$ Cumulative number of splits by initialization list. $7$ The current lowest level containing no split boxes. $8$–$100$ Reserved for future use.
8: $\mathbf{iuser}\left(*\right)$Integer array User Workspace
9: $\mathbf{ruser}\left(*\right)$Real (Kind=nag_wp) array User Workspace
10: $\mathbf{cpuser}$Type (c_ptr) User Workspace
iuser, ruser and cpuser are not used by e05kbf, but are passed directly to objfun and monit and may be used to pass information to these routines. If you do not need to reference cpuser, it should be initialized to c_null_ptr.
11: $\mathbf{ifail}$Integer Input/Output
On entry: ifail must be set to $0$, $-1$ or $1$ to set behaviour on detection of an error; these values have no effect when no error is detected.
A value of $0$ causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of $-1$ means that an error message is printed while a value of $1$ means that it is not.
If halting is not appropriate, the value $-1$ or $1$ is recommended. If message printing is undesirable, then the value $1$ is recommended. Otherwise, the value $-1$ is recommended since useful values can be provided in some output arguments even when ${\mathbf{ifail}}\ne {\mathbf{0}}$ on exit. When the value $-\mathbf{1}$ or $\mathbf{1}$ is used it is essential to test the value of ifail on exit.
On exit: ${\mathbf{ifail}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see Section 6).

## 6Error Indicators and Warnings

If on entry ${\mathbf{ifail}}=0$ or $-1$, explanatory error messages are output on the current error message unit (as defined by x04aaf).
Errors or warnings detected by the routine:
Note: in some cases e05kbf may return useful information.
${\mathbf{ifail}}=1$
The supplied handle does not define a valid handle to the data structure for the NAG optimization modelling suite. It has not been properly initialized or it has been corrupted.
${\mathbf{ifail}}=2$
The problem is already being solved.
This solver does not support the model defined in the handle.
${\mathbf{ifail}}=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 (${\mathbf{MCS Initialization Method}}=\mathrm{CUSTOM}$) 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.
${\mathbf{ifail}}=4$
On entry, ${\mathbf{nvar}}=⟨\mathit{\text{value}}⟩$, expected $\mathrm{value}=⟨\mathit{\text{value}}⟩$.
Constraint: nvar must match the current number of variables of the model in the handle.
${\mathbf{ifail}}=5$
On entry, user-supplied $\mathit{initpt}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{initpt}\left(\mathit{i}\right)\ge 1$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$.
For more details, please refer to Section 9.2.
On entry, user-supplied $\mathit{initpt}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$ and $\mathit{sdlist}=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{initpt}\left(\mathit{i}\right)\le \mathit{sdlist}$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$.
For more details, please refer to Section 9.2.
On entry, user-supplied $\mathit{list}\left(i,\mathit{j}\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$, $\mathit{j}=⟨\mathit{\text{value}}⟩$, and $\mathit{bl}\left(i\right)=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{list}\left(\mathit{i},\mathit{j}\right)\ge \mathit{bl}\left(\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$ and $\mathit{j}=1,2,\dots ,\mathit{numpts}\left(\mathit{i}\right)$.
For more details, please refer to Section 9.2.
On entry, user-supplied $\mathit{list}\left(i,\mathit{j}\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$, $\mathit{j}=⟨\mathit{\text{value}}⟩$, and $\mathit{bu}\left(i\right)=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{list}\left(i,\mathit{j}\right)\le \mathit{bu}\left(\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$ and $\mathit{j}=1,2,\dots ,\mathit{numpts}\left(\mathit{i}\right)$.
For more details, please refer to Section 9.2.
On entry, user-supplied $\mathit{numpts}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{numpts}\left(\mathit{i}\right)\ge 3$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$.
For more details, please refer to Section 9.2.
On entry, user-supplied $\mathit{numpts}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$ and $\mathit{sdlist}=⟨\mathit{\text{value}}⟩$.
Constraint: if ${\mathbf{x}}\left(i\right)$ is not fixed then $\mathit{numpts}\left(\mathit{i}\right)\le \mathit{sdlist}$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$.
For more details, please refer to Section 9.2.
On entry, user-supplied section $\mathit{list}\left(i,1:\mathit{numpts}\left(i\right)\right)$ contained $\mathit{ndist}$ distinct elements, and $\mathit{ndist}<\mathit{numpts}\left(i\right)$: $\mathit{ndist}=⟨\mathit{\text{value}}⟩$, $\mathit{numpts}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$.
For more details, please refer to Section 9.2.
On entry, user-supplied section $\mathit{list}\left(i,1:\mathit{numpts}\left(i\right)\right)$ was not in ascending order: $\mathit{numpts}\left(i\right)=⟨\mathit{\text{value}}⟩$, $i=⟨\mathit{\text{value}}⟩$.
For more details, please refer to Section 9.2.
The dimension of the array ‘MCS List’ is not a multiple of ${\mathbf{nvar}}$. For more details, please refer to Section 9.2.
${\mathbf{ifail}}=6$
There were $n=⟨\mathit{\text{value}}⟩$ variables and the optional parameter MCS Splits Limit $\mathit{smax}=⟨\mathit{\text{value}}⟩$.
Constraint: $\mathit{smax}>n+2$.
Use e04zmf to set compatible option values.
${\mathbf{ifail}}=7$
The dummy objfun routine was called but the problem requires these values. Please provide a proper objfun routine.
${\mathbf{ifail}}=8$
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. It is likely that your objective function is badly scaled: try rescaling it. Also, try relaxing the bounds or use a different initialization method. If the problem persists, please contact NAG quoting error code $⟨\mathit{\text{value}}⟩$.
${\mathbf{ifail}}=11$
The optional parameter MCS Initialization Method was set to CUSTOM but one of the arrays $\mathit{list}$, $\mathit{initpt}$ or $\mathit{numpts}$ was not passed correctly to the handle.
For more details, please refer to Section 9.2.
${\mathbf{ifail}}=20$
User-supplied monitoring routine requested termination.
User-supplied objective function requested termination.
${\mathbf{ifail}}=21$
The function evaluations limit was exceeded.
Approximately MCS Max Objective Calls function calls have been made without your chosen termination criterion being satisfied.
${\mathbf{ifail}}=24$
The division procedure completed but your target value could not be reached.
Despite every sub-box being processed MCS Splits Limit times, the target value you provided in MCS Target Objective Value could not be found to the tolerances given in MCS Target Objective Error and MCS Target Objective Safeguard. You could try reducing MCS Splits Limit or the objective tolerances.
${\mathbf{ifail}}=-99$
See Section 7 in the Introduction to the NAG Library FL Interface for further information.
${\mathbf{ifail}}=-399$
Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library FL Interface for further information.
${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.
See Section 9 in the Introduction to the NAG Library FL Interface for further information.

## 7Accuracy

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

## 8Parallelism and Performance

e05kbf is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
e05kbf makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the Users' Note for your implementation for any additional implementation-specific information.

### 9.1Description of the Printed Output

The solver can print information to give an overview of the problem and the progress of the computation. The output may be sent to two independent unit numbers which are set by optional parameters Print File and Monitoring File. Optional parameters Print Level, Print Options, Monitoring Level and Print Solution determine the exposed level of detail. This allows, for example, a detailed log file to be generated while the condensed information is displayed on the screen.
By default (${\mathbf{Print File}}=6$, ${\mathbf{Print Level}}=2$), four sections are printed to the standard output: a header, a list of options, an iteration log and a summary.
The header is a message indicating the start of the solver. It should look like:
``` -------------------------------------------------
E05KB, MCS method for bound constrained problems
-------------------------------------------------```
Optional parameters list
If ${\mathbf{Print Options}}=\mathrm{YES}$, a list of the optional parameters and their values is printed. The list shows all options of the solver, each displayed on one line. The line contains the option name, its current value and an indicator for how it was set. The options left at their defaults are noted by ‘d’ and the ones you have set are noted by ‘U’. Note that the output format is compatible with the file format expected by e04zpf. The output looks as follows:
```     Mcs Initialization Method     =       Simple Bounds     * d
Mcs Local Searches            =                  On     * d
Mcs Local Searches Limit      =                  50     * d
Mcs Local Searches Tolerance  =         2.22045E-16     * d
Mcs Max Objective Calls       =                   0     * d
Mcs Monitor Frequency         =                   0     * d
Mcs Repeatability             =                 Off     * d
Mcs Splits Limit              =                   0     * d
Mcs Static Limit              =                   0     * d
Mcs Target Objective Error    =         1.02648E-04     * d
Mcs Target Objective Safeguard=         1.05367E-08     * d
Mcs Target Objective Value    =        -1.7977E+308     * d
Mcs Print Frequency           =                   1     * d```
Problem statistics
If ${\mathbf{Print Level}}\ge 2$, statistics on the problem are printed, for example:
``` Problem Statistics
No of variables                  4
free (unconstrained)           1
bounded                        3
Objective function       Nonlinear```
Iteration log
If ${\mathbf{Print Level}}\ge 2$, the solver will print a summary line for each split. Each line shows the split number, the value of the objective function (obj), the number of sweeps (nsweeps), and the cumulative number of objective function evaluations (nf). The output looks as follows:
```  ----------------------------------------------
Splits |    obj        nsweeps    |    nf   |
----------------------------------------------
1 | -8.00E+00        1       |     8   |
2 | -9.00E+00        1       |     9   |
3 | -9.00E+00        1       |    10   |
4 | -9.00E+00        1       |    11   |```
Summary
Once the solver finishes, a summary is produced:
```Status: Converged, objective function is static

Value of the objective                   -9.00000E+00
Number of objective function evaluations           81
Number of sweeps                                   10```
Note that only the iterations that decrease the objective function are printed in the iteration log, meaning that objective evaluations are likely to happen between the last printed iteration and the convergence. This leads to a small difference between the last line of the iteration log and the final summary in terms of the number of function evaluations.
Optionally, if ${\mathbf{Stats Time}}=\mathrm{YES}$, the timings are printed:
```Timings
Total time spent in the solver            0.056
Time spent in the objective evaluation    0.012```
Additionally, if ${\mathbf{Print Solution}}=\mathrm{YES}$, the solution is printed along with the bounds:
``` Primal variables:
idx   Lower bound       Value       Upper bound
1   1.00000E+00    1.00000E+00    3.00000E+00
2  -2.00000E+00   -8.52342E-02    0.00000E+00
3       -inf       4.09304E-01         inf
4   1.00000E+00    1.00000E+00    3.00000E+00```

### 9.2Creating a Custom Initialization List

By default, the initialization list as described in Section 11.1 is computed internally by the solver using one of the methods defined by the optional parameter MCS Initialization Method. It is however possible to customize the initialization method by first setting MCS Initialization Method to $\mathrm{CUSTOM}$ and providing the following information to the handle before the solver call through the utility routines e04rxf and e04rwf:
Numpts(nvar): integer array
$\mathit{numpts}$ encodes the number of splitting points in each coordinate direction; e.g., for $i=1,\mathrm{...},{\mathbf{nvar}}$, $\mathrm{numpts}\left(i\right)$ is the number of splits intially created in the $i$th direction. The coordinates of the initial boxes must be provided in the argument $\mathit{list}\left({\mathbf{nvar}}*\mathit{sdlist}\right)$, where sdlist is defined implicitly when list is passed to the problem structure.
To provide the array numpts, the routine e04rwf must be called with the command string ${\mathbf{cmdstr}}=\text{'}\mathrm{MCS number points}\text{'}$.
Constraint:
• $\mathit{numpts}\left(\mathit{i}\right)\ge 3$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$
• ${\mathrm{max}}_{i}\mathit{numpts}\left(i\right)=\mathit{sdlist}$
Initpt(nvar): integer array
$\mathit{initpt}$ designates a point stored in $\mathit{list}$ that you wish e05kbf to consider as an ‘initial point’ for the purposes of the splitting procedure. Call this initial point ${\mathit{x}}^{0}$. The coordinates of ${\mathit{x}}^{0}$ correspond to a set of indices ${\mathit{j}}_{i}$, $i=1,\mathrm{...},{\mathbf{nvar}}$, such that ${\mathit{x}}_{i}^{0}$ is stored in $\mathit{list}\left(i+\left({j}_{i}-1\right)×{\mathbf{nvar}}\right)$. You must set $\mathit{initpt}\left(i\right)={j}_{i}$ for $i=1,\mathrm{...},{\mathbf{nvar}}$.
To provide the array initpt, the routine e04rwf must be called with the command string ${\mathbf{cmdstr}}=\text{'}\mathrm{MCS initial points}\text{'}$.
Constraint:
• $1\le \mathit{initpt}\left(\mathit{i}\right)\le \mathit{sdlist}={\mathrm{max}}_{\mathit{i}}\mathit{numpts}\left(\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,{\mathbf{nvar}}$
List(${\mathbf{nvar}}×\mathit{sdlist}$): real array
$\mathit{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 $i$ the split is done at the values $\mathit{list}\left(i+\left(j-1\right)×{\mathbf{nvar}}\right)$, $j=1,\mathrm{...},\mathit{numpts}\left(i\right)$, as well as at some adaptively chosen intermediate points. For all $i$, the values $\mathit{list}\left(i+\left(j-1\right)×{\mathbf{nvar}}\right)$, $j=1,\mathrm{...},\mathit{numpts}\left(i\right)$, must be in ascending order with each entry being distinct.
The dimension sdlist is defined implicitly as $\mathit{sdlist}=\frac{\mathit{dim}\left(\mathit{list}\right)}{{\mathbf{nvar}}}$ when list is passed to the handle.
To provide the array list, the routine e04rxf must be called with the command string ${\mathbf{cmdstr}}=\text{'}\mathrm{MCS list}\text{'}$.
Constraint:
• for all $j=1,\mathrm{...},\mathit{numpts}\left(\mathit{i}\right)$ and $i=1,\mathrm{...},{\mathbf{nvar}}$, $\mathit{bl}\left(\mathit{i}\right)\le \mathit{list}\left(\mathit{i}+\left(j-1\right)×{\mathbf{nvar}}\right)\le \mathit{bu}\left(\mathit{i}\right)$, where bl and bu are the bounds defined on the variables.

## 10Example

This example finds the global minimum of the ‘peaks’ function in two dimensions
 $F(x,y) = 3 (1-x) 2 exp(- x 2 - (y+1) 2 ) -10  ( x 5 - x 3 - y 5 ) exp(- x 2 - y 2 ) - 1 3 exp(- (x+1) 2 - y 2 )$
on the box $\left[-3,3\right]×\left[-3,3\right]$.
The function $F$ has several local minima and one global minimum in the given box. The global minimum is approximately located at $\left(0.23,-1.63\right)$, where the function value is approximately $-6.55$.

### 10.1Program Text

Program Text (e05kbfe.f90)

None.

### 10.3Program Results

Program Results (e05kbfe.r)

## 11Algorithmic Details

Here we summarise the main features of the MCS algorithm used in e05kbf, and we introduce some terminology used in the description of the subroutine and its arguments. We assume throughout that we will only do any work in coordinates $i$ in which ${x}_{i}$ is free to vary. The MCS algorithm is fully described in Huyer and Neumaier (1999).

### 11.1Initialization and Sweeps

Each sub-box is determined by a basepoint $\mathbf{x}$ and an opposite point $\mathbf{y}$. We denote such a sub-box by $B\left[\mathbf{x},\mathbf{y}\right]$. 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 $i$ for the first time (in the initialization procedure or later), the splitting is done at three or more user-defined values ${\left\{{x}_{i}^{j}\right\}}_{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$
and a vector $\mathbf{p}$ that, for each $i$, locates within ${\left\{{x}_{i}^{j}\right\}}_{j}$ the $i$th coordinate of an initial point ${\mathbf{x}}^{0}$; that is, if ${x}_{i}^{0}={x}_{i}^{j}$ for some $j=1,2,\dots ,{L}_{i}$, then ${p}_{i}=j$. A good guess for the global optimum can be used as ${\mathbf{x}}^{0}$.
The initialization points and the vectors $\mathbf{\ell }$ and $\mathbf{p}$ 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 e05kbf: if you provide them yourself, the arrays $\mathit{list}\left(i+\left(j-1\right)×{\mathbf{nvar}}\right)$ should contain ${x}_{i}^{j}$, $\mathit{numpts}\left(i\right)$ should contain ${L}_{i}$, and $\mathit{initpt}\left(i\right)$ should contain ${p}_{\mathit{i}}$, for $\mathit{i}=1,2,\dots ,n$ and $\mathit{j}=1,2,\dots ,{L}_{\mathit{i}}$; please refer to Section 9.2 for more information on how to define a custom initialization list. If you wish e05kbf 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 $\mathrm{subint}\left(x,y\right)$ (which is also used during the splitting procedure, and is described in Section 11.2). If infinite values persist, e05kbf exits with ${\mathbf{ifail}}={\mathbf{3}}$. There is also the option to generate an initialization list with the aid of linesearches (by setting ${\mathbf{MCS Initialization Method}}=SIMPLE OFF-BOUNDS$). 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 (${\mathbf{MCS Initialization Method}}=RANDOM$) generates a randomized list, so that independent multiple runs may be made if you suspect a global optimum has not been found. Depending on whether you set the optional parameter MCS Repeatability to $\mathrm{ON}$ or $\mathrm{OFF}$, the random state is initialized to give a repeatable or non-repeatable sequence. 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.
Given an initialization list (preset or otherwise), e05kbf evaluates $F$ at ${\mathbf{x}}^{0}$ and sets the initial estimate of the global minimum, ${\mathbf{x}}^{*}$, to ${\mathbf{x}}^{0}$. Then, for $i=1,2,\dots ,n$, the objective function $F$ is evaluated at ${L}_{i}-1$ points that agree with ${\mathbf{x}}^{*}$ in all but the $i$th coordinate. We obtain pairs $\left({\stackrel{^}{\mathbf{x}}}^{\mathit{j}},{f}_{i}^{\mathit{j}}\right)$, for $\mathit{j}=1,2,\dots ,{L}_{i}$, with: ${\mathbf{x}}^{*}={\stackrel{^}{\mathbf{x}}}^{{j}_{1}}$, say; with, for $j\ne {j}_{1}$,
 $x^kj = { xk* if k≠i; xkj otherwise;$
and with
 $fij = F ( x^ j ) .$
The point having the smallest function value is renamed ${\mathbf{x}}^{*}$ and the procedure is repeated with the next coordinate.
Once e05kbf 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\left[\mathbf{x},\mathbf{y}\right]=\left[\mathbf{\ell },\mathbf{u}\right]$, having basepoint $\mathbf{x}={\mathbf{x}}^{0}$. The opposite point $\mathbf{y}$ is a corner of $\left[\mathbf{\ell },\mathbf{u}\right]$ farthest away from $\mathbf{x}$, in some sense. The point $\mathbf{x}$ need not be a vertex of $\left[\mathbf{\ell },\mathbf{u}\right]$, and $\mathbf{y}$ is entitled to have infinite coordinates. We loop over each coordinate $i$, splitting the current box along coordinate $i$ into $2{L}_{i}-2$, $2{L}_{i}-1$ or $2{L}_{i}$ sub-intervals with exactly one of the ${\stackrel{^}{x}}_{i}^{j}$ as endpoints, depending on whether two, one or none of the ${\stackrel{^}{x}}_{i}^{j}$ are on the boundary. Thus, as well as splitting at ${\stackrel{^}{x}}_{i}^{\mathit{j}}$, for $\mathit{j}=1,2,\dots ,{L}_{i}$, we split at additional points ${z}_{i}^{\mathit{j}}$, for $\mathit{j}=2,3,\dots ,{L}_{i}$. These additional ${z}_{i}^{j}$ are such that
 $zij = x^ i j-1 + qm ( x^ i j - x^ i j-1 ) , j=2,…,Li ,$
where $q$ is the golden-section ratio $\left(\sqrt{5}-1\right)/2$, and the exponent $m$ takes the value $1$ or $2$, 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 ${\mathbf{x}}^{*}$ by changing ${x}_{i}^{*}$ to the ${x}_{i}^{j}$ that is a boundary point of the corresponding $i$th coordinate interval; this new basepoint, therefore, has function value ${f}_{i}^{j}$. The opposite point is derived from $\mathbf{y}$ by changing ${y}_{i}$ to the other end of that interval.
e05kbf can now rank the coordinates based on an estimated variability of $F$. For each $i$ we compute the union of the ranges of the quadratic interpolant through any three consecutive ${\stackrel{^}{x}}_{i}^{j}$, taking the difference between the upper and lower bounds obtained as a measure of the variability of $F$ in coordinate $i$. A vector $\mathbf{\pi }$ is populated in such a way that coordinate $i$ has the ${\pi }_{i}$th highest estimated variability. For tiebreaks, when the ${\mathbf{x}}^{*}$ obtained after splitting coordinate $i$ belongs to two sub-boxes, the one that contains the minimizer of the quadratic models is designated the current sub-box for coordinate $i+1$.
Boxes are assigned levels in the following manner. The root box is given level $1$. When a sub-box of level $s$ is split, the child with the smaller fraction of the golden-section split receives level $s+2$; all other children receive level $s+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,\dots ,{n}_{r}+2$, so it is meaningful to choose ${s}_{\mathrm{max}}$ much larger than ${n}_{r}$. Note that the internal structure of e05kbf demands that ${s}_{\mathrm{max}}$ be at least ${n}_{r}+3$.
Examples of initializations in two dimensions are given in Figure 1. In both cases the initial point is ${\mathbf{x}}^{0}=\left(\mathbf{\ell }+\mathbf{u}\right)/2$; on the left the initialization points are
 $x1 = ℓ , x2 = (ℓ+u) / 2 , x3 = u ,$
while on the right the points are
 $x1 = (5ℓ+u) / 6 , x2 = (ℓ+u) / 2 , x3 = (ℓ+5u) / 6 .$
In Figure 1, basepoints and levels after initialization are displayed. Note that these initialization lists correspond to ${\mathbf{MCS Initialization Method}}=SIMPLE BOUNDS$ and ${\mathbf{MCS Initialization Method}}=SIMPLE OFF-BOUNDS$, respectively.
After initialization, a series of sweeps through levels is begun. A sweep is defined by three steps:
1. (i)scan the list of non-split sub-boxes. Fill a record list $\mathbf{b}$ according to ${b}_{s}=0$ if there is no box at level $s$, and with ${b}_{s}$ pointing to a sub-box with the lowest function value among all sub-boxes with level $s$ otherwise, for $0;
2. (ii)the sub-box with label ${b}_{s}$ is a candidate for splitting. If the sub-box is not to be split, according to the rules described in Section 11.2, increase its level by $1$ and update ${b}_{s+1}$ if necessary. If the sub-box is split, mark it so, insert its children into the list of sub-boxes, and update $\mathbf{b}$ if any child with level ${s}^{\prime }$ yields a strict improvement of $F$ over those sub-boxes at level ${s}^{\prime }$;
3. (iii)increment $s$ by $1$. If $s={s}_{\mathrm{max}}$ then display the monitoring information and start a new sweep; else if ${b}_{s}=0$ then repeat this step; else display monitoring information and go to the previous step.
Clearly, each sweep ends after at most ${s}_{\mathrm{max}}-1$ visits of the third step.

### 11.2Splitting

Each sub-box is stored by e05kbf 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\left[\mathbf{x},\mathbf{y}\right]$ of level $s<{s}_{\mathrm{max}}$ is a candidate for splitting, as described in Section 11.1, we recover $\mathbf{x}$, $\mathbf{y}$, and the number, ${n}_{j}$, of times coordinate $j$ has been split in the history of $B$. Sub-box $B$ could be split in one of two ways.
1. (i)Splitting by rank
If $s>2{n}_{r}\left(\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}+1\right)$, the box is always split. The splitting index is set to a coordinate $i$ such that ${n}_{i}=\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}$.
2. (ii)Splitting by expected gain
If $s\le 2{n}_{r}\left(\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}+1\right)$, 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 $2{n}_{r}+1$ function values. If the expected gain is too small the sub-box is not split at all, and its level is increased by $1$.
Eventually, a sub-box that is not eligible for splitting by expected gain will reach level $2{n}_{r}\left(\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}+1\right)+1$ and then be split by rank, as long as ${s}_{\mathrm{max}}$ is large enough. As ${s}_{\mathrm{max}}\to \infty$, 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 $i$th coordinate interval $▯\left\{{x}_{i},{y}_{i}\right\}$, where we define $▯\left\{{x}_{i},{y}_{i}\right\}=\left[\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({x}_{i},{y}_{i}\right),\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({x}_{i},{y}_{i}\right)\right]$, for ${x}_{i}\in R$ and ${y}_{i}\in \overline{R}$, and where $\mathbf{x}$ 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 ${x}_{i}$. Moreover, if ${y}_{i}$ is large we want the new splitting value to not be too large, so we force it to belong to some smaller interval $▯\left\{{\xi }^{\prime },{\xi }^{\prime \prime }\right\}$, determined by
 $ξ′′ = subint (xi,yi) , ξ′ = xi + (ξ′′-xi) / 10 ,$
where the function $\mathrm{subint}$ 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.$

#### 11.2.1Splitting by rank

Consider a sub-box $B$ with level $s>2{n}_{r}\left(\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}+1\right)$. 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 $i$ such that ${n}_{i}=\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}$ for $B$, select the splitting index to be the coordinate with the lowest ${\pi }_{i}$ (and hence highest variability rank). ‘Splitting by rank’ refers to the ranking of the coordinates by ${n}_{i}$ and ${\pi }_{i}$.
If ${n}_{i}=0$, so that $B$ has never been split along coordinate $i$, the splitting is done according to the initialization list and the adaptively chosen golden-section split points, as described in Section 11.1. 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 $\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left\{s+2,{s}_{\mathrm{max}}\right\}$. All other children are given level $s+1$.
Otherwise, $B$ ranges between ${x}_{i}$ and ${y}_{i}$ in the $i$th coordinate direction. The splitting value is selected to be ${z}_{i}={x}_{i}+2\left(\mathrm{subint}\left({x}_{i},{y}_{i}\right)-{x}_{i}\right)/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 ${z}_{i}$ may not be optimal. Sub-box $B$ is split at ${z}_{i}$ and the golden-section split point, producing three parts and requiring only one additional function evaluation, at the point ${\mathbf{x}}^{\prime }$ obtained from $\mathbf{x}$ by changing the $i$th coordinate to ${z}_{i}$. The child with the smaller fraction of the golden-section split is given level $\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left\{s+2,{s}_{\mathrm{max}}\right\}$, while the other two parts are given level $s+1$. Basepoints are assigned as follows: the basepoint of the first child is taken to be $\mathbf{x}$, and the basepoint of the second and third children is the point ${\mathbf{x}}^{\prime }$. Opposite points are obtained by changing ${y}_{i}$ to the other end of the $i$th coordinate-interval of the corresponding child.

#### 11.2.2Splitting by expected gain

When a sub-box $B$ has level $s\le 2{n}_{r}\left(\mathrm{min}\phantom{\rule{0.25em}{0ex}}{n}_{j}+1\right)$, 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 $B$: whenever the $i$th coordinate was split in the history of $B$, we obtained values that can be used for the current quadratic interpolation in coordinate $i$.
We loop over $i$; for each coordinate we pursue the history of $B$ 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 $\mathbf{x}$. If the current coordinate has not yet been split we use the initialization list. Then we generate a local separable model $e\left(\mathbf{\xi }\right)$ for $F\left(\mathbf{\xi }\right)$ by interpolation at $\mathbf{x}$ and the $2{n}_{r}$ additional points just collected:
 $e(ξ) = F(x) + ∑ i=1 n ei (ξi) .$
We define the expected gain ${\stackrel{^}{e}}_{i}$ in function value when we evaluate at a new point obtained by changing coordinate $i$ in the basepoint, for each $i$, based on two cases:
1. (i)${n}_{i}=0$. We compute the expected gain as
 $e^i = min 1≤j≤ Li {fij} - f i pi .$
Again, we split according to the initialization list, with the new basepoints and opposite points being as before.
2. (ii)${n}_{i}>0$. Now, the $i$th component of our sub-box ranges from ${x}_{i}$ to ${y}_{i}$. Using the quadratic partial correction function
 $ei (ξi) = αi (ξi-xi) + βi (ξi-xi) 2$
we can approximate the maximal gain expected when changing ${x}_{i}$ only. We will choose the splitting value from $▯\left\{{\xi }^{\prime },{\xi }^{\prime \prime }\right\}$. We compute
 $e^i = min ξi ∈ ▯ {ξ′,ξ′′} ei (ξi)$
and call ${z}_{i}$ the minimizer in $▯\left\{{\xi }^{\prime },{\xi }^{\prime \prime }\right\}$.
If the expected best function value ${f}_{\mathrm{exp}}$ satisfies
 $fexp = F(x) + min 1≤i≤n e^i < fbest ,$ (1)
where ${f}_{\mathrm{best}}$ 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 ${\stackrel{^}{e}}_{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 ${z}_{i}$. The sub-box is split at ${z}_{i}$ as long as ${z}_{i}\ne {y}_{i}$, and at the golden-section split point; two or three children are produced. The larger fraction of the golden-section split receives level $s+1$, while the smaller fraction receives level $\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left\{s+2,{s}_{\mathrm{max}}\right\}$. If it is the case that ${z}_{i}\ne {y}_{i}$ and the third child is larger than the smaller of the two children from the golden-section split, the third child receives level $s+1$. Otherwise it is given the level $\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left\{s+2,{s}_{\mathrm{max}}\right\}$. The basepoint of the first child is set to $\mathbf{x}$, and the basepoint of the second (and third if it exists) is obtained by changing the $i$th coordinate of $\mathbf{x}$ to ${z}_{i}$. The opposite points are again derived by changing ${y}_{i}$ to the other end of the $i$th coordinate interval of $B$.
If equation (1) does not hold, we expect no improvement. We do not split, and we increase the level of $B$ by $1$.

### 11.3Local Search

The local optimization algorithm used by e05kbf 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 $F$. 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 $F$ 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 $r$ 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 $\mathit{loclim}$ visits to this part of the local search have occurred, where $\mathit{loclim}$ is the value of the optional parameter MCS 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 MCS Max Objective Calls). The final criterion checks if no improvement can be made to the function value, or whether the approximated gradient $\mathbf{g}$ is small, in the sense that
 $|g|T max(|x|,| x old |) < loctol (f0-f) .$
The vector ${\mathbf{x}}_{\mathrm{old}}$ is the best point at the start of the current loop in this iterative local-search procedure, the constant $\mathit{loctol}$ is the value of the optional parameter MCS Local Searches Tolerance, $f$ is the objective value at $\mathbf{x}$, and ${f}_{0}$ is the smallest function value found by the initialization procedure.
Next, e05kbf 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 $r$ 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 e05kbf uses a general nonlinear optimizer to minimize the models.

## 12Optional Parameters

Several optional parameters in e05kbf define choices in the problem specification or the algorithm logic. In order to reduce the number of formal arguments of e05kbf 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 optional parameters can be changed by calling e04zmf anytime between the initialization of the handle and the call to the solver. Modification of the optional parameters during intermediate monitoring stops is not allowed. Once the solver finishes, the optional parameters can be altered again for the next solve.
The option values may be retrieved by e04znf.
The following is a list of the optional parameters available. A full description of each optional parameter is provided in Section 12.1.

### 12.1Description 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:
• the keywords;
• a parameter value, where the letters $a$, $i$ and $r$ denote options that take character, integer and real values respectively;
• the default value, where the symbol $\epsilon$ is a generic notation for machine precision (see x02ajf).
All options accept the value $\mathrm{DEFAULT}$ to return single options to their default states.
Keywords and character values are case and white space insensitive.
 Defaults
This special keyword may be used to reset all optional parameters to their default values. Any value given with this keyword will be ignored.
 MCS Initialization Method $a$ Default $\text{}=\mathrm{SIMPLE BOUNDS}$
The solver provides several methods to produce the initial set of sub-boxes. When ${\mathbf{MCS Initialization Method}}=\mathrm{SIMPLE BOUNDS}$, $\mathrm{SIMPLE OFF-BOUNDS}$, $\mathrm{LINESEARCH}$ or $\mathrm{RANDOM}$, e05kbf will generate the sub-boxes internally. Please refer to Section 11.1 for more details about the differences between the automatic methods.
It is also possible to create a custom initialization list by setting ${\mathbf{MCS Initialization Method}}=\mathrm{CUSTOM}$. Section 9.2 describes how to provide the custom initialization list to the solver.
Constraint: ${\mathbf{MCS Initialization Method}}=\mathrm{SIMPLE BOUNDS}$, $\mathrm{SIMPLE OFF-BOUNDS}$, $\mathrm{LINESEARCH}$, $\mathrm{RANDOM}$ or $\mathrm{CUSTOM}$
 MCS Max Objective Calls $i$ Default $\text{}=100{n}_{r}^{2}$
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 MCS Max Objective Calls may be made.
Constraint: ${\mathbf{MCS Max Objective Calls}}\ge 1$
 MCS Monitor Frequency $i$ Default $\text{}=0$
This parameter specifies the frequency on which to call the monitor routine monit. If zero, the monitor routine will not be called.
Constraint: ${\mathbf{MCS Monitor Frequency}}\ge 0$.
 MCS Print Frequency $i$ Default $\text{}=1$
This parameter specifies the frequency with which to print information regarding each sweep to Print File and/or Monitoring File. By default, it will print information of every sweep.
Constraint: ${\mathbf{MCS Print Frequency}}\ge 1$.
 MCS Local Searches $a$ Default $\text{}=\mathrm{ON}$
If you want to try to accelerate convergence of e05kbf by starting local searches from candidate minima, you will require MCS Local Searches to be $\mathrm{ON}$.
Constraint: ${\mathbf{MCS Local Searches}}=\mathrm{ON}$ or $\mathrm{OFF}$
 MCS Local Searches Limit $i$ Default $\text{}=50$
This defines the maximal number of iterations to be used in the trust region loop of the local-search procedure.
Constraint: ${\mathbf{MCS Local Searches Limit}}>0$
 MCS Local Searches Tolerance $r$ Default $\text{}=2\epsilon$
The value of MCS Local Searches Tolerance is the multiplier used during local searches as a stopping criterion for when the approximated gradient is small, in the sense described in Section 11.3.
Constraint: ${\mathbf{MCS Local Searches Tolerance}}\ge 2\epsilon$
 MCS Repeatability $a$ Default $\text{}=\mathrm{OFF}$
For use with random initialization lists (${\mathbf{MCS Initialization Method}}=RANDOM$). When set to $\mathrm{ON}$, The initialization list will be consistent for every call of e05kbf.
Constraint: ${\mathbf{MCS Repeatability}}=\mathrm{ON}$ or $\mathrm{OFF}$
 MCS Splits Limit $i$ Default $\text{}=⌊d\left({n}_{r}+2\right)/3⌋$
Along with the initialization list, this defines a limit on the number of times the root box will be split along any single coordinate direction. If MCS Local Searches is $\mathrm{OFF}$ you may find the default value to be too small.
Constraint: ${\mathbf{MCS Splits Limit}}>{n}_{r}+2$
 MCS Static Limit $i$ Default $\text{}=3{n}_{r}$
As the default termination criterion, computation stops when the best function value is static for MCS Static Limit sweeps through levels. This parameter is ignored if you have specified a target value to reach in MCS Target Objective Value.
Constraint: ${\mathbf{MCS Static Limit}}>0$
 MCS Target Objective Error $r$ Default $\text{}={\epsilon }^{\frac{1}{4}}$
If you have given a target objective value to reach in MCS Target Objective Value, this optional parameter sets your desired relative error between F(x) and MCS Target Objective Value, as described in Section 7. See also the description of the optional parameter MCS Target Objective Safeguard.
Constraint: ${\mathbf{MCS Target Objective Error}}\ge 2\epsilon$
 MCS Target Objective Safeguard $r$ Default $\text{}={\epsilon }^{\frac{1}{2}}$
If you have given a target objective value to reach in MCS Target Objective Value (the value of the optional parameter MCS Target Objective Value), this optional parameter sets your desired safeguarded termination tolerance, for when MCS Target Objective Value is close to zero.
Constraint: ${\mathbf{MCS Target Objective Safeguard}}\ge 2\epsilon$
 MCS Target Objective Value $r$
This parameter may be set if you wish e05kbf to use a specific value as the target function value to reach during the optimization. Setting MCS Target Objective Value overrides the default termination criterion determined by the optional parameter MCS Target Objective Value.
 Infinite Bound Size $r$ Default $\text{}={r}_{\mathrm{max}}^{\frac{1}{4}}$
This defines the ‘infinite’ bound $\mathit{bigbnd}$ in the definition of the problem constraints. Any upper bound greater than or equal to $\mathit{bigbnd}$ will be regarded as $+\infty$ (and similarly any lower bound less than or equal to $-\mathit{bigbnd}$ will be regarded as $-\infty$). Note that a modification of this optional parameter does not influence constraints which have already been defined; only the constraints formulated after the change will be affected.
Constraint: ${\mathbf{Infinite Bound Size}}\ge 1000$.
 Monitoring File $i$ Default $\text{}=-1$
If $i\ge 0$, the unit number for the secondary (monitoring) output. If ${\mathbf{Monitoring File}}=-1$, no secondary output is provided. The information output to this unit is controlled by Monitoring Level.
Constraint: ${\mathbf{Monitoring File}}\ge -1$.
 Monitoring Level $i$ Default $=4$
This parameter sets the amount of information detail that will be printed by the solver to the secondary output. The meaning of the levels is the same as with Print Level.
Constraint: $0\le {\mathbf{Monitoring Level}}\le 5$.
 Print File $i$ Default $=\text{advisory message unit number}$
If $i\ge 0$, the unit number for the primary output of the solver. If ${\mathbf{Print File}}=-1$, the primary output is completely turned off independently of other settings. The default value is the advisory message unit number as defined by x04abf at the time of the optional parameters initialization, e.g., at the initialization of the handle. The information output to this unit is controlled by Print Level.
Constraint: ${\mathbf{Print File}}\ge -1$.
 Print Level $i$ Default $=2$
This parameter defines how detailed information should be printed by the solver to the primary and secondary output.
$\mathbit{i}$ Output
$0$ No output from the solver.
$1$ The Header and Summary.
$2$, $3$, $4$, $5$ Additionally, the Iteration log.
Constraint: $0\le {\mathbf{Print Level}}\le 5$.
 Print Options $a$ Default $=\mathrm{YES}$
If ${\mathbf{Print Options}}=\mathrm{YES}$, a listing of optional parameters will be printed to the primary output and is always printed to the secondary output.
Constraint: ${\mathbf{Print Options}}=\mathrm{YES}$ or $\mathrm{NO}$.
 Print Solution $a$ Default $=\mathrm{NO}$
If ${\mathbf{Print Solution}}=\mathrm{YES}$, the solution will be printed to the primary and secondary output.
Constraint: ${\mathbf{Print Solution}}=\mathrm{YES}$ or $\mathrm{NO}$.
 Stats Time $a$ Default $=\mathrm{NO}$
This parameter turns on timings of various parts of the algorithm to give a better overview of where most of the time is spent. This might be helpful for a choice of different solving approaches. It is possible to choose between CPU and wall clock time. Choice $\mathrm{YES}$ is equivalent to $\mathrm{WALL CLOCK}$.
Constraint: ${\mathbf{Stats Time}}=\mathrm{YES}$, $\mathrm{NO}$, $\mathrm{CPU}$ or $\mathrm{WALL CLOCK}$.
 Task $a$ Default $=\mathrm{MINIMIZE}$
This parameter specifies the required direction of the optimization. If , the objective function (if set) is ignored and the algorithm stops as soon as a feasible point is found with respect to the given tolerance.
Constraint: ${\mathbf{Task}}=\mathrm{MINIMIZE}$, $\mathrm{MAXIMIZE}$ or $\mathrm{FEASIBLE POINT}$.


text:
#### Tagged as

, 4 Feb 2014 CPOL
 Rate this:
You can subclass Task - but it is tricky. Here's how you do it by "tying the knot" - using lambdas and variable capture to implement lazy evaluation

## Introduction

Perhaps you have some calculation or other operation that you want to run asynchronously, in your C# program, and you want to use a Task1. The expected way to do this is to just provide the task factory with the code you want to run, frequently in the form of a lambda:

    Task t = Task.Factory.StartNew( () => MyLongAsynchronousCalculation() );

That works well if the lambda is pretty much self-contained, or uses resources from the object that's creating the asynchronous activity.

But perhaps the asynchronous activity is somewhat complicated, and has an ongoing state, and, in other words, is the kind of thing you'd like to encapsulate in a class. There are two basic Object Oriented ways to approach this: Use subclassing, or use composition.

    I wouldn't personally extend Task<T>, I'd compose it instead.3

As soon as you try subclassing Task you discover the major problem: Task takes the code it is to run as an Action delegate parameter to its constructor and it doesn't provide any setter for that delegate that you can use after the Task is constructed—even if you're not starting the Task right away. In the constructor you can't make that delegate refer to the Task subclass instance you're creating.

    class D : Task {
public D() : base(??) { }
public void Run() { … }
}
…
D d = new D();

What do you put where the question marks are? You can't simply refer to D's method Run: It is not a static member so you need an object reference—and you don't have one. You can't use the this keyword there: it isn't allowed by the language. There's no good choice, so you're stuck?

## Using a value before it is computed: Tying the Knot

In a functional language all values are immutable. If the value has multiple fields they can't be modified. So how do you construct a data structure that refers to itself? This comes up in (for example) cyclic lists, or graph structures that aren't DAGs.

Consider, for example, representing rational numbers in the range [0..1) as a linked list of base 10 digits. The rational $$\frac{1}{8} = 0.125$$ is finite, so that's easy. But what about $$\frac{1}{7} = 0.\overline{142857}$$, which is a non-terminating repeating decimal fraction?9

In an imperative language it isn't hard since you can just clobber a field after it is created.

    var OneEighth = new SinglyLinkedList<int> { 1, 2, 5 };
Console.WriteLine("1/8 = " + string.Join(",", OneEighth.Take(30).Select(e => e.ToString())));
Console: 1/8 = 1,2,5

var OneSeventh = new SinglyLinkedList<int> { 1, 4, 2, 8, 5, 7 };
OneSeventh.Next.Next.Next.Next.Next.Next = OneSeventh;
Console.WriteLine("1/7 = " + string.Join(",", OneSeventh.Take(30).Select(e => e.ToString())));
Console: 1/7 = 1,4,2,8,5,7,1,4,2,8,5,7,1,4,2,8,5,7,1,4,2,8,5,7,1,4,2,8,5,7
    oneEighth = 1 : 2 : 5

oneSeventh = let x = 1 : 4 : 2 : 8 : 5 : 7 : x
in  x

Here, the name x refers to the list under construction and is also used to construct the tail of the list. It works because x refers to a memory location that isn't going to be referenced until some code uses the variable oneSeventh and traverses past the 6th element of the list. (Note the difference between a value and a variable: The variable is the location that can hold a value.

That is tying the knot!

## Tying the Knot in C#: closures and variable capture

Given that C# is not a lazy evaluation language12 how is tying the knot to be implemented?

The combination works like this:

1. Create a variable to provide a binding location for some value, but don't provide the value.
2. Create a lambda expression that closes over that variable, and returns the value in the variable. Since creating the lambda doesn't dereference the variable to get the value, it is fine.
3. Create the data structure, passing in the lambda to the creation routine, which will save it somewhere but not invoke it yet.
4. Store the data structure into the closed-over-variable.
5. Evaluate the lambda, which dereferences the variable to get the data structure, and does something with it (like store it in some field internal to itself.

To make this concrete, suppose we have a class T where the constructor takes an Action, called A, which it stores away in a readonly field, so that there is no way to (re-)set A after the instance of T is constructed. And then it has a method T.M which is called sometime later, after construction is finished, and invokes the Action A.

(And further suppose that we can't change T to provide a setter for A or anything else to fix up this situation.)

Now suppose that we subclass T with a class D and we want the Action A to run a method on itself, an instance of D. Normally our Action A would look something like this:

class D : T {
public D(Action a) : base(a) { }
…
Action f = () => this.Foo();
…
}

But that won't work because we don't have this yet, and actually we can't get a reference to our new instance D until after its base constructor, a member of T, runs and the D constructor returns.

The way to solve this is to use variable capture in creating the lambda expression:

    …
D d = null;			// 1: Create a variable to hold a D.
Action g = () => d.Foo();	// 2: Capture the variable in our lambda.
d = new D(g);			// 3: Create a new D, passing in our lambda,
//    and store it in the variable we captured.
d.M();				// 4: Execute the method that is going to
//    run our Action and call d.Foo().
…

The real trick is going to be to get an API that can accept an Action or Func instead of the object—an instance of some specific type—it is expecting. But fortunately, in our problem of deriving from Task, that problem is solved because Task takes two parameters: An Action and an arbitrary object—and we'll pass our delaying lambda in as the arbitrary object!

## Back to the problem: Subclassing the Task class

Back to the Task that's left to us: How do we subclass Task?

To answer the first question: The Action, as well as the derived instance, will be created in a factory method. The constructor of our derived class will have protected access so a developer can't create one directly.

To answer the second question: It would be easy enough to have our derived class have a constructor that, in addition to its other arguments for initializing itself, took the Action and immediately passed it into its base class:

    protected D(Action a) : base(a) { … }

And then this article would be over. But I prefer to not have to write this code more than once. I would like to provide a generic abstract class to do all the work. It will derive from Task and all my various subclasses will derive from it.

As soon as I do this, however, I run into a problem: my factory method will be generic in my most-derived subclass. That's so that, as a factory, it can return an instance of that subclass (instead of a superclass—its own type—that would need to be cast to the actual subclass). But if it is generic in the type that it is creating (and returning) then it can only create instances of types that have a zero-argument constructor (due to the way the generic constraint new() works). And with that zero-argument constructor, how can the Action be passed in?

The answer is somewhat unsatisfactory: It will be passed in via a static field that the constructor can reference. And that's our last problem to solve: Ensure that the factory method is serialized so that it is safe to set the static field and then immediately create a new object that refers to that field in its constructor, so that if the factory method is called on two different threads simultaneously there are no race conditions that would lead to one of the new instances getting the Action for the other instance.

Really, this is a complication that is only necessary if you want to have a resusable abstract base class to own all the code that handles passing the Action in to the Task.

Anyway, the code for the abstract generic class DeriveFromTaskBase is in the zip archive associated with this article, so I'll only comment on highlights here.

### The public API of DeriveFromTaskBase

The public API of DeriveFromTaskBase consists of the factory method Create that creates an instance and starts it, and an abstract method Run that must be overriden to provide the subclass-specific computation that is the entire purpose for subclassing Task.

The Create factory method takes an optional Action<T>, called beforeStartInitializer, which is run just before the instance is started. Its purpose is to provide a chance to initialize the instance and make up for the fact that there is only a zero-argument constructor. The Action<T> you provide will be given the instance itself and can set properties or run methods on that instance. (Remember that when you create the Action with a lambda expression you can capture any values you need at that time.) If you also, or alternatively, have things you can do at construction time (that, necessarily, can't rely on any outside inputs), you can (optionally) override the method Constructor and do that initialization.

    ///
/// Abstract base class for classes that want to derive from Task.
///
{
#region Public interface
public static T Create<T>(Action<T> beforeStartInitializer = null)
{
…
}

public virtual void Constructor() { }

public abstract void Run();
#endregion

### The construction of the derived instance

The constructor is fairly simple: referring to a static field that holds the indirection Action it simply passes that Action to the base Task and then calls the (optional) Constructor method.

        private static Action thisDeferred;

{
Constructor();
}

### The factory which "ties the knot" and creates the derived instance

Interestingly, there are two knots to be tied!

The factory method grabs a lock to ensure serialization. Then it provides a location for the to-be-created instance and provides a location for the true Action. It ties the first knot by creating the indirection-Action that invokes the true Action by capturing its location. It then performs new T() to finally create the derived instance you're really looking for. And it ties the second knot by creating the true Action that captures the location of the new instance. After all that it starts the instance—and the Task calls its start Action which calls the true Action which calls the instance's Run method, and at last the Task is going!

        private static object createLock = new object();
private static Action thisDeferred;

private static T Create<T>(Action<T> beforeStartInitializer = null)
{
T t = null;

lock (createLock)
{
Action thisDeferredInner = null;

thisDeferred = () => thisDeferredInner();

t = new T();

thisDeferredInner = () =>
{
if (null != beforeStartInitializer)
beforeStartInitializer(t);
t.Run();
};
}
t.Start();
return t;
}

(By now you shouldn't need comments to understand the above code…but no worries: There are comments in the sources that are in the zip.)

## Article Summary

Is it worth it? Well, that depends. For the particular case in hand, utilizing a Task to run complicated code that needs ongoing state, if you were starting from scratch then it would be easiest to write a class that, intead of deriving from Task, simply owns a Task instance. That is, use composition.

Then, also, in the general case (moving away from Tasks), knowing how to tie the knot as a technique in general, and the particular way you specifically do it in C#, with lambdas and variable capture, may prove useful to you (as well as interesting.) I hope so.

Oh, one more thing. Now that you understand the previous section you should be easily able to understand the diagram at the top of this article. What's that you say? You can't? Well, it's not your fault. Fact is, it isn't a very good diagram; it's just the best I could do. Please feel free to provide me with a better diagram (let me know in a comment) and I'll be glad to put it up to replace the one I came up with and give you the credit for it.

## Article Revision History

• 01-FEB-2014: Original article.

## Footnotes

2Jon Skeet's C# In Depth is an excellent book, as is Jon Skeet's blog.

3 See Jon Skeet's answer on SO.

6 A full description of Tying the Knot. The canonical algorithm for tying the knot is repmin, which is a one pass algorithm for building a tree which is the same shape as the given tree except each leaf has the minimum leaf value from the original tree.

12 Lazy evaluation (wikipedia) means that an expression is not evaluated when it is bound to a variable, but only when it is used. Nearly all languages use strict evaluation where expressions are fully evaluated when they are bound to a variable (or procedure argument). In fact, this is why some logging libraries sometimes go out of their way to use non-language mechanisms, like preprocessor macros in C/C++, to provide efficient processing of formatted messages and their arguments: They are trying to improve performance and reduce the overhead of logging by avoiding the evaluation of message arguments unless the message's log/trace level is high enough that the log message will actually be written to some sink.

    We can solve any problem by introducing an extra level of indirection..

14 There's some question about whether C# closures are "true" closures or not. There was quite a discussion on Wikipedia, and also at the programming language blog Lambda The Ultimate, with respect to the Wikipedia article on closures. Some people will accept nothing less than closures which "bind return" (that is, allow for call-with-current-continuation). Other people believe that the closures in C# 3.0 and Javascript are as close to "true" closures as makes no difference.

## Share

Architect Bakin's Bits
United States
I've been programming for 35 years and I'm still learning new things. My major interests recently are programming models and methods for concurrency, applications of functional programming, and doing interesting things with Mathematica.

 First Prev Next
 My vote of 3 Paulo Zemek 11-Mar-14 9:11
 [ThreadStatic] and more. Paulo Zemek 11-Mar-14 9:11
 What was the problem, again? John Brett 5-Feb-14 1:15
 My vote of 5 Klaus Luedenscheidt 4-Feb-14 20:06
 Very nice article ! Volynsky Alex 4-Feb-14 10:19
 Alternative solutions Henrik Jonsson 4-Feb-14 9:42
 Last Visit: 31-Dec-99 19:00     Last Update: 31-Jan-15 14:27 Refresh 1

format: plain


text:
### Rigid Body

A ball of the radius r rotates with the angular velocity $\omega$ around the horizontal axis passing through the center of mass. The ball is carefully put on a horizontal surface with the friction coefficient $\mu$. Find $v_{cm}(t)$.


text:
Data determined by PrimitiveRecoverySchemes at a single grid point. More...

#include <PrimitiveRecoveryData.hpp>

## Public Attributes

double rest_mass_density

double lorentz_factor

double pressure

double rho_h_w_squared

## Detailed Description

Data determined by PrimitiveRecoverySchemes at a single grid point.

rho_h_w_squared is $$\rho h W^2$$ where $$\rho$$ is the rest mass density, $$h$$ is the specific enthalpy, and $$W$$ is the Lorentz factor.

The documentation for this struct was generated from the following file:


text:
# How do you write the partial fraction decomposition of the rational expression (x^2 + 1)/(x(x + 1)(x^2 + 2))?

Dec 14, 2015

$\frac{1}{2 x} - \frac{2}{3 \left(x + 1\right)} + \frac{x + 2}{6 \left({x}^{2} + 2\right)}$

#### Explanation:

$\frac{{x}^{2} + 1}{x \left(x + 1\right) \left({x}^{2} + 2\right)} = \frac{A}{x} + \frac{B}{x + 1} + \frac{C x + D}{{x}^{2} + 2}$

${x}^{2} + 1 = A \left(x + 1\right) \left({x}^{2} + 2\right) + B \left(x\right) \left({x}^{2} + 2\right) + \left(C x + D\right) \left(x\right) \left(x + 1\right)$

${x}^{2} + 1 = A {x}^{3} + A {x}^{2} + 2 A x + 2 A + B {x}^{3} + 2 B x + C {x}^{3} + C {x}^{2} + D {x}^{2} + D x$

${x}^{2} + 1 = {x}^{3} \left(A + B + C\right) + {x}^{2} \left(A + C + D\right) + x \left(2 A + 2 B + D\right) + 1 \left(2 A\right)$

Thus, $\left\{\begin{matrix}A + B + C = 0 \\ A + C + D = 1 \\ 2 A + 2 B + D = 0 \\ 2 A = 1\end{matrix}\right.$

Solve this to see that $\left\{\begin{matrix}A = \frac{1}{2} \\ B = - \frac{2}{3} \\ C = \frac{1}{6} \\ D = \frac{1}{3}\end{matrix}\right.$

Therefore,

$\frac{{x}^{2} + 1}{x \left(x + 1\right) \left({x}^{2} + 2\right)} = \frac{1}{2 x} - \frac{2}{3 \left(x + 1\right)} + \frac{x + 2}{6 \left({x}^{2} + 2\right)}$


text:
# Scraping Instagram for Hashtag data

Over the past few months, I've been actively using python and I have made a few scripts to scrape #hashtag data from Instagram.

It all started with some basic script I had made early 2017 and I have been adding and modifying it ever since. Over the last few months, I made progress in my own skill of Python, successfully adding things like user agent and proxy rotation.

Now that I have a tool that does exactly what I want, I'm looking to:

• Optimize code structure (it's really copying and pasting mostly) and removing 'crappy' code.

Therefore I'm hoping SO can help me analyze my code and suggest optimizations.

My script does the following:

• It analyzes hashtags from the input file (hashtags.txt)
• It then scrapes data from Instagram (like post count, average engagement,...)
• This data is then stored in a .csv. Which is being processed again afterward to remove duplicates.

I also included user agent randomization and proxy rotation.

However, I feel like my code is far from optimal and when I want to add additional things (like catching HTTP errors, retrying on proxy timeouts,...) I'm just adding more levels of indentation so I'm pretty sure there are other options there!

Any help or feedback to optimize my code below is GREATLY appreciated!

    # This script is written for personal research and is not endorsed by Instagram.
# Use at your own risk!
# -*- coding: utf-8 -*-
import csv
import requests
import urllib.request
import json
import re
import random
import time
from fake_useragent import UserAgent
from random import randint
from time import sleep

ua = UserAgent(cache=False)
ts = time.gmtime()
timestamp = time.strftime("%d-%m-%Y %H-%M", ts)

fieldnames = ['Hashtag','Active Days Ago','Post Count','AVG. Likes','MAX. Likes','MIN. Likes','AVG. Comments','Hashtag URL','Post Ready Tag']
return fieldnames

with open(filename, 'w', newline='') as f_out:
return

with open(t_file) as f:
return keyword_list

with open(p_file) as f:
return proxy_list

#file
data_filename = 'Hashtag Scrape ' + timestamp + '.csv'
KEYWORD_FILE = './hashtags.txt'
DATA_FILE = './' + data_filename
PROXY_FILE = './proxies.txt'

low = input("Please enter minimal delay time (in seconds): ")
low_random = int(low)
high = input("Please enter maximal delay time (in seconds): ")
high_random = int(high)

#get the data
for keyword in keywords:
import urllib, json
if len(proxies)!=0:
proxy_ip = random.choice(proxies)
proxy_support = urllib.request.ProxyHandler({'https':proxy_ip})
opener = urllib.request.build_opener(proxy_support)
urllib.request.install_opener(opener)
prepare_url = urllib.request.Request(
'https://www.instagram.com/explore/tags/' + urllib.parse.quote_plus(keyword) + '/?__a=1',
'User-Agent': ua.random
}
)
url = urllib.request.urlopen(prepare_url)
post_info = {}
response = json.load(url) #response is the JSON dump of the url.

#defining some script helpers
x = len(response['graphql']['hashtag']['edge_hashtag_to_top_posts']['edges'])
i = avg_post_likes = 0
likes_value = []

#Getting the general tag data
hashtag_name = response['graphql']['hashtag']['name']
post_count = response['graphql']['hashtag']['edge_hashtag_to_media']['count']
hashtag_url = 'https://www.instagram.com/explore/tags/' + keyword
top_posts = response['graphql']['hashtag']['edge_hashtag_to_top_posts']['edges']

#calculate the active days ago
most_recent_post = response['graphql']['hashtag']['edge_hashtag_to_media']['edges'][0]['node']['taken_at_timestamp']
import datetime
from dateutil import relativedelta
post_datetime = datetime.datetime.fromtimestamp(most_recent_post).strftime('%Y-%m-%d %H:%M:%S')
post_cleandate = datetime.datetime.fromtimestamp(most_recent_post).strftime('%Y-%m-%d')
from datetime import datetime, date
most_recent_clean = datetime.strptime(post_cleandate, '%Y-%m-%d')
today = datetime.strptime(str(date.today()),'%Y-%m-%d')
posted_days_ago = relativedelta.relativedelta(today, most_recent_clean).days

while i <=x-1:
#Getting data from top posts
top_post_likes = response['graphql']['hashtag']['edge_hashtag_to_top_posts']['edges'][i]['node']['edge_liked_by']
post_like = response['graphql']['hashtag']['edge_hashtag_to_top_posts']['edges'][i]['node']['edge_liked_by']['count']
post_comment = response['graphql']['hashtag']['edge_hashtag_to_top_posts']['edges'][i]['node']['edge_media_to_comment']['count']
likes_value.append(post_like)
i += 1
print('Writing ' + keyword + ' to output file')
with open(data_filename, 'a', newline='',  encoding='utf-8') as data_out:
post_info["Hashtag"] = hashtag_name
post_info["Active Days Ago"] = posted_days_ago
post_info["Post Count"] = post_count
post_info["AVG. Likes"] = round(sum(likes_value)/len(likes_value),2)
post_info["MAX. Likes"] = max(likes_value)
post_info["MIN. Likes"] = min(likes_value)
post_info["Hashtag URL"] = hashtag_url
csv_writer.writerow(post_info)

#Randomly pause script based on input values
sleep(randint(low_random,high_random))
#cleaning up the file:
destination = data_filename[:-4] + '_unique.csv'
data = open(data_filename, 'r',encoding='utf-8')
target = open(destination, 'w',encoding='utf-8')
# Let the user know you are starting, in case you are de-dupping a huge file
print("\nRemoving duplicates from %r" % data_filename)

# Initialize variables and counters
unique_lines = set()
source_lines = 0
duplicate_lines = 0

# Loop through data, write uniques to output file, skip duplicates.
for line in data:
source_lines += 1
# Strip out the junk for an easy set check, also saves memory
line_to_check = line.strip('\r\n')
if line_to_check in unique_lines: # Skip if line is already in set
duplicate_lines += 1
continue
else: # Write if new and append stripped line to list of seen lines
target.write(line)
# Be nice and close out the files
target.close()
data.close()
import os
os.remove(data_filename)
os.rename(destination, data_filename)
print("SUCCESS: Removed %d duplicate line(s) from file with %d line(s)." % \
(duplicate_lines, source_lines))
print("Wrote output to %r\n" % data_filename)
print("\n" + 'ALL DONE !!!! ')

For those interested, this is how the output file looks:

• Why aren't you using their API? instagram.com/developer – Reinderien Jan 4 '19 at 18:01
• Because I don't need to (yet), and don't really want to either :). I can get this data without using tokens and login credentials. So that's my preferred approach. – ThomasSt Jan 4 '19 at 18:11
• I coded something as similar as yours: codereview.stackexchange.com/questions/210613/… If you're interested in more web scraping. – austingae Jan 4 '19 at 19:07
• Please do not update the code in your question to incorporate feedback from answers, doing so goes against the Question + Answer style of Code Review. This is not a forum where you should keep the most updated version in your question. Please see what you may and may not do after receiving answers. – Zeta Jan 5 '19 at 8:26

This function:

def get_csv_header(top_numb):
fieldnames = ['Hashtag','Active Days Ago','Post Count','AVG. Likes','MAX. Likes','MIN. Likes','AVG. Comments','Hashtag URL','Post Ready Tag']
return fieldnames

has a few issues. top_numb is unused, so delete it. You can both construct and return the list in the same statement, but due to its length I suggest that you add some linebreaks in that list. Finally: per Python 3 docs, fieldnames must be a sequence but needn't be a list - so make this a tuple () and not a list [] because the data are immutable.

Otherwise:

## Remove redundant returns

i.e. the no-op return seen in write_csv_header.

## Make a main function

...for all of your global code, for a couple of reasons - to clean up the global namespace, and to make your code callable as a library for other applications.

## Use f-strings

...for strings like this:

data_filename = 'Hashtag Scrape ' + timestamp + '.csv'

that can be:

data_filename = f'Hashtag Scrape {timestamp}.csv'

## Write more subroutines

The bulk of your logic within the main for keyword in keywords loop is quite long. Break this up into several subroutines for legibility and maintainability.

## Use requests

You're calling into urllib.request.Request, but there's usually no good reason to do this. Use requests instead, which is better in nearly every way.

## Apply a linter

This will catch non-PEP8 whitespace (or lack thereof) such as that seen in this statement:

if len(proxies)!=0:

## Imports at the top

In the middle of your source, we see:

import datetime
from dateutil import relativedelta
post_datetime = datetime.datetime.fromtimestamp(most_recent_post).strftime('%Y-%m-%d %H:%M:%S')
post_cleandate = datetime.datetime.fromtimestamp(most_recent_post).strftime('%Y-%m-%d')
from datetime import datetime, date

It's usually considered better practice to do all of your imports at the top of the source file.

## Don't declare indices that you don't use

This loop:

i = avg_post_likes = 0
while i <=x-1:
# ...
i += 1

should be

for _ in range(x):
# ...

You also need a better name for x.

## Use dict.update

This code:

        post_info["Hashtag"] = hashtag_name
post_info["Active Days Ago"] = posted_days_ago
post_info["Post Count"] = post_count
post_info["AVG. Likes"] = round(sum(likes_value)/len(likes_value),2)
post_info["MAX. Likes"] = max(likes_value)
post_info["MIN. Likes"] = min(likes_value)
post_info["Hashtag URL"] = hashtag_url

can be greatly simplified by use of update:

post_info.update({
'Hashtag': hashtag_name,
'Active Days Ago': posted_days_ago,
# ...

## Use context management

You were doing so well elsewhere in the file! But then we see this:

data = open(data_filename, 'r',encoding='utf-8')
target = open(destination, 'w',encoding='utf-8')

Those should also use with. You can keep the indentation from getting out-of-control by writing more subroutines.