Multi-phase flows

A multi-phase flow concatenates several single-phase flows end-to-end, with switching times between phases and optional jump functions applied at each switch. The integration is exact: each phase is solved independently, and the output of one phase becomes the initial condition of the next.

Concatenation operator `*`

Use the * operator to chain two flows at a switching time:

# flow1 on [t0, 1.0], then flow2 on [1.0, tf]
mpf = flow1 * (1.0, flow2)

MultiPhaseStateFlow
  phases: 2
  systems: VectorFieldSystem, VectorFieldSystem
  integrators: SciML, SciML
  switching_times: Real[1.0]
  jumps: Any[nothing]

The result is a MultiPhaseStateFlow. You can chain more phases:

mpf3 = flow1 * (0.5, flow2) * (1.0, flow3)
MultiPhase.n_phases(mpf3)

Switching times must be strictly increasing. If they are not, a PreconditionError is thrown.

Phases with different dynamics

Each phase may wrap a different vector field — the dynamics, system type, and integrator are all free per phase. The only requirement is that all phases share the same TimeDependence and VariableDependence traits (both autonomous or both non-autonomous, both fixed or both variable).

# flow1: f(x) = -x   |   flow_het: f(x) = -2x  — different functions, same TD/VD
mpf_het = flow1 * (1.0, flow_het)
MultiPhase.n_phases(mpf_het)

This enables the typical use cases of multi-phase optimal control:

Bang-bang trajectories: flow_plus * (t_switch, flow_minus) with distinct controlled dynamics per arc.
Multi-regime systems: phases with different physical models (e.g. free flight / contact phase).
Different solvers per phase: stiff phases can use an implicit integrator while smooth phases use an explicit one.

Mixing state flows with Hamiltonian flows is not allowed — attempting it raises a PreconditionError.

Jumps at switching times

A jump function $g$ is applied to the state at the switching time before starting the next phase. Pass it as the second element of the switching tuple:

jump = x -> 2.0 .* x   # double the state at the switch

mpf_jump = flow1 * (1.0, jump, flow2)

MultiPhaseStateFlow
  phases: 2
  systems: VectorFieldSystem, VectorFieldSystem
  integrators: SciML, SciML
  switching_times: Real[1.0]
  jumps: Any[Main.var"#2#3"()]

Pass nothing (or use the two-element tuple) for a continuous switch with no jump.

Calling a multi-phase flow

Multi-phase flows share the same call interface as single-phase flows:

x0 = [1.0, 0.0]

# Point integration: final state only
xf = mpf(0.0, x0, 2.0)

2-element Vector{Float64}:
 0.13533528334857758
 0.0

# Trajectory integration: full time history (phases merged)
sol = mpf((0.0, 2.0), x0)

VectorFieldTrajectory
  result: SciMLIntegrationResult

Inspecting a multi-phase flow

MultiPhase.n_phases(mpf)                 # number of phases
MultiPhase.get_flow(mpf, 1)              # flow for phase 1
MultiPhase.get_switching_time(mpf, 1)    # switching time after phase 1
MultiPhase.get_jump(mpf, 1)              # jump at that switching time (nothing here)

MultiPhase.get_flows(mpf)                # tuple of all phase flows
MultiPhase.get_switching_times(mpf)      # all switching times
MultiPhase.get_jumps(mpf)                # all jump functions

1-element Vector{Any}:
 nothing

Hamiltonian multi-phase flows

The same operators work for Hamiltonian flows:

hvf_f(x, p) = (p, -x)
hvf = Data.HamiltonianVectorField(hvf_f)
hflow1 = Flows.Flow(hvf; reltol=1e-10)
hflow2 = Flows.Flow(hvf; reltol=1e-10)

hmpf = hflow1 * (1.0, hflow2)
typeof(hmpf)

MultiPhaseHamiltonianFlow{Autonomous, Fixed, Tuple{HamiltonianFlow{Autonomous, Fixed, HamiltonianVectorFieldSystem{typeof(hvf_f), Autonomous, Fixed, OutOfPlace}, SciML{StrategyOptions{@NamedTuple{internalnorm::OptionValue{typeof(real_norm)}, alg::OptionValue{Tsit5{typeof(trivial_limiter!), typeof(trivial_limiter!), Serial}}, reltol::OptionValue{Float64}, save_everystep::OptionValue{Symbol}, abstol::OptionValue{Float64}, save_start::OptionValue{Symbol}, dense::OptionValue{Symbol}}}, Dict{Symbol, Any}, Dict{Symbol, Any}}}, HamiltonianFlow{Autonomous, Fixed, HamiltonianVectorFieldSystem{typeof(hvf_f), Autonomous, Fixed, OutOfPlace}, SciML{StrategyOptions{@NamedTuple{internalnorm::OptionValue{typeof(real_norm)}, alg::OptionValue{Tsit5{typeof(trivial_limiter!), typeof(trivial_limiter!), Serial}}, reltol::OptionValue{Float64}, save_everystep::OptionValue{Symbol}, abstol::OptionValue{Float64}, save_start::OptionValue{Symbol}, dense::OptionValue{Symbol}}}, Dict{Symbol, Any}, Dict{Symbol, Any}}}}, Vector{Real}, Vector{Any}} (alias for CTFlows.MultiPhase.MultiPhaseFlow{CTBase.Traits.Autonomous, CTBase.Traits.Fixed, CTBase.Traits.HamiltonianDynamics, Tuple{CTFlows.Flows.Flow{CTBase.Traits.Autonomous, CTBase.Traits.Fixed, CTBase.Traits.HamiltonianDynamics, CTFlows.Systems.HamiltonianVectorFieldSystem{typeof(Main.hvf_f), CTBase.Traits.Autonomous, CTBase.Traits.Fixed, CTBase.Traits.OutOfPlace}, CTFlows.Integrators.SciML{CTBase.Strategies.StrategyOptions{@NamedTuple{internalnorm::CTBase.Options.OptionValue{typeof(CTFlows.Common.real_norm)}, alg::CTBase.Options.OptionValue{Tsit5{typeof(OrdinaryDiffEqCore.trivial_limiter!), typeof(OrdinaryDiffEqCore.trivial_limiter!), FastBroadcast.Serial}}, reltol::CTBase.Options.OptionValue{Float64}, save_everystep::CTBase.Options.OptionValue{Symbol}, abstol::CTBase.Options.OptionValue{Float64}, save_start::CTBase.Options.OptionValue{Symbol}, dense::CTBase.Options.OptionValue{Symbol}}}, Dict{Symbol, Any}, Dict{Symbol, Any}}}, CTFlows.Flows.Flow{CTBase.Traits.Autonomous, CTBase.Traits.Fixed, CTBase.Traits.HamiltonianDynamics, CTFlows.Systems.HamiltonianVectorFieldSystem{typeof(Main.hvf_f), CTBase.Traits.Autonomous, CTBase.Traits.Fixed, CTBase.Traits.OutOfPlace}, CTFlows.Integrators.SciML{CTBase.Strategies.StrategyOptions{@NamedTuple{internalnorm::CTBase.Options.OptionValue{typeof(CTFlows.Common.real_norm)}, alg::CTBase.Options.OptionValue{Tsit5{typeof(OrdinaryDiffEqCore.trivial_limiter!), typeof(OrdinaryDiffEqCore.trivial_limiter!), FastBroadcast.Serial}}, reltol::CTBase.Options.OptionValue{Float64}, save_everystep::CTBase.Options.OptionValue{Symbol}, abstol::CTBase.Options.OptionValue{Float64}, save_start::CTBase.Options.OptionValue{Symbol}, dense::CTBase.Options.OptionValue{Symbol}}}, Dict{Symbol, Any}, Dict{Symbol, Any}}}}, Array{Real, 1}, Array{Any, 1}})

x0, p0 = [1.0, 0.0], [0.0, 1.0]
xf, pf = hmpf(0.0, x0, p0, 2.0)
(xf, pf)

(-0.41614683654602175, 0.909297426625848)

Design notes

Concatenation is an associative binary operation on flows: (f1 * (t, f2)) * (s, f3).
Phases may have different system and integrator types (S, I): each phase can wrap a distinct function and use its own solver. The flows field is a heterogeneous tuple, so the compiler specializes the integration loop per phase combination. The required uniformity across phases is TD, VD, and the dynamics family (state vs Hamiltonian), enforced at construction time.
The trajectory integration merges phase results via Integrators.merge, which concatenates the time grids and result vectors.