Subproblem Solvers¶
This module provides the machinery to solve 1- and N-dimensional trust-region subproblems.
Functions Summary
|
Derivative of the secular equation |
|
Computes the derivative of the solution \(s(\lambda)\) with respect to lambda, where \(s\) is the subproblem solution according to |
|
|
|
Computes the quadratic form \(x^TQx + x^Tp\) |
|
Secular equation |
|
Computes the solution \(s(\lambda)\) as subproblem solution according to |
|
Solves the special case of a one-dimensional subproblem |
|
This function exactly solves the n-dimensional subproblem. |
Functions
- fides.subproblem.dsecular(lam, w, eigvals, eigvecs, delta)[source]¶
Derivative of the secular equation
\(\phi(\lambda) = \frac{1}{||s||} - \frac{1}{\Delta}\)
with respect to \(\lambda\)
- Parameters
lam (
float
) – \(\lambda\)w (
numpy.ndarray
) – precomputed eigenvector coefficients for -geigvals (
numpy.ndarray
) – precomputed eigenvalues of Beigvecs (
numpy.ndarray
) – precomputed eigenvectors of Bdelta (
float
) – trust region radius \(\Delta\)
- Returns
\(\frac{\partial \phi(\lambda)}{\partial \lambda}\)
- fides.subproblem.dslam(lam, w, eigvals, eigvecs)[source]¶
Computes the derivative of the solution \(s(\lambda)\) with respect to lambda, where \(s\) is the subproblem solution according to
\(-(B + \lambda I)s = g\)
- Parameters
lam (
float
) – \(\lambda\)w (
numpy.ndarray
) – precomputed eigenvector coefficients for -geigvals (
numpy.ndarray
) – precomputed eigenvalues of Beigvecs (
numpy.ndarray
) – precomputed eigenvectors of B
- Returns
\(\frac{\partial s(\lambda)}{\partial \lambda}\)
- fides.subproblem.quadratic_form(Q, p, x)[source]¶
Computes the quadratic form \(x^TQx + x^Tp\)
- Parameters
Q (
numpy.ndarray
) – Matrixp (
numpy.ndarray
) – Vectorx (
numpy.ndarray
) – Input
- Return type
- Returns
Value of form
- fides.subproblem.secular(lam, w, eigvals, eigvecs, delta)[source]¶
Secular equation
\(\phi(\lambda) = \frac{1}{||s||} - \frac{1}{\Delta}\)
Subproblem solutions are given by the roots of this equation
- Parameters
lam (
float
) – \(\lambda\)w (
numpy.ndarray
) – precomputed eigenvector coefficients for -geigvals (
numpy.ndarray
) – precomputed eigenvalues of Beigvecs (
numpy.ndarray
) – precomputed eigenvectors of Bdelta (
float
) – trust region radius \(\Delta\)
- Returns
\(\phi(\lambda)\)
- fides.subproblem.slam(lam, w, eigvals, eigvecs)[source]¶
Computes the solution \(s(\lambda)\) as subproblem solution according to
\(-(B + \lambda I)s = g\)
- Parameters
lam (
float
) – \(\lambda\)w (
numpy.ndarray
) – precomputed eigenvector coefficients for -geigvals (
numpy.ndarray
) – precomputed eigenvalues of Beigvecs (
numpy.ndarray
) – precomputed eigenvectors of B
- Return type
- Returns
\(s(\lambda)\)
- fides.subproblem.solve_1d_trust_region_subproblem(B, g, s, delta, s0)[source]¶
Solves the special case of a one-dimensional subproblem
- Parameters
B (
numpy.ndarray
) – Hessian of the quadratic subproblemg (
numpy.ndarray
) – Gradient of the quadratic subproblems (
numpy.ndarray
) – Vector defining the one-dimensional search directiondelta (
float
) – Norm boundary for the solution of the quadratic subproblems0 (
numpy.ndarray
) – reference point from where search is started, also counts towards norm of step
- Return type
- Returns
Proposed step-length
- fides.subproblem.solve_nd_trust_region_subproblem(B, g, delta, logger=None)[source]¶
This function exactly solves the n-dimensional subproblem.
\(argmin_s\{s^T B s + s^T g = 0: ||s|| <= \Delta, s \in \mathbb{ R}^n\}\)
The solution is characterized by the equation \(-(B + \lambda I)s = g\). If B is positive definite, the solution can be obtained by \(\lambda = 0\) if \(Bs = -g\) satisfies \(||s|| <= \Delta\). If B is indefinite or \(Bs = -g\) satisfies \(||s|| > \Delta\) and an appropriate \(\lambda\) has to be identified via 1D rootfinding of the secular equation
\(\phi(\lambda) = \frac{1}{||s(\lambda)||} - \frac{1}{\Delta} = 0\)
with \(s(\lambda)\) computed according to an eigenvalue decomposition of B. The eigenvalue decomposition, although being more expensive than a cholesky decomposition, has the advantage that eigenvectors are invariant to changes in \(\lambda\) and eigenvalues are linear in \(\lambda\), so factorization only has to be performed once. We perform the linesearch via Newton’s algorithm and Brent-Q as fallback. The hard case is treated separately and serves as general fallback.
- Parameters
B (
numpy.ndarray
) – Hessian of the quadratic subproblemg (
numpy.ndarray
) – Gradient of the quadratic subproblemdelta (
float
) – Norm boundary for the solution of the quadratic subproblemlogger (
typing.Optional
[logging.Logger
]) – Logger instance to be used for logging
- Return type
- Returns
s: Selected step, step_type: Type of solution that was obtained