Home > manopt > tools > hessianextreme.m

hessianextreme

PURPOSE ^

Compute an extreme eigenvector / eigenvalue of the Hessian of a problem.

SYNOPSIS ^

function [y, lambda, info] = hessianextreme(problem, x, side, y0, options, storedb, key)

DESCRIPTION ^

 Compute an extreme eigenvector / eigenvalue of the Hessian of a problem.

 [u, lambda, info] = hessianextreme(problem, x)
 [u, lambda, info] = hessianextreme(problem, x, side)
 [u, lambda, info] = hessianextreme(problem, x, side, u0)
 [u, lambda, info] = hessianextreme(problem, x, side, u0, options)
 [u, lambda, info] = hessianextreme(problem, x, side, u0, options, storedb)
 [u, lambda, info] = hessianextreme(problem, x, side, u0, options, storedb, key)
 
 (For side, u0 and options, pass [] to omit any.)

 Given a Manopt problem structure and a point x on the manifold problem.M,
 this function computes a tangent vector u at x of unit norm such that the
 Hessian quadratic form is minimized or maximized:

    minimize or maximize <u, Hess f(x)[u]> such that <u, u> = 1,

 where <.,.> is the Riemannian metric on the tangent space at x. Choose
 between minimizing and maximizing by setting side = 'min' or 'max', with
 'min' being the default. The value attained is returned as lambda, and
 is the minimal or maximal eigenvalue of the Hessian (actually, the last
 value attained when the solver stopped). This is a real number since the
 Hessian is a symmetric operator.

 If u0 is specified, it should be a unit-norm tangent vector at x. It is
 then used as initial guess to solve the above problem. Pass [] to omit.

 The options structure, if provided, will be passed along to manoptsolve.
 As such, you may choose which solver to use to solve the above
 optimization problem by setting options.solver. See manoptsolve's help.
 The other options will be passed along to the chosen solver too.
 Pass [] to omit.

 Often times, it is only necessary to compute a vector u such that the
 quadratic form is negative, if that is at all possible. To do so, set the
 following stopping criterion: options.tolcost = -1e-10; (for example)
 and side = 'min'. The solver will return as soon as the quadratic cost
 defined above drops below the set value (or sooner if another stopping
 criterion triggers first.)

 storedb is a StoreDB object, key is the StoreDB key to point x.

 info is the info struct-array returned by the solver.

 See also: hessianspectrum manoptsolve tangentspherefactory

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [y, lambda, info] = hessianextreme(problem, x, side, y0, options, storedb, key)
0002 % Compute an extreme eigenvector / eigenvalue of the Hessian of a problem.
0003 %
0004 % [u, lambda, info] = hessianextreme(problem, x)
0005 % [u, lambda, info] = hessianextreme(problem, x, side)
0006 % [u, lambda, info] = hessianextreme(problem, x, side, u0)
0007 % [u, lambda, info] = hessianextreme(problem, x, side, u0, options)
0008 % [u, lambda, info] = hessianextreme(problem, x, side, u0, options, storedb)
0009 % [u, lambda, info] = hessianextreme(problem, x, side, u0, options, storedb, key)
0010 %
0011 % (For side, u0 and options, pass [] to omit any.)
0012 %
0013 % Given a Manopt problem structure and a point x on the manifold problem.M,
0014 % this function computes a tangent vector u at x of unit norm such that the
0015 % Hessian quadratic form is minimized or maximized:
0016 %
0017 %    minimize or maximize <u, Hess f(x)[u]> such that <u, u> = 1,
0018 %
0019 % where <.,.> is the Riemannian metric on the tangent space at x. Choose
0020 % between minimizing and maximizing by setting side = 'min' or 'max', with
0021 % 'min' being the default. The value attained is returned as lambda, and
0022 % is the minimal or maximal eigenvalue of the Hessian (actually, the last
0023 % value attained when the solver stopped). This is a real number since the
0024 % Hessian is a symmetric operator.
0025 %
0026 % If u0 is specified, it should be a unit-norm tangent vector at x. It is
0027 % then used as initial guess to solve the above problem. Pass [] to omit.
0028 %
0029 % The options structure, if provided, will be passed along to manoptsolve.
0030 % As such, you may choose which solver to use to solve the above
0031 % optimization problem by setting options.solver. See manoptsolve's help.
0032 % The other options will be passed along to the chosen solver too.
0033 % Pass [] to omit.
0034 %
0035 % Often times, it is only necessary to compute a vector u such that the
0036 % quadratic form is negative, if that is at all possible. To do so, set the
0037 % following stopping criterion: options.tolcost = -1e-10; (for example)
0038 % and side = 'min'. The solver will return as soon as the quadratic cost
0039 % defined above drops below the set value (or sooner if another stopping
0040 % criterion triggers first.)
0041 %
0042 % storedb is a StoreDB object, key is the StoreDB key to point x.
0043 %
0044 % info is the info struct-array returned by the solver.
0045 %
0046 % See also: hessianspectrum manoptsolve tangentspherefactory
0047 
0048 % This file is part of Manopt: www.manopt.org.
0049 % Original author: Nicolas Boumal, Aug. 13, 2014.
0050 % Contributors:
0051 % Change log:
0052 %
0053 %   April 3, 2015 (NB):
0054 %       Works with the new StoreDB class system.
0055 %
0056 %   May 7, 2015 (NB):
0057 %       Default solver options: verbosity = 0 and defaults to trustregions.
0058 %
0059 %   Nov 27, 2015 (NB):
0060 %       The function now also returns the info struct-array.
0061 
0062     
0063     % By default, minimize
0064     if ~exist('side', 'var') || isempty(side)
0065         side = 'min';
0066     end
0067     
0068     % If no initial guess was specified, prepare the empty one.
0069     if ~exist('y0', 'var')
0070         y0 = [];
0071     end
0072 
0073     % Merge default solver options with potential user-specified options.
0074     % Set local defaults here
0075     localdefaults.verbosity = 0;
0076     localdefaults.solver = @trustregions;
0077     if ~exist('options', 'var') || isempty(options)
0078         options = struct();
0079     end
0080     options = mergeOptions(localdefaults, options);
0081 
0082     % Allow omission of the key, and even of storedb.
0083     if ~exist('key', 'var')
0084         if ~exist('storedb', 'var')
0085             storedb = StoreDB();
0086         end
0087         key = storedb.getNewKey();
0088     end
0089     
0090     % Convert the side into a sign.
0091     % Since Manopt minimizes, 'min' asks for no sign change.
0092     switch lower(side)
0093         case 'min'
0094             sign = +1;
0095         case 'max'
0096             sign = -1;
0097         otherwise
0098             error('The side should be either ''min'' or ''max''.');
0099     end
0100 
0101     % We define a manifold that is actually the unit sphere on the tangent
0102     % space to problem.M at x. A generalization would be to consider
0103     % Stiefel or Grassmann on the tangent space, but this would require
0104     % manipulating collections of tangent vectors, which in full generality
0105     % may be more complex (from a programming point of view).
0106     % Points are represented as tangent vectors of unit norm.
0107     % Tangent vectors are represented as tangent vectors orthogonal to the
0108     % root point, with respect to the Riemannian metric on the tangent
0109     % space.
0110     
0111     % M is the original manifold. x is a point on M.
0112     M = problem.M;
0113     
0114     % N is the manifold we build. y will be a point on N, thus also a
0115     % tangent vector to M at x. This is a typical Riemannian submanifold of
0116     % a Euclidean space, hence it is easy to describe in terms of the tools
0117     % available for M.
0118     N = tangentspherefactory(M, x);
0119     
0120     % It is usually a good idea to force a gradient computation to make
0121     % sure precomputable things are precomputed.
0122     if canGetGradient(problem)
0123         [unused1, unused2] = getCostGrad(problem, x, storedb, key); %#ok
0124     end
0125     
0126     % This is the star operator of this party.
0127     hessian = @(y) getHessian(problem, x, y, storedb, key);
0128     
0129     % Start a Manopt problem structure for the quadratic optimization
0130     % problem on the sphere N.
0131     new_problem.M = N;
0132     
0133     % Define the cost function, its gradient and its Hessian.
0134 
0135     new_problem.cost = @cost;
0136     function [f, store] = cost(y, store)
0137         store = prepare(y, store);
0138         f = sign*store.f;
0139     end
0140 
0141     new_problem.grad = @grad;
0142     function [g, store] = grad(y, store)
0143         store = prepare(y, store);
0144         g = N.lincomb(y, sign*2, store.Hy, sign*(-2)*store.f, y);
0145     end
0146 
0147     new_problem.hess = @hess;
0148     function [h, store] = hess(y, ydot, store)
0149         store = prepare(y, store);
0150         Hydot = hessian(ydot);
0151         h = N.lincomb(y, sign*2, Hydot, sign*(-2)*store.f, ydot);
0152         h = N.proj(y, h);
0153     end
0154 
0155     % This helper makes sure we do not duplicate Hessian computations.
0156     function store = prepare(y, store)
0157         if ~isfield(store, 'ready')
0158             Hy = hessian(y);
0159             store.f = M.inner(x, y, Hy);
0160             store.Hy = Hy;
0161             store.ready = true;
0162         end
0163     end
0164     
0165     % Call a Manopt solver to solve the quadratic optimization problem on
0166     % the abstract sphere N.
0167     [y, lambda, info] = manoptsolve(new_problem, y0, options);
0168     lambda = sign*lambda;
0169 
0170 end

Generated on Fri 30-Sep-2022 13:18:25 by m2html © 2005