GLTR#
purpose#
The gltr package uses a Krylov-subspace iteration to find an approximation of the global minimizer of quadratic objective function within an ellipsoidal region; this is commonly known as the trust-region subproblem. The aim is to minimize the quadratic objective function
See Section 4 of $GALAHAD/doc/gltr.pdf for additional details.
method#
The required solution
The method is iterative. Starting with the vector
If the solution to (1) lies interior to the constraint, the required
solution
A cheaper alternative is to use the Steihuag-Toint strategy, which is simply
to stop at the first boundary point encountered along the piecewise
linear path generated by the conjugate-gradient iterates. Note that if
reference#
The method is described in detail in
N. I. M. Gould, S. Lucidi, M. Roma and Ph. L. Toint, ``Solving the trust-region subproblem using the Lanczos method’’. SIAM Journal on Optimization 9(2) (1999), 504-525.
functions#
- gltr.initialize()#
Set default option values and initialize private data
Returns:
- optionsdict
- dictionary containing default control options:
- errorint
error and warning diagnostics occur on stream error.
- outint
general output occurs on stream out.
- problemint
unit to write problem data into file problem_file.
- print_levelint
the level of output required is specified by print_level. Possible values are
<=0
gives no output,
1
gives a one-line summary for every iteration.
2
gives a summary of the inner iteration for each iteration.
>=3
gives increasingly verbose (debugging) output.
- itmaxint
the maximum number of iterations allowed (-ve = no bound).
- Lanczos_itmaxint
the maximum number of iterations allowed once the boundary has been encountered (-ve = no bound).
- extra_vectorsint
the number of extra work vectors of length n used.
- ritz_printout_deviceint
the unit number for writing debug Ritz values.
- stop_relativefloat
the iteration stops successfully when the gradient in the
norm is smaller than max( stop_relative
* norm initial gradient,stop_absolute
).- stop_absolutefloat
see stop_relative.
- fraction_optfloat
an estimate of the solution that gives at least
fraction_opt
times the optimal objective value will be found.- f_minfloat
the iteration stops if the objective-function value is lower than f_min.
- rminvr_zerofloat
the smallest value that the square of the M norm of the gradient of the the objective may be before it is considered to be zero.
- f_0float
the constant term,
, in the objective function. - unitmbool
is
the identity matrix ?. - steihaug_tointbool
should the iteration stop when the Trust-region is first encountered ?.
- boundarybool
is the solution thought to lie on the constraint boundary ?.
- equality_problembool
is the solution required to lie on the constraint boundary?.
- space_criticalbool
if
space_critical
is True, every effort will be made to use as little space as possible. This may result in longer computation time.- deallocate_error_fatalbool
if
deallocate_error_fatal
is True, any array/pointer deallocation error will terminate execution. Otherwise, computation will continue.- print_ritz_valuesbool
should the Ritz values be written to the debug stream?.
- ritz_file_namestr
name of debug file containing the Ritz values.
- prefixstr
all output lines will be prefixed by the string contained in quotes within
prefix
, e.g. ‘word’ (note the qutoes) will result in the prefix word.
- gltr.load_options(options=None)#
Import control options into internal storage prior to solution.
Parameters:
- optionsdict, optional
dictionary of control options (see
gltr.initialize
).
- gltr.solve_problem(status, n, radius, r, v)#
Find the global moinimizer of the quadratic objective function
within a trust-region of radius . Parameters:
- statusint
holds the entry status. Possible values are
1
an initial entry with r set to
.
4
a restart entry with
unchanged but a smaller radius .
other
the value returned from the previous call, see Returns below.
- nint
holds the number of variables.
- radiusfloat
holds the strictly positive trust-region radius,
. - rndarray(n)
holds the values of the linear term
in the objective function when initial or return status = 1, 4 or 5 (see below). - vndarray(n)
holds the result vector when return status = 2 or 3 (see below).
Returns:
- statusint
holds the exit status. Possible values are
0
the solution has been found, no further reentry is required
2
the inverse of
must be applied to the vector returned in v, the result placed in v, and the function recalled with status set to 2. This will only occur if options[‘unitm’] is False.
3
the product of
with the vector returned in v must be formed, the result placed in v, and the function recalled with status set to 3.
5
the iteration must be restarted by setting r to
, and the function recalled with status set to 5.
<0
an error occurred, see
status
ingltr.information
for further details.- xndarray(n)
holds the values of the approximate minimizer
. - rndarray(n)
holds the values of the gradient
at the current . - vndarray(n)
holds the return vector when return status = 2 or 3 (see above).
- [optional] gltr.information()
Provide optional output information
Returns:
- informdict
- dictionary containing output information:
- statusint
return status. Possible values are:
0
The run was successful.
-1
An allocation error occurred. A message indicating the offending array is written on unit options[‘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 options[‘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.
-3
The restriction n > 0 or
has been violated. -15
appears to be indefinite. -18
The iteration limit has been exceeded.
-30
The trust-region has been encountered in Steihaug-Toint mode.
-31
The function value is smaller than options[‘f_min’].
- alloc_statusint
the status of the last attempted allocation/deallocation.
- bad_allocstr
the name of the array for which an allocation/deallocation error occurred.
- iterint
the total number of iterations required.
- iter_pass2int
the total number of pass-2 iterations required if the solution lies on the trust-region boundary.
- objfloat
the value of the quadratic function.
- multiplierfloat
the Lagrange multiplier corresponding to the trust-region constraint.
- mnormxfloat
the
-norm of , . - pivfloat
the latest pivot in the Cholesky factorization of the Lanczos tridiagona.
- curvfloat
the most negative cuurvature encountered.
- rayleighfloat
the current Rayleigh quotient.
- leftmostfloat
an estimate of the leftmost generalized eigenvalue of the pencil
. - negative_curvaturebool
was negative curvature encountered ?.
- hard_casebool
did the hard case occur ?.
- gltr.terminate()#
Deallocate all internal private storage.
example code#
from galahad import gltr
import numpy as np
np.set_printoptions(precision=4,suppress=True,floatmode='fixed')
print("\n** python test: gltr")
# set problem data
n = 100
v = np.empty(n)
hv = np.empty(n)
g = np.empty(n)
g[:] = 1.0
# set default options
options = gltr.initialize()
# set some non-default options
options['unitm'] = False
options['print_level'] = 0
#print("options:", options)
# load non-default options
gltr.load_options(options)
# first problem
print(" solve problem 1")
status = 1
radius = 1.0
r = g
while status > 0:
status, x, r, v = gltr.solve_problem(status, n, radius, r, v)
if status == 2: # precondition v -> v/2
for i in range(n):
v[i] = v[i] / 2.0;
elif status == 3: # form v -> H v
hv[0] = 2.0 * v[0] + v[1]
for i in range(1,n-1):
hv[i] = v[i-1] + 2.0 * v[i] + v[i+1]
hv[n-1] = v[n-2] + 2.0 * v[n-1]
for i in range(n):
v[i] = hv[i]
elif status == 5: # reset r
r = g
# get information
inform = gltr.information()
print(" optimal f: %.4f" % inform['obj'])
#print(" x:",x)
# second problem, same data with smaller radius
print("\n solve problem 2 with smaller radius")
status = 4
radius = 0.1
r = g
while status > 0:
status, x, r, v = gltr.solve_problem(status, n, radius, r, v)
if status == 2: # precondition v -> v/2
for i in range(n):
v[i] = v[i] / 2.0;
elif status == 3: # form v -> H v
hv[0] = 2.0 * v[0] + v[1]
for i in range(1,n-1):
hv[i] = v[i-1] + 2.0 * v[i] + v[i+1]
hv[n-1] = v[n-2] + 2.0 * v[n-1]
for i in range(n):
v[i] = hv[i]
elif status == 5: # reset r
r = g
# get information
inform = gltr.information()
print(" optimal f: %.4f" % inform['obj'])
#print(" x:",x)
print('** gltr exit status:', inform['status'])
# deallocate internal data
gltr.terminate()
This example code is available in $GALAHAD/src/gltr/Python/test_gltr.py .