henry2004y · henry2004y · Apr 16, 2026 · Apr 16, 2026 · Apr 16, 2026 · Apr 16, 2026
diff --git a/README.md b/README.md
@@ -56,7 +56,7 @@ sol = TestParticle.solve(prob; dt, savestepinterval)[1]
 
 Besides the standard Boris method, we also support various versions of Boris solvers including:
 - **Multistep Boris**: fixed time step with `n > 1` substeps.
-- **Adaptive Boris**: uses `AdaptiveBoris()` for automatic time step selection based on local gyroperiod.
+- **Adaptive Boris**: uses `Boris(safety = 0.1)` to automatically adjust time step based on local gyroperiod.
 
 For plotting with Makie,
 

diff --git a/benchmark/benchmarks.jl b/benchmark/benchmarks.jl
@@ -168,7 +168,7 @@ SUITE["trace"]["numerical field"]["Hyper Boris (n=2, N=4)"] = @benchmarkable TP.
 SUITE["trace"]["numerical field"]["Hyper Boris (n=2, N=6)"] = @benchmarkable TP.solve(
     $prob_boris; dt = 1 / 7, savestepinterval = 10, n = 2, N = 6
 )
-alg_adaptive = AdaptiveBoris(safety = 0.1)
+alg_adaptive = Boris(safety = 0.1)
 SUITE["trace"]["numerical field"]["Adaptive Boris"] = @benchmarkable TP.solve(
     $prob_boris, $alg_adaptive
 )

diff --git a/docs/examples/features/demo_adaptive_accuracy.jl b/docs/examples/features/demo_adaptive_accuracy.jl
@@ -1,7 +1,7 @@
 # # Adaptive Solver Accuracy
 #
 # This example compares the accuracy and energy conservation of solvers that
-# use **changing time steps**: the native `AdaptiveBoris` solver and the
+# use **changing time steps**: the native `Boris(safety = η)` solver and the
 # adaptive ODE solvers `Tsit5` and `Vern6`.
 #
 # Unlike the fixed-step tests in [Solver Accuracy Analysis](@ref) and
@@ -91,16 +91,16 @@ prob_ode1 = ODEProblem(trace_normalized!, u0_1, tspan1, param1)
 
 adaptive_solvers_1 = [
     (
-        "AdaptiveBoris (η=1/16)",
-        AdaptiveBoris(; safety = 1 / 16),
+        "Boris (η=1/16)",
+        Boris(; safety = 1 / 16),
     ),
     (
-        "AdaptiveBoris (η=1/8)",
-        AdaptiveBoris(; safety = 1 / 8),
+        "Boris (η=1/8)",
+        Boris(; safety = 1 / 8),
     ),
     (
-        "AdaptiveBoris (η=1/4)",
-        AdaptiveBoris(; safety = 1 / 4),
+        "Boris (η=1/4)",
+        Boris(; safety = 1 / 4),
     ),
 ]
 
@@ -134,7 +134,7 @@ ax1b = Axis(
     f1[1, 2];
     xlabel = "Time [Gyroperiod]",
     ylabel = L"\Delta t / T_\mathrm{gyro}",
-    title = "AdaptiveBoris: Time Step History",
+    title = "Adaptive Boris: Time Step History",
 )
 
 for (i, (name, alg)) in enumerate(adaptive_solvers_1)
@@ -201,7 +201,7 @@ plot_table(results1) #hide
 #
 # With E = [0, 0.5, 0.1] and B = [0, 0, 1], an exact velocity
 # solution exists (Section 6, Zenitani & Kato 2025). We sweep the
-# `safety` parameter η (representing dt / T_gyro) for `AdaptiveBoris`
+# `safety` parameter η (representing dt / T_gyro) for `Boris(safety = η)`
 # and the tolerance for the ODE solvers, then plot maximum velocity
 # error vs. number of time steps (computational cost).
 
@@ -245,7 +245,7 @@ boris_steps = Int[]
 boris_errors = Float64[]
 
 for s in safety_values
-    alg = AdaptiveBoris(; safety = s)
+    alg = Boris(; safety = s)
     sol = TP.solve(prob_tp2, alg)[1]
     push!(boris_steps, length(sol.t))
     push!(boris_errors, max_velocity_error(sol))
@@ -288,7 +288,7 @@ ax2 = Axis(
 
 scatterlines!(
     ax2, boris_steps, boris_errors;
-    label = "AdaptiveBoris",
+    label = "Adaptive Boris",
     marker = :circle, linewidth = 2,
 )
 
@@ -335,12 +335,12 @@ prob_ode3 = ODEProblem(trace_normalized!, u0_2, tspan3, param3)
 
 adaptive_solvers_3 = [
     (
-        "AdaptiveBoris (η=0.1)",
-        AdaptiveBoris(; safety = 0.1),
+        "Boris (η=0.1)",
+        Boris(; safety = 0.1),
     ),
     (
-        "AdaptiveBoris (η=0.2)",
-        AdaptiveBoris(; safety = 0.2),
+        "Boris (η=0.2)",
+        Boris(; safety = 0.2),
     ),
 ]
 

diff --git a/docs/examples/features/demo_boris.jl b/docs/examples/features/demo_boris.jl
@@ -99,7 +99,7 @@ sol_boris_2 = TP.solve(prob; dt, n = 2)[1];
 sol_boris_4 = TP.solve(prob; dt, n = 4)[1];
 sol_boris_hyper = TP.solve(prob; dt, n = 2, N = 4)[1];
 
-alg_adaptive = AdaptiveBoris(safety = 0.1)
+alg_adaptive = Boris(safety = 0.1)
 sol_boris_adaptive = TP.solve(prob, alg_adaptive)[1];
 
 # Let's compare against the default ODE solver `Tsit5` from DifferentialEquations.jl, in both fixed time step mode and adaptive mode:
@@ -148,7 +148,7 @@ sol_boris_4 = TP.solve(prob_boris; dt, savestepinterval = 36, n = 4)[1];
 sol_boris_hyper = TP.solve(prob_boris; dt, savestepinterval = 36, n = 2, N = 4)[1];
 sol_boris_adaptive = TP.solve(
     prob_boris,
-    AdaptiveBoris(safety = 0.1)
+    Boris(safety = 0.1)
 )[1]
 sol1 = solve(prob, Tsit5(); adaptive = false, dt, dense = false, saveat = dt);
 sol2 = solve(prob, Tsit5());

diff --git a/docs/src/tutorial/advanced_boris.md b/docs/src/tutorial/advanced_boris.md
@@ -74,7 +74,7 @@ Combining both $n > 1$ and $N > 2$ ensures ultra-high stability tracking over dr
 
 ## 3. Adaptive Boris Method
 
-The `AdaptiveBoris` solver adjusts the time step $\Delta t$ dynamically based on the local cyclotron frequency $\Omega_c = |q B / m|$. This is particularly useful in systems with strong magnetic field gradients, such as magnetic mirrors or planetary magnetospheres, where the required resolution varies significantly along the particle's trajectory.
+The `Boris(safety = η)` solver adjusts the time step $\Delta t$ dynamically based on the local cyclotron frequency $\Omega_c = |q B / m|$. This is particularly useful in systems with strong magnetic field gradients, such as magnetic mirrors or planetary magnetospheres, where the required resolution varies significantly along the particle's trajectory.
 
 The time step is determined by:
 ```math
@@ -113,10 +113,10 @@ sol = TestParticle.solve(prob; dt, n=2, N=4)
 
 ### Adaptive Boris
 
-You can use the adaptive solver by passing an `AdaptiveBoris` object as the second argument to `solve`.
+You can use the adaptive solver by passing `Boris(safety = η)` as the second argument to `solve`.
 
 ```julia
 # Adaptive Boris with safety factor 0.05 (20 steps per period)
-alg = AdaptiveBoris(safety=0.05)
+alg = Boris(safety=0.05)
 sol = TestParticle.solve(prob, alg)[1]
 ```
diff --git a/lib/OrdinaryDiffEqBoris/Project.toml b/lib/OrdinaryDiffEqBoris/Project.toml
@@ -0,0 +1,27 @@
+name = "OrdinaryDiffEqBoris"
+uuid = "1e8b2ed2-ea1c-4b6e-b3f5-67c29d6bd197"
+version = "0.1.0"
+
+[deps]
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+MuladdMacro = "46d2c3a1-f734-5fdb-9937-b9b9aeba4221"
+OrdinaryDiffEqCore = "bbf590c4-e513-4bbe-9b18-05decba2e5d8"
+RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
+Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
+StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
+
+[compat]
+MuladdMacro = "0.2"
+OrdinaryDiffEqCore = "3"
+RecursiveArrayTools = "3.24, 4"
+Reexport = "1.2.2"
+SciMLBase = "2, 3"
+StaticArrays = "1"
+julia = "1.10"
+
+[extras]
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+
+[targets]
+test = ["Test"]
diff --git a/lib/OrdinaryDiffEqBoris/src/OrdinaryDiffEqBoris.jl b/lib/OrdinaryDiffEqBoris/src/OrdinaryDiffEqBoris.jl
@@ -0,0 +1,21 @@
+module OrdinaryDiffEqBoris
+
+using Reexport
+@reexport using SciMLBase
+import OrdinaryDiffEqCore: OrdinaryDiffEqAlgorithm, OrdinaryDiffEqAdaptiveAlgorithm,
+    OrdinaryDiffEqMutableCache, OrdinaryDiffEqConstantCache,
+    alg_order, alg_cache, isfsal, initialize!, perform_step!
+using RecursiveArrayTools
+using StaticArrays
+using MuladdMacro
+using LinearAlgebra
+
+include("algorithms.jl")
+include("alg_utils.jl")
+include("boris_caches.jl")
+include("boris_perform_step.jl")
+
+export Boris
+export MultistepBoris2, MultistepBoris4, MultistepBoris6
+
+end
diff --git a/lib/OrdinaryDiffEqBoris/src/alg_utils.jl b/lib/OrdinaryDiffEqBoris/src/alg_utils.jl
@@ -0,0 +1,5 @@
+alg_order(alg::Boris) = 2
+isfsal(alg::Boris) = false
+
+alg_order(alg::MultistepBoris{N}) where {N} = 2
+isfsal(alg::MultistepBoris{N}) where {N} = false
diff --git a/lib/OrdinaryDiffEqBoris/src/algorithms.jl b/lib/OrdinaryDiffEqBoris/src/algorithms.jl
@@ -0,0 +1,45 @@
+"""
+    Boris()
+
+The standard Boris method for particle pushing in electric and magnetic fields.
+
+This solver expects a problem where `p` is structured as `(q2m, m, E, B, ...)`, which matches the signature used by TestParticle.jl's `TraceProblem`.
+"""
+struct Boris{T} <: OrdinaryDiffEqAlgorithm
+    safety::T
+end
+Boris(; safety = 0.0) = Boris(safety)
+
+"""
+    MultistepBoris{N}(; n=1)
+
+The Multistep/Hyper Boris method of order `N`.
+`n` specifies the number of subcycles.
+`N` specifies the gyrophase correction order. `N=2` corresponds to the Multicycle solver, while `N=4` or `N=6` are the Hyper Boris solvers.
+"""
+struct MultistepBoris{N, T} <: OrdinaryDiffEqAlgorithm
+    n::Int
+    safety::T
+end
+MultistepBoris{N}(; n::Int = 1, safety = 0.0) where {N} = MultistepBoris{N, typeof(safety)}(n, safety)
+
+"""
+    MultistepBoris2(; n=1)
+
+The Multicycle Boris method (MultistepBoris with N=2).
+"""
+const MultistepBoris2 = MultistepBoris{2}
+
+"""
+    MultistepBoris4(; n=1)
+
+The 4th order Hyper Boris method (MultistepBoris with N=4).
+"""
+const MultistepBoris4 = MultistepBoris{4}
+
+"""
+    MultistepBoris6(; n=1)
+
+The 6th order Hyper Boris method (MultistepBoris with N=6).
+"""
+const MultistepBoris6 = MultistepBoris{6}
diff --git a/lib/OrdinaryDiffEqBoris/src/boris_caches.jl b/lib/OrdinaryDiffEqBoris/src/boris_caches.jl
@@ -0,0 +1,49 @@
+struct BorisConstantCache <: OrdinaryDiffEqConstantCache end
+
+struct BorisCache{uType, rateType} <: OrdinaryDiffEqMutableCache
+    u::uType
+    uprev::uType
+    tmp::uType
+    k::rateType
+end
+
+function alg_cache(
+        alg::Boris, u, rate_prototype, ::Type{uEltypeNoUnits},
+        ::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t, dt, reltol, p, calck,
+        ::Val{false}, args...; kwargs...
+    ) where {uEltypeNoUnits, uBottomEltypeNoUnits, tTypeNoUnits}
+    return BorisConstantCache()
+end
+
+function alg_cache(
+        alg::Boris, u, rate_prototype, ::Type{uEltypeNoUnits},
+        ::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t, dt, reltol, p, calck,
+        ::Val{true}, args...; kwargs...
+    ) where {uEltypeNoUnits, uBottomEltypeNoUnits, tTypeNoUnits}
+    return BorisCache(u, uprev, similar(u), similar(rate_prototype))
+end
+
+struct MultistepBorisConstantCache <: OrdinaryDiffEqConstantCache end
+
+struct MultistepBorisCache{uType, rateType} <: OrdinaryDiffEqMutableCache
+    u::uType
+    uprev::uType
+    tmp::uType
+    k::rateType
+end
+
+function alg_cache(
+        alg::MultistepBoris{N}, u, rate_prototype, ::Type{uEltypeNoUnits},
+        ::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t, dt, reltol, p, calck,
+        ::Val{false}, args...; kwargs...
+    ) where {uEltypeNoUnits, uBottomEltypeNoUnits, tTypeNoUnits, N}
+    return MultistepBorisConstantCache()
+end
+
+function alg_cache(
+        alg::MultistepBoris{N}, u, rate_prototype, ::Type{uEltypeNoUnits},
+        ::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t, dt, reltol, p, calck,
+        ::Val{true}, args...; kwargs...
+    ) where {uEltypeNoUnits, uBottomEltypeNoUnits, tTypeNoUnits, N}
+    return MultistepBorisCache(u, uprev, similar(u), similar(rate_prototype))
+end