0001 function M = spherecomplexfactory(n, m)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038 if ~exist('m', 'var')
0039 m = 1;
0040 end
0041
0042 if m == 1
0043 M.name = @() sprintf('Complex sphere S^%d', n-1);
0044 else
0045 M.name = @() sprintf('Unit F-norm %dx%d complex matrices', n, m);
0046 end
0047
0048 M.dim = @() 2*(n*m)-1;
0049
0050 M.inner = @(x, d1, d2) real(d1(:)'*d2(:));
0051
0052 M.norm = @(x, d) norm(d, 'fro');
0053
0054 M.dist = @(x, y) real(2*asin(.5*norm(x - y, 'fro')));
0055
0056 M.typicaldist = @() pi;
0057
0058 M.proj = @(x, d) reshape(d(:) - x(:)*(real(x(:)'*d(:))), n, m);
0059
0060
0061
0062 M.egrad2rgrad = M.proj;
0063
0064 M.ehess2rhess = @ehess2rhess;
0065 function rhess = ehess2rhess(x, egrad, ehess, u)
0066 rhess = M.proj(x, ehess) - real((x(:)'*egrad(:)))*u;
0067 end
0068
0069 M.tangent = M.proj;
0070
0071 M.exp = @exponential;
0072
0073 M.retr = @retraction;
0074
0075 M.log = @logarithm;
0076 function v = logarithm(x1, x2)
0077 v = M.proj(x1, x2 - x1);
0078 di = M.dist(x1, x2);
0079
0080 if di > 1e-6
0081 nv = norm(v, 'fro');
0082 v = v * (di / nv);
0083 end
0084 end
0085
0086 M.hash = @(x) ['z' hashmd5([real(x(:)) ; imag(x(:))])];
0087
0088 M.rand = @() random(n, m);
0089
0090 M.randvec = @(x) randomvec(n, m, x);
0091
0092 M.lincomb = @matrixlincomb;
0093
0094 M.zerovec = @(x) zeros(n, m);
0095
0096 M.transp = @(x1, x2, d) M.proj(x2, d);
0097
0098 M.pairmean = @pairmean;
0099 function y = pairmean(x1, x2)
0100 y = x1+x2;
0101 y = y / norm(y, 'fro');
0102 end
0103
0104 mn = m*n;
0105 M.vec = @(x, u_mat) [real(u_mat(:)) ; imag(u_mat(:))];
0106 M.mat = @(x, u_vec) reshape(u_vec(1:mn), m, n) + 1i*reshape(u_vec((mn+1):end), m, n);
0107 M.vecmatareisometries = @() true;
0108
0109 end
0110
0111
0112 function y = exponential(x, d, t)
0113
0114 if nargin == 2
0115
0116 td = d;
0117 else
0118 td = t*d;
0119 end
0120
0121 nrm_td = norm(td, 'fro');
0122
0123 if nrm_td > 0
0124 y = x*cos(nrm_td) + td*(sin(nrm_td)/nrm_td);
0125 else
0126 y = x;
0127 end
0128
0129 end
0130
0131
0132 function y = retraction(x, d, t)
0133
0134 if nargin == 2
0135 t = 1;
0136 end
0137
0138 y = x+t*d;
0139 y = y/norm(y, 'fro');
0140
0141 end
0142
0143
0144 function x = random(n, m)
0145
0146 x = randn(n, m) + 1i*randn(n, m);
0147 x = x/norm(x, 'fro');
0148
0149 end
0150
0151
0152 function d = randomvec(n, m, x)
0153
0154 d = randn(n, m) + 1i*randn(n, m);
0155 d = reshape(d(:) - x(:)*(real(x(:)'*d(:))), n, m);
0156 d = d / norm(d, 'fro');
0157
0158 end