GALAHAD UGO package#

purpose#

The ugo package aims to find the global minimizer of a univariate twice-continuously differentiable function $f(x)$ of a single variable over the finite interval $x^l <= x <= x^u$. Function and derivative values are provided via a subroutine call. Second derivatives may be used to advantage if they are available.

See Section 4 of $GALAHAD/doc/ugo.pdf for additional details.

method#

The algorithm starts by splitting the interval $[x^l,x^u]$ into a specified number of subintervals $[x_i,x_{i+1}]$ of equal length, and evaluating $f$ and its derivatives at each $x_i$. A surrogate (approximating) lower bound function is constructed on each subinterval using the function and derivative values at each end, and an estimate of the first- and second-derivative Lipschitz constant. This surrogate is minimized, the true objective evaluated at the best predicted point, and the corresponding interval split again at this point. Any interval whose surrogate lower bound value exceeds an evaluated actual value is discarded. The method continues until only one interval of a maximum permitted width remains.

reference#

Many ingredients in the algorithm are based on the paper

D. Lera and Ya. D. Sergeyev, “Acceleration of univariate global optimization algorithms working with Lipschitz functions and Lipschitz first derivatives” SIAM J. Optimization 23(1), (2013) 508–529

but adapted to use second derivatives.

introduction to function calls#

To solve a given problem, functions from the ugo package must be called in the following order:

ugo_initialize - provide default control parameters and set up initial data structures
ugo_read_specfile (optional) - override control values by reading replacement values from a file
ugo_import - set up problem data structures and fixed values
ugo_reset_control (optional) - possibly change control parameters if a sequence of problems are being solved
solve the problem by calling one of
- ugo_solve_direct - solve using function calls to evaluate function and derivative values, or
- ugo_solve_reverse - solve returning to the calling program to obtain function and derivative values
ugo_information (optional) - recover information about the solution and solution process
ugo_terminate - deallocate data structures

See the examples section for illustrations of use.

parametric real type T and integer type INT#

Below, the symbol T refers to a parametric real type that may be Float32 (single precision), Float64 (double precision) or, if supported, Float128 (quadruple precision). The symbol INT refers to a parametric integer type that may be Int32 (32-bit integer) or Int64 (64-bit integer).

callable functions#

    function ugo_initialize(T, INT, data, control, status)

Set default control values and initialize private data

Parameters:

data

holds private internal data

control

is a structure containing control information (see ugo_control_type)

status

is a scalar variable of type INT, that gives the exit status from the package. Possible values are (currently):

0
The initialization was successful.

    function ugo_read_specfile(T, INT, control, specfile)

Read the content of a specification file, and assign values associated with given keywords to the corresponding control parameters. An in-depth discussion of specification files is available, and a detailed list of keywords with associated default values is provided in $GALAHAD/src/ugo/UGO.template. See also Table 2.1 in the Fortran documentation provided in $GALAHAD/doc/ugo.pdf for a list of how these keywords relate to the components of the control structure.

Parameters:

control	is a structure containing control information (see ugo_control_type)
specfile	is a one-dimensional array of type Vararg{Cchar} that must give the name of the specification file

    function ugo_import(T, INT, control, data, status, x_l, x_u)

Import problem data into internal storage prior to solution.

Parameters:

control	is a structure whose members provide control parameters for the remaining procedures (see ugo_control_type)
data	holds private internal data
status	is a scalar variable of type INT, that gives the exit status from the package. Possible values are: 1 The import was successful, and the package is ready for the solve phase -1 An allocation error occurred. A message indicating the offending array is written on unit control.error, and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -2 A deallocation error occurred. A message indicating the offending array is written on unit control.error and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively.
x_l	is a scalar variable of type T, that holds the value $x^l$ of the lower bound on the optimization variable $x$.
x_u	is a scalar variable of type T, that holds the value $x^u$ of the upper bound on the optimization variable $x$.

    function ugo_reset_control(T, INT, control, data, status)

Reset control parameters after import if required.

Parameters:

control

is a structure whose members provide control parameters for the remaining procedures (see ugo_control_type)

data

holds private internal data

status

is a scalar variable of type INT, that gives the exit status from the package. Possible values are:

1
The import was successful, and the package is ready for the solve phase

    function ugo_solve_direct(T, INT, data, userdata, status, x, f, g, h, eval_fgh)

Find an approximation to the global minimizer of a given univariate function with a Lipschitz gradient in an interval.

This version is for the case where all function/derivative information is available by function calls.

Parameters:

data	holds private internal data
userdata	is a structure that allows data to be passed into the function and derivative evaluation programs (see below).
status	is a scalar variable of type INT, that gives the entry and exit status from the package. On initial entry, status must be set to 1. Possible exit values are: 0 The run was successful -1 An allocation error occurred. A message indicating the offending array is written on unit control.error, and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -2 A deallocation error occurred. A message indicating the offending array is written on unit control.error and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -7 The objective function appears to be unbounded from below -18 Too many iterations have been performed. This may happen if control.maxit is too small, but may also be symptomatic of a badly scaled problem. -19 The CPU time limit has been reached. This may happen if control.cpu_time_limit is too small, but may also be symptomatic of a badly scaled problem. -40 The user has forced termination of solver by removing the file named control.alive_file from unit unit control.alive_unit.
x	is a scalar variable of type T, that holds the value of the approximate global minimizer $x$ after a successful (status = 0) call.
f	is a scalar variable of type T, that holds the the value of the objective function $f(x)$ at the approximate global minimizer $x$ after a successful (status = 0) call.
g	is a scalar variable of type T, that holds the the value of the gradient of the objective function $f^{\prime}(x)$ at the approximate global minimizer $x$ after a successful (status = 0) call.
h	is a scalar variable of type T, that holds the the value of the second derivative of the objective function $f^{\prime\prime}(x)$ at the approximate global minimizer $x$ after a successful (status = 0) call.
eval_fgh	is a user-provided function that must have the following signature: function eval_fgh(x, f, g, h, userdata) The value of the objective function $f(x)$ and its first derivative $f^{\prime}(x)$ evaluated at x= $x$ must be assigned to f and g respectively, and the function return value set to 0. In addition, if control.second_derivatives_available has been set to true, when calling ugo_import, the user must also assign the value of the second derivative $f^{\prime\prime}(x)$ in h; it need not be assigned otherwise. If the evaluation is impossible at x, return should be set to a nonzero value.

    function ugo_solve_reverse(T, INT, data, status, eval_status, x, f, g, h)

Find an approximation to the global minimizer of a given univariate function with a Lipschitz gradient in an interval.

This version is for the case where function/derivative information is only available by returning to the calling procedure.

Parameters:

data	holds private internal data
status	is a scalar variable of type INT, that gives the entry and exit status from the package. On initial entry, status must be set to 1. Possible exit values are: 0 The run was successful -1 An allocation error occurred. A message indicating the offending array is written on unit control.error, and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -2 A deallocation error occurred. A message indicating the offending array is written on unit control.error and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -7 The objective function appears to be unbounded from below -18 Too many iterations have been performed. This may happen if control.maxit is too small, but may also be symptomatic of a badly scaled problem. -19 The CPU time limit has been reached. This may happen if control.cpu_time_limit is too small, but may also be symptomatic of a badly scaled problem. -40 The user has forced termination of solver by removing the file named control.alive_file from unit unit control.alive_unit. 3 The user should compute the objective function value $f(x)$ and its first derivative $f^{\prime}(x)$, and then re-enter the function. The required values should be set in f and g respectively, and eval_status (below) should be set to 0. If the user is unable to evaluate $f(x)$ or $f^{\prime}(x)$ - for instance, if the function or its first derivative are undefined at x - the user need not set f or g, but should then set eval_status to a non-zero value. This value can only occur when control.second_derivatives_available = false. 4 The user should compute the objective function value $f(x)$ and its first two derivatives $f^{\prime}(x)$ and $f^{\prime\prime}(x)$ at x= $x$, and then re-enter the function. The required values should be set in f, g and h respectively, and eval_status (below) should be set to 0. If the user is unable to evaluate $f(x)$, $f^{\prime}(x)$ or $f^{\prime\prime}(x)$ - for instance, if the function or its derivatives are undefined at x - the user need not set f, g or h, but should then set eval_status to a non-zero value. This value can only occur when control.second_derivatives_available = true.
eval_status	is a scalar variable of type INT, that is used to indicate if objective function and its derivatives can be provided (see above).
x	is a scalar variable of type T, that holds the next value of $x$ at which the user is required to evaluate the objective (and its derivatives) when status > 0, or the value of the approximate global minimizer when status = 0
f	is a scalar variable of type T, that must be set by the user to hold the value of $f(x)$ if required by status > 0 (see above), and will return the value of the approximate global minimum when status = 0
g	is a scalar variable of type T, that must be set by the user to hold the value of $f^{\prime}(x)$ if required by status > 0 (see above), and will return the value of the first derivative of $f$ at the approximate global minimizer when status = 0
h	is a scalar variable of type T, that must be set by the user to hold the value of $f^{\prime\prime}(x)$ if required by status > 0 (see above), and will return the value of the second derivative of $f$ at the approximate global minimizer when status = 0

    function ugo_information(T, INT, data, inform, status)

Provides output information

Parameters:

data

holds private internal data

inform

is a structure containing output information (see ugo_inform_type)

status

is a scalar variable of type INT, that gives the exit status from the package. Possible values are (currently):

0
The values were recorded successfully

    function ugo_terminate(T, INT, data, control, inform)

Deallocate all internal private storage

Parameters:

data	holds private internal data
control	is a structure containing control information (see ugo_control_type)
inform	is a structure containing output information (see ugo_inform_type)

available structures#

ugo_control_type structure#

    struct ugo_control_type{T,INT}
      error::INT
      out::INT
      print_level::INT
      start_print::INT
      stop_print::INT
      print_gap::INT
      maxit::INT
      initial_points::INT
      storage_increment::INT
      buffer::INT
      lipschitz_estimate_used::INT
      next_interval_selection::INT
      refine_with_newton::INT
      alive_unit::INT
      alive_file::NTuple{31,Cchar}
      stop_length::T
      small_g_for_newton::T
      small_g::T
      obj_sufficient::T
      global_lipschitz_constant::T
      reliability_parameter::T
      lipschitz_lower_bound::T
      cpu_time_limit::T
      clock_time_limit::T
      second_derivative_available::Bool
      space_critical::Bool
      deallocate_error_fatal::Bool
      prefix::NTuple{31,Cchar}

detailed documentation#

control derived type as a Julia structure

components#

INT error

error and warning diagnostics occur on stream error

INT out

general output occurs on stream out

INT print_level

the level of output required. Possible values are:

$\leq$ 0 no output,
1 a one-line summary for every improvement
2 a summary of each iteration
$\geq$ 3 increasingly verbose (debugging) output

INT start_print

any printing will start on this iteration

INT stop_print

any printing will stop on this iteration

INT print_gap

the number of iterations between printing

INT maxit

the maximum number of iterations allowed

INT initial_points

the number of initial (uniformly-spaced) evaluation points (<2 reset to 2)

INT storage_increment

incremenets of storage allocated (less that 1000 will be reset to 1000)

INT buffer

unit for any out-of-core writing when expanding arrays

INT lipschitz_estimate_used

what sort of Lipschitz constant estimate will be used:

1 = global contant provided
2 = global contant estimated
3 = local costants estimated

INT next_interval_selection

how is the next interval for examination chosen:

1 = traditional
2 = local_improvement

INT refine_with_newton

try refine_with_newton Newton steps from the vacinity of the global minimizer to try to improve the estimate

INT alive_unit

removal of the file alive_file from unit alive_unit terminates execution

NTuple{31,Cchar} alive_file

see alive_unit

T stop_length

overall convergence tolerances. The iteration will terminate when the step is less than .stop_length

T small_g_for_newton

if the absolute value of the gradient is smaller than small_g_for_newton, the next evaluation point may be at a Newton estimate of a local minimizer

T small_g

if the absolute value of the gradient at the end of the interval search is smaller than small_g, no Newton serach is necessary

T obj_sufficient

stop if the objective function is smaller than a specified value

T global_lipschitz_constant

the global Lipschitz constant for the gradient (-ve means unknown)

T reliability_parameter

the reliability parameter that is used to boost insufficiently large estimates of the Lipschitz constant (-ve means that default values will be chosen depending on whether second derivatives are provided or not)

T lipschitz_lower_bound

a lower bound on the Lipscitz constant for the gradient (not zero unless the function is constant)

T cpu_time_limit

the maximum CPU time allowed (-ve means infinite)

T clock_time_limit

the maximum elapsed clock time allowed (-ve means infinite)

Bool second_derivative_available

if .second_derivative_available is true, the user must provide them when requested. The package is generally more effective if second derivatives are available.

Bool space_critical

if .space_critical is true, every effort will be made to use as little space as possible. This may result in longer computation time

Bool deallocate_error_fatal

if .deallocate_error_fatal is true, any array/pointer deallocation error will terminate execution. Otherwise, computation will continue

NTuple{31,Cchar} prefix

all output lines will be prefixed by .prefix(2:LEN(TRIM(.prefix))-1) where .prefix contains the required string enclosed in quotes, e.g. “string” or ‘string’

ugo_time_type structure#

    struct ugo_time_type{T}
      total::Float32
      clock_total::T

detailed documentation#

time derived type as a Julia structure

components#

Float32 total

the total CPU time spent in the package

T clock_total

the total clock time spent in the package

ugo_inform_type structure#

    struct ugo_inform_type{T,INT}
      status::INT
      eval_status::INT
      alloc_status::INT
      bad_alloc::NTuple{81,Cchar}
      iter::INT
      f_eval::INT
      g_eval::INT
      h_eval::INT
      time::ugo_time_type{T}

detailed documentation#

inform derived type as a Julia structure

components#

INT status

return status. See UGO_solve for details

INT eval_status

evaluation status for reverse communication interface

INT alloc_status

the status of the last attempted allocation/deallocation

NTuple{81,Cchar) bad_alloc

the name of the array for which an allocation/deallocation error occurred

INT iter

the total number of iterations performed

INT f_eval

the total number of evaluations of the objective function

INT g_eval

the total number of evaluations of the gradient of the objective function

INT h_eval

the total number of evaluations of the Hessian of the objective function

ugo_time_type{T} ugo_time_type time

timings (see above)

Index

example calls#

This is an example of how to use the package to minimize a univariate function; the code is available in $GALAHAD/src/ugo/Julia/test_ugo.jl .

# test_ugo.jl
# Simple code to test the Julia interface to UGO

using GALAHAD
using Test
using Printf
using Accessors
using Quadmath

function test_ugo(::Type{T}, ::Type{INT}) where {T,INT}
  # Test problem objective
  function objf(x::T)
    a = 10.0
    res = x * x * cos(a * x)
    return res
  end

  # Test problem first derivative
  function gradf(x::T)
    a = 10.0
    res = -a * x * x * sin(a * x) + 2.0 * x * cos(a * x)
    return res
  end

  # Test problem second derivative
  function hessf(x::T)
    a = 10.0
    res = -a * a * x * x * cos(a * x) - 4.0 * a * x * sin(a * x) + 2.0 * cos(a * x)
    return res
  end

  # Derived types
  data = Ref{Ptr{Cvoid}}()
  control = Ref{ugo_control_type{T,INT}}()
  inform = Ref{ugo_inform_type{T,INT}}()

  # Initialize UGO
  status = Ref{INT}()
  eval_status = Ref{INT}()
  ugo_initialize(T, INT, data, control, status)

  # Set user-defined control options
  @reset control[].print_level = INT(1)

  # Test problem bounds
  x_l = Ref{T}(-1.0)
  x_u = Ref{T}(2.0)

  # Test problem objective, gradient, Hessian values
  x = Ref{T}(0.0)
  f = Ref{T}(objf(x[]))
  g = Ref{T}(gradf(x[]))
  h = Ref{T}(hessf(x[]))

  # import problem data
  ugo_import(T, INT, control, data, status, x_l, x_u)

  # Set for initial entry
  status[] = 1

  # Solve the problem: min f(x), x_l ≤ x ≤ x_u
  terminated = false
  while !terminated
    # Call UGO_solve
    ugo_solve_reverse(T, INT, data, status, eval_status, x, f, g, h)

    # Evaluate f(x) and its derivatives as required
    if (status[] ≥ 2)  # need objective
      f[] = objf(x[])
      if (status[] ≥ 3)  # need first derivative
        g[] = gradf(x[])
        if (status[] ≥ 4) # need second derivative
          h[] = hessf(x[])
        end
      end
    else  # the solution has been found (or an error has occured)
      terminated = true
    end
  end

  # Record solution information
  ugo_information(T, INT, data, inform, status)

  if inform[].status == 0
    @printf("%i evaluations. Optimal objective value = %5.2f status = %1i\n",
            inform[].f_eval, f[], inform[].status)
  else
    @printf("UGO_solve exit status = %1i\n", inform[].status)
  end

  # Delete internal workspace
  ugo_terminate(T, INT, data, control, inform)

  return 0
end

for (T, INT, libgalahad) in ((Float32 , Int32, GALAHAD.libgalahad_single      ),
                             (Float32 , Int64, GALAHAD.libgalahad_single_64   ),
                             (Float64 , Int32, GALAHAD.libgalahad_double      ),
                             (Float64 , Int64, GALAHAD.libgalahad_double_64   ),
                             (Float128, Int32, GALAHAD.libgalahad_quadruple   ),
                             (Float128, Int64, GALAHAD.libgalahad_quadruple_64))
  if isfile(libgalahad)
    @testset "UGO -- $T -- $INT" begin
      @test test_ugo(T, INT) == 0
    end
  end
end