Home > manopt > manifolds > symfixedrank > symfixedrankYYcomplexfactory.m

# symfixedrankYYcomplexfactory

## PURPOSE Manifold of n x n complex Hermitian pos. semidefinite matrices of rank k.

## SYNOPSIS function M = symfixedrankYYcomplexfactory(n, k)

## DESCRIPTION ``` Manifold of n x n complex Hermitian pos. semidefinite matrices of rank k.

function M = symfixedrankYYcomplexfactory(n, k)

Manifold of n-by-n complex Hermitian positive semidefinite matrices of
fixed rank k. This follows the quotient geometry described
in Sarod Yatawatta's 2013 paper:
"Radio interferometric calibration using a Riemannian manifold", ICASSP.

A point X on the manifold M is parameterized as YY^*, where
Y is a complex matrix of size nxk. For any point Y on the manifold M,
given any kxk complex unitary matrix U, we say Y*U  is equivalent to Y,
i.e., YY^* does not change. Therefore, M is the set of equivalence
classes and is a Riemannian quotient manifold C^{nk}/SU(k).
The metric is the usual real-trace inner product, that is,
it is the usual metric for the complex plane identified with R^2.

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.

The geometry is based on the following papers (and references therein).
Please cite the Manopt paper as well as the research papers:

@INPROCEEDINGS{Yatawatta2013A,
author={Yatawatta, S.},
booktitle={Acoustics, Speech and Signal Processing (ICASSP), 2013 IEEE International Conference on},
title={Radio interferometric calibration using a {R}iemannian manifold},
year={2013},
month={May},
pages={3866--3870},
doi={10.1109/ICASSP.2013.6638382},
ISSN={1520-6149},
}

@article{Yatawatta2013B,
author = {Yatawatta, S.},
title = {On the interpolation of calibration solutions obtained in radio interferometry},
volume = {428},
number = {1},
pages = {828--833},
year = {2013},
doi = {10.1093/mnras/sts069},
journal = {Monthly Notices of the Royal Astronomical Society}
}

## CROSS-REFERENCE INFORMATION This function calls:
• 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 = symfixedrankYYcomplexfactory(n, k)
0002 % Manifold of n x n complex Hermitian pos. semidefinite matrices of rank k.
0003 %
0004 % function M = symfixedrankYYcomplexfactory(n, k)
0005 %
0006 % Manifold of n-by-n complex Hermitian positive semidefinite matrices of
0007 % fixed rank k. This follows the quotient geometry described
0008 % in Sarod Yatawatta's 2013 paper:
0009 % "Radio interferometric calibration using a Riemannian manifold", ICASSP.
0010 %
0012 %
0013 % A point X on the manifold M is parameterized as YY^*, where
0014 % Y is a complex matrix of size nxk. For any point Y on the manifold M,
0015 % given any kxk complex unitary matrix U, we say Y*U  is equivalent to Y,
0016 % i.e., YY^* does not change. Therefore, M is the set of equivalence
0017 % classes and is a Riemannian quotient manifold C^{nk}/SU(k).
0018 % The metric is the usual real-trace inner product, that is,
0019 % it is the usual metric for the complex plane identified with R^2.
0020 %
0021 % Notice that this manifold is not complete: if optimization leads Y to be
0022 % rank-deficient, the geometry will break down. Hence, this geometry should
0023 % only be used if it is expected that the points of interest will have rank
0024 % exactly k. Reduce k if that is not the case.
0025 %
0026 % The geometry is based on the following papers (and references therein).
0027 % Please cite the Manopt paper as well as the research papers:
0028 %
0029 % @INPROCEEDINGS{Yatawatta2013A,
0030 %  author={Yatawatta, S.},
0031 %  booktitle={Acoustics, Speech and Signal Processing (ICASSP), 2013 IEEE International Conference on},
0032 %  title={Radio interferometric calibration using a {R}iemannian manifold},
0033 %  year={2013},
0034 %  month={May},
0035 %  pages={3866--3870},
0036 %  doi={10.1109/ICASSP.2013.6638382},
0037 %  ISSN={1520-6149},
0038 % }
0039 %
0040 % @article{Yatawatta2013B,
0041 %  author = {Yatawatta, S.},
0042 %  title = {On the interpolation of calibration solutions obtained in radio interferometry},
0043 %  volume = {428},
0044 %  number = {1},
0045 %  pages = {828--833},
0046 %  year = {2013},
0047 %  doi = {10.1093/mnras/sts069},
0048 %  journal = {Monthly Notices of the Royal Astronomical Society}
0049 % }
0050 %
0052
0053
0054 % This file is part of Manopt: www.manopt.org.
0055 % Original author: Sarod Yatawatta, June 29, 2015.
0056 % Contributors: Bamdev Mishra.
0057 % Change log:
0058 %
0059 %   June 28, 2016 (NB):
0060 %       Metric scaled down by factor 2 to match the metric used in
0061 %       euclideancomplexfactory.
0062 %
0063 %   Apr. 18, 2018 (NB):
0064 %       Removed lyap depedency.
0065 %
0066 %   Sep.  6, 2018 (NB):
0067 %       Removed M.exp() as it was not implemented.
0068
0069     M.name = @() sprintf('YY'' quotient manifold of Hermitian %dx%d complex matrices of rank %d.', n, n, k);
0070
0071     M.dim = @() 2*k*n - k*k; % SY: dim of ambient space (2*k*n) - dim of kxk unitary matrix  (k^2).
0072
0073     % Euclidean metric on the total space.
0074     % BM: equivalent to real(trace(eta'*zeta)), but more efficient.
0075     M.inner = @(Y, eta, zeta) real(eta(:)'*zeta(:));
0076
0077     M.norm = @(Y, eta) sqrt(M.inner(Y, eta, eta));
0078
0079     % Find unitary U to minimize ||Y - Z*U||,
0080     % i.e., the Procrustes problem, with svd(Y'*Z).
0081     M.dist = @(Y, Z) distance;
0082     function distval = distance(Y, Z)
0083         [u, ignore, v] = svd(Z'*Y); %#ok<ASGLU>
0084         E = Y - Z*u*v'; % SY: checked.
0085         distval = real(E(:)'*E(:));
0086     end
0087
0088     M.typicaldist = @() 10*k; % BM: To do.
0089
0090     M.proj = @projection;
0091     function etaproj = projection(Y, eta)
0092         % Projection onto the horizontal space
0093         xx = Y'*Y;
0094         rr = Y'*eta - eta'*Y;
0095         Omega = lyapunov_symmetric(xx, rr);
0096         etaproj = eta - Y*Omega;
0097     end
0098
0099     M.tangent = M.proj;
0100     M.tangent2ambient = @(Y, eta) eta;
0101
0102     M.retr = @retraction;
0103     function Ynew = retraction(Y, eta, t)
0104         if nargin < 3
0105             t = 1.0;
0106         end
0107         Ynew = Y + t*eta;
0108     end
0109
0110
0112     M.ehess2rhess = @(Y, egrad, ehess, U) M.proj(Y, ehess);
0113
0114
0115     % Notice that the hash of two equivalent points will be different...
0116     M.hash = @(Y) ['z' hashmd5([real(Y(:)); imag(Y(:))])];
0117
0118     M.rand = @random;
0119     function Y = random()
0120         Y = randn(n, k) + 1i*randn(n,k);
0121     end
0122
0123     M.randvec = @randomvec;
0124     function eta = randomvec(Y)
0125         eta = randn(n, k) + 1i*randn(n,k);
0126         eta = projection(Y, eta);
0127         nrm = M.norm(Y, eta);
0128         eta = eta / nrm;
0129     end
0130
0131     M.lincomb = @matrixlincomb;
0132
0133     M.zerovec = @(Y) zeros(n, k);
0134
0135     M.transp = @(Y1, Y2, d) projection(Y2, d);
0136
0137     M.vec = @(Y, u_mat) [real(u_mat(:)); imag(u_mat(:))];
0138     M.mat = @(Y, u_vec) reshape(u_vec(1 : n*k), [n, k]) + 1i*reshape(u_vec(n*k + 1: end), [n, k]);
0139     M.vecmatareisometries = @() true;
0140
0141 end```

Generated on Fri 30-Sep-2022 13:18:25 by m2html © 2005