CAPI trlib_tri_factor¶

Functions¶

trlib_int_t trlib_tri_factor_min(trlib_int_t nirblk, trlib_int_t *irblk, trlib_flt_t *diag, trlib_flt_t *offdiag, trlib_flt_t *neglin, trlib_flt_t radius, trlib_int_t itmax, trlib_flt_t tol_rel, trlib_flt_t tol_newton_tiny, trlib_int_t pos_def, trlib_int_t equality, trlib_int_t *warm0, trlib_flt_t *lam0, trlib_int_t *warm, trlib_flt_t *lam, trlib_int_t *warm_leftmost, trlib_int_t *ileftmost, trlib_flt_t *leftmost, trlib_int_t *warm_fac0, trlib_flt_t *diag_fac0, trlib_flt_t *offdiag_fac0, trlib_int_t *warm_fac, trlib_flt_t *diag_fac, trlib_flt_t *offdiag_fac, trlib_flt_t *sol0, trlib_flt_t *sol, trlib_flt_t *ones, trlib_flt_t *fwork, trlib_int_t refine, trlib_int_t verbose, trlib_int_t unicode, char *prefix, FILE *fout, trlib_int_t *timing, trlib_flt_t *obj, trlib_int_t *iter_newton, trlib_int_t *sub_fail)¶

Solves tridiagonal trust region subproblem

Computes minimizer to

\(\min \frac 12 \langle h, T h \rangle + \langle g, h \rangle\) subject to the trust region constraint \(\Vert h \Vert \le r\), where \(T \in \mathbb R^{n \times n}\) is symmetric tridiagonal.

Let \(T = \begin{pmatrix} T_1 & & \\ & \ddots & \\ & & T_\ell \end{pmatrix}\) be composed into irreducible blocks \(T_i\).

The minimizer is a global minimizer (modulo floating point).

The algorithm is the Moré-Sorensen-Method as decribed as Algorithm 5.2 in [Gould1999].

Convergence

Exit with success is reported in several cases:

interior solution: the stationary point \(T h = - g\) is suitable, iterative refinement is used in the solution of this system.

hard case: the smallest eigenvalue \(-\theta\) was found to be degenerate and \(h = v + \alpha w\) is returned with \(v\) solution of \((T + \theta I) v = -g\), \(w\) eigenvector corresponding to \(- \theta\) and \(\alpha\) chosen to satisfy the trust region constraint and being a minimizer.

boundary and hard case: convergence in Newton iteration is reported if \(\Vert h(\lambda) \Vert - r \le \texttt{tol}\_\texttt{rel} \, r\) with \((T+\lambda I) h(\lambda) = -g\), Newton breakdown reported if \(\vert d\lambda \vert \le \texttt{macheps} \vert \lambda \vert\).

Parameters:

nirblk (trlib_int_t, input) – number of irreducible blocks \(\ell\), ensure \(\ell > 0\)
irblk (trlib_int_t, input) –
pointer to indices of irreducible blocks, length nirblk+1:
- irblk[i] is the start index of block \(i\) in diag and offdiag
- irblk[i+1] - 1 is the stop index of block \(i\)
- irblk[i+1] - irred[i] the dimension \(n_\ell\) of block \(i\)
- irblk[nirred] the dimension of \(T\)
diag (trlib_flt_t, input) – pointer to array holding diagonal of \(T\), length irblk[nirblk]
offdiag (trlib_flt_t, input) – pointer to array holding offdiagonal of \(T\), length irblk[nirblk]
neglin (trlib_flt_t, input) – pointer to array holding \(-g\), length n
radius (trlib_flt_t, input) – trust region constraint radius \(r\)
itmax (trlib_int_t, input) – maximum number of Newton iterations
tol_rel (trlib_flt_t, input) – relative stopping tolerance for residual in Newton iteration for lam, good default may be \(\texttt{macheps}\) (TRLIB_EPS)
tol_newton_tiny (trlib_flt_t, input) – stopping tolerance for step in Newton iteration for lam, good default may be \(10 \texttt{macheps}^{.75}\)
pos_def (trlib_int_t, input) – set 1 if you know \(T\) to be positive definite, otherwise 0
equality (trlib_int_t, input) – set 1 if you want to enfore trust region constraint as equality, otherwise 0
warm0 (trlib_int_t, input/output) – set 1 if you provide a valid value in lam0, otherwise 0
lam0 (trlib_flt_t, input/output) – Lagrange multiplier such that \(T_0+ \mathtt{lam0} I\) positive definite and \((T_0+ \mathtt{lam0} I) \mathtt{sol0} = \mathtt{neglin}\) - on entry: estimate suitable for warmstart - on exit: computed multiplier
warm (trlib_int_t, input/output) – set 1 if you provide a valid value in lam, otherwise 0
lam (trlib_flt_t, input/output) – Lagrange multiplier such that \(T+ \mathtt{lam} I\) positive definite and \((T+ \mathtt{lam} I) \mathtt{sol} = \mathtt{neglin}\) - on entry: estimate suitable for warmstart - on exit: computed multiplier
warm_leftmost (trlib_int_t, input/output) –
- on entry: set 1 if you provide a valid value in leftmost section in fwork
- on exit: 1 if leftmost section in fwork holds valid exit value, otherwise 0
ileftmost (trlib_int_t, input/output) – index to block with smallest leftmost eigenvalue
leftmost (trlib_flt_t, input/output) – array holding leftmost eigenvalues of blocks
warm_fac0 (trlib_int_t, input/output) – set 1 if you provide a valid factoricization in diag_fac0, offdiag_fac0
diag_fac0 (trlib_flt_t, input/output) – pointer to array holding diagonal of Cholesky factorization of \(T_0 + \texttt{lam0} I\), length irblk[1] - on entry: factorization corresponding to provided estimation lam0 on entry - on exit: factorization corresponding to computed multiplier lam0
offdiag_fac0 (trlib_flt_t, input/output) – pointer to array holding offdiagonal of Cholesky factorization of \(T_0 + \texttt{lam0} I\), length irblk[1]-1 - on entry: factorization corresponding to provided estimation lam0 on entry - on exit: factorization corresponding to computed multiplier lam0
warm_fac (trlib_int_t, input/output) – set 1 if you provide a valid factoricization in diag_fac, offdiag_fac
diag_fac (trlib_flt_t, input/output) – pointer to array of length n that holds the following: - let \(j = \texttt{ileftmost}\) and \(\theta = - \texttt{leftmost[ileftmost]}\) - on position \(0, \ldots, \texttt{irblk[1]}\): diagonal of factorization of \(T_0 + \theta I\) - other positions have to be allocated
offdiag_fac (trlib_flt_t, input/output) – pointer to array of length n-1 that holds the following: - let \(j = \texttt{ileftmost}\) and \(\theta = - \texttt{leftmost[ileftmost]}\) - on position \(0, \ldots, \texttt{irblk[1]}-1\): offdiagonal of factorization of \(T_0 + \theta I\) - other positions have to be allocated
sol0 (trlib_flt_t, input/output) – pointer to array holding solution, length irblk[1]
sol (trlib_flt_t, input/output) – pointer to array holding solution, length n
ones (trlib_flt_t, input) – array with every value 1.0, length n
fwork (trlib_flt_t, input/output) –
floating point workspace, must be allocated memory on input of size trlib_tri_factor_memory_size() (call with argument n) and can be discarded on output, memory layout:

start end (excl) description

0 n auxiliary vector

n 2 n holds diagonal of \(T + \lambda I\)

2 n 4 n workspace for iterative refinement
refine (trlib_int_t, input) – set to 1 if iterative refinement should be used on solving linear systems, otherwise to 0
verbose (trlib_int_t, input) – determines the verbosity level of output that is written to fout
unicode (trlib_int_t, input) – set to 1 if fout can handle unicode, otherwise to 0
prefix (char, input) – string that is printed before iteration output
fout (FILE, input) – output stream
timing (trlib_int_t, input/output) –
gives timing details, provide allocated zero initialized memory of length trlib_tri_timing_size()

block description

0 total duration

1 timing of linear algebra calls

2 timing from trlib_leftmost_irreducible()

3 timing from trlib_eigen_inverse()
obj (trlib_flt_t, output) – objective function value at solution point
iter_newton (trlib_int_t, output) – number of Newton iterations
sub_fail (trlib_int_t, output) – status code of subroutine if failure occured in subroutines called

Returns:

status

TRLIB_TTR_CONV_BOUND success with solution on boundary
TRLIB_TTR_CONV_INTERIOR success with interior solution
TRLIB_TTR_HARD success, but hard case encountered and solution may be approximate
TRLIB_TTR_NEWTON_BREAK most likely success with accurate result; premature end of Newton iteration due to tiny step
TRLIB_TTR_HARD_INIT_LAM hard case encountered without being able to find suitable initial \(\lambda\) for Newton iteration, returned approximate stationary point that maybe suboptimal
TRLIB_TTR_ITMAX iteration limit exceeded
TRLIB_TTR_FAIL_FACTOR failure on matrix factorization
TRLIB_TTR_FAIL_LINSOLVE failure on backsolve
TRLIB_TTR_FAIL_EIG failure on eigenvalue computation. status code of trlib_eigen_inverse() in sub_fail
TRLIB_TTR_FAIL_LM failure on leftmost eigenvalue computation. status code of trlib_leftmost_irreducible() in sub_fail

Return type:

trlib_int_t

trlib_int_t trlib_tri_factor_regularized_umin(trlib_int_t n, trlib_flt_t *diag, trlib_flt_t *offdiag, trlib_flt_t *neglin, trlib_flt_t lam, trlib_flt_t *sol, trlib_flt_t *ones, trlib_flt_t *fwork, trlib_int_t refine, trlib_int_t verbose, trlib_int_t unicode, char *prefix, FILE *fout, trlib_int_t *timing, trlib_flt_t *norm_sol, trlib_int_t *sub_fail)¶

Computes minimizer of regularized unconstrained problem

Computes minimizer of

\(\min \frac 12 \langle h, (T + \lambda I) h \rangle + \langle g, h \rangle\), where \(T \in \mathbb R^{n \times n}\) is symmetric tridiagonal and \(\lambda\) such that \(T + \lambda I\) is spd.

Parameters:

n (trlib_int_t, input) – dimension, ensure \(n > 0\)
diag (trlib_flt_t, input) – pointer to array holding diagonal of \(T\), length n
offdiag (trlib_flt_t, input) – pointer to array holding offdiagonal of \(T\), length n-1
neglin (trlib_flt_t, input) – pointer to array holding \(-g\), length n
lam (trlib_flt_t, input) – regularization parameter \(\lambda\)
sol (trlib_flt_t, input/output) – pointer to array holding solution, length n
ones (trlib_flt_t, input) – array with every value 1.0, length n

fwork (trlib_flt_t, input/output) –

floating point workspace, must be allocated memory on input of size trlib_tri_factor_memory_size() (call with argument n) and can be discarded on output, memory layout:

start	end (excl)	description
0	`n`	holds diagonal of \(T + \lambda I\)
`n`	2 `n`	holds diagonal of factorization \(T + \lambda I\)
2 `n`	3 `n`	holds offdiagonal of factorization \(T + \lambda I\)
3 `n`	5 `n`	workspace for iterative refinement

refine (trlib_int_t, input) – set to 1 if iterative refinement should be used on solving linear systems, otherwise to 0
verbose (trlib_int_t, input) – determines the verbosity level of output that is written to fout
unicode (trlib_int_t, input) – set to 1 if fout can handle unicode, otherwise to 0
prefix (char, input) – string that is printed before iteration output
fout (FILE, input) – output stream
timing (trlib_int_t, input/output) –
gives timing details, provide allocated zero initialized memory of length trlib_tri_timing_size()

block description

0 total duration

1 timing of linear algebra calls
norm_sol (trlib_flt_t, output) – norm of solution fector
sub_fail (trlib_int_t, output) – status code of subroutine if failure occured in subroutines called

Returns:

status

TRLIB_TTR_CONV_INTERIOR success with interior solution
TRLIB_TTR_FAIL_FACTOR failure on matrix factorization
TRLIB_TTR_FAIL_LINSOLVE failure on backsolve

Return type:

trlib_int_t

trlib_int_t trlib_tri_factor_get_regularization(trlib_int_t n, trlib_flt_t *diag, trlib_flt_t *offdiag, trlib_flt_t *neglin, trlib_flt_t *lam, trlib_flt_t sigma, trlib_flt_t sigma_l, trlib_flt_t sigma_u, trlib_flt_t *sol, trlib_flt_t *ones, trlib_flt_t *fwork, trlib_int_t refine, trlib_int_t verbose, trlib_int_t unicode, char *prefix, FILE *fout, trlib_int_t *timing, trlib_flt_t *norm_sol, trlib_int_t *sub_fail)¶

Computes regularization parameter needed in trace

Let \(s(\lambda)\) be solution of \((T + \lambda I) s(\lambda) + g = 0\), where \(T \in \mathbb R^{n \times n}\) is symmetric tridiagonal and \(\lambda\) such that \(T + \lambda I\) is spd.

Then find \(\lambda\) with \(\sigma_{\text l} \le \frac{\lambda}{s(\lambda)} \le \sigma_{\text u}\).

Parameters:

n (trlib_int_t, input) – dimension, ensure \(n > 0\)
diag (trlib_flt_t, input) – pointer to array holding diagonal of \(T\), length n
offdiag (trlib_flt_t, input) – pointer to array holding offdiagonal of \(T\), length n-1
neglin (trlib_flt_t, input) – pointer to array holding \(-g\), length n
lam (trlib_flt_t, input/output) – regularization parameter \(\lambda\), on input initial guess, on exit desired parameter
sigma (trlib_flt_t, input) – value that is used in root finding, e.g. \(\frac 12 (\sigma_{\text l} + \sigma_{\text u})\)
sigma_l (trlib_flt_t, input) – lower bound
sigma_u (trlib_flt_t, input) – upper bound
sol (trlib_flt_t, input/output) – pointer to array holding solution, length n
ones (trlib_flt_t, input) – array with every value 1.0, length n

fwork (trlib_flt_t, input/output) –

floating point workspace, must be allocated memory on input of size trlib_tri_factor_memory_size() (call with argument n) and can be discarded on output, memory layout:

start	end (excl)	description
0	`n`	holds diagonal of \(T + \lambda I\)
`n`	2 `n`	holds diagonal of factorization \(T + \lambda I\)
2 `n`	3 `n`	holds offdiagonal of factorization \(T + \lambda I\)
3 `n`	5 `n`	workspace for iterative refinement
5 `n`	6 `n`	auxiliary vector

refine (trlib_int_t, input) – set to 1 if iterative refinement should be used on solving linear systems, otherwise to 0
verbose (trlib_int_t, input) – determines the verbosity level of output that is written to fout
unicode (trlib_int_t, input) – set to 1 if fout can handle unicode, otherwise to 0
prefix (char, input) – string that is printed before iteration output
fout (FILE, input) – output stream
timing (trlib_int_t, input/output) –
gives timing details, provide allocated zero initialized memory of length trlib_tri_timing_size()

block description

0 total duration

1 timing of linear algebra calls
norm_sol (trlib_flt_t, output) – norm of solution fector
sub_fail (trlib_int_t, output) – status code of subroutine if failure occured in subroutines called

Returns:

status

TRLIB_TTR_CONV_INTERIOR success with interior solution
TRLIB_TTR_FAIL_FACTOR failure on matrix factorization
TRLIB_TTR_FAIL_LINSOLVE failure on backsolve

Return type:

trlib_int_t

trlib_int_t trlib_tri_factor_regularize_posdef(trlib_int_t n, trlib_flt_t *diag, trlib_flt_t *offdiag, trlib_flt_t tol_away, trlib_flt_t security_step, trlib_flt_t *regdiag)¶

Compute diagonal regularization to make tridiagonal matrix positive definite

Parameters:	n (`trlib_int_t`, input) – dimension, ensure \(n > 0\) diag (`trlib_flt_t`, input) – pointer to array holding diagonal of \(T\), length `n` offdiag (`trlib_flt_t`, input) – pointer to array holding offdiagonal of \(T\), length `n-1` tol_away (`trlib_flt_t`, input) – tolerance that diagonal entries in factorization should be away from zero, relative to previous entry. Good default \(10^{-12}\). security_step (`trlib_flt_t`, input) – factor greater `1.0` that defines a margin to get away from zero in the step taken. Good default `2.0`. regdiag (`trlib_flt_t`, input/output) – pointer to array holding regularization term, length `n`
Returns:	`0`
Return type:	`trlib_int_t`

trlib_int_t trlib_tri_factor_memory_size(trlib_int_t n)¶

Gives information on memory that has to be allocated for trlib_tri_factor_min()

Parameters:	n (`trlib_int_t`, input) – dimension, ensure \(n > 0\) fwork_size (`trlib_flt_t`, output) – size of floating point workspace fwork that has to be allocated for `trlib_tri_factor_min()`
Returns:	`0`
Return type:	`trlib_int_t`

trlib_int_t trlib_tri_timing_size(void)¶: size that has to be allocated for timing in trlib_tri_factor_min()

Definitions¶

TRLIB_TTR_CONV_BOUND¶: 0

TRLIB_TTR_CONV_INTERIOR¶: 1

TRLIB_TTR_HARD¶: 2

TRLIB_TTR_NEWTON_BREAK¶: 3

TRLIB_TTR_HARD_INIT_LAM¶: 4

TRLIB_TTR_ITMAX¶: -1

TRLIB_TTR_FAIL_FACTOR¶: -2

TRLIB_TTR_FAIL_LINSOLVE¶: -3

TRLIB_TTR_FAIL_EIG¶: -4

TRLIB_TTR_FAIL_LM¶: -5

TRLIB_TTR_FAIL_HARD¶: -10