Integrators

What are Integrators?

Integrators discretize continuous-time dynamics into constraints for the NLP solver. They implement the relationship:

\[x_{k+1} = \Phi(x_k, u_k, \Delta t_k)\]

where Φ approximates the continuous evolution ẋ = f(x, u, t).

using DirectTrajOpt
using NamedTrajectories
using LinearAlgebra

BilinearIntegrator

Overview

Used for control-linear (bilinear) dynamics:

\[\dot{x} = (G_0 + \sum_i u_i G_i) x\]

where:

  • G₀ is the drift term (dynamics with no control)
  • Gᵢ are the drive terms (how controls affect the system)
  • uᵢ are the control inputs

How it Works

Uses the matrix exponential for exact integration:

\[x_{k+1} = \exp(\Delta t \cdot G(u_k)) x_k\]

where G(u) = G₀ + Σᵢ uᵢ Gᵢ.

Example: Simple 2D System

N = 50
traj = NamedTrajectory(
    (x = randn(2, N), u = randn(1, N), Δt = fill(0.1, N));
    timestep = :Δt,
    controls = :u,
    initial = (x = [1.0, 0.0],),
    final = (x = [0.0, 1.0],),
)
N = 50, (x = 1:2, u = 3:3, → Δt = 4:4)

Define drift (natural dynamics) and drives (control terms)

G_drift = [-0.1 1.0; -1.0 -0.1]     # Damped oscillator
G_drives = [[0.0 1.0; 1.0 0.0]]     # Symmetric control coupling
1-element Vector{Matrix{Float64}}:
 [0.0 1.0; 1.0 0.0]

Create generator function

G = u -> G_drift + sum(u .* G_drives)
#2 (generic function with 1 method)

Create integrator

integrator = BilinearIntegrator(G, :x, :u, traj)
BilinearIntegrator{DirectTrajOpt.Integrators.var"#7#8"{Main.var"#2#3"}}(DirectTrajOpt.Integrators.var"#7#8"{Main.var"#2#3"}(Main.var"#2#3"()), :x, :u, 2, 6, 98)

Multiple Drives Example

traj_multi = NamedTrajectory(
    (x = randn(3, N), u = randn(2, N), Δt = fill(0.1, N));
    timestep = :Δt,
    controls = :u,
)

G_drift_3d = [
    0.0 1.0 0.0;
    -1.0 0.0 0.0;
    0.0 0.0 -0.1
]

G_drives_3d = [
    [1.0 0.0 0.0; 0.0 0.0 0.0; 0.0 0.0 0.0],  # Drive 1
    [0.0 0.0 0.0; 0.0 1.0 0.0; 0.0 0.0 0.0],   # Drive 2
]

G_multi = u -> G_drift_3d + sum(u .* G_drives_3d)

integrator_multi = BilinearIntegrator(G_multi, :x, :u, traj_multi)
BilinearIntegrator{DirectTrajOpt.Integrators.var"#7#8"{Main.var"#5#6"}}(DirectTrajOpt.Integrators.var"#7#8"{Main.var"#5#6"}(Main.var"#5#6"()), :x, :u, 3, 9, 147)

When to Use BilinearIntegrator

✓ Quantum systems (Hamiltonian evolution) ✓ Rotating systems (attitude dynamics) ✓ Systems linear in controls ✓ When you want exact integration (no discretization error)

TimeDependentBilinearIntegrator

Overview

For time-varying bilinear dynamics:

\[\dot{x} = (G_0(t) + \sum_i u_i(t) G_i(t)) x\]

The generator function now depends on both control and time.

Example: Periodic Disturbance

traj_td = NamedTrajectory(
    (
        x = randn(2, N),
        u = randn(1, N),
        t = collect(range(0, 5, N)),  # time variable
        Δt = fill(0.1, N),
    );
    timestep = :Δt,
    controls = :u,
)
N = 50, (x = 1:2, u = 3:3, t = 4:4, → Δt = 5:5)

Time-dependent generator

G_td = (u, t) -> [-0.1 + 0.5*sin(t) 1.0; -1.0 -0.1] + u[1] * [0.0 1.0; 1.0 0.0]

integrator_td = TimeDependentBilinearIntegrator(G_td, :x, :u, :t, traj_td)
TimeDependentBilinearIntegrator{DirectTrajOpt.Integrators.var"#21#22"{@NamedTuple{}, SciMLBase.ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, DirectTrajOpt.Integrators.var"#f!#20"{Main.var"#8#9"}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing, Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Nothing, @NamedTuple{}}, SciMLBase.StandardODEProblem}}}(DirectTrajOpt.Integrators.var"#21#22"{@NamedTuple{}, SciMLBase.ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, DirectTrajOpt.Integrators.var"#f!#20"{Main.var"#8#9"}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing, Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Nothing, @NamedTuple{}}, SciMLBase.StandardODEProblem}}(NamedTuple(), SciMLBase.ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, DirectTrajOpt.Integrators.var"#f!#20"{Main.var"#8#9"}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing, Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Nothing, @NamedTuple{}}, SciMLBase.StandardODEProblem}(SciMLBase.ODEFunction{true, SciMLBase.AutoSpecialize, DirectTrajOpt.Integrators.var"#f!#20"{Main.var"#8#9"}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing, Nothing, Nothing}(DirectTrajOpt.Integrators.var"#f!#20"{Main.var"#8#9"}(Main.var"#8#9"(), Core.Box(DirectTrajOpt.Integrators.var"#18#19"{Int64}(1))), LinearAlgebra.UniformScaling{Bool}(true), nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, SciMLBase.DEFAULT_OBSERVED, nothing, nothing, nothing, nothing), [0.0, 0.0], (0.0, 1.0), [0.0, 0.0, 1.0, 0.0], Base.Pairs{Symbol, Union{}, Nothing, @NamedTuple{}}(), SciMLBase.StandardODEProblem())), :x, :u, :t, 1, 2, 1, 98)

When to Use TimeDependentBilinearIntegrator

✓ Time-varying Hamiltonians ✓ Systems with periodic forcing ✓ Carrier wave modulation (e.g., rotating frame transformations)

DerivativeIntegrator

Overview

Enforces derivative relationships between trajectory components:

\[\frac{d(\text{var})}{dt} = \text{deriv}\]

This is used for smoothness or when controls are derivatives of other variables.

Example: Smooth Controls

traj_smooth = NamedTrajectory(
    (
        x = randn(2, N),
        u = randn(2, N),
        du = zeros(2, N),   # control derivative
        Δt = fill(0.1, N),
    );
    timestep = :Δt,
    controls = :u,
    initial = (u = [0.0, 0.0],),
    final = (u = [0.0, 0.0],),
)
N = 50, (x = 1:2, u = 3:4, du = 5:6, → Δt = 7:7)

Enforce du/dt = du

deriv_integrator = DerivativeIntegrator(:u, :du, traj_smooth)
DerivativeIntegrator{DirectTrajOpt.Integrators.var"#2#3"}(DirectTrajOpt.Integrators.var"#2#3"(), :u, :du, 2, 7, 98)

Now you can penalize du to get smooth controls: obj = QuadraticRegularizer(:u, trajsmooth, 1e-2) obj += QuadraticRegularizer(:du, trajsmooth, 1e-1) # Smoothness penalty

Multiple Derivative Orders

traj_smooth2 = NamedTrajectory(
    (
        x = randn(2, N),
        u = randn(1, N),
        du = zeros(1, N),
        ddu = zeros(1, N),
        Δt = fill(0.1, N),
    );
    timestep = :Δt,
    controls = :u,
)
N = 50, (x = 1:2, u = 3:3, du = 4:4, ddu = 5:5, → Δt = 6:6)

Chain derivatives: d(u)/dt = du, d(du)/dt = ddu

deriv_u = DerivativeIntegrator(:u, :du, traj_smooth2)
deriv_du = DerivativeIntegrator(:du, :ddu, traj_smooth2)
DerivativeIntegrator{DirectTrajOpt.Integrators.var"#2#3"}(DirectTrajOpt.Integrators.var"#2#3"(), :du, :ddu, 1, 4, 49)

When to Use DerivativeIntegrator

✓ Enforce smooth, implementable controls ✓ Acceleration limits (when control is jerk) ✓ Tracking derivative information

Combining Multiple Integrators

You can use multiple integrators simultaneously:

traj_combined = NamedTrajectory(
    (x = randn(2, N), u = randn(2, N), du = zeros(2, N), Δt = fill(0.1, N));
    timestep = :Δt,
    controls = :u,
    initial = (x = [0.0, 0.0], u = [0.0, 0.0]),
    final = (u = [0.0, 0.0],),
)
N = 50, (x = 1:2, u = 3:4, du = 5:6, → Δt = 7:7)

Create problem with multiple integrators

G_drift_simple = [-0.1 1.0; -1.0 -0.1]
G_drives_simple = [[0.0 1.0; 1.0 0.0], [1.0 0.0; 0.0 1.0]]
G_simple = u -> G_drift_simple + sum(u .* G_drives_simple)

obj = QuadraticRegularizer(:u, traj_combined, 1e-2)
obj += QuadraticRegularizer(:du, traj_combined, 1e-1)

integrators_combined = [
    BilinearIntegrator(G_simple, :x, :u, traj_combined),
    DerivativeIntegrator(:u, :du, traj_combined),
]

prob = DirectTrajOptProblem(traj_combined, obj, integrators_combined)
DirectTrajOptProblem
   timesteps            = 50
   duration             = 4.899999999999999
   variable names       = (:x, :u, :du, :Δt)
   knot point dimension = 7

Integration Methods Comparison

IntegratorDynamics TypeAccuracyUse Case
BilinearIntegratorControl-linearExactQuantum, rotation
TimeDependentBilinearIntegratorTime-varying control-linearExactModulated systems
DerivativeIntegratorDerivative relationExactSmoothness

Custom Integrators

You can implement custom integrators by subtyping AbstractIntegrator and defining the constraint function. See the Advanced Topics section for details.

Interface Requirements

struct MyIntegrator <: AbstractIntegrator
    # ... fields ...
end

# Implement constraint evaluation
function (int::MyIntegrator)(δ, zₖ, zₖ₊₁, k)
    # Compute constraint: δ = xₖ₊₁ - Φ(xₖ, uₖ, Δtₖ)
    # where Φ is your integration scheme
end

Best Practices

Initialization

  • Start with good initial guesses for states and controls
  • For smooth control problems, initialize derivatives to zero
  • Use linear interpolation for states between boundary conditions

Performance

  • Matrix exponential (BilinearIntegrator) is efficient for small systems (n < 20)
  • For large systems, consider sparse representations
  • DerivativeIntegrator is cheap (just finite differences)

Numerical Stability

  • Keep time steps reasonable (not too large)
  • For stiff systems, smaller time steps help
  • BilinearIntegrator handles stiff systems well

Common Patterns

Pattern 1: Basic Bilinear Problem

G_basic = u -> [-0.1 1.0; -1.0 -0.1] + u[1] * [0.0 1.0; 1.0 0.0]
#14 (generic function with 1 method)

integrator = BilinearIntegrator(G_basic, :x, :u, traj)

Pattern 2: Smooth Control Problem

integrators = [ BilinearIntegrator(G, :x, :u, traj), DerivativeIntegrator(:u, :du, traj) ]

Pattern 3: Time-Dependent with Smoothness

integrators = [ TimeDependentBilinearIntegrator(G_td, :x, :u, :t, traj), DerivativeIntegrator(:u, :du, traj) ]

Summary

Key Takeaways:

  1. Integrators convert continuous dynamics to discrete constraints
  2. BilinearIntegrator is the workhorse for control-linear systems
  3. DerivativeIntegrator adds smoothness
  4. You can combine multiple integrators
  5. Good initialization helps convergence

Next Steps

  • Objectives: Learn how to define cost functions
  • Constraints: Add bounds and path constraints
  • Tutorials: See integrators in complete examples

This page was generated using Literate.jl.