Home > manopt > tools > smallestinconvexhull.m

smallestinconvexhull

PURPOSE ^

Computes a minimal norm convex combination of given tangent vectors in Manopt.

SYNOPSIS ^

function [u_norm, coeffs, u] = smallestinconvexhull(M, x, U)

DESCRIPTION ^

 Computes a minimal norm convex combination of given tangent vectors in Manopt.

 function [u_norm, coeffs, u] = smallestinconvexhull(M, x, U)

 M is a manifold as returned by a Manopt factory.
 x is a point on this manifold.
 U is a cell containing N tangent vectors U{1} to U{N} at x.
 
 This function computes u, a tangent vector at x contained in the convex
 hull spanned by the N vectors U{i}, with minimal norm (according to the
 Riemannian metric on M). This is obtained by solving a convex quadratic
 program involving the Gram matrix of the given tangent vectors.
 The quadratic program is solved using Matlab's built-in quadprog,
 which requires the optimization toolbox. If this toolbox is not
 available, consider replacing with CVX for example.


 u_norm is the norm of the smallest vector u.
 coeffs is a vector of length N with entries in [0, 1] summing to 1.
 u is the sought vector: u = coeffs(1)*U{1} + ... + coeffs(N)*U{N}.

 Nicolas Boumal, Feb. 19, 2013
 Modified April 6, 2016 to work with Manopt.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [u_norm, coeffs, u] = smallestinconvexhull(M, x, U)
0002 % Computes a minimal norm convex combination of given tangent vectors in Manopt.
0003 %
0004 % function [u_norm, coeffs, u] = smallestinconvexhull(M, x, U)
0005 %
0006 % M is a manifold as returned by a Manopt factory.
0007 % x is a point on this manifold.
0008 % U is a cell containing N tangent vectors U{1} to U{N} at x.
0009 %
0010 % This function computes u, a tangent vector at x contained in the convex
0011 % hull spanned by the N vectors U{i}, with minimal norm (according to the
0012 % Riemannian metric on M). This is obtained by solving a convex quadratic
0013 % program involving the Gram matrix of the given tangent vectors.
0014 % The quadratic program is solved using Matlab's built-in quadprog,
0015 % which requires the optimization toolbox. If this toolbox is not
0016 % available, consider replacing with CVX for example.
0017 %
0018 %
0019 % u_norm is the norm of the smallest vector u.
0020 % coeffs is a vector of length N with entries in [0, 1] summing to 1.
0021 % u is the sought vector: u = coeffs(1)*U{1} + ... + coeffs(N)*U{N}.
0022 %
0023 % Nicolas Boumal, Feb. 19, 2013
0024 % Modified April 6, 2016 to work with Manopt.
0025 
0026 % This file is part of Manopt: www.manopt.org.
0027 % Original author: Nicolas Boumal, June 28, 2016.
0028 % Contributors:
0029 % Change log:
0030 %
0031 %   June 28, 2016 (NB):
0032 %       Adapted for Manopt from original code by same author (Feb. 19, 2013)
0033 
0034 % Example code: pick a manifold, a point, and a collection of tangent
0035 % vectors at that point, then get the smallest vector in the convex hull
0036 % of those:
0037 %
0038 % M = spherefactory(5);
0039 % x = M.rand();
0040 % N = 3;
0041 % U = cell(N,1);
0042 % for k = 1 : N, U{k} = M.randvec(x); end
0043 % [u_norm, coeffs, u] = smallestinconvexhull(M, x, U)
0044 
0045     % We simply need to solve the following quadratic program:
0046     % minimize ||u||^2 such that u = sum_i s_i U_i, 0 <= s_i <= 1
0047     %                            and sum_i s_i = 1
0048     %
0049     % This is equivalent to solving:
0050     %  min s'*G*s s.t. 0 <= s <= 1, s'*ones = 1, with G(i, j) = <U_i, U_j> (Gram matrix)
0051     % Then our solution is s_1 U_1 + ... + s_N U_N.
0052     
0053     
0054     % Compute the Gram matrix of the given tangent vectors
0055     N = numel(U);
0056     G = grammatrix(M, x, U);
0057     
0058     % Solve the quadratic program.
0059     % If the optimization toolbox is not available, consider replacing with
0060     % CVX.
0061     
0062     opts = optimset('Display', 'off');
0063     [s_opt, cost_opt] ...
0064           = quadprog(G, zeros(N, 1),     ...  % objective (squared norm)
0065                      [], [],             ...  % inequalities (none)
0066                      ones(1, N), 1,      ...  % equality (sum to 1)
0067                      zeros(N, 1),        ...  % lower bounds (s_i >= 0)
0068                      ones(N, 1),         ...  % upper bounds (s_i <= 1)
0069                      [],                 ...  % we do not specify an initial guess
0070                      opts);
0071 
0072     % Norm of the smallest tangent vector in the convex hull:
0073     u_norm = real(sqrt(2*cost_opt));
0074 
0075     % Keep track of optimal coefficients
0076     coeffs = s_opt;
0077     
0078     % If required, construct the vector explicitly.
0079     if nargout >= 3
0080         u = lincomb(M, x, U, coeffs);
0081     end
0082                  
0083 end

Generated on Sat 12-Nov-2016 14:11:22 by m2html © 2005