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

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:
• norm NORM Norm of a TT/MPS tensor.
• norm NORM Norm of a TT/MPS block-mu tensor.
• hashmd5 Computes the MD5 hash of input data.
• lyapunov_symmetric Solves AX + XA = C when A = A', as a pseudo-inverse.
• matrixlincomb Linear combination function for tangent vectors represented as matrices.
This function is called by:

## 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".
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".
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
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
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