Home > manopt > tools > hessianmatrix.m

# hessianmatrix

## PURPOSE Computes a matrix which represents the Hessian in some tangent basis.

## SYNOPSIS function [H, basis] = hessianmatrix(problem, x, basis)

## DESCRIPTION ``` Computes a matrix which represents the Hessian in some tangent basis.

[H, basis] = hessianmatrix(problem, x)
[H, basis] = hessianmatrix(problem, x, basis)

problem is a Manopt problem structure with a manifold and cost function.
x is a point on the manifold problem.M.
basis (optional) is an orthonormal basis for the tangent space to the
manifold at x. If no basis is supplied, one will be generated at random.
If the basis spans only a subspace of the tangent space at x,
then the returned matrix represents the Hessian restricted to that subspace.

H is an n-by-n symmetric matrix (with n the number of vectors in the basis)
such that H(i, j) is the inner product between basis{i}
and Hess(basis{j}), with respect to the metric on the tangent space to
problem.M at x, where Hess(basis{j}) is the vector obtained after
applying the Hessian at x to basis{j}.

For optimization, it is usually not useful to compute the Hessian matrix,
as this quickly becomes expensive. This tool is provided mostly for
exploration and debugging rather than to be used algorithmically in
solvers. To access the spectrum of the Hessian, it may be more practical
to call hessianextreme or hessianspectrum. This should coincide with eig(H).

Example of equivalence:

Hu = getHessian(problem, x, u)

is equivalent to (but much faster than):

B = tangentorthobasis(M, x);
H = hessianmatrix(problem, x, B);
u_vec = tangent2vec(M, x, B, u);
Hu_vec = H*u_vec;
Hu = lincomb(M, x, B, Hu_vec);

Note that there will be some error due to numerical round-off.

## CROSS-REFERENCE INFORMATION This function calls:
• StoreDB
• canGetApproxHessian Checks whether an approximate Hessian can be computed for this problem.
• canGetHessian Checks whether the Hessian can be computed for a problem structure.
• getHessian Computes the Hessian of the cost function at x along d.
• tangentorthobasis Returns an orthonormal basis of tangent vectors in the Manopt framework.
This function is called by:

## SOURCE CODE ```0001 function [H, basis] = hessianmatrix(problem, x, basis)
0002 % Computes a matrix which represents the Hessian in some tangent basis.
0003 %
0004 % [H, basis] = hessianmatrix(problem, x)
0005 % [H, basis] = hessianmatrix(problem, x, basis)
0006 %
0007 % problem is a Manopt problem structure with a manifold and cost function.
0008 % x is a point on the manifold problem.M.
0009 % basis (optional) is an orthonormal basis for the tangent space to the
0010 % manifold at x. If no basis is supplied, one will be generated at random.
0011 % If the basis spans only a subspace of the tangent space at x,
0012 % then the returned matrix represents the Hessian restricted to that subspace.
0013 %
0014 % H is an n-by-n symmetric matrix (with n the number of vectors in the basis)
0015 % such that H(i, j) is the inner product between basis{i}
0016 % and Hess(basis{j}), with respect to the metric on the tangent space to
0017 % problem.M at x, where Hess(basis{j}) is the vector obtained after
0018 % applying the Hessian at x to basis{j}.
0019 %
0020 % For optimization, it is usually not useful to compute the Hessian matrix,
0021 % as this quickly becomes expensive. This tool is provided mostly for
0022 % exploration and debugging rather than to be used algorithmically in
0023 % solvers. To access the spectrum of the Hessian, it may be more practical
0024 % to call hessianextreme or hessianspectrum. This should coincide with eig(H).
0025 %
0026 %
0027 % Example of equivalence:
0028 %
0029 %     Hu = getHessian(problem, x, u)
0030 %
0031 % is equivalent to (but much faster than):
0032 %
0033 %     B = tangentorthobasis(M, x);
0034 %     H = hessianmatrix(problem, x, B);
0035 %     u_vec = tangent2vec(M, x, B, u);
0036 %     Hu_vec = H*u_vec;
0037 %     Hu = lincomb(M, x, B, Hu_vec);
0038 %
0039 % Note that there will be some error due to numerical round-off.
0040 %
0041 %
0043
0044 % This file is part of Manopt: www.manopt.org.
0045 % Original author: Nicolas Boumal, July 14, 2016.
0046 % Contributors:
0047 % Change log:
0048
0049 % TODO: refactor using operator2matrix
0050
0051
0052     % No warning if an approximate Hessian is available, as then the user
0053     % is presumably aware of what they are doing.
0054     if ~canGetHessian(problem) && ~canGetApproxHessian(problem)
0055         warning('manopt:hessianmatrix:nohessian', ...
0056                 ['The Hessian appears to be unavailable.\n' ...
0057                  'Will try to use an approximate Hessian instead.\n'...
0058                  'Since this approximation may not be linear or '...
0059                  'symmetric,\nthe computation might fail and the '...
0060                  'results (if any)\nmight make no sense.']);
0061     end
0062
0063
0064     % Unless an orthonormal basis for the tangent space at x is provided,
0065     % pick a random one.
0066     if ~exist('basis', 'var') || isempty(basis)
0067         n = problem.M.dim();
0068         basis = tangentorthobasis(problem.M, x, n);
0069     else
0070         n = numel(basis);
0071     end
0072
0073     % Create a store database and get a key for x
0074     storedb = StoreDB(1);
0075     key = storedb.getNewKey();
0076
0077     % Apply the Hessian at x to each basis vector
0078     Hbasis = cell(n, 1);
0079     for k = 1 : numel(Hbasis)
0080         Hbasis{k} = getHessian(problem, x, basis{k}, storedb, key);
0081     end
0082
0083     % H is the matrix which contains the inner products of
0084     % the ((basis vectors)) with the ((Hessian applied to basis vectors)).
0085     H = zeros(n);
0086     for i = 1 : n
0087         H(i, i) = problem.M.inner(x, basis{i}, Hbasis{i});
0088         for j = (i+1) : n
0089             H(i, j) = problem.M.inner(x, basis{i}, Hbasis{j});
0090             H(j, i) = H(i, j);
0091         end
0092     end
0093
0094 end```

Generated on Sun 05-Sep-2021 17:57:00 by m2html © 2005