Fix more RAT-v4 iteration semantics misses

Three more spots that iterate or index `AbstractVectorOfArray`
with pre-RAT-v4 semantics, surfaced by CI on the previous commit:

1. `lib/BoundaryValueDiffEqMIRK/src/adaptivity.jl:39,93,148` —
   `[maximum(abs, d) for d in errors]` where `cache.errors` is
   a `VectorOfArray`. Under RAT v4 this iterates **scalars**,
   producing a `Vector{Float64}` of length `prod(size(errors))`
   instead of one entry per timestep, which then fails
   `ŝ ./= mesh_dt` with `DimensionMismatch`. Fix: iterate
   `errors.u`. (Lines 204–205 already had this fix.)
2. `lib/BoundaryValueDiffEqFIRK/src/adaptivity.jl:263` —
   identical pattern with `defect`.
3. `lib/BoundaryValueDiffEqShooting/src/multiple_shooting.jl:531-532`
   — `length(first(u0_))` / `similar(first(u0_), …)` where `u0_`
   is the `VectorOfArray` returned by `__initial_guess_on_mesh`.
   Under RAT v4 `first(u0_)` is the first scalar, so
   `similar(::Float64, ::Int)` throws `MethodError`. Fix:
   `first(u0_.u)`.

Bumps MIRK / FIRK / Shooting patch versions to 1.16.1.

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>

Author: ChrisRackauckas

lib/BoundaryValueDiffEqFIRK/Project.toml

@@ -1,6 +1,6 @@
 name = "BoundaryValueDiffEqFIRK"
 uuid = "85d9eb09-370e-4000-bb32-543851f73618"
-version = "1.16.0"
+version = "1.16.1"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"

lib/BoundaryValueDiffEqFIRK/src/adaptivity.jl

@@ -260,7 +260,7 @@ Generate new mesh based on the defect.
 
     info = ReturnCode.Success
 
-    ŝ = [maximum(abs, d) for d in defect]  # Broadcasting breaks GPU Compilation
+    ŝ = [maximum(abs, d) for d in defect.u]  # Broadcasting breaks GPU Compilation
     ŝ .= (ŝ ./ abstol) .^ (T(1) / (order + 1))
     r₁ = maximum(ŝ)
     r₂ = sum(ŝ)

lib/BoundaryValueDiffEqMIRK/Project.toml

@@ -1,6 +1,6 @@
 name = "BoundaryValueDiffEqMIRK"
 uuid = "1a22d4ce-7765-49ea-b6f2-13c8438986a6"
-version = "1.16.0"
+version = "1.16.1"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"

lib/BoundaryValueDiffEqMIRK/src/adaptivity.jl

@@ -36,7 +36,7 @@ Generate new mesh based on the defect or the global error.
 
     info = ReturnCode.Success
 
-    ŝ = [maximum(abs, d) for d in errors]  # Broadcasting breaks GPU Compilation
+    ŝ = [maximum(abs, d) for d in errors.u]  # Broadcasting breaks GPU Compilation
     ŝ .= (ŝ ./ abstol) .^ (T(1) / (order + 1))
     r₁ = maximum(ŝ)
     r₂ = sum(ŝ)
@@ -90,7 +90,7 @@ end
 
     info = ReturnCode.Success
 
-    ŝ = [maximum(abs, d) for d in errors]
+    ŝ = [maximum(abs, d) for d in errors.u]
     ŝ .= (ŝ ./ abstol) .^ (T(1) / order)
     r₁ = maximum(ŝ)
     r₂ = sum(ŝ)
@@ -145,7 +145,7 @@ end
 
     info = ReturnCode.Success
 
-    ŝ = [maximum(abs, d) for d in errors]
+    ŝ = [maximum(abs, d) for d in errors.u]
     ŝ .= (ŝ ./ abstol) .^ (T(1) / (order + 1))
     r₁ = maximum(ŝ)
     r₂ = sum(ŝ)

lib/BoundaryValueDiffEqShooting/Project.toml

@@ -1,6 +1,6 @@
 name = "BoundaryValueDiffEqShooting"
 uuid = "ed55bfe0-3725-4db6-871e-a1dc9f42a757"
-version = "1.16.0"
+version = "1.16.1"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"

lib/BoundaryValueDiffEqShooting/src/multiple_shooting.jl

@@ -528,8 +528,8 @@ end
     nodes .= range(tspan[1], tspan[2]; length = nshoots + 1)
     u0_ = __initial_guess_on_mesh(u0, nodes, p)
 
-    N = length(first(u0_))
-    u_at_nodes = similar(first(u0_), (nshoots + 1) * N)
+    N = length(first(u0_.u))
+    u_at_nodes = similar(first(u0_.u), (nshoots + 1) * N)
     recursive_flatten!(u_at_nodes, u0_.u)
 
     return u_at_nodes