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