Grassmannian manifold

Manifolds.Grassmann β€” Type
Grassmann{n,k,𝔽} <: AbstractEmbeddedManifold{𝔽,DefaultIsometricEmbeddingType}

The Grassmann manifold $\operatorname{Gr}(n,k)$ consists of all subspaces spanned by $k$ linear independent vectors $𝔽^n$, where $𝔽 ∈ \{ℝ, β„‚\}$ is either the real- (or complex-) valued vectors. This yields all $k$-dimensional subspaces of $ℝ^n$ for the real-valued case and all $2k$-dimensional subspaces of $β„‚^n$ for the second.

The manifold can be represented as

\[\operatorname{Gr}(n,k) := \bigl\{ \operatorname{span}(p) : p ∈ 𝔽^{n Γ— k}, p^\mathrm{H}p = I_k\},\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transpose or Hermitian and $I_k$ is the $k Γ— k$ identity matrix. This means, that the columns of $p$ form an unitary basis of the subspace, that is a point on $\operatorname{Gr}(n,k)$, and hence the subspace can actually be represented by a whole equivalence class of representers. Another interpretation is, that

\[\operatorname{Gr}(n,k) = \operatorname{St}(n,k) / \operatorname{O}(k),\]

i.e the Grassmann manifold is the quotient of the Stiefel manifold and the orthogonal group $\operatorname{O}(k)$ of orthogonal $k Γ— k$ matrices.

The tangent space at a point (subspace) $x$ is given by

\[T_x\mathrm{Gr}(n,k) = \bigl\{ X ∈ 𝔽^{n Γ— k} : X^{\mathrm{H}}p + p^{\mathrm{H}}X = 0_{k} \bigr\},\]

where $0_k$ is the $k Γ— k$ zero matrix and $\overline{\cdot}$ the complex conjugate.

Note that a point $p ∈ \operatorname{Gr}(n,k)$ might be represented by different matrices (i.e. matrices with unitary column vectors that span the same subspace). Different representations of $p$ also lead to different representation matrices for the tangent space $T_p\mathrm{Gr}(n,k)$

The manifold is named after Hermann G. Graßmann (1809-1877).

Constructor

Grassmann(n,k,field=ℝ)

Generate the Grassmann manifold $\operatorname{Gr}(n,k)$, where the real-valued case field = ℝ is the default.

source
Base.exp β€” Method
exp(M::Grassmann, p, X)

Compute the exponential map on the Grassmann M$= \mathrm{Gr}(n,k)$ starting in p with tangent vector (direction) X. Let $X = USV$ denote the SVD decomposition of $X$. Then the exponential map is written using

\[z = p V\cos(S)V^\mathrm{H} + U\sin(S)V^\mathrm{H},\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian and the cosine and sine are applied element wise to the diagonal entries of $S$. A final QR decomposition $z=QR$ is performed for numerical stability reasons, yielding the result as

\[\exp_p X = Q.\]
source
Base.log β€” Method
log(M::Grassmann, p, q)

Compute the logarithmic map on the Grassmann M$ = \mathcal M=\mathrm{Gr}(n,k)$, i.e. the tangent vector X whose corresponding geodesic starting from p reaches q after time 1 on M. The formula reads

\[\log_p q = V\cdot \operatorname{atan}(S) \cdot U^\mathrm{H},\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian. The matrices $U$ and $V$ are the unitary matrices, and $S$ is the diagonal matrix containing the singular values of the SVD-decomposition

\[USV = (q^\mathrm{H}p)^{-1} ( q^\mathrm{H} - q^\mathrm{H}pp^\mathrm{H}).\]

In this formula the $\operatorname{atan}$ is meant elementwise.

source
Manifolds.uniform_distribution β€” Method
uniform_distribution(M::Grassmann{n,k,ℝ}, p)

Uniform distribution on given (real-valued) Grassmann M. Specifically, this is the normalized Haar measure on M. Generated points will be of similar type as p.

The implementation is based on Section 2.5.1 in [Chikuse2003]; see also Theorem 2.2.2(iii) in [Chikuse2003].

source
ManifoldsBase.check_manifold_point β€” Method
check_manifold_point(M::Grassmann{n,k,𝔽}, p)

Check whether p is representing a point on the Grassmann M, i.e. its a n-by-k matrix of unitary column vectors and of correct eltype with respect to 𝔽.

source
ManifoldsBase.check_tangent_vector β€” Method
check_tangent_vector(M::Grassmann{n,k,𝔽}, p, X; check_base_point = true, kwargs...)

Check whether X is a tangent vector in the tangent space of p on the Grassmann M, i.e. that X is of size and type as well as that

\[ p^{\mathrm{H}}X + \overline{X^{\mathrm{H}}p} = 0_k,\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transpose or Hermitian, $\overline{\cdot}$ the (elementwise) complex conjugate, and $0_k$ the $k Γ— k$ zero matrix. The optional parameter check_base_point indicates, whether to call check_manifold_point for p.

source
ManifoldsBase.distance β€” Method
distance(M::Grassmann, p, q)

Compute the Riemannian distance on Grassmann manifold M$= \mathrm{Gr}(n,k)$.

Let $USV = p^\mathrm{H}q$ denote the SVD decomposition of $p^\mathrm{H}q$, where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian. Then the distance is given by

\[d_{\mathrm{Gr}(n,k)}(p,q) = \operatorname{norm}(\operatorname{Re}(b)).\]

where

\[b_{i}=\begin{cases} 0 & \text{if} \; S_i β‰₯ 1\\ \arccos(S_i) & \, \text{if} \; S_i<1. \end{cases}\]
source
ManifoldsBase.inner β€” Method
inner(M::Grassmann, p, X, Y)

Compute the inner product for two tangent vectors X, Y from the tangent space of p on the Grassmann manifold M. The formula reads

\[g_p(X,Y) = \operatorname{tr}(X^{\mathrm{H}}Y),\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

source
ManifoldsBase.inverse_retract β€” Method
inverse_retract(M::Grassmann, p, q, ::PolarInverseRetraction)

Compute the inverse retraction for the PolarRetraction, on the Grassmann manifold M, i.e.,

\[\operatorname{retr}_p^{-1}q = q*(p^\mathrm{H}q)^{-1} - p,\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

source
ManifoldsBase.inverse_retract β€” Method
inverse_retract(M, p, q, ::QRInverseRetraction)

Compute the inverse retraction for the QRRetraction, on the Grassmann manifold M, i.e.,

\[\operatorname{retr}_p^{-1}q = q(p^\mathrm{H}q)^{-1} - p,\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

source
ManifoldsBase.project β€” Method
project(M::Grassmann, p, X)

Project the n-by-k X onto the tangent space of p on the Grassmann M, which is computed by

\[\operatorname{proj_p}(X) = X - pp^{\mathrm{H}}X,\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

source
ManifoldsBase.retract β€” Method
retract(M::Grassmann, p, X, ::PolarRetraction)

Compute the SVD-based retraction PolarRetraction on the Grassmann M. With $USV = p + X$ the retraction reads

\[\operatorname{retr}_p X = UV^\mathrm{H},\]

where $\cdot^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

source
ManifoldsBase.retract β€” Method
retract(M::Grassmann, p, X, ::QRRetraction )

Compute the QR-based retraction QRRetraction on the Grassmann M. With $QR = p + X$ the retraction reads

\[\operatorname{retr}_p X = QD,\]

where D is a $m Γ— n$ matrix with

\[D = \operatorname{diag}( \operatorname{sgn}(R_{ii}+0,5)_{i=1}^n ).\]
source
ManifoldsBase.vector_transport_to β€” Method
vector_transport_to(M::Grassmann,p,X,q,::ProjectionTransport)

compute the projection based transport on the Grassmann M by interpreting X from the tangent space at p as a point in the embedding and projecting it onto the tangent space at q.

source