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
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