Home > manopt > manifolds > symfixedrank > symfixedrankYYfactory.m

symfixedrankYYfactory

PURPOSE ^

Manifold of n-by-n symmetric positive semidefinite matrices of rank k.

SYNOPSIS ^

function M = symfixedrankYYfactory(n, k)

DESCRIPTION ^

 Manifold of n-by-n symmetric positive semidefinite matrices of rank k.

 function M = symfixedrankYYfactory(n, k)

 A point X on the manifold is parameterized as YY^T where Y is a matrix of
 size nxk. As such, X is symmetric, positive semidefinite. We restrict to
 full-rank Y's, such that X has rank exactly k. The point X is numerically
 represented by Y (this is more efficient than working with X, which may
 be big). Tangent vectors are represented as matrices of the same size as
 Y, call them Ydot, so that Xdot = Y Ydot' + Ydot Y'. The metric is the
 canonical Euclidean metric on Y.
 
 Since for any orthogonal Q of size k, it holds that (YQ)(YQ)' = YY',
 we "group" all matrices of the form YQ in an equivalence class. The set
 of equivalence classes is a Riemannian quotient manifold, implemented
 here.

 Notice that this manifold is not complete: if optimization leads Y to be
 rank deficient, the geometry will break down. Hence, this geometry should
 only be used if it is expected that the points of interest will have rank
 exactly k. Reduce k if that is not the case.
 
 An alternative, complete, geometry for positive semidefinite matrices of
 rank k is described in Bonnabel and Sepulchre 2009, "Riemannian Metric
 and Geometric Mean for Positive Semidefinite Matrices of Fixed Rank",
 SIAM Journal on Matrix Analysis and Applications.


 The geometry here implemented is the simplest case of the 2010 paper:
 M. Journee, P.-A. Absil, F. Bach and R. Sepulchre,
 "Low-Rank Optimization on the Cone of Positive Semidefinite Matrices".
 Paper link: http://www.di.ens.fr/~fbach/journee2010_sdp.pdf
 The expressions for the distance and the logarithm come from the 2018
 preprint:
 Estelle Massart, P.-A. Absil, 
 "Quotient geometry with simple geodesics for the manifold of fixed-rank
 positive-semidefinite matrices".
 Paper link: https://sites.uclouvain.be/absil/2018-06/quotient_tech_report.pdf
 
 
 Please cite the Manopt paper as well as the research papers:
     @Article{journee2010low,
       Title   = {Low-rank optimization on the cone of positive semidefinite matrices},
       Author  = {Journ{\'e}e, M. and Bach, F. and Absil, P.-A. and Sepulchre, R.},
       Journal = {SIAM Journal on Optimization},
       Year    = {2010},
       Number  = {5},
       Pages   = {2327--2351},
       Volume  = {20},
       Doi     = {10.1137/080731359}
     }

     @TechReport{massart2018quotient,
       author        = {Massart, Estelle and Absil, P.-A.},
       title         = {Quotient geometry with simple geodesics for the manifold of fixed-rank positive-semidefinite matrices},
       institution   = {UCLouvain},
       year          = {2018},
       number        = {UCL-INMA-2018.06-v2},
       note          = {Preprint: \url{http://sites.uclouvain.be/absil/2018.06}},
     }

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function M = symfixedrankYYfactory(n, k)
0002 % Manifold of n-by-n symmetric positive semidefinite matrices of rank k.
0003 %
0004 % function M = symfixedrankYYfactory(n, k)
0005 %
0006 % A point X on the manifold is parameterized as YY^T where Y is a matrix of
0007 % size nxk. As such, X is symmetric, positive semidefinite. We restrict to
0008 % full-rank Y's, such that X has rank exactly k. The point X is numerically
0009 % represented by Y (this is more efficient than working with X, which may
0010 % be big). Tangent vectors are represented as matrices of the same size as
0011 % Y, call them Ydot, so that Xdot = Y Ydot' + Ydot Y'. The metric is the
0012 % canonical Euclidean metric on Y.
0013 %
0014 % Since for any orthogonal Q of size k, it holds that (YQ)(YQ)' = YY',
0015 % we "group" all matrices of the form YQ in an equivalence class. The set
0016 % of equivalence classes is a Riemannian quotient manifold, implemented
0017 % here.
0018 %
0019 % Notice that this manifold is not complete: if optimization leads Y to be
0020 % rank deficient, the geometry will break down. Hence, this geometry should
0021 % only be used if it is expected that the points of interest will have rank
0022 % exactly k. Reduce k if that is not the case.
0023 %
0024 % An alternative, complete, geometry for positive semidefinite matrices of
0025 % rank k is described in Bonnabel and Sepulchre 2009, "Riemannian Metric
0026 % and Geometric Mean for Positive Semidefinite Matrices of Fixed Rank",
0027 % SIAM Journal on Matrix Analysis and Applications.
0028 %
0029 %
0030 % The geometry here implemented is the simplest case of the 2010 paper:
0031 % M. Journee, P.-A. Absil, F. Bach and R. Sepulchre,
0032 % "Low-Rank Optimization on the Cone of Positive Semidefinite Matrices".
0033 % Paper link: http://www.di.ens.fr/~fbach/journee2010_sdp.pdf
0034 % The expressions for the distance and the logarithm come from the 2018
0035 % preprint:
0036 % Estelle Massart, P.-A. Absil,
0037 % "Quotient geometry with simple geodesics for the manifold of fixed-rank
0038 % positive-semidefinite matrices".
0039 % Paper link: https://sites.uclouvain.be/absil/2018-06/quotient_tech_report.pdf
0040 %
0041 %
0042 % Please cite the Manopt paper as well as the research papers:
0043 %     @Article{journee2010low,
0044 %       Title   = {Low-rank optimization on the cone of positive semidefinite matrices},
0045 %       Author  = {Journ{\'e}e, M. and Bach, F. and Absil, P.-A. and Sepulchre, R.},
0046 %       Journal = {SIAM Journal on Optimization},
0047 %       Year    = {2010},
0048 %       Number  = {5},
0049 %       Pages   = {2327--2351},
0050 %       Volume  = {20},
0051 %       Doi     = {10.1137/080731359}
0052 %     }
0053 %
0054 %     @TechReport{massart2018quotient,
0055 %       author        = {Massart, Estelle and Absil, P.-A.},
0056 %       title         = {Quotient geometry with simple geodesics for the manifold of fixed-rank positive-semidefinite matrices},
0057 %       institution   = {UCLouvain},
0058 %       year          = {2018},
0059 %       number        = {UCL-INMA-2018.06-v2},
0060 %       note          = {Preprint: \url{http://sites.uclouvain.be/absil/2018.06}},
0061 %     }
0062 %
0063 
0064 % See also: elliptopefactory spectrahedronfactory symfixedrankYYcomplexfactory
0065 
0066 % This file is part of Manopt: www.manopt.org.
0067 % Original author: Bamdev Mishra, Dec. 30, 2012.
0068 % Contributors: Estelle Massart
0069 % Change log:
0070 %
0071 %   July 10, 2013 (NB):
0072 %       Added vec, mat, tangent, tangent2ambient ;
0073 %       Correction for the dimension of the manifold.
0074 %
0075 %   Apr.  2, 2015 (NB):
0076 %       Replaced trace(A'*B) by A(:)'*B(:) (equivalent but faster).
0077 %
0078 %   Apr. 17, 2018 (NB):
0079 %       Removed dependence on lyap.
0080 %
0081 %   Sep.  6, 2018 (NB):
0082 %       Removed M.exp() as it was not implemented.
0083 %
0084 %   June 7, 2019  (EM):
0085 %       Added M.dist, M.exp, M.log and M.invretr.
0086 
0087     M.name = @() sprintf('YY'' quotient manifold of %dx%d psd matrices of rank %d', n, k);
0088 
0089     M.dim = @() k*n - k*(k-1)/2;
0090 
0091     % Euclidean metric on the total space
0092     M.inner = @(Y, eta, zeta) eta(:)'*zeta(:);
0093 
0094     M.norm = @(Y, eta) sqrt(M.inner(Y, eta, eta));
0095 
0096     M.dist = @(Y, Z) norm(logarithm(Y,Z),'fro');
0097 
0098     M.typicaldist = @() 10*k;
0099 
0100     M.proj = @projection;
0101     function etaproj = projection(Y, eta)
0102         % Projection onto the horizontal space
0103         YtY = Y'*Y;
0104         SS = YtY;
0105         AS = Y'*eta - eta'*Y;
0106         % Omega = lyap(SS, -AS);
0107         Omega = lyapunov_symmetric(SS, AS);
0108         etaproj = eta - Y*Omega;
0109     end
0110 
0111     M.tangent = M.proj;
0112     M.tangent2ambient = @(Y, eta) eta;
0113 
0114     M.exp = @exponential;
0115     function Ynew = exponential(Y, eta, t)
0116         if nargin < 3
0117             t = 1.0;
0118         end
0119         Ynew = Y + t*eta;
0120     end
0121 
0122     M.retr = M.exp;
0123     
0124     M.log = @logarithm;
0125     function eta = logarithm(Y, Z)
0126         YtZ = Y'*Z;
0127         [U, ~, V] = svd(YtZ);
0128         Qt = V*U';
0129         eta = Z*Qt - Y;
0130     end
0131 
0132     M.invretr = M.log;
0133 
0134     M.egrad2rgrad = @(Y, eta) eta;
0135     M.ehess2rhess = @(Y, egrad, ehess, U) M.proj(Y, ehess);
0136 
0137     % Notice that the hash of two equivalent points will be different...
0138     M.hash = @(Y) ['z' hashmd5(Y(:))];
0139 
0140     M.rand = @random;
0141     function Y = random()
0142         Y = randn(n, k);
0143     end
0144 
0145     M.randvec = @randomvec;
0146     function eta = randomvec(Y)
0147         eta = randn(n, k);
0148         eta = projection(Y, eta);
0149         nrm = M.norm(Y, eta);
0150         eta = eta / nrm;
0151     end
0152 
0153     M.lincomb = @matrixlincomb;
0154 
0155     M.zerovec = @(Y) zeros(n, k);
0156 
0157     M.transp = @(Y1, Y2, d) projection(Y2, d);
0158         
0159     M.vec = @(Y, u_mat) u_mat(:);
0160     M.mat = @(Y, u_vec) reshape(u_vec, [n, k]);
0161     M.vecmatareisometries = @() true;
0162 
0163 end

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