michakraus · Oct 31, 2026 · Oct 30, 2026 · Oct 30, 2026 · Oct 31, 2026 · Oct 31, 2026
diff --git a/CHANGELOG.md b/CHANGELOG.md
diff --git a/KNOWN_ISSUES.md b/KNOWN_ISSUES.md
@@ -76,24 +76,24 @@ What is known to be wrong in GeometricMachineLearning and is not fixed.
   asserts the fields with `hasfield`. The fields stay because code outside this package may need
   them; removing them changes the type's parameters, its layout and those tests.
 
-### K41 · `CrossModalCycleLoss(nn, x_a, x_b)` with a mixed-component `NamedTuple` raises a `MethodError`, not the documented `ArgumentError`.
+### K41 · `CrossModalLoss(nn, x_a, x_b)` with a mixed-component `NamedTuple` raises a `MethodError`, not the documented `ArgumentError`.
 
-- **location:** `src/loss/multimodal_losses.jl:323`
+- **location:** `src/loss/multimodal_losses.jl:382`
 - **kind:** defect
 - **found:** #370
 - **evidence:**
-  The three `(nn, x_a, x_b)` methods (lines 323, 329 and 335) annotate `x_a` and `x_b` with
+  The three `(nn, x_a, x_b)` methods (lines 382, 388 and 394) annotate `x_a` and `x_b` with
   `ArrayOrNamedTuple`, and AbstractNeuralNetworks defines that as
   `Union{AbstractArray{T}, NamedTupleOfArrays{T}}`, with one element type `T` for every component.
   A `NamedTuple` whose `q` and `p` differ in element type matches no `T`, so dispatch fails before
-  `_check_eltype` (line 204) runs. In a fresh `julia --startup-file=no --project=.` process:
+  `_check_eltype` (line 251) runs. In a fresh `julia --startup-file=no --project=.` process:
 
   ```julia
   T = typeof((q = zeros(Float32, 2, 3), p = zeros(Float64, 2, 3)))
   T <: AbstractNeuralNetworks.ArrayOrNamedTuple   # false
   ```
 
-  The pair form `CrossModalCycleLoss(nn, (x_a, x_b))` takes `Tuple{Any, Any}`, reaches
+  The pair form `CrossModalLoss(nn, (x_a, x_b))` takes `Tuple{Any, Any}`, reaches
   `_check_eltype`, and is tested. The docstring says that another element type raises an
   `ArgumentError`. Closing this entry means widening the `x_a` and `x_b` annotations of the three
   methods, for example to `Union{AbstractArray, NamedTuple}`. That is a design call about what the

diff --git a/docs/src/layers/stiefel_layer.md b/docs/src/layers/stiefel_layer.md
@@ -1,9 +1,14 @@
 # Bilinear and Low-Rank Stiefel Layers
 
 `BilinearStiefelLayer`, `TransposeBilinearStiefelLayer`, `LowRankStiefelLayer`,
-`DecodeLowRankStiefelLayer`, and `FastDLRATraceLoss` provide canonical velocity-parity
-Stiefel frontends for two-sided phase-space grid fields ``F \in \mathbb{R}^{N_x \times N_v}``
-and dynamical low-rank approximation (DLRA) states ``(X, S, V)``.
+`DecodeLowRankStiefelLayer`, and `FastDLRATraceLoss` provide velocity-parity Stiefel
+frontends for two-sided phase-space grid fields ``F \in \mathbb{R}^{N_x \times N_v}`` and
+dynamical low-rank approximation (DLRA) states ``(X, S, V)``.
+
+Their ``(q, p)`` is canonical only for the structure ``\mathbb{J}`` that they declare on the
+parity coordinates. A frontend that is linear in ``f`` is not a Poisson map for the Vlasov
+bracket: with ``L^2``-orthonormal even and odd bases, the bracket linearised at ``f_0`` gives a
+``\{Q, P\}`` that is not the identity and can be singular.
 
 ## Bilinear Grid and Spline Stiefel Layers
 

diff --git a/docs/src/layers/symmetric_particle_layers.md b/docs/src/layers/symmetric_particle_layers.md
@@ -11,7 +11,7 @@ P_\sigma^T \mathbb{J}_{2Md} P_\sigma = \mathbb{J}_{2Md}, \qquad
 P_\sigma = \begin{bmatrix} \Pi_\sigma \otimes I_d & \mathbb{O} \\ \mathbb{O} & \Pi_\sigma \otimes I_d \end{bmatrix}.
 ```
 
-A standard [`PSDLayer`](@ref), [`GradientLayerQ`](@ref), or [`GradientLayerP`](@ref) breaks ``S_M`` symmetry unless the particle state is canonicalized or processed by ``S_M``-equivariant symplectic layers. `GeometricMachineLearning` provides three particle frontends and the composite architecture [`SymmetricSymplecticAutoencoder`](@ref). The sorting and slot frontends are ``S_M``-invariant; the pairwise shears are ``S_M``-equivariant.
+A standard [`PSDLayer`](@ref), [`GradientLayerQ`](@ref), or [`GradientLayerP`](@ref) breaks ``S_M`` symmetry unless the particle state is canonicalized or processed by ``S_M``-equivariant symplectic layers. `GeometricMachineLearning` provides three particle frontends and the composite architecture [`SymmetricSymplecticAutoencoder`](@ref). The sorting and slot frontends are invariant under a relabelling of particles with distinct coordinates along the sorting axis; the pairwise shears are ``S_M``-equivariant. Sorting is discontinuous in ``p`` where two particles cross along the sorting axis, so the Jacobian identity of a sorting frontend holds only away from those crossings. The sorting and slot frontends suit systems that keep their particle order, not a collisionless plasma, whose sorted snapshots are not low-rank.
 
 ## Canonical Weyl-Chamber Sorting
 
@@ -36,7 +36,7 @@ InversePairwiseSymplecticLayerP
 
 ## Canonical Slot Reduction
 
-[`CanonicalSlotLayer`](@ref) sorts the particles along coordinate row 1 of ``q`` with [`CanonicalSortLayer`](@ref), separates the center-of-mass mode ``(Q_{\mathrm{cm}}, P_{\mathrm{cm}})`` from the ``M - 1`` zero-sum relative modes via the orthogonal Helmert / Jacobi transformation, sorts the ``M - 1`` relative modes canonically by descending squared phase-space norm ``\|q_{\mathrm{rel}, k}\|_2^2 + \|p_{\mathrm{rel}, k}\|_2^2``, and projects onto ``N_{\mathrm{slots}} - 1`` relative slots via a learned Stiefel matrix ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)``. Its analytic decoder lift is [`InverseCanonicalSlotLayer`](@ref).
+[`CanonicalSlotLayer`](@ref) sorts the particles along coordinate row 1 of ``q`` with [`CanonicalSortLayer`](@ref), separates the center-of-mass mode ``(Q_{\mathrm{cm}}, P_{\mathrm{cm}})`` from the ``M - 1`` zero-sum relative modes via the orthogonal Helmert / Jacobi transformation, sorts the ``M - 1`` relative modes by descending squared phase-space norm ``\|q_{\mathrm{rel}, k}\|_2^2 + \|p_{\mathrm{rel}, k}\|_2^2``, and projects onto ``N_{\mathrm{slots}} - 1`` relative slots via a learned Stiefel matrix ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)``. Its analytic decoder lift is [`InverseCanonicalSlotLayer`](@ref). The particles are already sorted, so the second sort removes no symmetry, and it makes the output discontinuous where two relative-mode norms cross. The layer is a linear projection of sorted states for a fixed particle count ``M``; it takes no variable particle count and has no Löwdin-orthogonalised affinity slots.
 
 ```@docs
 CanonicalSlotLayer
@@ -45,7 +45,7 @@ InverseCanonicalSlotLayer
 
 ## Symmetric Symplectic Autoencoder
 
-[`SymmetricSymplecticAutoencoder`](@ref) combines a chosen particle frontend (`:sort`, `:pairwise`, or `:slot`) with [`WeightTiedSymplecticAutoencoder`](@ref). Only `:sort` and `:slot` give a permutation-invariant latent state; `:pairwise` is equivariant only in its frontend, because the core is not equivariant.
+[`SymmetricSymplecticAutoencoder`](@ref) combines a chosen particle frontend (`:sort`, `:pairwise`, or `:slot`) with [`WeightTiedSymplecticAutoencoder`](@ref). Only `:sort` and `:slot` give a latent state that is invariant under a relabelling of particles with distinct coordinates along the sorting axis; `:pairwise` is equivariant only in its frontend, because the core is not equivariant.
 
 ```@docs
 SymmetricSymplecticAutoencoder

diff --git a/docs/src/reduced_order_modeling/losses.md b/docs/src/reduced_order_modeling/losses.md
@@ -48,7 +48,7 @@ AutoEncoderLoss
 ReducedLoss
 CharbonnierRelativeLoss
 IsometricRolloutLoss
-CrossModalCycleLoss
+CrossModalLoss
 GeometricMachineLearning.ParametricLoss
 GeometricMachineLearning.apply_parametric
 projection_error

diff --git a/src/GeometricMachineLearning.jl b/src/GeometricMachineLearning.jl
@@ -276,7 +276,7 @@ export LowRankStiefelLayer, DecodeLowRankStiefelLayer
 export FastDLRATraceLoss
 export SymmetricSymplecticAutoencoder
 export MultiModalSymplecticAutoencoder, encode, decode, translate
-export CharbonnierRelativeLoss, IsometricRolloutLoss, CrossModalCycleLoss
+export CharbonnierRelativeLoss, IsometricRolloutLoss, CrossModalLoss
 
 include("map_to_cpu.jl")
 

diff --git a/src/architectures/multimodal_symplectic_autoencoder.jl b/src/architectures/multimodal_symplectic_autoencoder.jl
@@ -29,7 +29,8 @@ ChainRulesCore.@non_differentiable _check_modality(::Any, ::Any)
 
 Make an instance of a multi-modal symplectic autoencoder mapping a non-empty `NamedTuple` of modality frontends into a shared canonical latent space ``\xi = (\tilde{q}, \tilde{p}) \in \mathbb{R}^{2n}`` via `core`.
 
-Each entry of `frontends` maps a modality symbol (such as `:particles`, `:spline`, `:grid`, or `:dlra`) to a symplectic or Poisson frontend layer whose canonical output dimension equals `core.full_dim = 2N`:
+Each entry of `frontends` maps a modality symbol (such as `:particles`, `:spline`, `:grid`, or `:dlra`) to a frontend layer whose output dimension equals `core.full_dim = 2N`.
+The latent state is canonical for the structure that each frontend declares; for the two parity frontends, [`BilinearStiefelLayer`](@ref) and [`LowRankStiefelLayer`](@ref), that structure is the ``\mathbb{J}_{2N}`` that they declare on the parity coordinates, not the Vlasov bracket. The frontends are:
 - [`CanonicalSortLayer`](@ref), [`CanonicalSlotLayer`](@ref), [`PairwiseSymplecticLayerQ`](@ref), [`PairwiseSymplecticLayerP`](@ref), or a `Chain` of particle layers with output dimension ``2N``;
 - [`BilinearStiefelLayer`](@ref)`(Nx, Nv, Mx, Mv)` with ``2 M_x M_v = 2N`` (with or without `metric_factors` for spline vs grid modalities);
 - [`LowRankStiefelLayer`](@ref)`(Nx, Nv, Mx, Mv; r_out = ...)` with ``2 M_x M_v = 2N``.

diff --git a/src/architectures/symmetric_symplectic_autoencoder.jl b/src/architectures/symmetric_symplectic_autoencoder.jl
@@ -6,6 +6,7 @@ Make an instance of a symplectic autoencoder with a particle frontend for `M` pa
 `SymmetricSymplecticAutoencoder` combines a particle frontend (`frontend = :sort`, `:pairwise`, or `:slot`) with the weight-tied symplectic core [`WeightTiedSymplecticAutoencoder`](@ref) so that [`encoder`](@ref)`(nn)` and [`decoder`](@ref)`(nn)` share a single parameter tree.
 
 With `frontend = :sort` or `frontend = :slot`, the latent state is invariant under a relabelling of particles with distinct coordinates along the sorting axis (`sort_axis` for `:sort`, coordinate 1 for `:slot`). With `frontend = :pairwise`, only the frontend is permutation-equivariant: the core is not, so the latent state is not invariant.
+The sort of `:sort` and `:slot` is discontinuous in `p` where two particles cross along the sorting axis, so the Jacobian identity of the frontend holds only away from those crossings. `:sort` and `:slot` suit systems that keep their particle order, not a collisionless plasma, whose sorted snapshots are not low-rank.
 
 The latent right-inversion ``\mathrm{encoder}(nn)(\mathrm{decoder}(nn)(\xi)) = \xi`` holds to machine precision for every ``\xi \in \mathbb{R}^{2n}`` with `frontend = :pairwise`. With `frontend = :sort` and `frontend = :slot`, it holds for ``\xi`` whose decoded state lies in the canonical chamber, where the sorting permutations on the decoded state are the identity; the decoder does not map into this chamber in general.
 

diff --git a/src/layers/bilinear_stiefel_layer.jl b/src/layers/bilinear_stiefel_layer.jl
@@ -30,6 +30,12 @@ When `metric_factors = (Rx, Rv)` is supplied with Cholesky or square-root mass f
 `Rx` of size `(Nx, Nx)` and `Rv` of size `(Nv ÷ 2, Nv ÷ 2)` or `(Nv, Nv)`, the field is
 mapped into the Euclidean metric by `(Rx, Rv)` prior to Stiefel projection.
 
+The coordinates ``(q, p)`` are canonical only for the structure ``\mathbb{J}_{2n}`` that the
+layer and [`LowRankStiefelLayer`](@ref) declare on the parity coordinates. A frontend that is
+linear in ``F`` is not a Poisson map for the Vlasov bracket: with ``L^2``-orthonormal even and
+odd bases, the bracket linearised at ``f_0`` gives a ``\{Q, P\}`` that is not the identity
+and can be singular.
+
 See also [`TransposeBilinearStiefelLayer`](@ref), [`LowRankStiefelLayer`](@ref), and
 [`inverse_layer`](@ref).
 """

diff --git a/src/layers/canonical_slot_layer.jl b/src/layers/canonical_slot_layer.jl
@@ -41,15 +41,20 @@ end
     CanonicalSlotLayer{N, NOut, M, NSlots}()
     CanonicalSlotLayer(M::Integer, d::Integer, NSlots::Integer)
 
-Make an instance of a canonical slot reduction layer mapping `M` particles in ``d`` spatial dimensions to `NSlots` collective canonical slots (`M >= 2`, `1 <= NSlots <= M`).
+Make an instance of a canonical slot reduction layer mapping `M` particles in ``d`` spatial dimensions to `NSlots` canonical slots (`M >= 2`, `1 <= NSlots <= M`).
 
 The first two type parameters are the phase-space dimensions `N = 2Md` of the particles and `NOut = 2 NSlots d` of the slots, which `input_dimension` and `output_dimension` return.
 The layer first sorts the particles along coordinate row 1 of `q` with [`CanonicalSortLayer`](@ref), so that its output is invariant under a relabelling of particles with distinct first coordinates. It then applies the orthogonal Helmert / Jacobi transformation separating the center-of-mass mode
 
 ```math
 (Q_{\mathrm{cm}}, P_{\mathrm{cm}}) = \frac{1}{\sqrt{M}}\sum_{i=1}^M (q_i, p_i)
 ```
-from the ``M - 1`` zero-sum relative modes ``(q_{\mathrm{rel}}, p_{\mathrm{rel}}) \in \mathbb{R}^{d \times (M - 1)}``, sorts the ``M - 1`` relative modes canonically by descending squared phase-space norm ``\|q_{\mathrm{rel}, k}\|_2^2 + \|p_{\mathrm{rel}, k}\|_2^2`` (`sortperm` with `rev = true, alg = MergeSort`, applying the identical permutation ``\sigma_{\mathrm{rel}} \in S_{M - 1}`` to both ``q_{\mathrm{rel}}`` and ``p_{\mathrm{rel}}``), and projects the sorted relative modes onto ``N_{\mathrm{slots}} - 1`` relative slots via a learned [`StiefelManifold`](@extref GeometricOptimizers GeometricOptimizers.StiefelManifold) parameter ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)`` while keeping ``(Q_{\mathrm{cm}}, P_{\mathrm{cm}})`` intact as slot 1.
+from the ``M - 1`` zero-sum relative modes ``(q_{\mathrm{rel}}, p_{\mathrm{rel}}) \in \mathbb{R}^{d \times (M - 1)}``, sorts the ``M - 1`` relative modes by descending squared phase-space norm ``\|q_{\mathrm{rel}, k}\|_2^2 + \|p_{\mathrm{rel}, k}\|_2^2`` (`sortperm` with `rev = true, alg = MergeSort`, applying the identical permutation ``\sigma_{\mathrm{rel}} \in S_{M - 1}`` to both ``q_{\mathrm{rel}}`` and ``p_{\mathrm{rel}}``), and projects the sorted relative modes onto ``N_{\mathrm{slots}} - 1`` relative slots via a learned [`StiefelManifold`](@extref GeometricOptimizers GeometricOptimizers.StiefelManifold) parameter ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)`` while keeping ``(Q_{\mathrm{cm}}, P_{\mathrm{cm}})`` intact as slot 1.
+The particles are already sorted, so this second sort removes no symmetry, and it makes the output discontinuous where two relative-mode norms cross.
+
+The layer is a linear projection of sorted states for a fixed particle count `M`; it takes no variable particle count and has no Löwdin-orthogonalised affinity slots.
+The sort of the particles is discontinuous in `p` where two particles cross along coordinate row 1, so the identities of the layer's Jacobian hold only away from those crossings.
+The layer suits systems that keep their particle order, not a collisionless plasma, whose sorted snapshots are not low-rank.
 """
 struct CanonicalSlotLayer{N, NOut, M, NSlots} <: AbstractExplicitLayer{N, NOut}
     function CanonicalSlotLayer{N, NOut, M, NSlots}() where {N, NOut, M, NSlots}
@@ -181,7 +186,7 @@ end
     InverseCanonicalSlotLayer{N, NOut, M, NSlots}()
     InverseCanonicalSlotLayer(M::Integer, d::Integer, NSlots::Integer)
 
-Analytic symplectic decoder lift (`inverse_layer(CanonicalSlotLayer(M, d, NSlots))`) mapping `NSlots` collective canonical slots back to `M` particles via ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)`` and the inverse Helmert matrix.
+Analytic symplectic decoder lift (`inverse_layer(CanonicalSlotLayer(M, d, NSlots))`) mapping `NSlots` canonical slots back to `M` particles via ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)`` and the inverse Helmert matrix.
 
 It takes the type parameters of [`CanonicalSlotLayer`](@ref), so its input dimension is `NOut = 2 NSlots d` and its output dimension is `N = 2Md`.
 For any ``\Phi_{\mathrm{rel}} \in \mathrm{St}(N_{\mathrm{slots}} - 1, M - 1)``, the decoder lift satisfies ``(\nabla \Psi_{\mathrm{slot}}^{\mathrm{dec}})^T \mathbb{J}_{2Md} \nabla \Psi_{\mathrm{slot}}^{\mathrm{dec}} = \mathbb{J}_{2 N_{\mathrm{slots}} d}`` and preserves ``(Q_{\mathrm{cm}}, P_{\mathrm{cm}})``.

diff --git a/src/layers/canonical_sort_layer.jl b/src/layers/canonical_sort_layer.jl
@@ -34,6 +34,9 @@ Make an instance of a cotangent-lifted canonical Weyl-chamber particle sorting l
 The first type parameter is the phase-space dimension `N = 2Md`, which `input_dimension` and `output_dimension` return.
 Given a `(q, p)` `NamedTuple` with `q, p` of shape `(d, M)` (single snapshot) or `(d, M, B)` (batch of ``B`` snapshots), `CanonicalSortLayer{N, M, Axis}()` sorts the `M` particles along coordinate row `Axis` of `q` using a stable merge sort (`MergeSort`), applies the **identical** permutation ``\sigma_b \in S_M`` to both `q` and `p` in each batch sample ``b``, and returns `(q = q_sorted, p = p_sorted, perm = σ)`.
 On a flattened phase-space vector of length ``2Md`` or snapshot matrix of shape ``(2Md, B)``, it reshapes each half to `(d, M)` or `(d, M, B)` and returns the sorted flattened state.
+
+The sort is discontinuous in `p` where two particles cross along the sorting axis, because the permutation changes there, so the Jacobian identity ``(\nabla\Psi)^T \mathbb{J}_{2Md} \nabla\Psi = \mathbb{J}_{2Md}`` holds only away from those crossings.
+Sorting suits systems that keep their particle order, not a collisionless plasma, whose sorted snapshots are not low-rank.
 """
 struct CanonicalSortLayer{N, M, Axis} <: AbstractExplicitLayer{N, N}
     function CanonicalSortLayer{N, M, Axis}() where {N, M, Axis}

diff --git a/src/layers/low_rank_stiefel_layer.jl b/src/layers/low_rank_stiefel_layer.jl
@@ -24,6 +24,12 @@ which produces the exact same canonical coordinates `(q = vec(Q_red), p = vec(P_
 orthogonal gauge transformations ``(X, S, V) \mapsto (X O_x, O_x^T S O_v, V O_v)`` for all
 ``O_x, O_v \in O(r_{\mathrm{in}})``.
 
+The coordinates ``(q, p)`` are canonical only for the structure ``\mathbb{J}_{2n}`` that the layer and
+[`BilinearStiefelLayer`](@ref) declare on the parity coordinates. A frontend that is linear in
+``f = X S V^T`` is not a Poisson map for the Vlasov bracket: with ``L^2``-orthonormal even and
+odd bases, the bracket linearised at ``f_0`` gives a ``\{Q, P\}`` that is not the identity
+and can be singular.
+
 See also [`DecodeLowRankStiefelLayer`](@ref), [`BilinearStiefelLayer`](@ref),
 [`FastDLRATraceLoss`](@ref), and [`inverse_layer`](@ref).
 """