diff --git a/.github/workflows/Documenter.yml b/.github/workflows/Documenter.yml
index ed5e92c95..bccff7933 100644
--- a/.github/workflows/Documenter.yml
+++ b/.github/workflows/Documenter.yml
@@ -19,10 +19,8 @@ jobs:
       - uses: actions/checkout@v7
       - name: Install packages for generating tikz images
         run: |
-          sudo apt-get install imagemagick
-          sudo apt-get install poppler-utils
-          sudo apt-get install texlive-xetex
-          sudo apt-get install texlive-science
+          sudo apt-get update
+          sudo apt-get install -y imagemagick poppler-utils texlive-xetex texlive-science
       # `@latest` is not a moving alias — it is a literal tag, and upstream stopped moving it
       # in November 2024. Pin the major instead.
       - uses: julia-actions/setup-julia@v3
diff --git a/.github/workflows/Latex.yml b/.github/workflows/Latex.yml
index 1e72c716b..9a1ce6c05 100644
--- a/.github/workflows/Latex.yml
+++ b/.github/workflows/Latex.yml
@@ -17,8 +17,8 @@ jobs:
       - uses: actions/checkout@v7
       - name: install imagemagick, poppler-utils and texlive
         run: |
-            sudo apt-get install imagemagick
-            sudo apt-get install poppler-utils
+            sudo apt-get update
+            sudo apt-get install -y imagemagick poppler-utils
             sudo apt-get install texlive-xetex
             sudo apt-get install texlive-science
             sudo apt-get install latex-cjk-all
diff --git a/CHANGELOG.md b/CHANGELOG.md
index fb2cc7d11..6c9faae0b 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -13,45 +13,12 @@ breaking release).
 > alongside the work. Where a release removed exported names the list is given; where it is a
 > reconstruction of intent, it says so.
 
-## [Unreleased] — 0.6.1
+## [Unreleased] — 0.8.0
 
-### Changed
-
-- **Julia 1.11 is the minimum**, up from the 1.10 LTS, in step with the rest of this family of
-  packages. Nothing here was accommodating 1.10; the two places that name it —
-  `test/training_parameters.jl` and `test/optimizers/optimizer_convergence_tests/psd_optim.jl` — record
-  where a measurement was taken and how a defect was found, and both stay as the history they are.
-
-- **`GeometricOptimizers = "0.6"`.** 0.6.0 is what takes a whole `NetworkParameters` through the
-  optimizer, and it is breaking: four `outer!`/`_mul!` methods, one `update!(::BFGSState, …)` and
-  `add!` on a parameter set are gone, and two `l2norm` methods moved upstream to `GeometricBase`. None
-  of them had a caller here — the `add!` arm was deleted on this side in 0.7 for the same reason — so
-  the bump is a compat entry and no code. The full suite is green against it.
-
-- **The fifteen sites spelling `Union{NamedTuple, NetworkParameters}` inline now name it.** It is
-  `NeuralNetworkParameters.ParameterSet`, added upstream in 0.2.2: the same union, in the package that
-  owns the type. `SymbolicNeuralNetworks` had it as `EquationSet`, `AbstractNeuralNetworks` spelled it
-  out eight times and `GeometricOptimizers` sixteen; all four name one thing now. Definitionally
-  identical, so no dispatch changes anywhere. Compat is `NeuralNetworkParameters = "0.2.2"`.
-
-### Fixed
-
-- **The comment on `_make_optimizer_cache`'s branch ordering had been wrong for a release.** It said
-  `NetworkParameters` is not one of the types `GeometricOptimizers.OptimizerSolution` unions, and
-  therefore that asking the structural question before the capability question was invisible. 0.5.0 made
-  the container a member, so `_use_go_cache` is true at the *root* now and the ordering is load-bearing:
-  without it a whole network would get one cache instead of one per layer, silently, with
-  `_leaf_optim_step!` handed the entire tree. The code was already right — the ordering was put in
-  ahead of the change it anticipated — and only the comment needed correcting. Found in the review of
-  [GeometricOptimizers #68](https://github.com/JuliaGNI/GeometricOptimizers.jl/pull/68).
-
-  It also records what the better end state is, now that `GeometricOptimizers` 0.6.0 takes a whole
-  container as one solution: one cache for the network, which means giving `_GMLGradient` a
-  `NetworkParameters` method and dropping `_tree_optim_step!`. That is a change of behaviour — one
-  `GlobalSection` tree and one `Q` across every layer instead of one per layer — so it wants its own
-  release rather than being folded in here.
-
-## [0.7.0]
+> [!NOTE]
+> Not released. This was written as `[0.7.0]` before 0.6.1 and 0.7.0 were cut, and the number was
+> taken by the release below; the work it describes -- `apply_toNT`, `_eltype`, `map_to_cpu` -- is
+> still to do.
 
 **The traversal of a parameter set now belongs to the package that owns the parameters, and the
 traversal of a `NamedTuple` belongs to `Base`.** 0.6.0 handed the HDF5 walk over to
@@ -215,6 +182,135 @@ so it cannot coexist with `GeometricOptimizers` 0.5.
   docs build — here, one release behind. The old one-liner also called `save` unqualified, and
   `DocInventories` does not export it.
 
+
+## [0.7.0]
+
+**A layer is wrapped at the `GeometricOptimizers` boundary, and a whole set of parameters is a
+`NetworkParameters`.** Breaking: this release tracks `NeuralNetworkParameters` 0.3.0,
+`AbstractNeuralNetworks` 0.8.0, `SymbolicNeuralNetworks` 0.8.0 and `GeometricOptimizers` 0.7.0.
+
+### Changed
+
+- **`GeometricOptimizers = "0.7"`.** 0.6.0 is what takes a whole `NetworkParameters` through the
+  optimizer, and it was breaking: four `outer!`/`_mul!` methods, one `update!(::BFGSState, …)` and
+  `add!` on a parameter set are gone, and two `l2norm` methods moved upstream to `GeometricBase`. None
+  of them had a caller here — the `add!` arm was deleted on this side in 0.7 for the same reason.
+
+  0.7.0 goes further and **does** need code here: a whole set of parameters reaches that package only
+  as a `NetworkParameters` now, and a bare `NamedTuple` is turned away. This package hands it one
+  **layer** at a time, so it cannot convert by wrapping at the root; see the entry below. The full
+  suite is green against it, `test/reduced_system.jl` aside — that one needs `GeometricIntegrators`,
+  which still pins `SimpleSolvers` 0.12 and is the ecosystem-wide blocker.
+
+- **`ParameterSet` is gone from the signatures here.** `NeuralNetworkParameters` 0.3.0 removes it, so
+  `Optimizer` and `optimize_for_one_epoch!` take a `NetworkParameters`, `_get_contents` has a method
+  per shape, and the loss functors drop the annotation entirely — it dispatched nothing there, since
+  the model and the input/output types settle every method, and both a container and the bare
+  `NamedTuple` a reverse pass produces reach them.
+
+  `_tree_optim_step!`'s test for a section tree is now `_is_section_tree`, a named predicate rather than
+  an `isa` against a union. It answers one question — descend, or hand this whole section to the layer
+  — and naming it is what keeps the container arm from looking like breadth for its own sake.
+
+- **A *flat* set of parameters gets one optimizer cache again, and the branch order is what decides it.**
+  `_make_optimizer_cache` asks `_is_layer` — "is this one cache's worth?" — *before* asking whether the
+  argument is a tree to descend into. A flat set answers yes to both, and the tree reading is wrong: it
+  gives every individual weight its own cache, and for a manifold weight that is not merely wasteful,
+  because a bare `Manifold` then reaches `_GMLGradient` where a whole layer should have, which is
+  ambiguous against `GeometricOptimizers`' own `Gradient` functor. `_is_layer` accepts a container as
+  well as a bare `NamedTuple` for this reason.
+
+- **A layer is wrapped at the `GeometricOptimizers` boundary, and the flatness rule is written out
+  here.** `_as_go_solution` wraps a layer in a `NetworkParameters` before handing it to
+  `OptimizerCache`, `OptimizerState` or `update!`, and `_leaf_optim_step!` uses the wrapped form
+  throughout. The wrap **shares the leaf arrays**, so the step still writes through to the network's
+  own weights and nothing is copied back — which is what makes this a boundary change and not a
+  change of behaviour.
+
+  What did have to be restated is the test that separates a *layer* from a *subtree*. Upstream made it
+  by dispatch: `ArrayNamedTuple{T}` bounded a `NamedTuple`'s values by `AbstractArray{T}`, so a
+  `NamedTuple` of branches was not one and got a cache per layer rather than one for the whole of it.
+  That alias is gone — it was an alias for `Base.NamedTuple`, so every method on it was a method on a
+  `Base` type — so `_is_layer` says it here, where it is this package's rule rather than a property of
+  whichever types upstream happens to accept this month.
+
+  **One behaviour follows from the wrap rather than from the old bound.** A layer whose weights do not
+  share an element type is now one cache with a promoted `T`, where the `Vararg` bound admitted no `T`
+  at all and the layer fell through to one `GMLEuclideanState` per weight — i.e. no manifold handling.
+  `_adapt_method_to_T` reads the promotion, so the cache is built for the type the layer actually has.
+
+  `_GMLGradient` gains its `NetworkParameters` method as part of this, and `_gml_rgrad` walks with
+  `mapparameters` instead of `Base.map`: it recurses on the branches, so `_gml_rgrad` is only ever
+  called on leaves, and it rebuilds in the shape of its first argument — a container, which is what
+  `GeometricOptimizers._copyto!(gradient_array(cache), ·)` has a method for.
+
+### Fixed
+
+- **`_GMLGradient` on a bare `Manifold` was ambiguous, and the branch ordering was the only thing
+  preventing it.** `_GMLGradient{T}` is a `GeometricOptimizers.Gradient{T}` and `Manifold{T}` is an
+  `AbstractMatrix{T}`, so `(::_GMLGradient{T})(::AbstractArray{T})` here and upstream's
+  `(::Gradient{T})(::Manifold{T})` are neither of them more specific on the intersection —
+  `MethodError: … is ambiguous`, several frames into a solve. Asking `_is_layer` before descending
+  keeps a bare `Manifold` from arriving, and that ordering is a correctness fix in its own right, but
+  a routing decision must not be the only thing standing between a caller and an ambiguity.
+
+  There is now a method for it, and the answer it gives is this package's: the Euclidean gradient is
+  already computed and carried in `dp`, so the projection is `rgrad(x, dp)` rather than upstream's
+  `rgrad(x, reshape(grad(vec(x)), …))`, which would evaluate an inner gradient this functor does not
+  have. `test/optimizers/gml_gradient_dispatch.jl` pins it, along with the shapes the optimizer
+  actually produces.
+
+- **An empty set was one cache's worth.** `all` over an empty collection is vacuously true, so
+  `_is_layer(NamedTuple())` was `true` and `OptimizerCache` would have been asked for the element type
+  of a set with no leaves to promote. An empty set is *zero* caches' worth and descends.
+
+- `_get_contents` takes a method per shape in its wrapped forms too, rather than a `Tuple` and an
+  `AbstractVector` over a union of the two. Same rule as everywhere else in this release.
+
+## [0.6.1]
+
+### Changed
+
+- **Julia 1.11 is the minimum**, up from the 1.10 LTS, in step with the rest of this family of
+  packages. Nothing here was accommodating 1.10; the two places that name it —
+  `test/training_parameters.jl` and `test/optimizers/optimizer_convergence_tests/psd_optim.jl` — record
+  where a measurement was taken and how a defect was found, and both stay as the history they are.
+
+- **`GeometricOptimizers = "0.6"`.** 0.6.0 is what takes a whole `NetworkParameters` through the
+  optimizer, and it is breaking: four `outer!`/`_mul!` methods, one `update!(::BFGSState, …)` and
+  `add!` on a parameter set are gone, and two `l2norm` methods moved upstream to `GeometricBase`. None
+  of them had a caller here — the `add!` arm was deleted on this side in 0.7 for the same reason — so
+  the bump is a compat entry and no code. The full suite is green against it.
+
+- **The fifteen sites spelling `Union{NamedTuple, NetworkParameters}` inline now name it.** It is
+  `NeuralNetworkParameters.ParameterSet`, added upstream in 0.2.2: the same union, in the package that
+  owns the type. `SymbolicNeuralNetworks` had it as `EquationSet`, `AbstractNeuralNetworks` spelled it
+  out eight times and `GeometricOptimizers` sixteen; all four name one thing now. Definitionally
+  identical, so no dispatch changes anywhere. Compat is `NeuralNetworkParameters = "0.2.2"`.
+- **The fifteen sites spelling `Union{NamedTuple, NetworkParameters}` inline now name it.** It is
+  `NeuralNetworkParameters.ParameterSet`, added upstream in 0.2.2: the same union, in the package that
+  owns the type. `SymbolicNeuralNetworks` had it as `EquationSet`, `AbstractNeuralNetworks` spelled it
+  out eight times and `GeometricOptimizers` sixteen; all four name one thing now. Definitionally
+  identical, so no dispatch changes anywhere. Compat is `NeuralNetworkParameters = "0.2.2"`.
+
+### Fixed
+
+- **The comment on `_make_optimizer_cache`'s branch ordering had been wrong for a release.** It said
+  `NetworkParameters` is not one of the types `GeometricOptimizers.OptimizerSolution` unions, and
+  therefore that asking the structural question before the capability question was invisible. 0.5.0 made
+  the container a member, so `_use_go_cache` is true at the *root* now and the ordering is load-bearing:
+  without it a whole network would get one cache instead of one per layer, silently, with
+  `_leaf_optim_step!` handed the entire tree. The code was already right — the ordering was put in
+  ahead of the change it anticipated — and only the comment needed correcting. Found in the review of
+  [GeometricOptimizers #68](https://github.com/JuliaGNI/GeometricOptimizers.jl/pull/68).
+
+  It also records what the better end state is, now that `GeometricOptimizers` takes a whole container
+  as one solution: one cache for the network. `_GMLGradient` has its `NetworkParameters` method as of
+  this release, but that is one wrapped **layer** and not the tree — the container branch still splits
+  the network first. Dropping `_tree_optim_step!` and handing the whole tree to one cache is a change
+  of behaviour — one `GlobalSection` tree and one `Q` across every layer instead of one per layer — so
+  it still wants its own release rather than being folded in here.
+
 ## [0.6.0] — 2026-08-24
 
 **The parameter container moves out to [NeuralNetworkParameters.jl][nnp].**
diff --git a/Project.toml b/Project.toml
index b65825b43..e2f285467 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,6 +1,6 @@
 name = "GeometricMachineLearning"
 uuid = "194d25b2-d3f5-49f0-af24-c124f4aa80cc"
-version = "0.6.1"
+version = "0.7.0"
 authors = ["Michael Kraus <michael.kraus@ipp.mpg.de>"]
 
 [deps]
@@ -31,7 +31,7 @@ HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
 HDF5Ext = "HDF5"
 
 [compat]
-AbstractNeuralNetworks = "0.7"
+AbstractNeuralNetworks = "0.8"
 ChainRulesCore = "1"
 ChainRulesTestUtils = "1"
 Distances = "0.10"
@@ -39,15 +39,15 @@ ForwardDiff = "1"
 GeometricBase = "0.14"
 GeometricEquations = "0.21"
 GeometricIntegrators = "0.18.2"
-GeometricOptimizers = "0.6"
+GeometricOptimizers = "0.7"
 GeometricSolutions = "0.6"
 HDF5 = "0.16, 0.17"
 KernelAbstractions = "0.9"
 NNlib = "0.8, 0.9"
-NeuralNetworkParameters = "0.2.2"
+NeuralNetworkParameters = "0.3"
 ProgressMeter = "1"
 SafeTestsets = "0.1"
-SymbolicNeuralNetworks = "0.7"
+SymbolicNeuralNetworks = "0.8"
 Symbolics = "7"
 TimerOutputs = "0.5, 1"
 Zygote = "0.7"
diff --git a/codecov.yml b/codecov.yml
new file mode 100644
index 000000000..3555747ff
--- /dev/null
+++ b/codecov.yml
@@ -0,0 +1,18 @@
+# Without this file Codecov's defaults apply: `target: auto` with `threshold: 0%`, which fails a pull
+# request for *any* decrease. That is sharper than it reads, because Codecov counts a *modified* line
+# as part of the diff: changing a signature's annotation, or splitting one into a method per shape,
+# re-attributes lines that were already uncovered to whoever touched them.
+#
+# 1% on both, matching `GeometricOptimizers` and `SymbolicNeuralNetworks`. `informational` is
+# deliberately *not* set: the check should still be able to fail, and on #255 it did — six genuinely
+# untested methods, which were tested rather than waved through.
+coverage:
+  status:
+    project:
+      default:
+        target: auto
+        threshold: 1%
+    patch:
+      default:
+        target: auto
+        threshold: 1%
diff --git a/docs/src/optimizers/optimizer.md b/docs/src/optimizers/optimizer.md
index 708b5b558..ca0b095ee 100644
--- a/docs/src/optimizers/optimizer.md
+++ b/docs/src/optimizers/optimizer.md
@@ -26,4 +26,5 @@ opt = Optimizer(nn; AdamOptimizerWithDecay(n_epochs, Float32)...)
 Optimizer
 optimize_for_one_epoch!
 optimization_step!
+GeometricMachineLearning._as_go_solution
 ```
diff --git a/src/GeometricMachineLearning.jl b/src/GeometricMachineLearning.jl
index 9dd1ea4d9..2248e94fa 100644
--- a/src/GeometricMachineLearning.jl
+++ b/src/GeometricMachineLearning.jl
@@ -11,7 +11,7 @@ using AbstractNeuralNetworks
 # `parameter_eltype`, which promotes over the leaves of a set. A plain `NamedTuple` of layers is
 # walked with `Base.map`, which needs nothing from anybody.
 import NeuralNetworkParameters: NetworkParameters
-using NeuralNetworkParameters: mapstorage, parameter_eltype, ParameterSet
+using NeuralNetworkParameters: mapparameters, mapstorage, parameter_eltype
 using ChainRulesCore
 # `sqeuclidean` is the default distance of every `TrainingMethod` in `src/training_method/`.
 using Distances
diff --git a/src/data_loader/optimize.jl b/src/data_loader/optimize.jl
index 91ece127b..188733ee1 100644
--- a/src/data_loader/optimize.jl
+++ b/src/data_loader/optimize.jl
@@ -48,7 +48,7 @@ number_of_batches(dl, batch)
 ```
 """
 function optimize_for_one_epoch!(   opt::Optimizer, 
-                                    model, ps::ParameterSet, 
+                                    model, ps::NetworkParameters, 
                                     dl::DataLoader{T}, 
                                     batch::Batch, 
                                     loss::NetworkLoss, 
@@ -58,7 +58,7 @@ end
 
 function optimize_for_one_epoch!(   opt::Optimizer, 
                                     model, 
-                                    ps::ParameterSet, 
+                                    ps::NetworkParameters, 
                                     dl::DataLoader{T}, 
                                     batch::Batch, 
                                     _pullback::AbstractPullback, 
diff --git a/src/loss/hnn_loss.jl b/src/loss/hnn_loss.jl
index 3b575fade..9ddbada94 100644
--- a/src/loss/hnn_loss.jl
+++ b/src/loss/hnn_loss.jl
@@ -29,13 +29,13 @@ end
 AbstractNeuralNetworks.NetworkLoss(arch::HamiltonianArchitecture) = HNNLoss(arch)
 
 function (loss::HNNLoss)(::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet,
+        ps,
         input::QPTOAT,
         output::QPTOAT)
     loss(ps, input, output)
 end
 
-function (loss::HNNLoss)(ps::ParameterSet,
+function (loss::HNNLoss)(ps,
         input::QPTOAT,
         output::QPTOAT)
     norm(loss.hvf(input, ps) - output) / norm(output)
diff --git a/src/loss/losses.jl b/src/loss/losses.jl
index e8769609c..d19b52ce7 100644
--- a/src/loss/losses.jl
+++ b/src/loss/losses.jl
@@ -73,7 +73,7 @@ function crop_array_for_transformer_loss(nn_output::AT,
 end
 
 function (loss::TransformerLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet, input::AT,
+        ps, input::AT,
         output::AT) where {T, AT <: AbstractArray{T, 3}}
     input_dim, input_seq_length = size(input)
     output_dim, output_prediction_window = size(output)
@@ -88,13 +88,13 @@ function (loss::TransformerLoss)(model::Union{Chain, AbstractExplicitLayer},
 end
 
 function (loss::TransformerLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet, input::AT,
+        ps, input::AT,
         output::AT) where {T, AT <: AbstractArray{T, 2}}
     loss(model, ps, reshape(input, size(input)..., 1), reshape(output, size(output)..., 1))
 end
 
 function (loss::TransformerLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet,
+        ps,
         input::T, output::T) where {T <: QPT}
     loss(model, ps, vcat(input.q, input.p), vcat(output.q, output.p))
 end
@@ -120,7 +120,7 @@ end
 struct ClassificationTransformerLoss <: NetworkLoss end
 
 function (loss::ClassificationTransformerLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet,
+        ps,
         input::AbstractArray, output::AbstractArray)
     predicted_output_uncropped = model(input, ps)
     # predicted_output_cropped = crop_array_for_transformer_loss(predicted_output_uncropped, output)
@@ -179,12 +179,12 @@ function (loss::AutoEncoderLoss)(nn::NeuralNetwork, input::QPTOAT)
 end
 
 function (loss::AutoEncoderLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet, input::QPTOAT)
+        ps, input::QPTOAT)
     loss(model, ps, input, input)
 end
 
 function (loss::AutoEncoderLoss)(model::Union{Chain, AbstractExplicitLayer},
-        ps::ParameterSet, input::QPTOAT, output::QPTOAT)
+        ps, input::QPTOAT, output::QPTOAT)
     FeedForwardLoss()(model, ps, input, output)
 end
 
diff --git a/src/optimizers/optimizer.jl b/src/optimizers/optimizer.jl
index f4a651d7c..1c73d7557 100644
--- a/src/optimizers/optimizer.jl
+++ b/src/optimizers/optimizer.jl
@@ -11,12 +11,30 @@ end
 
 _gml_rgrad(x::Manifold, dp) = rgrad(x, dp)
 _gml_rgrad(x, dp) = dp
-_gml_rgrad(x::NamedTuple, dp::NamedTuple) = map(_gml_rgrad, x, dp)
+# `mapparameters` and not `map`: it recurses on the branches, so `_gml_rgrad` is only ever called on
+# leaves, and it rebuilds in the shape of its *first* argument -- a container, which is what
+# `GeometricOptimizers._copyto!(gradient_array(cache), ·)` has a method for. It also normalises its
+# trailing arguments, so `dp` may stay the plain `NamedTuple` the pullback produced.
+_gml_rgrad(x::NetworkParameters, dp) = mapparameters(_gml_rgrad, x, dp)
 
-(g::_GMLGradient{T})(x::GeometricOptimizers.ArrayNamedTuple{T}) where {T} =
-    _gml_rgrad(x, g.dp)
+(g::_GMLGradient{T})(x::NetworkParameters{T}) where {T} = _gml_rgrad(x, g.dp)
 (g::_GMLGradient{T})(x::AbstractArray{T}) where {T} = g.dp
 
+# `Manifold` and not `AbstractArray` alone, because the two are otherwise **ambiguous** and not merely
+# overlapping. `_GMLGradient{T}` is a `GeometricOptimizers.Gradient{T}`, `Manifold{T}` is an
+# `AbstractMatrix{T}`, and `GeometricOptimizers` has `(::Gradient{T})(::Manifold{T})`; so on a bare
+# manifold neither `(::_GMLGradient{T})(::AbstractArray{T})` above nor that one is more specific, and
+# the call is a run-time `MethodError: ... is ambiguous`. This method is more specific than both.
+#
+# The answer it gives is this package's: the gradient is already computed, so the projection is
+# `rgrad(x, dp)` and not upstream's `rgrad(x, reshape(grad(vec(x)), ...))`, which would evaluate an
+# inner gradient this functor does not have. `_gml_rgrad` is where that lives.
+#
+# `_make_optimizer_cache` below asks `_is_layer` first and so keeps a bare `Manifold` from arriving
+# here at all. That ordering is a correctness fix in its own right -- one cache per layer rather than
+# per weight -- but it must not be the *only* thing standing between a caller and an ambiguity.
+(g::_GMLGradient{T})(x::Manifold{T}) where {T} = _gml_rgrad(x, g.dp)
+
 # State for Euclidean (non-manifold) parameters.
 mutable struct GMLEuclideanState{T, AT<:AbstractArray{T}}
     iterations::Int
@@ -44,48 +62,92 @@ _adapt_method_to_T(method::GeometricOptimizers.MomentumMethod, ::Type{T}) where
     GeometricOptimizers.MomentumMethod(T(method.α))
 _adapt_method_to_T(method, ::Type) = method
 
-_use_go_cache(method, x) =
-    _is_go_native_method(method) && x isa GeometricOptimizers.OptimizerSolution
+# Whether one `GeometricOptimizers` cache covers `x` whole.
+#
+# A **layer** is a flat `NamedTuple` of arrays, and it is exactly what one cache is for. A `NamedTuple`
+# whose values are branches is a *subtree* and gets one cache per layer instead, which is why the test
+# is flatness and not merely "a `NamedTuple`".
+#
+# The rule lives here rather than upstream because a whole set of parameters reaches
+# `GeometricOptimizers` only as a `NetworkParameters`, and a container says nothing about how deep the
+# tree beneath it is. Splitting a network into layers is this package's decision, so this package makes
+# it; [`_as_go_solution`](@ref) then does the wrap.
+#
+# A layer whose weights do not share an element type is still one cache, with `T` the promotion over
+# them. `_adapt_method_to_T` reads that promotion, so the cache is built for the element type the
+# layer actually has rather than for one its weights are required to agree on.
+# `!isempty` and not `all` alone: `all` over an empty collection is vacuously true, so without it an
+# empty set would be one cache's worth -- and `OptimizerCache` would then be asked for the element
+# type of a set that has no leaves to promote. An empty set is zero caches' worth, so it descends and
+# yields an empty tree of caches, which is what the step below iterates over zero times.
+_is_layer(x::NamedTuple) = _all_leaves(values(x))
+_is_layer(x::NetworkParameters) = _all_leaves(values(x))
+_is_layer(_) = false
+
+_all_leaves(vs) = !isempty(vs) && all(v -> v isa AbstractArray, vs)
+
+_use_go_cache(method, x) = _is_go_native_method(method) &&
+    (x isa GeometricOptimizers.OptimizerSolution || _is_layer(x))
+
+"""
+    _as_go_solution(x)
+
+`x` in the shape `GeometricOptimizers` takes a solution in.
+
+A layer given as a bare `NamedTuple` is wrapped; anything already in the right shape — a container, a
+bare `Manifold`, an `AbstractArray` — is passed through. The wrap **shares the leaf arrays**, so an
+in-place optimizer step writes through to the network's own weights and nothing has to be copied back.
+"""
+_as_go_solution(x::NetworkParameters) = x
+_as_go_solution(x) = _is_layer(x) ? NetworkParameters(x) : x
 
 # A `NetworkParameters` is always a tree of layers to descend into, never a single
 # `GeometricOptimizers` leaf, so its branch comes *first* — ahead of `_use_go_cache`.
 #
-# That ordering is **load-bearing**, and was not always. It used to be invisible, because
-# `NetworkParameters` was not one of the types `GeometricOptimizers.OptimizerSolution` unions, so
-# `_use_go_cache` was false at the root anyway and control reached the container branch either way.
-# `GeometricOptimizers` 0.5.0 added the container to that union, so `_use_go_cache` is now true at the
-# *root*: without the test above, a whole network would get one cache instead of one per layer --
-# silently, and with `_leaf_optim_step!` handed the entire tree. Asking the structural question before
-# the capability question is what makes the shape of the cache depend on the shape of the parameters
-# rather than on which types upstream happens to accept this month.
+# That ordering is **load-bearing**. `NetworkParameters` *is* one of the types
+# `GeometricOptimizers.OptimizerSolution` unions, so `_use_go_cache` is true at the root: without the
+# test above, a whole network would get one cache instead of one per layer -- silently, and with
+# `_leaf_optim_step!` handed the entire tree. Asking the structural question before the capability
+# question is what makes the shape of the cache depend on the shape of the parameters rather than on
+# which types upstream happens to accept.
+#
+# One cache for the whole network is the better end state, and it is **not** what this does.
+# `_GMLGradient` does take a `NetworkParameters`, but that is one wrapped *layer* and not the tree,
+# because the container branch above splits the network first. Getting to one cache means dropping
+# `_tree_optim_step!` and handing the whole tree over, which is a change of behaviour -- one
+# `GlobalSection` tree and one `Q` across every layer -- so it wants its own release. Nothing here
+# depends on it.
 #
-# One cache for the whole network is the better end state and is what
-# `GeometricOptimizers` 0.6.0 makes possible; getting there means giving `_GMLGradient` a
-# `NetworkParameters` method and dropping `_tree_optim_step!`, which is a change of behaviour -- one
-# `GlobalSection` tree and one `Q` across every layer -- and so its own release.
+# **"Is this one cache's worth?" is asked before "is this a tree?"**, and a *flat* set answers yes to
+# both. `_is_layer` is that question: a set whose values are all leaves is a layer, whether it arrives
+# wrapped or as a branch, and one `GeometricOptimizers` cache is exactly what a layer is for. Asking
+# the tree question first would descend into such a set and give every individual weight its own cache
+# -- which for a manifold weight is not merely wasteful but wrong, since a bare `Manifold` then reaches
+# the gradient functor where a whole layer should have.
 #
-# The `NamedTuple` branch stays *after* `_use_go_cache`, and that asymmetry is the point: a layer is a
-# `NamedTuple` of arrays, which is exactly what one `GeometricOptimizers` cache is for. Hoisting it
-# too would descend into the individual weights and give each its own `GMLEuclideanState`.
+# The tree branch comes second and covers both carriers, because a network is a tree of layers whether
+# it is wrapped or not.
 function _make_optimizer_cache(method, x)
-    if x isa NetworkParameters
+    if _is_go_native_method(method) && _is_layer(x)
+        GeometricOptimizers.OptimizerCache(_adapt_method_to_T(method, parameter_eltype(x)),
+                                          _as_go_solution(x))
+    elseif x isa NetworkParameters || x isa NamedTuple
         NamedTuple{keys(x)}(Tuple(_make_optimizer_cache(method, x[k]) for k in keys(x)))
     elseif _use_go_cache(method, x)
         GeometricOptimizers.OptimizerCache(_adapt_method_to_T(method, parameter_eltype(x)), x)
-    elseif x isa NamedTuple
-        NamedTuple{keys(x)}(Tuple(_make_optimizer_cache(method, x[k]) for k in keys(x)))
     else
         GMLEuclideanState(x)
     end
 end
 
 function _make_optimizer_state(method, x)
-    if x isa NetworkParameters
+    if _is_go_native_method(method) && _is_layer(x)
+        GeometricOptimizers.OptimizerState(_adapt_method_to_T(method, parameter_eltype(x)),
+                                          _as_go_solution(x))
+    elseif x isa NetworkParameters || x isa NamedTuple
         NamedTuple{keys(x)}(Tuple(_make_optimizer_state(method, x[k]) for k in keys(x)))
     elseif _use_go_cache(method, x)
         GeometricOptimizers.OptimizerState(_adapt_method_to_T(method, parameter_eltype(x)), x)
-    elseif x isa NamedTuple
-        NamedTuple{keys(x)}(Tuple(_make_optimizer_state(method, x[k]) for k in keys(x)))
     else
         GMLEuclideanState(x)
     end
@@ -148,7 +210,7 @@ function Optimizer(method::GeometricOptimizers.OptimizerMethod, nn::NeuralNetwor
 end
 
 function Optimizer(method::GeometricOptimizers.OptimizerMethod,
-        ps::ParameterSet;
+        ps::NetworkParameters;
         retraction = GeometricOptimizers.cayley,
         step_size = _default_step_size(method))
     Optimizer(method, _make_optimizer_cache(method, ps), _make_optimizer_state(method, ps),
@@ -158,7 +220,7 @@ end
 # The keyword form, so that the `(algorithm, linesearch)` pairing `GeometricOptimizers` returns from
 # `AdamOptimizerWithDecay` splats in unchanged. `linesearch` is the step size under the name
 # upstream gives it; a `Static` carries its own `α`, which is then the fixed learning rate.
-function Optimizer(nn_or_ps::Union{NeuralNetwork, NamedTuple, NetworkParameters};
+function Optimizer(nn_or_ps::Union{NeuralNetwork, NetworkParameters};
         algorithm::GeometricOptimizers.OptimizerMethod,
         linesearch = nothing,
         retraction = GeometricOptimizers.cayley,
@@ -221,15 +283,20 @@ function _go_update_leaf!(cache, state, local_grad,
                                 GeometricOptimizers.NoHessian{T}(), ps_leaf)
 end
 
-# GO-managed leaf step (manifolds, vectors, ArrayNamedTuples)
+# GO-managed leaf step (manifolds, vectors, whole layers)
+#
+# `ps` and not `ps_leaf` from here on: a layer is wrapped for the same reason its cache is, and the
+# wrap shares the arrays, so `_copyto!(ps, solution(cache))` below writes into the weights the caller
+# holds. A bare `Manifold` or `AbstractArray` passes through untouched.
 function _leaf_optim_step!(cache::GeometricOptimizers.OptimizerCache,
         state::GeometricOptimizers.OptimizerState,
         dp_leaf, ps_leaf, λY_leaf, method, retraction, step_size)
     T = parameter_eltype(ps_leaf)
+    ps = _as_go_solution(ps_leaf)
     local_grad = _GMLGradient{T, typeof(dp_leaf)}(dp_leaf)
     adapted = _adapt_method_to_T(method, T)
     state.iterations += 1
-    _go_update_leaf!(cache, state, local_grad, adapted, ps_leaf)
+    _go_update_leaf!(cache, state, local_grad, adapted, ps)
     GeometricOptimizers._rmul!(GeometricOptimizers.direction(cache), step_size)
     GeometricOptimizers.update_section!(GeometricOptimizers.section(cache),
                                          GeometricOptimizers.section(state),
@@ -237,7 +304,7 @@ function _leaf_optim_step!(cache::GeometricOptimizers.OptimizerCache,
                                          retraction)
     GeometricOptimizers._copyto!(GeometricOptimizers.solution(cache),
                                   GeometricOptimizers.section(cache))
-    GeometricOptimizers._copyto!(ps_leaf, GeometricOptimizers.solution(cache))
+    GeometricOptimizers._copyto!(ps, GeometricOptimizers.solution(cache))
     GeometricOptimizers._copyto!(λY_leaf, GeometricOptimizers.section(cache))
     # `section(cache)` is `update_section!(section(state), direction, retraction)`, so copying it is
     # the same thing as retracting a second time -- and a retraction on a manifold is `O(N³)` where
@@ -267,6 +334,16 @@ function _leaf_optim_step!(cache::GMLEuclideanState, state::GMLEuclideanState,
     nothing
 end
 
+# Whether `λY` is a tree of sections to descend into, or a single `GlobalSection` standing in for a
+# whole subtree. Named rather than written as an `isa` union: it is one question, asked once, and the
+# answer decides whether a layer gets its own section or the whole tree. A section tree is a plain
+# `NamedTuple` today, because `GeometricOptimizers`' `GlobalSection(::NetworkParameters)` unwraps the
+# container; the container arm is here so that a change there is a `MethodError` and not a silent
+# hand-off of the whole tree to every layer.
+_is_section_tree(::NamedTuple) = true
+_is_section_tree(::NetworkParameters) = true
+_is_section_tree(_) = false
+
 # Recursive dispatcher over the *cache* tree.
 #
 # Deliberately hand-written rather than `NeuralNetworkParameters.foreachparameters`, which walks a
@@ -285,17 +362,17 @@ end
 #
 # The `λY` test names both container types for the same reason `_make_optimizer_cache` asks the
 # structural question first: a section tree is something to descend into, whichever type carries it.
-# Today only a `NamedTuple` ever arrives, because the `GlobalSection(::NetworkParameters)` method at
-# the top of this file unwraps the container before `GeometricOptimizers` sees it. That method is
-# GML's to give back once `GeometricOptimizers` depends on `NeuralNetworkParameters` — and its
-# replacement there is expected to return a *container* of sections, at which point `isa NamedTuple`
-# alone would be false and every layer would be handed the whole tree instead of its own section.
+# Only a `NamedTuple` ever arrives: `GeometricOptimizers`' own
+# `GlobalSection(ps::NetworkParameters) = GlobalSection(params(ps))` unwraps the container, so a
+# section tree is a plain `NamedTuple` whether the parameters are wrapped or not. Should that ever
+# return a *container* of sections instead, `isa NamedTuple` alone would be false here and every layer
+# would be handed the whole tree instead of its own section.
 function _tree_optim_step!(caches, states, dp, ps, λY, method, retraction, step_size)
     if caches isa NamedTuple
         for k in keys(caches)
             dp_k = dp[k]
             dp_k === nothing && continue
-            λY_k = λY isa ParameterSet ? λY[k] : λY
+            λY_k = _is_section_tree(λY) ? λY[k] : λY
             _tree_optim_step!(caches[k], states[k], dp_k, ps[k], λY_k,
                               method, retraction, step_size)
         end
diff --git a/src/pullbacks/zygote_pullback.jl b/src/pullbacks/zygote_pullback.jl
index 177df48bf..b92933f6e 100644
--- a/src/pullbacks/zygote_pullback.jl
+++ b/src/pullbacks/zygote_pullback.jl
@@ -32,6 +32,14 @@ end
 (_pullback::ZygotePullback)(ps, model, input_nt_output_nt::Tuple{<:QPTOAT, <:QPTOAT})::Tuple = Zygote.pullback(
     ps -> _pullback.loss(model, ps, input_nt_output_nt...), ps)
 
+# Both shapes, because a pullback produces both. `NeuralNetworkParameters` rewraps the tangent of a
+# `NetworkParameters` into one, while a reverse pass seeded with a bare `NamedTuple` hands one back.
+# A method per shape here as everywhere else in this release, and in the wrapped forms too: a union
+# over `NamedTuple` reads as breadth for its own sake, where two methods say which shape arrived.
+#
+# The comment goes *above* the docstring and not between it and the first method: a docstring attaches
+# to whatever follows it immediately, so a comment in between orphans it -- `docs/check_references.jl`
+# reports that as "documented in nowhere".
 """
     _get_contents(returned_pullback)
 
@@ -39,9 +47,16 @@ Unwrap the single element `Zygote` may wrap a pullback result in.
 
 Together with [`_get_params`](@ref) this makes up [`_processing`](@ref).
 """
-_get_contents(nt::ParameterSet) = nt
-_get_contents(nt::Tuple{<:ParameterSet}) = nt[1]
-function _get_contents(nt::AbstractVector{<:ParameterSet})
+_get_contents(nt::NetworkParameters) = nt
+_get_contents(nt::NamedTuple) = nt
+
+_get_contents(nt::Tuple{<:NetworkParameters}) = nt[1]
+_get_contents(nt::Tuple{<:NamedTuple}) = nt[1]
+
+_get_contents(nt::AbstractVector{<:NetworkParameters}) = _only_set(nt)
+_get_contents(nt::AbstractVector{<:NamedTuple}) = _only_set(nt)
+
+function _only_set(nt)
     length(nt) == 1 || throw(ArgumentError(
         "the pullback returned $(length(nt)) parameter sets, expected one."))
     nt[1]
diff --git a/test/optimizers/gml_gradient_dispatch.jl b/test/optimizers/gml_gradient_dispatch.jl
new file mode 100644
index 000000000..1d75f7eb1
--- /dev/null
+++ b/test/optimizers/gml_gradient_dispatch.jl
@@ -0,0 +1,114 @@
+# What `_GMLGradient` is called on, and by which method.
+#
+# `_GMLGradient` is a `GeometricOptimizers.Gradient`, so its functor methods share a method table with
+# upstream's. Two of them overlap on a shape this package actually produces, and the overlap is an
+# *ambiguity* rather than a choice — which is a run-time error several frames into a solve, not a
+# signature a reader would notice. So it is asserted here rather than left to whichever call path
+# happens to reach it.
+
+using GeometricMachineLearning
+using GeometricOptimizers: StiefelManifold, GrassmannManifold, rgrad
+using NeuralNetworkParameters: NetworkParameters
+using Test
+
+GML = GeometricMachineLearning
+
+# A bare `Manifold`. `_GMLGradient{T} <: GeometricOptimizers.Gradient{T}` and `Manifold{T} <:
+# AbstractMatrix{T}`, so `(::_GMLGradient{T})(::AbstractArray{T})` and upstream's
+# `(::Gradient{T})(::Manifold{T})` are ambiguous here: neither is more specific. Without a method of
+# its own this call is `MethodError: ... is ambiguous`.
+#
+# The value is this package's answer and not upstream's. Upstream's body is
+# `rgrad(x, reshape(grad(vec(x)), size(x)...))`, which evaluates an *inner* gradient — and a
+# `_GMLGradient` has none, since the whole point of it is that the Euclidean gradient was computed
+# already and is being carried in `dp`.
+@testset "a bare Manifold reaches this package's method, not upstream's" begin
+    for M in (StiefelManifold, GrassmannManifold)
+        Y  = rand(M, 4, 2)
+        dp = randn(4, 2)
+        g  = GML._GMLGradient{Float64, typeof(dp)}(dp)
+
+        @test g(Y) == rgrad(Y, dp)
+        @test which(g, Tuple{typeof(Y)}).module === GeometricMachineLearning
+    end
+end
+
+# The two shapes the optimizer actually hands it, so that the method above is not the only one
+# covered and a later narrowing of either shows up here.
+@testset "a wrapped layer and a plain array" begin
+    dp = randn(3, 3)
+    g  = GML._GMLGradient{Float64, typeof(dp)}(dp)
+    # a plain array leaf: the gradient is already Euclidean, so it passes through
+    @test g(randn(3, 3)) === dp
+
+    ps = NetworkParameters((weight = rand(StiefelManifold, 4, 2),))
+    dps = (weight = randn(4, 2),)
+    gps = GML._GMLGradient{Float64, typeof(dps)}(dps)
+    out = gps(ps)
+    # rebuilt in the shape of its *first* argument -- a container, which is what
+    # `GeometricOptimizers._copyto!` has a method for
+    @test out isa NetworkParameters
+    @test out.weight == rgrad(ps.weight, dps.weight)
+end
+
+# `_is_layer` is the rule for how a network is split into caches, so its boundary cases are stated
+# rather than inherited. An empty set is *zero* caches' worth and descends; `all` over an empty
+# collection would otherwise make it vacuously one.
+@testset "_is_layer at the boundary" begin
+    @test GML._is_layer((W = randn(2, 2), b = randn(2)))
+    @test GML._is_layer(NetworkParameters((W = randn(2, 2), b = randn(2))))
+    # a set whose values are branches is a subtree, and gets one cache per layer
+    @test !GML._is_layer((L1 = (W = randn(2, 2),), L2 = (W = randn(2, 2),)))
+    @test !GML._is_layer(NamedTuple())
+    @test !GML._is_layer(NetworkParameters(NamedTuple()))
+    @test !GML._is_layer(randn(2, 2))
+
+    # a layer whose weights do not share an element type is still one cache, with `T` the promotion
+    @test GML._is_layer((W = randn(Float32, 2, 2), b = randn(Float64, 2)))
+end
+
+# `_as_go_solution` is the wrap itself, and the property that makes the wrap free: it shares the leaf
+# arrays, so the optimizer step writes through to the network's own weights and nothing is copied
+# back. That is what makes this a boundary change rather than a change of behaviour, and it is the one
+# claim of this release that a reader cannot check by reading a signature.
+@testset "_as_go_solution wraps a layer and shares its arrays" begin
+    W, b = randn(2, 2), randn(2)
+    layer = (W = W, b = b)
+
+    wrapped = GML._as_go_solution(layer)
+    @test wrapped isa NetworkParameters
+    @test wrapped.W === W          # shared, not copied
+    @test wrapped.b === b
+
+    # already in the right shape: passed through untouched
+    ps = NetworkParameters(layer)
+    @test GML._as_go_solution(ps) === ps
+    Y = rand(StiefelManifold, 4, 2)
+    @test GML._as_go_solution(Y) === Y
+    v = randn(3)
+    @test GML._as_go_solution(v) === v
+
+    # a subtree is not a layer, so it is left alone for the tree branch to descend
+    tree = (L1 = (W = randn(2, 2),), L2 = (W = randn(2, 2),))
+    @test GML._as_go_solution(tree) === tree
+end
+
+# `_get_contents` unwraps whatever shape the pullback handed back. Six methods, because a reverse pass
+# produces a container or the bare `NamedTuple` it wraps, and `Zygote` may hand either back bare, in a
+# one-tuple, or in a one-element vector. Each is exercised here: splitting a union into a method per
+# shape means each method needs its own case, or the split trades one tested path for six untested
+# ones.
+@testset "_get_contents unwraps every shape a pullback returns" begin
+    nt = (L1 = (W = randn(2, 2), b = randn(2)),)
+    ps = NetworkParameters(nt)
+
+    for x in (nt, ps)
+        @test GML._get_contents(x) === x
+        @test GML._get_contents((x,)) === x
+        @test GML._get_contents([x]) === x
+    end
+
+    # more than one set is a caller error and says so, rather than silently taking the first
+    @test_throws ArgumentError GML._get_contents([ps, ps])
+    @test_throws ArgumentError GML._get_contents([nt, nt])
+end
diff --git a/test/optimizers/gradient_optimizer.jl b/test/optimizers/gradient_optimizer.jl
index 7ff70bbed..c79bd455c 100644
--- a/test/optimizers/gradient_optimizer.jl
+++ b/test/optimizers/gradient_optimizer.jl
@@ -13,7 +13,7 @@ function gradient_optimizer_euclidean(N; n_steps = 20, step_size = 1e-3)
     B = inv(rand(N, N))
     loss(ps) = norm(ps.A - B) ^ 2
     A = randn(N, N)
-    ps = (A = A,)
+    ps = NetworkParameters((A = A,))
     loss1 = loss(ps)
     o = Optimizer(GradientMethod(), ps; step_size = step_size)
     for _ in 1:n_steps
@@ -33,7 +33,7 @@ function gradient_optimizer_stiefel(N, n; n_steps = 20, step_size = 1e-3)
     B  = YB * YB'
     loss(ps) = norm(ps.Y * ps.Y' - B) ^ 2
     Y  = rand(StiefelManifold, N, n)
-    ps = (Y = Y,)
+    ps = NetworkParameters((Y = Y,))
     loss1 = loss(ps)
     o = Optimizer(GradientMethod(), ps; step_size = step_size)
     for _ in 1:n_steps
diff --git a/test/optimizers/momentum_optimizer.jl b/test/optimizers/momentum_optimizer.jl
index 4a9322355..0351140a1 100644
--- a/test/optimizers/momentum_optimizer.jl
+++ b/test/optimizers/momentum_optimizer.jl
@@ -15,7 +15,7 @@ For a *constant* gradient ``g`` the momentum ``p^{(t)} = \alpha{}p^{(t-1)} + g``
 growing one that never saturates, which is what this test rules out.
 """
 function momentum_matches_geometric_sum(; α = 0.9, η = 0.1, n_steps = 10)
-    ps = (A = zeros(2, 2),)
+    ps = NetworkParameters((A = zeros(2, 2),))
     o = Optimizer(MomentumMethod(α), ps; step_size = η)
     λY = GlobalSection(ps)
     previous = 0.0
@@ -36,7 +36,7 @@ function momentum_optimizer_stiefel(N, n; n_steps = 50, step_size = 1e-3, α = 0
     YB = rand(StiefelManifold, N, n)
     B = YB * YB'
     loss(ps) = norm(ps.Y * ps.Y' - B)^2
-    ps = (Y = rand(StiefelManifold, N, n),)
+    ps = NetworkParameters((Y = rand(StiefelManifold, N, n),))
     loss1 = loss(ps)
     o = Optimizer(MomentumMethod(α), ps; step_size = step_size)
     λY = GlobalSection(ps)
diff --git a/test/optimizers/optimizer_convergence_tests/adam_with_learning_rate_decay.jl b/test/optimizers/optimizer_convergence_tests/adam_with_learning_rate_decay.jl
index b9b53cb99..4f5648d83 100644
--- a/test/optimizers/optimizer_convergence_tests/adam_with_learning_rate_decay.jl
+++ b/test/optimizers/optimizer_convergence_tests/adam_with_learning_rate_decay.jl
@@ -81,7 +81,7 @@ schedule, which `test/adam_optimizer_with_decay.jl` upstream asserts does not ha
 """
 function schedule_starts_at_one(; n_epochs = 100, η₁ = 1e-2, η₂ = 1e-6)
     method = AdamOptimizerWithDecay(n_epochs, Float64; η₁ = η₁, η₂ = η₂)
-    o = Optimizer((weight = zeros(2, 2),); method...)
+    o = Optimizer(NetworkParameters((weight = zeros(2, 2),)); method...)
 
     γ = exp(log(η₂ / η₁) / n_epochs)
     @test o.step_size isa DecayingStatic
@@ -94,7 +94,7 @@ function schedule_starts_at_one(; n_epochs = 100, η₁ = 1e-2, η₂ = 1e-6)
 
     # and the same thing through the public entry point: one step of a Euclidean parameter with a
     # gradient of `1` moves it by `α₁ / (√1 + δ) ≈ α₁`, so the distance travelled reports the α used
-    ps = (weight = zeros(2, 2),)
+    ps = NetworkParameters((weight = zeros(2, 2),))
     opt = Optimizer(ps; method...)
     optimization_step!(opt, GlobalSection(ps), ps, (weight = ones(2, 2),))
     @test opt.iterations == 1
diff --git a/test/runtests.jl b/test/runtests.jl
index 252118a94..a4c5e9686 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -101,6 +101,9 @@ end
 @safetestset "Optimizers with structured (non-manifold) weights                               " begin
     include("optimizers/structured_array_parameters.jl")
 end
+@safetestset "_GMLGradient dispatch                                                           " begin
+    include("optimizers/gml_gradient_dispatch.jl")
+end
 @safetestset "map_to_cpu                                                                      " begin
     include("map_to_cpu_tests.jl")
 end
diff --git a/test/transformer_related/multi_head_attention_stiefel_optim_cache.jl b/test/transformer_related/multi_head_attention_stiefel_optim_cache.jl
index 932ca2ee9..be54d7566 100644
--- a/test/transformer_related/multi_head_attention_stiefel_optim_cache.jl
+++ b/test/transformer_related/multi_head_attention_stiefel_optim_cache.jl
@@ -6,6 +6,13 @@ using GeometricMachineLearning, Test
 import GeometricOptimizers
 using GeometricOptimizers: GradientCache, MomentumCache, AdamCache
 import Random, LinearAlgebra
+# The cache and the state hold a `NetworkParameters` per layer: a layer is wrapped at the
+# `GeometricOptimizers` boundary, and the wrap shares the leaf arrays. `mapparameters` and not `map`,
+# because it recurses on the branches and so reaches the leaves whichever shape it is handed.
+#
+# Each helper below has two methods, because two shapes arrive: a wrapped layer is a
+# `NetworkParameters` and the tree above it is a plain `NamedTuple`.
+using NeuralNetworkParameters: NetworkParameters, mapparameters
 
 Random.seed!(1234)
 
@@ -16,6 +23,7 @@ The direction δ and moments m₁, m₂ should be zero-initialised.
 """
 _check_slahm_zero(A::StiefelLieAlgHorMatrix{T}, tol) where T =
     (@test typeof(A) <: StiefelLieAlgHorMatrix; @test LinearAlgebra.norm(A) < tol)
+_check_slahm_zero(A::NetworkParameters, tol) = foreach(v -> _check_slahm_zero(v, tol), values(A))
 _check_slahm_zero(A::NamedTuple, tol) = foreach(v -> _check_slahm_zero(v, tol), values(A))
 
 function check_adam_cache(C::GeometricOptimizers.OptimizerCache{T}, tol=T(10) * eps(T)) where T
@@ -23,7 +31,8 @@ function check_adam_cache(C::GeometricOptimizers.OptimizerCache{T}, tol=T(10) *
     _check_slahm_zero(C.δ, tol)
     _check_slahm_zero(C.m₁, tol)
 end
-check_adam_cache(B::NamedTuple) = map(check_adam_cache, B)
+check_adam_cache(B::NetworkParameters) = mapparameters(check_adam_cache, B)
+check_adam_cache(B::NamedTuple) = mapparameters(check_adam_cache, B)
 
 @doc raw"""
 This checks if the momentum cache was set up in the right way.
@@ -34,7 +43,8 @@ function check_momentum_cache(C::GeometricOptimizers.OptimizerCache{T}, tol=T(10
     @test C isa MomentumCache
     _check_slahm_zero(C.δ, tol)
 end
-check_momentum_cache(B::NamedTuple) = map(check_momentum_cache, B)
+check_momentum_cache(B::NetworkParameters) = mapparameters(check_momentum_cache, B)
+check_momentum_cache(B::NamedTuple) = mapparameters(check_momentum_cache, B)
 
 @doc raw"""
 This checks if the gradient cache was set up in the right way.
@@ -44,7 +54,8 @@ function check_gradient_cache(C::GeometricOptimizers.OptimizerCache{T}) where T
     @test C isa GradientCache
     @test hasproperty(C, :δ)
 end
-check_gradient_cache(B::NamedTuple) = map(check_gradient_cache, B)
+check_gradient_cache(B::NetworkParameters) = mapparameters(check_gradient_cache, B)
+check_gradient_cache(B::NamedTuple) = mapparameters(check_gradient_cache, B)
 
 @doc raw"""
 This checks if all the caches are set up in the right way for the `MultiHeadAttention` layer with Stiefel weights.
diff --git a/test/transformer_related/multi_head_attention_stiefel_retraction.jl b/test/transformer_related/multi_head_attention_stiefel_retraction.jl
index f113a3ce0..c7c55966b 100644
--- a/test/transformer_related/multi_head_attention_stiefel_retraction.jl
+++ b/test/transformer_related/multi_head_attention_stiefel_retraction.jl
@@ -5,6 +5,13 @@ using GeometricMachineLearning: cayley
 # so they are named qualified rather than through a re-export.
 using GeometricOptimizers: MomentumCache
 import Random, Test, LinearAlgebra, KernelAbstractions
+# The cache and the state hold a `NetworkParameters` per layer: a layer is wrapped at the
+# `GeometricOptimizers` boundary, and the wrap shares the leaf arrays. `mapparameters` and not `map`,
+# because it recurses on the branches and so reaches the leaves whichever shape it is handed.
+#
+# Each helper below has two methods, because two shapes arrive: a wrapped layer is a
+# `NetworkParameters` and the tree above it is a plain `NamedTuple`.
+using NeuralNetworkParameters: NetworkParameters, mapparameters
 
 Random.seed!(1234)
 
@@ -16,7 +23,8 @@ function check_retraction_geodesic(A::AbstractMatrix{T}, tol=eps(T)) where T
     @test typeof(A_retracted) <: StiefelManifold
     @test LinearAlgebra.norm(A_retracted - StiefelProjection(A_retracted)) < tol
 end
-check_retraction_geodesic(cache::NamedTuple) = map(check_retraction_geodesic, cache)
+check_retraction_geodesic(cache::NetworkParameters) = mapparameters(check_retraction_geodesic, cache)
+check_retraction_geodesic(cache::NamedTuple) = mapparameters(check_retraction_geodesic, cache)
 check_retraction_geodesic(B::MomentumCache) = check_retraction_geodesic(B.δ)
 
 @doc raw"""
@@ -27,7 +35,8 @@ function check_retraction_cayley(A::AbstractMatrix{T}, tol=eps(T)) where T
     @test typeof(A_retracted) <: StiefelManifold
     @test LinearAlgebra.norm(A_retracted - StiefelProjection(A_retracted)) < tol
 end
-check_retraction_cayley(cache::NamedTuple) = map(check_retraction_cayley, cache)
+check_retraction_cayley(cache::NetworkParameters) = mapparameters(check_retraction_cayley, cache)
+check_retraction_cayley(cache::NamedTuple) = mapparameters(check_retraction_cayley, cache)
 check_retraction_cayley(B::MomentumCache) = check_retraction_cayley(B.δ)
 
 @doc raw"""
diff --git a/test/transformer_related/multi_head_attention_stiefel_setup.jl b/test/transformer_related/multi_head_attention_stiefel_setup.jl
index 52e8a19a6..20f38d2c7 100644
--- a/test/transformer_related/multi_head_attention_stiefel_setup.jl
+++ b/test/transformer_related/multi_head_attention_stiefel_setup.jl
@@ -2,6 +2,13 @@ using GeometricMachineLearning, Test
 # see the note in `multi_head_attention_stiefel_retraction.jl`
 using GeometricOptimizers: MomentumCache
 import Random, Test, LinearAlgebra, KernelAbstractions
+# The cache and the state hold a `NetworkParameters` per layer: a layer is wrapped at the
+# `GeometricOptimizers` boundary, and the wrap shares the leaf arrays. `mapparameters` and not `map`,
+# because it recurses on the branches and so reaches the leaves whichever shape it is handed.
+#
+# Each helper below has two methods, because two shapes arrive: a wrapped layer is a
+# `NetworkParameters` and the tree above it is a plain `NamedTuple`.
+using NeuralNetworkParameters: NetworkParameters, mapparameters
 
 Random.seed!(1234)
 
@@ -12,7 +19,8 @@ function check_setup(A::AbstractMatrix{T}, tol=T(10)*eps(T)) where T
     @test typeof(A) <: StiefelManifold
     @test check(A) < tol
 end
-check_setup(ps::NamedTuple) = map(check_setup, ps)
+check_setup(ps::NetworkParameters) = mapparameters(check_setup, ps)
+check_setup(ps::NamedTuple) = mapparameters(check_setup, ps)
 check_setup(ps::NetworkParameters) = check_setup(GeometricMachineLearning.params(ps))
 
 @doc raw"""
@@ -22,7 +30,8 @@ function check_grad_setup(B::AbstractMatrix{T}, tol=T(10)*eps(T)) where T
     @test typeof(B) <: StiefelLieAlgHorMatrix
     @test LinearAlgebra.norm(B) < tol
 end
-check_grad_setup(gx::NamedTuple) = map(check_grad_setup, gx)
+check_grad_setup(gx::NetworkParameters) = mapparameters(check_grad_setup, gx)
+check_grad_setup(gx::NamedTuple) = mapparameters(check_grad_setup, gx)
 check_grad_setup(B::MomentumCache) = check_grad_setup(B.δ)
 
 @doc raw"""