Home > manopt > manifolds > complexcircle > complexcirclefactory.m

complexcirclefactory

PURPOSE ^

Returns a manifold struct to optimize over unit-modulus complex numbers.

SYNOPSIS ^

function M = complexcirclefactory(n)

DESCRIPTION ^

 Returns a manifold struct to optimize over unit-modulus complex numbers.

 function M = complexcirclefactory()
 function M = complexcirclefactory(n)

 Description of vectors z in C^n (complex) such that each component z(i)
 has unit modulus. The manifold structure is the Riemannian submanifold
 structure from the embedding space R^2 x ... x R^2, i.e., the complex
 circle is identified with the unit circle in the real plane.

 By default, n = 1.

 See also spherecomplexfactory

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function M = complexcirclefactory(n)
0002 % Returns a manifold struct to optimize over unit-modulus complex numbers.
0003 %
0004 % function M = complexcirclefactory()
0005 % function M = complexcirclefactory(n)
0006 %
0007 % Description of vectors z in C^n (complex) such that each component z(i)
0008 % has unit modulus. The manifold structure is the Riemannian submanifold
0009 % structure from the embedding space R^2 x ... x R^2, i.e., the complex
0010 % circle is identified with the unit circle in the real plane.
0011 %
0012 % By default, n = 1.
0013 %
0014 % See also spherecomplexfactory
0015 
0016 % This file is part of Manopt: www.manopt.org.
0017 % Original author: Nicolas Boumal, Dec. 30, 2012.
0018 % Contributors:
0019 % Change log:
0020 %
0021 %   July 7, 2014 (NB): Added ehess2rhess function.
0022 %
0023 %   Sep. 3, 2014 (NB): Correction to the dist function (extract real part).
0024 %
0025 %   April 13, 2015 (NB): Fixed logarithm.
0026 %
0027 %   Oct. 8, 2016 (NB)
0028 %       Code for exponential was simplified to only treat the zero vector
0029 %       as a particular case.
0030 %
0031 %   July 20, 2017 (NB)
0032 %       The distance function is now even more accurate. Improved logarithm
0033 %       accordingly.
0034     
0035     if ~exist('n', 'var')
0036         n = 1;
0037     end
0038 
0039     M.name = @() sprintf('Complex circle (S^1)^%d', n);
0040     
0041     M.dim = @() n;
0042     
0043     M.inner = @(z, v, w) real(v'*w);
0044     
0045     M.norm = @(x, v) norm(v);
0046     
0047     M.dist = @(x, y) norm(real(2*asin(.5*abs(x - y))));
0048     
0049     M.typicaldist = @() pi*sqrt(n);
0050     
0051     M.proj = @(z, u) u - real( conj(u) .* z ) .* z;    
0052     
0053     M.tangent = M.proj;
0054     
0055     % For Riemannian submanifolds, converting a Euclidean gradient into a
0056     % Riemannian gradient amounts to an orthogonal projection.
0057     M.egrad2rgrad = M.proj;
0058     
0059     M.ehess2rhess = @ehess2rhess;
0060     function rhess = ehess2rhess(z, egrad, ehess, zdot)
0061         rhess = M.proj(z, ehess - real(z.*conj(egrad)).*zdot);
0062     end
0063     
0064     M.exp = @exponential;
0065     function y = exponential(z, v, t)
0066         
0067         if nargin == 2
0068             % t = 1;
0069             tv = v;
0070         else
0071             tv = t*v;
0072         end
0073 
0074         y = zeros(n, 1);
0075 
0076         nrm_tv = abs(tv);
0077         
0078         % We need to be careful for zero steps.
0079         mask = (nrm_tv > 0);
0080         y(mask) = z(mask).*cos(nrm_tv(mask)) + ...
0081                   tv(mask).*(sin(nrm_tv(mask))./nrm_tv(mask));
0082         y(~mask) = z(~mask);
0083         
0084     end
0085     
0086     M.retr = @retraction;
0087     function y = retraction(z, v, t)
0088         if nargin == 2
0089             % t = 1;
0090             tv = v;
0091         else
0092             tv = t*v;
0093         end
0094         y = sign(z+tv);
0095     end
0096 
0097     M.log = @logarithm;
0098     function v = logarithm(x1, x2)
0099         v = M.proj(x1, x2 - x1);
0100         di = real(2*asin(.5*abs(x1 - x2)));
0101         nv = abs(v);
0102         factors = di ./ nv;
0103         factors(di <= 1e-10) = 1;
0104         v = v .* factors;
0105     end
0106     
0107     M.hash = @(z) ['z' hashmd5( [real(z(:)) ; imag(z(:))] ) ];
0108     
0109     M.rand = @random;
0110     function z = random()
0111         z = sign(randn(n, 1) + 1i*randn(n, 1));
0112     end
0113     
0114     M.randvec = @randomvec;
0115     function v = randomvec(z)
0116         % i*z(k) is a basis vector of the tangent vector to the k-th circle
0117         v = randn(n, 1) .* (1i*z);
0118         v = v / norm(v);
0119     end
0120     
0121     M.lincomb = @matrixlincomb;
0122     
0123     M.zerovec = @(x) zeros(n, 1);
0124     
0125     M.transp = @(x1, x2, d) M.proj(x2, d);
0126     
0127     M.pairmean = @pairmean;
0128     function z = pairmean(z1, z2)
0129         z = sign(z1+z2);
0130     end
0131 
0132     M.vec = @(x, u_mat) [real(u_mat) ; imag(u_mat)];
0133     M.mat = @(x, u_vec) u_vec(1:n) + 1i*u_vec((n+1):end);
0134     M.vecmatareisometries = @() true;
0135 
0136 end

Generated on Fri 08-Sep-2017 12:43:19 by m2html © 2005