Trust Region Steps

This module provides the machinery to compute different trust-region( -reflective) step proposals and select among them based on to their performance according to the quadratic approximation of the objective function

Classes Summary

GradientStep(x, sg, hess, scaling, …)

This class provides the machinery to compute a gradient step.

Step(x, sg, hess, scaling, g_dscaling, …)

Base class for the computation of a proposal step

TRStep2D(x, sg, hess, scaling, g_dscaling, …)

This class provides the machinery to compute an approximate solution of the trust region subproblem according to a 2D subproblem

TRStepFull(x, sg, hess, scaling, g_dscaling, …)

This class provides the machinery to compute an exact solution of the trust region subproblem.

TRStepReflected(x, sg, hess, scaling, …)

This class provides the machinery to compute a reflected step based on trust region subproblem solution thathit the boundaries.

Functions Summary

normalize(v)

Inplace normalization of a vector

trust_region_reflective(x, g, hess, scaling, …)

Compute a step according to the solution of the trust-region subproblem.

Classes

class fides.trust_region.GradientStep(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]

This class provides the machinery to compute a gradient step.

__init__(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]
Parameters
  • x – Reference point

  • sg – Gradient in rescaled coordinates

  • hess – Hessian in unscaled coordinates

  • scaling – Matrix that defines scaling transformation

  • g_dscaling – Unscaled gradient multiplied by derivative of scaling transformation

  • delta – Trust region Radius in scaled coordinates

  • theta – Stepback parameter that controls how close steps are allowed to get to the boundary

  • ub – Upper boundary

  • lb – Lower boundary

class fides.trust_region.Step(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]

Base class for the computation of a proposal step

Variables
  • x – current state of optimization variables

  • s – proposed step

  • sc – coefficients in the 1D/2D subspace that defines the affine transformed step ss: ss = subspace * sc

  • ss – affine transformed step: s = scaling * ss

  • og_s – s without step back

  • og_sc – st without step back

  • og_ss – ss without step back

  • sg – rescaled gradient scaling * g

  • hess – hessian of the objective function at x

  • g_dscaling – diag(g) * dscaling

  • delta – trust region radius in the transformed space defined by scaling matrix

  • theta – controls step back, fraction of step to take if full step would reach breakpoint

  • lb – lower boundaries for x

  • ub – upper boundaries for x

  • minbr – maximal fraction of step s that can be taken to reach first breakpoint

  • ipt – index of x that specifies the variable that will hit the breakpoint if a step minbr * s is taken

  • qpval0 – value to the quadratic subproblem at x

  • qpval – value of the quadratic subproblem for the proposed step

  • shess – matrix of the full quadratic problem

  • cg – projection of the g_hat to the subspace

  • chess – projection of the B to the subspace

__init__(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]
Parameters
  • x (numpy.ndarray) – Reference point

  • sg (numpy.ndarray) – Gradient in rescaled coordinates

  • hess (numpy.ndarray) – Hessian in unscaled coordinates

  • scaling (scipy.sparse.csc.csc_matrix) – Matrix that defines scaling transformation

  • g_dscaling (scipy.sparse.csc.csc_matrix) – Unscaled gradient multiplied by derivative of scaling transformation

  • delta (float) – Trust region Radius in scaled coordinates

  • theta (float) – Stepback parameter that controls how close steps are allowed to get to the boundary

  • ub (numpy.ndarray) – Upper boundary

  • lb (numpy.ndarray) – Lower boundary

calculate()[source]

Calculates step and the expected objective function value according to the quadratic approximation

compute_step()[source]

Compute the step as solution to the trust region subproblem. Special code is used for the special case 1-dimensional subspace case

Return type

None

reduce_to_subspace()[source]

This function projects the matrix shess and the vector sg to the subspace

Return type

None

step_back()[source]

This function truncates the step based on the distance of the current point to the boundary.

class fides.trust_region.TRStep2D(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb, subspace)[source]

This class provides the machinery to compute an approximate solution of the trust region subproblem according to a 2D subproblem

__init__(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb, subspace)[source]
Parameters

subspace – Precomputed subspace

class fides.trust_region.TRStepFull(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]

This class provides the machinery to compute an exact solution of the trust region subproblem.

__init__(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb)[source]
Parameters
  • x – Reference point

  • sg – Gradient in rescaled coordinates

  • hess – Hessian in unscaled coordinates

  • scaling – Matrix that defines scaling transformation

  • g_dscaling – Unscaled gradient multiplied by derivative of scaling transformation

  • delta – Trust region Radius in scaled coordinates

  • theta – Stepback parameter that controls how close steps are allowed to get to the boundary

  • ub – Upper boundary

  • lb – Lower boundary

class fides.trust_region.TRStepReflected(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb, tr_step)[source]

This class provides the machinery to compute a reflected step based on trust region subproblem solution thathit the boundaries.

__init__(x, sg, hess, scaling, g_dscaling, delta, theta, ub, lb, tr_step)[source]
Parameters

tr_step – Trust-region step that is reflected

Functions

fides.trust_region.normalize(v)[source]

Inplace normalization of a vector

Parameters

v (numpy.ndarray) – vector to be normalized

Return type

None

fides.trust_region.trust_region_reflective(x, g, hess, scaling, tr_subspace, delta, dv, theta, lb, ub, subspace_dim)[source]

Compute a step according to the solution of the trust-region subproblem. If step-back is necessary, gradient and reflected trust region step are also evaluated in terms of their performance according to the local quadratic approximation

Parameters
  • x (numpy.ndarray) – Current values of the optimization variables

  • g (numpy.ndarray) – Objective function gradient at x

  • hess (numpy.ndarray) – (Approximate) objective function Hessian at x

  • scaling (scipy.sparse.csc.csc_matrix) – Scaling transformation according to distance to boundary

  • tr_subspace (numpy.ndarray) – Precomputed subspace from previous iteration for reuse if proposed step was not accepted

  • delta (float) – Trust region radius, note that this applies after scaling transformation

  • dv (numpy.ndarray) – derivative of scaling transformation

  • theta (float) – parameter regulating stepback

  • lb (numpy.ndarray) – lower optimization variable boundaries

  • ub (numpy.ndarray) – upper optimization variable boundaries

  • subspace_dim (fides.constants.SubSpaceDim) – Subspace dimension in which the subproblem will be solved. Larger subspaces require more compute time but can yield higher quality step proposals.

Return type

typing.Tuple[numpy.ndarray, numpy.ndarray, float, numpy.ndarray, str]

Returns

s: proposed step, ss: rescaled proposed step, qpval: expected function value according to local quadratic approximation, subspace: computed subspace for reuse if proposed step is not accepted, steptype: type of step that was selected for proposal