3.2.5. horton.correlatedwfn.trustregionopt – Optimization of Newton step

Module contains different trust radius optimization algorithms

class horton.correlatedwfn.trustregionopt.Dogleg(pn, g, H, radius)

Powell’s single dogleg steps approximate the trust region step over the two dimensional subspace spanned by the Cauchy step and any Newton step. Depending on the size of current step length, the step is adjusted yielding one of the following steps:

  • if p(newton) < radius: Do full Newton step (not calculated here)
  • elif p(cauchy) > radius: Do Cauchy step
  • else: Do dogleg step

** Arguements **

pn
The full Newton step
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance
radius
The current trustradius

** Returns **

step:final step, a OneIndex instance
stepNorm:Euclidian norm of the step
__call__()

The dogleg algorithm.

__init__(pn, g, H, radius)

Powell’s single dogleg steps approximate the trust region step over the two dimensional subspace spanned by the Cauchy step and any Newton step. Depending on the size of current step length, the step is adjusted yielding one of the following steps:

  • if p(newton) < radius: Do full Newton step (not calculated here)
  • elif p(cauchy) > radius: Do Cauchy step
  • else: Do dogleg step

** Arguements **

pn
The full Newton step
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance
radius
The current trustradius

** Returns **

step:final step, a OneIndex instance
stepNorm:Euclidian norm of the step
class horton.correlatedwfn.trustregionopt.DoubleDogleg(pn, g, H, radius)

Double dogleg steps approximate the trust region step

  • if ||p(g)|| >= radius: Do step in g direction
  • elif c*||p(n)|| >= radius: Do ddl step
  • elif ||p(n)|| >= radius: Do scaled Newton step
  • else: Do Newton step

** Arguements **

pn
The full Newton step
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance
radius
The current trustradius

** Returns **

step:final step, a OneIndex instance
stepNorm:Euclidian norm of the step
__call__()

The double dogleg algorithm.

__init__(pn, g, H, radius)

Double dogleg steps approximate the trust region step

  • if ||p(g)|| >= radius: Do step in g direction
  • elif c*||p(n)|| >= radius: Do ddl step
  • elif ||p(n)|| >= radius: Do scaled Newton step
  • else: Do Newton step

** Arguements **

pn
The full Newton step
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance
radius
The current trustradius

** Returns **

step:final step, a OneIndex instance
stepNorm:Euclidian norm of the step
class horton.correlatedwfn.trustregionopt.TruncatedCG(lf, g, H, radius, **kwargs)

Solve the quadratic trust-region subproblem

minimize < g, s > + 1/2 < s, Hs > subject to < s, s > <= radius

by means of the truncated conjugate gradient algorithm as described in Trond Steihaug, SIAM Journal on Numerical Analysis (1983), 20, 626-637.

The algorithm stops as soon as the preconditioned norm of the gradient falls under

max( abstol, reltol * g0 )

where g0 is the preconditioned norm of the initial gradient (or the Euclidian norm if no preconditioner is given), or as soon as the iterates cross the boundary of the trust region.

** Arguements **

lf
A linalg factory instance
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance

** Returns **

step:final step
niter:number of iterations
stepNorm:Euclidian norm of the step
__call__(**kwargs)

Solve the trust-region subproblem.

** Keywords **

s0:initial guess (default: [0,0,…,0]),
abstol:absolute stopping tolerance (default: 1.0e-8),
reltol:relative stopping tolerance (default: 1.0e-6),
maxiter:maximum number of iterations (default: 2n),
prec:a user-defined preconditioner.
__init__(lf, g, H, radius, **kwargs)

Solve the quadratic trust-region subproblem

minimize < g, s > + 1/2 < s, Hs > subject to < s, s > <= radius

by means of the truncated conjugate gradient algorithm as described in Trond Steihaug, SIAM Journal on Numerical Analysis (1983), 20, 626-637.

The algorithm stops as soon as the preconditioned norm of the gradient falls under

max( abstol, reltol * g0 )

where g0 is the preconditioned norm of the initial gradient (or the Euclidian norm if no preconditioner is given), or as soon as the iterates cross the boundary of the trust region.

** Arguements **

lf
A linalg factory instance
g
The gradient, a OneIndex instance
H
The Hessian, a OneIndex instance

** Returns **

step:final step
niter:number of iterations
stepNorm:Euclidian norm of the step
to_boundary(s, p, radius, ss=None)

Given vectors s and p and a trust-region radius radius > 0, return the positive scalar sigma such that

|| s + sigma * p || = radius

in Euclidian norm. If known, supply optional argument ss whose value should be the squared Euclidian norm of s.