Home > manopt > manifolds > positive > positivefactory.m

positivefactory

PURPOSE

Manifold of m-by-n matrices with positive entries; scale invariant metric

SYNOPSIS

function M = positivefactory(m, n)

DESCRIPTION

``` Manifold of m-by-n matrices with positive entries; scale invariant metric

function M = positivefactory(m)
function M = positivefactory(m, n)

A point X on the manifold M is represented as a matrix X of size mxn with
all individual entries real, strictly positive. By default, n = 1.

A tangent vector at X is represented as a matrix of the same size as X.
Entries of tangent vectors are free (in particular, not necessarily
positive.)

The Riemannian metric for each individual entry is the bi-invariant
metric for positive scalars, as a particular case of the bi-invariant
metric for positive definite matrices studied in Chapter 6 of the book

"Positive definite matrices" by Rajendra Bhatia,
Princeton University Press, 2007.

The Riemannian structure of M is obtained as the Cartesian product of the
geometry for mxn positive real numbers.

It should be stressed that matrices with one or more zero entries do not
belong to this manifold: they appear to be infinitely far away as a
result of the metric scaling like X.^(-1). Thus, if the solutions of an
optimization problem have entries equal to zero, these solutions are not
attainable on the manifold, which is likely to create serious numerical
issues. This geometry is best used when the solutions of the optimization
problem are indeed entry-wise positive, yet may have very different
scales (with some entries being very small, and some entries being very
large, relatively.)

CROSS-REFERENCE INFORMATION

This function calls:
• hashmd5 Computes the MD5 hash of input data.
• matrixlincomb Linear combination function for tangent vectors represented as matrices.
This function is called by:

SOURCE CODE

```0001 function M = positivefactory(m, n)
0002 % Manifold of m-by-n matrices with positive entries; scale invariant metric
0003 %
0004 % function M = positivefactory(m)
0005 % function M = positivefactory(m, n)
0006 %
0007 % A point X on the manifold M is represented as a matrix X of size mxn with
0008 % all individual entries real, strictly positive. By default, n = 1.
0009 %
0010 % A tangent vector at X is represented as a matrix of the same size as X.
0011 % Entries of tangent vectors are free (in particular, not necessarily
0012 % positive.)
0013 %
0014 % The Riemannian metric for each individual entry is the bi-invariant
0015 % metric for positive scalars, as a particular case of the bi-invariant
0016 % metric for positive definite matrices studied in Chapter 6 of the book
0017 %
0018 %    "Positive definite matrices" by Rajendra Bhatia,
0019 %    Princeton University Press, 2007.
0020 %
0021 % The Riemannian structure of M is obtained as the Cartesian product of the
0022 % geometry for mxn positive real numbers.
0023 %
0024 % It should be stressed that matrices with one or more zero entries do not
0025 % belong to this manifold: they appear to be infinitely far away as a
0026 % result of the metric scaling like X.^(-1). Thus, if the solutions of an
0027 % optimization problem have entries equal to zero, these solutions are not
0028 % attainable on the manifold, which is likely to create serious numerical
0029 % issues. This geometry is best used when the solutions of the optimization
0030 % problem are indeed entry-wise positive, yet may have very different
0031 % scales (with some entries being very small, and some entries being very
0032 % large, relatively.)
0033 %
0035
0036 % This file is part of Manopt: www.manopt.org.
0037 % Original author: Bamdev Mishra, Dec 03, 2017.
0038
0039     if ~exist('n', 'var') || isempty(n)
0040         n = 1;
0041     end
0042
0043     M.name = @() sprintf('Element-wise positive %dx%d matrices', m, n);
0044
0045     M.dim = @() m*n;
0046
0047     % The metric is the scale invariant metric for scalars.
0048     M.inner = @myinner;
0049     function innerproduct = myinner(X, eta, zeta)
0050         innerproduct = (eta(:)./X(:))'*(zeta(:)./X(:));
0051     end
0052
0053     M.norm = @(X, eta) sqrt(myinner(X, eta, eta));
0054
0055     M.dist = @(X, Y) norm(log(Y./X), 'fro');
0056
0057     M.typicaldist = @() sqrt(m*n);
0058
0061         eta = X.*(eta).*X;
0062     end
0063
0064     M.ehess2rhess = @ehess2rhess;
0065     function Hess = ehess2rhess(X, egrad, ehess, eta)
0066         % Directional derivatives of the Riemannian gradient
0067         Hess = X.*(ehess).*X + 2*(eta.*(egrad).*X);
0068
0069         % Correction factor for the non-constant metric
0070         Hess = Hess - (eta.*(egrad).*X);
0071     end
0072
0073     % Since this manifold is an open subset of R^(nxm), the tangent space
0074     % at any X on M is all of R^(nxm).
0075     M.proj = @(X, eta) eta;
0076
0077     M.tangent = M.proj;
0078     M.tangent2ambient = @(X, eta) eta;
0079
0080     M.retr = @exponential;
0081
0082     M.exp = @exponential;
0083     function Y = exponential(X, eta, t)
0084         if nargin < 3
0085             t = 1.0;
0086         end
0087         % It is unclear whether this is the numerically most stable way to
0088         % implement this operation. If you run into trouble with this
0089         % factory, please get in touch on the forum.
0090         Y = (X.*(exp((t*eta)./X)));
0091     end
0092
0093     M.log = @logarithm;
0094     function H = logarithm(X, Y)
0095         % Same comment about numerical stability as for exp.
0096         H = (X.*(log(Y./X)));
0097     end
0098
0099     M.hash = @(X) ['z' hashmd5(X(:))];
0100
0101     % Generate a random element-wise positive matrix following a
0102     % certain distribution. The particular choice of a distribution is of
0103     % course arbitrary, and specific applications might require different
0104     % ones.
0105     M.rand = @random;
0106     function X = random()
0107         X = exp(randn(m, n));
0108     end
0109
0110     % Generate a uniformly random unit-norm tangent vector at X.
0111     M.randvec = @randomvec;
0112     function eta = randomvec(X)
0113         eta = randn(m, n).*X;
0114         nrm = M.norm(X, eta);
0115         eta = eta / nrm;
0116     end
0117
0118     M.lincomb = @matrixlincomb;
0119
0120     M.zerovec = @(X) zeros(m, n);
0121
0122
0123     M.transp = @(X1, X2, eta) eta;
0124
0125     % For reference, a proper vector transport is given here, following
0126     % work by Sra and Hosseini: "Conic geometric optimisation on the
0127     % manifold of positive definite matrices".
0128     % This is not used by default. To force the use of this transport,
0129     % execute "M.transp = M.paralleltransp;" on your M returned by the
0130     % present factory.
0131     M.paralleltransp = @parallel_transport;
0132     function zeta = parallel_transport(X, Y, eta)
0133         E = sqrt(Y./X);
0134         zeta = E.*eta.*E;
0135     end
0136
0137     % vec and mat are not isometries, because of the unusual inner metric.
0138     M.vec = @(X, U) U(:);
0139     M.mat = @(X, u) reshape(u, m, n);
0140     M.vecmatareisometries = @() true;
0141
0142 end```

Generated on Mon 10-Sep-2018 11:48:06 by m2html © 2005