Edit on GitHub

Introduction

Purpose

Given real $n$ by $n$ symmetric matrices $H$ and $M$ (with $M$ positive definite), a real $n$ vector $c$ and scalars $\Delta>0$ and $f_0$, this package finds an ** approximate minimizer of the quadratic objective function $\frac{1}{2} x^T H x+c^T x + f_0$, where the vector $x$ is required to satisfy the constraint $\|x\|_M \leq\Delta$**, and where the $M$-norm of $x$ is $\|x\|_M = \sqrt{x^T M x}$. This problem commonly occurs as a trust-region subproblem in nonlinear optimization calculations. The method may be suitable for large $n$ as no factorization of $H$ is required. Reverse communication is used to obtain matrix-vector products of the form $H z$ and $M^{-1} z$.

The package may also be used to solve the related problem in which $x$ is instead required to satisfy the equality constraint $\|x\|_M = \Delta$.

Authors

N. I. M. Gould, STFC-Rutherford Appleton Laboratory, England.

C interface, additionally J. Fowkes, STFC-Rutherford Appleton Laboratory.

Julia interface, additionally A. Montoison and D. Orban, Polytechnique Montréal.

Originally released

April 1997, C interface December 2021.

Terminology

Method

The required solution $x$ necessarily satisfies the optimality condition $H x + \lambda M x + c = 0$, where $\lambda \geq 0$is a Lagrange multiplier corresponding to the constraint $\|x\|_M\leq\Delta$. In addition, the matrix $H + \lambda M$ will be positive definite.

The method is iterative. Startingwith the vector $M^{-1} c$, a matrix of Lanczos vectors is built one column at a time so that the $k$-th column is generated during iteration $k$. These columns span a so-called Krylov space. The resulting $n$ by $k$ matrix $Q_k $ has the property that $Q_k^T H Q_k^{}=T_k^{}$, where $T_k$ is tridiagonal. An approximation to the required solution may then be expressed formally as $x_{k+1}=Q_k y_k,$ \n x{k+1}=Qk yk, \n where yk $ solves the “tridiagonal” subproblem of minimizing $(1) \;\;\; \frac{1}{2} y^T T_k y+ \|c\|_{M^{-1}} e_1^T y\;\mbox{subject to the constraint}\; \|y\|_2\leq\Delta,$ \n (1)1/2 y^T Tk y+ ||c||M^{-1} e1^T y subject to the constraint \|y\|2leqDelta, \n and where $e_1$ is the first unit vector.

If the solution to (1) lies interior to the constraint, the required solution $x_{k+1}$ may simply be found as the $k$-th (preconditioned) conjugate-gradient iterate. This solution can be obtained without the need to access the whole matrix $Q_k$. These conjugate-gradient iterates increase in $M$-norm, and thus once one of them exceeds $\Delta$ in $M$-norm, the solution must occur on the constraint boundary. Thereafter, the solution to (1) is less easy to obtain, but an efficient inner iteration to solve (1) is nonetheless achievable because $T_k $ is tridiagonal. It is possible to observe the optimality measure $\|H x+\lambda M x+c\|_{M^{-1}}$ without computing $x_{k+1}$, and thus without needing $Q_k $. Once this measure is sufficiently small, a second pass is required to obtain the estimate $x_{k+1} $ from $y_k $. As this second pass is an additional expense, a record is kept of the optimal objective function values for each value of $k$, and the second pass is only performed so far as to ensure a given fraction of the final optimal objective value. Large savings may be made in the second pass by choosing the required fraction to be significantly smaller than one.

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 $H$ is significantly indefinite, this strategy often produces a far from optimal point, but is effective when $H$ is positive definite or almost

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.

Call order

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

  • gltr_initialize - provide default control parameters and set up initial data structures
  • gltr_read_specfile (optional) - override control values by reading replacement values from a file
  • gltr_import_control - import control parameters prior to

solution

  • gltr_solve_problem - solve the problem by reverse

communication, a sequence of calls are made under control of a status parameter, each exit either asks the user to provide additional informaton and to re-enter, or reports that either the solution has been found or that an error has occurred

  • gltr_information (optional) - recover information about the solution and solution process
  • gltr_terminate - deallocate data structures

See Section~\ref{examples} for an example of use. See the <a href="examples.html">examples tab</a> for an illustration of use. See the examples section for an illustration of use.