Grassmannian manifold

Grassmann{𝔽, T} <: AbstractDecoratorManifold{𝔽}

The Grassmann manifold $\mathrm{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

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

where $⋅^{\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

\[\mathrm{Gr}(n,k) = \mathrm{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. Note that it doesn't matter whether we start from the Euclidean or canonical metric on the Stiefel manifold, the resulting quotient metric on Grassmann is the same.

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

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

where $0_k$ is the $k×k$ zero matrix. Note that this is a subspace (called the horizontal space) of the tangent space of the Stiefel manifold, whose tangent space requires $p^{\mathrm{H}}X = -p^{\mathrm{H}}X$, which is here fulfilled for all $X$ automatically, even for any “rotated basis” $pQ$, $Q ∈ \operatorname{O}(k)$

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).

Points can also be represented as orthogonal projectors, where

\[\operatorname{Gr}(n,k) := \bigl\{ p \in \mathbb 𝔽^{n×n} : p = p^{\mathrm{H}, p^2 = p, \operatorname{rank}(p) = k\},\]

with tangent space

\[T_p\mathrm{Gr}(n,k) = \bigl\{ X ∈ \mathbb 𝔽^{n×n} : X=X^{\mathrm{H}} \text{ and } X = pX+Xp \bigr\},\]

see also ProjectorPoint and ProjectorTangentVector.

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

A good overview can be found in[BZA20].

Constructor

Grassmann(n, k, field=ℝ, parameter::Symbol=:type)

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

GrassmannAtlas()

The standard atlas of the real Grassmann manifold $\mathrm{Gr}(n,k)$. Its $\binom{n}{k}$ charts are indexed by ordered tuples i of k row indices. The chart indexed by i contains the subspaces whose corresponding k-by-k minor in the rows i is invertible.

For a Stiefel representative p, let pᵢ denote this minor and let j be the complementary row indices. The coordinate map is

\[\varphi_i([p]) = \operatorname{vec}\left((p p_i^{-1})_j\right).\]

Its inverse inserts the reshaped coordinates into the rows j, inserts the identity into the rows i, and orthonormalizes the resulting matrix.

The construction is similar to the one presented in [UM14].

convert(::Type{ProjectorPoint}, p::AbstractMatrix)

Convert a point p on Stiefel that also represents a point (i.e. subspace) on Grassmann to a projector representation of said subspace, i.e. compute the canonical_project! for

\[ π^{\mathrm{SG}}(p) = pp^{\mathrm{T}}.\]

convert(::Type{ProjectorPoint}, ::StiefelPoint)

Convert a point p on Stiefel that also represents a point (i.e. subspace) on Grassmann to a projector representation of said subspace, i.e. compute the canonical_project! for

\[ π^{\mathrm{SG}}(p) = pp^{\mathrm{T}}.\]

isapprox(M::Grassmann, p, X, Y; kwargs...)

Test whether two tangent vectors X, Y at p on the Grassmann manifold M are (approximately) the same, that is whether $\lVert X - Y \rVert_p \approx 0$ with the classical isapprox. The kwargs... are passed on to this; atol defaults to sqrt(max_eps(X, Y)).

isapprox(M::Grassmann, p, q; kwargs...)

Test whether two points p, q on the Grassmann manifold M are (approximately) the same. Since this is a quotient manifold, the test is performed by checking their distance if they are not equal, that is $d_{\mathcal M}(p,q) \approx 0$ with the classical isapprox. The kwargs... are passed on to this; atol defaults to sqrt(max_eps(p, q)).

affine_connection!(M::Grassmann, Zc, A::GrassmannAtlas, i, a, Xc, Yc)

Store the Levi-Civita covariant derivative of Yc in direction Xc in Zc, using the standard GrassmannAtlas. Writing Z = reshape(a, n-k, k), U = reshape(Xc, n-k, k), V = reshape(Yc, n-k, k), and G = I + ZᵀZ, the coordinate expression is

\[\nabla_U V = -UG^{-1}Z^\mathsf{T}V - VG^{-1}Z^\mathsf{T}U.\]

det_local_metric(M::Grassmann, A::GrassmannAtlas, i::AbstractVector, a)

Return the determinant of the local metric in the standard GrassmannAtlas at coordinates a in chart i.

inverse_chart_injectivity_radius(M::Grassmann, A::GrassmannAtlas, i)

Return the injectivity radius of an affine chart in the standard GrassmannAtlas, which is infinite.

change_metric(M::Grassmann, ::EuclideanMetric, p, X)

Change X to the corresponding vector with respect to the metric of the Grassmann M, which is just the identity, since the manifold is isometrically embedded.

change_representer(M::Grassmann, ::EuclideanMetric, p, X)

Change X to the corresponding representer of a cotangent vector at p. Since the Grassmann manifold M is isometrically embedded, this is the identity.

default_retraction_method(M::Grassmann)
default_retraction_method(M::Grassmann, ::Type{StiefelPoint})
default_retraction_method(M::Grassmann, ::Type{ProjectorPoint})

Return ExponentialRetraction as the default on the Grassmann manifold for both representations.

default_vector_transport_method(M::Grassmann)

Return the default vector transport method for the Grassmann manifold, which is ParallelTransport().

get_total_space(::Grassmann)

Return the total space of the Grassmann manifold, which is the corresponding Stiefel manifold, independent of whether the points are represented already in the total space or as ProjectorPoints.

injectivity_radius(M::Grassmann)
injectivity_radius(M::Grassmann, p)

Return the injectivity radius on the Grassmann M, which is $\frac{π}{2}$.

inner(M::Grassmann, A::GrassmannAtlas, i::AbstractVector, a, Xc, Yc)

Compute the Riemannian inner product of coordinate vectors Xc and Yc at coordinates a in the chart i of the standard GrassmannAtlas. Writing Z = reshape(a, n-k, k), U = reshape(Xc, n-k, k), V = reshape(Yc, n-k, k), and G = I + ZᵀZ, the formula is

\[g_Z(U,V) = \operatorname{tr}\left(G^{-1} U^\mathsf{T} \left(I - ZG^{-1}Z^\mathsf{T}\right)V\right).\]

is_flat(M::Grassmann)

Return true if Grassmann M is one-dimensional.

manifold_dimension(M::Grassmann)

Return the dimension of the Grassmann(n,k,𝔽) manifold M, i.e.

\[\dim \operatorname{Gr}(n,k) = k(n-k) \dim_ℝ 𝔽,\]

where $\dim_ℝ 𝔽$ is the real_dimension of 𝔽.

mean(
    M::Grassmann,
    x::AbstractVector,
    [w::AbstractWeights,]
    method = GeodesicInterpolationWithinRadius(π/4);
    kwargs...,
)

Compute the Riemannian mean of x using GeodesicInterpolationWithinRadius.

StiefelPoint <: AbstractManifoldPoint

A point on a Stiefel manifold. This point is mainly used for representing points on the Grassmann where this is also the default representation and hence equivalent to using AbstractMatrices thereon. they can also used be used as points on Stiefel.

StiefelTangentVector <: AbstractTangentVector

A tangent vector on the Grassmann manifold represented by a tangent vector from the tangent space of a corresponding point from the Stiefel manifold, see StiefelPoint. This is the default representation so is can be used interchangeably with just abstract matrices.

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 $⋅^{\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.\]

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⋅ \operatorname{atan}(S) ⋅ U^\mathrm{H},\]

where $⋅^{\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.

rand(M::Grassmann; σ::Real=1.0, vector_at=nothing)

When vector_at is nothing, return a random point p on Grassmann manifold M by generating a random (Gaussian) matrix with standard deviation σ in matching size, which is orthonormal.

When vector_at is not nothing, return a (Gaussian) random vector from the tangent space $T_p\mathrm{Gr}(n,k)$ with mean zero and standard deviation σ by projecting a random Matrix onto the tangent space at vector_at.

riemannian_Hessian(M::Grassmann, p, G, H, X)

The Riemannian Hessian can be computed by adopting Eq. (6.6) [Ngu23], where we use for the EuclideanMetric $α_0=α_1=1$ in their formula. Let $\nabla f(p)$ denote the Euclidean gradient G, $\nabla^2 f(p)[X]$ the Euclidean Hessian H. Then the formula reads

\[ \operatorname{Hess}f(p)[X] = \operatorname{proj}_{T_p\mathcal M}\Bigl( ∇^2f(p)[X] - X p^{\mathrm{H}}∇f(p) \Bigr).\]

Compared to Eq. (5.6) also the metric conversion simplifies to the identity.

get_chart_index(M::Grassmann, A::GrassmannAtlas, p)
get_chart_index(M::Grassmann, A::GrassmannAtlas, i, a)

Return a chart index suitable for a point p or coordinates a in chart i. For a point, the index of the row minor with the largest absolute determinant is returned, breaking ties lexicographically. For coordinates, the represented point is reconstructed and the same selection rule is applied.

get_coordinates_induced_basis!(M::Grassmann, c, p, X, B::InducedBasis{<:Any, <:Any, <:GrassmannAtlas})

Store in c the coordinates of a tangent vector X at p with respect to the basis induced by the standard GrassmannAtlas.

get_parameters!(M::Grassmann, a, A::GrassmannAtlas, i, p)

Store the standard affine coordinates of p in a for the chart indexed by the row tuple i. The selected row minor of p must be invertible.

get_point!(M::Grassmann, p, A::GrassmannAtlas, i::AbstractVector, a)

Store in p the point represented by affine coordinates a in the chart of the standard GrassmannAtlas indexed by i.

get_vector_induced_basis!(M::Grassmann, X, p, c, B::InducedBasis{ℝ, TangentSpaceType, <:GrassmannAtlas})

Store in X the tangent vector at p represented by coordinates c with respect to the basis induced by the standard GrassmannAtlas.

distance(M::Grassmann, p, q)

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

The distance is the arc length

\[d_{\mathrm{Gr}(n,k)}(p,q) = \operatorname{norm}(\log_p(q)) = \Bigl(\sum_{i=1}^k θ_i^2\Bigr)^{1/2}\]

of the principal angles $θ_i$ between the spans of p and q, see Section 4.3 of [EAS98]. The angles are $θ_i = \arccos(σ_i)$ for the singular values $σ_i$ of $p^{\mathrm{H}}q$; the small ones are recomputed with the arc sine, see Algorithm 3.2 of [KA02].

get_coordinates(M::Grassmann{ℝ}, p, X, B::DefaultOrthonormalBasis)

Given a point p on the Grassmann manifold M in Stiefel representation, i.e. $p ∈ ℝ^{n×k}$ compute the coordinates $c ∈ ℝ^{k(n-k)}$ representing the tangent vector X with respect to the DefaultOrthonormalBasis.

The tangent space is characterized by $p^{\mathrm{T}}X = 0_{k×k}$, where $0_{k×k}$ is the $k×k$ matrix containing zeros A default orthonormal basis can is constructed as follows: Since $p$ has full column rank $k$, the null space of $p^{\mathrm{T}}$ is of dimension $n-k$. Let $v_1,…v_{n-k}$ be an orthonormal basis of that null space, then

\[X_1 = \bigl(v_1, 0_k,…0_k), X_2 = \bigl(v_2 0_k,…,0_k), …, X_{n-k} = \bigl(v_{n-k} 0_k,…0_k), X_{n-k1} = \bigl(0_k, v_1, 0_k, …, 0_k), …, X_{(n-k)k} = (0_k,\ldots,0_k,v_{n-k}),\]

where $0_k$ denotes the $k$-dimensional zero vector. Let $V = (v_1,…,v_{n-k}) ∈ ℝ^{n×(n-k)}$ be the matrix of the basis vectors

For a tangent vector $X$ we know that every column X_i,i=1,…,klies in the span ofv_1,…v_{n-k}by looking at every column of the matrix equationp^{\mathrm{T}}X = 0_{k×k}`.

Hence we can compute the coordinates $c$ in $k$ “blocks” $C_1,…C_k ∈ ℝ^{n-k}$ by solving the linear systems

\[VC_i = X_i, \qquad i=1,…,k.\]

get_vector(M::Grassmann{ℝ}, p, c, B::DefaultOrthonormalBasis)

Given a point p on the Grassmann manifold M in Stiefel representation, i.e. $p ∈ ℝ^{n×k}$ reconstruct a tangent vector with respect to the DefaultOrthonormalBasis given coefficients $c ∈ ℝ^{k(n-k)}$.

The tangent space is characterized by $p^{\mathrm{T}}X = 0_{k×k}$, where $0_{k×k}$ is the $k×k$ matrix containing zeros A default orthonormal basis can is constructed as follows: Since $p$ has full column rank $k$, the null space of $p^{\mathrm{T}}$ is of dimension $n-k$. Let $v_1,…v_{n-k}$ be an orthonormal basis of that null space, then

\[X_1 = \bigl(v_1, 0_k,…0_k), X_2 = \bigl(v_2 0_k,…,0_k), …, X_{n-k} = \bigl(v_{n-k} 0_k,…0_k), X_{n-k1} = \bigl(0_k, v_1, 0_k, …, 0_k), …, X_{(n-k)k} = (0_k,\ldots,0_k,v_{n-k}),\]

where $0_k$ denotes the $k$-dimensional zero vector. Let $V = (v_1,…,v_{n-k}) ∈ ℝ^{n×(n-k)}$ be the matrix of the basis vectors.

Reconstructing a tangent vector $X$from a vector of coefficientsc_1,…,c_{k(n-k)}can be done based on the ONBv_1,…v_{n-k}. Theith columnC_iofXis given by using “blocks”C_1,…C_k ∈ ℝ^{n-k}` of the coordinates and compute

\[X_i = VC_i = \sum_{j=1}^{n-k} c_{i(n-k) + j}v_j,\qquad i=1,…k\]

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 $⋅^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

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 $⋅^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

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 $⋅^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

is_vector(M::Grassmann, p, X; atol, rtol, kwargs...)

Check whether X is a tangent vector at p on the Grassmann M, that is whether $p^{\mathrm{H}}X = 0$ holds up to max(atol, rtol * norm(X)). The relative tolerance rtol refers to the size of X; its default is the one of isapprox.

parallel_transport_direction(M::Grassmann, p, X, Y)

Compute the parallel transport of $X \in T_p\mathcal M$ along the geodesic starting in direction $\dot γ (0) = Y$.

Let $Y = USV$ denote the SVD decomposition of $Y$. Then the parallel transport is given by the formula according to Equation (8.5) (p. 171) [AMS08] as

\[\mathcal P_{p,Y} X = -pV \sin(S)U^{\mathrm{T}}X + U\cos(S)U^{\mathrm{T}}X + (I-UU^{\mathrm{T}})X\]

where the sine and cosine applied to the diagonal matrix $S$ are meant to be elementwise

parallel_transport_to(M::Grassmann, p, X, q)

Compute the parallel transport of $X ∈ T_p\mathcal M$ along the geodesic connecting $p$ to $q$.

This method uses the logarithmic map and the parallel transport in that direction.

project(M::Grassmann, p)

Project p from the embedding onto the Grassmann M, i.e. compute q as the polar decomposition of $p$ such that $q^{\mathrm{H}}q$ is the identity, where $⋅^{\mathrm{H}}$ denotes the Hermitian, i.e. complex conjugate transposed.

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 $⋅^{\mathrm{H}}$ denotes the complex conjugate transposed or Hermitian.

representation_size(M::Grassmann)

Return the representation size or matrix dimension of a point on the Grassmann M, i.e. $(n,k)$ for both the real-valued and the complex value case.

retract(M::Grassmann, p, X, ::PolarRetraction)

Compute the SVD-based retraction PolarRetraction on the Grassmann M. With $USV^\mathrm{H} = p + X$ the retraction reads

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

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

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 $k×k$ matrix with

\[D = \operatorname{diag}\left( \operatorname{sgn}\left(R_{ii}+\frac{1}{2}\right)_{i=1}^k \right).\]

riemann_tensor(::Grassmann{ℝ}, p, X, Y, Z)

Compute the value of Riemann tensor on the real Grassmann manifold. The formula reads [Ren11]

\[R(X,Y)Z = (XY^\mathrm{T} - YX^\mathrm{T})Z + Z(Y^\mathrm{T}X - X^\mathrm{T}Y).\]

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.

zero_vector(M::Grassmann, p)

Return the zero tangent vector from the tangent space at p on the Grassmann M, which is given by a zero matrix the same size as p.

ProjectorPoint <: AbstractManifoldPoint

A type to represent points on a manifold Grassmann that are orthogonal projectors, i.e. a matrix $p ∈ \mathbb F^{n,n}$ projecting onto a $k$-dimensional subspace.

ProjectorTangentVector <: AbstractTangentVector

A type to represent tangent vectors to points on a Grassmann manifold that are orthogonal projectors.

exp(M::Grassmann, p::ProjectorPoint, X::ProjectorTangentVector)

Compute the exponential map on the Grassmann as

\[ \exp_p X = \operatorname{Exp}([X,p])p\operatorname{Exp}(-[X,p]),\]

where $\operatorname{Exp}$ denotes the matrix exponential and $[A,B] = A B - B A$ denotes the matrix commutator.

For details, see Proposition 3.2 in [BZA20].

norm(M::Grassmann, p::ProjectorPoint, X::ProjectorTangentVector)

Compute the norm of a ProjectorTangentVector X at p on the Grassmann manifold M.

It is the norm of the inner product, see Section 3.1, p. 15, of [BZA24],

\[\lVert X \rVert_p = \frac{1}{\sqrt{2}}\lVert X \rVert_{\mathrm{F}},\]

where $\lVert ⋅ \rVert_{\mathrm{F}}$ denotes the Frobenius norm.

canonical_project!(M::Grassmann, q::ProjectorPoint, p)

Compute the canonical projection $π(p)$ from the Stiefel manifold onto the Grassmann manifold when represented as ProjectorPoint, i.e.

\[ π^{\mathrm{SG}}(p) = pp^{\mathrm{T}}\]

check_point(::Grassmann, p::ProjectorPoint; kwargs...)

Check whether an orthogonal projector is a point from the Grassmann(n,k) manifold, i.e. the ProjectorPoint $p ∈ \mathbb F^{n×n}$, $\mathbb F ∈ \{\mathbb R, \mathbb C\}$ has to fulfill $p^{\mathrm{H}} = p$, $p^2=p$, and $\operatorname{rank} p = k$. Here $⋅^{\mathrm{H}}$ denotes the conjugate transpose, which is the transpose for real matrices. For $\mathbb F = \mathbb C$ these are the conditions in Section 1.2, p. 9, of [HM94].

check_size(M::Grassmann, p::ProjectorPoint; kwargs...)

Check that the ProjectorPoint is of correct size, i.e. from $\mathbb F^{n×n}$

check_vector(::Grassmann, p::ProjectorPoint, X::ProjectorTangentVector; kwargs...)

Check whether the ProjectorTangentVector X is from the tangent space $T_p\operatorname{Gr}(n,k)$ at the ProjectorPoint p on the Grassmann manifold $\operatorname{Gr}(n,k)$. This means that X has to be Hermitian, $X^{\mathrm{H}} = X$, and that

\[Xp + pX = X\]

must hold, where the kwargs can be used to check both $X^{\mathrm{H}} = X$ and this equality up to a certain tolerance.

diff_canonical_project(M::Grassmann, p::StiefelPoint, X::StiefelTangentVector)
diff_canonical_project!(M::Grassmann, Y::ProjectorTangentVector, p, X)

Compute the differential of canonical projection $π(p)$ from the Stiefel manifold onto the Grassmann manifold when represented as ProjectorPoint, i.e.

\[ Dπ^{\mathrm{SG}}(p)[X] = Xp^{\mathrm{T}} + pX^{\mathrm{T}}\]

distance(M::Grassmann, p::ProjectorPoint, q::ProjectorPoint)

Compute the Riemannian distance between two ProjectorPoints p and q on the Grassmann manifold M.

The projectors are written as $p = UU^{\mathrm{H}}$ and $q = VV^{\mathrm{H}}$, where the columns of $U, V ∈ 𝔽^{n×k}$ are the eigenvectors of p and q for their $k$ largest eigenvalues, which are one, see Section 2.1, p. 6, of [BZA24]. The result is the distance distance(M, U, V) of these Stiefel representatives,

\[d_{\mathrm{Gr}(n,k)}(p,q) = \Bigl(\sum_{i=1}^k θ_i^2\Bigr)^{1/2},\]

where $θ_i$ are the principal angles between the spans of $U$ and $V$, see Section 5.1, equation (5.3), p. 29, of [BZA24].

get_embedding(M::Grassmann, ::Type{<:ProjectorPoint})

Return the embedding of the ProjectorPoint representation of the Grassmann manifold, i.e. the Euclidean space $\mathbb F^{n×n}$.

horizontal_lift(N::Stiefel{n,k}, q, X::ProjectorTangentVector)

Compute the horizontal lift of X from the tangent space at $p=π(q)$ on the Grassmann manifold, i.e.

\[Y = Xq ∈ T_q\mathrm{St}(n,k)\]

inner(M::Grassmann, p::ProjectorPoint, X::ProjectorTangentVector, Y::ProjectorTangentVector)

Compute the inner product of two ProjectorTangentVectors X and Y at p on the Grassmann manifold M.

It is half the Frobenius inner product, see Section 3.1, equation (3.2), p. 14, of [BZA24],

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

parallel_transport_direction(
    M::Grassmann,
    p::ProjectorPoint,
    X::ProjectorTangentVector,
    d::ProjectorTangentVector
)

Compute the parallel transport of X from the tangent space at p into direction d, i.e. to $q=\exp_pd$. The formula is given in Proposition 3.5 of [BZA20] as

\[\mathcal{P}_{q ← p}(X) = \operatorname{Exp}([d,p])X\operatorname{Exp}(-[d,p]),\]

where $\operatorname{Exp}$ denotes the matrix exponential and $[A,B] = A B - B A$ denotes the matrix commutator.

representation_size(M::Grassmann, p::ProjectorPoint)

Return the representation size or matrix dimension of a point on the Grassmann M when using ProjectorPoints, i.e. $(n,n)$.