Constraints & duals

[ocots]

Constraints & duals — design note (roadmap v4, §3 then §2)

Date: July 8, 2026. Status: decisions validated (July 8, 2026) — A4 for Part A, B1 for Part B, the label-or-function keyword API, plus the follow-up remarks recorded in B.4 (keyword composition for NonFixed, and the :total/:partial semantics with frozen μ in :partial).

Scope: the two blocked items of roadmap-v4.md:

• Part A (§3) — how CTModels.Solutions.Solution should represent the presence/absence of dual variables (multipliers), the accessor contract, active-set bookkeeping, and display.
• Part B (§2) — how constrained flows are built in CTFlows: are constraints and multipliers separate objects or bundled, who computes the multiplier, and what the augmented state looks like.

Part A is treated first because its outcome (the storage question) gates Part B, per the v4 dependency graph. It turns out the gate opens easily: the recommended Part B design needs no new storage at all (multipliers are inputs to the flow, not integrated state), so Part A reduces to cleaning up what CTModels already has.

Sources consulted: current CTModels.Solutions (solution_types.jl, dual_model.jl, build_solution.jl, show.jl), the pre-refactor CTFlows in save/src/optimal_control_problem_utils.jl, the current CTFlows architecture (Systems, Flows, MultiPhase), CTBase.Data.ControlLaw, and the Handbook (types-traits-interfaces).

---

0. What exists today (facts, not opinions)

CTModels already stores duals. Solution has a dual::DualModelType field with DualModelType <: AbstractDualModel as a type parameter. The single concrete DualModel has 8 fields, each Union{Function,Nothing} / Union{ctVector,Nothing}:

Field	Type	Filled by
path_constraints_dual	t -> μ(t) or nothing	CTDirect (NLP duals)
boundary_constraints_dual	vector or nothing	CTDirect
state/control_constraints_{lb,ub}_dual	t -> ν(t) or nothing	CTDirect
variable_constraints_{lb,ub}_dual	vector or nothing	CTDirect

build_solution defaults every dual kwarg to nothing; CTFlows' _build_ocp_solution never passes any, so flow-built solutions carry an all-nothing DualModel. The public query is dual(sol, model, label) — per-label lookup across all constraint kinds, throwing IncorrectArgument for unknown labels. show branches at runtime on each field being nothing.

The old CTFlows treated multipliers as user inputs. Flow(ocp, u, g, μ) built H(t,x,p,v) = p·f + s p⁰ f⁰ + μ(t,x,p,v)'·g(t,x,u,v) (makeH overloads in save/). Separate wrapper types StateConstraint / MixedConstraint / Multiplier, paired positionally at the call site. The flow never computed a multiplier and never returned one; junction conditions (costate jumps) were the shooting method's business, supported by flow concatenation with jumps — which the current MultiPhase module still has (jumps field, jump functions at switching times).

CTModels already uses sentinel component types for absence: EmptyTimeGridModel <: AbstractTimeGridModel, EmptyDefinition <: AbstractDefinition. Presence/absence of a component is encoded in the component's type, not in a Solution-level flag.

CTModels can hand out a declared constraint by label: constraint(model, label) functors exist (Models/constraint_functors.jl), and ConstraintsModel stores labeled tuples (lb, f!, ub, labels) for path/boundary constraints and labeled boxes.

---

Part A — Solution with and without duals

A.1 The options

A1 — Status quo: one Solution, per-field nothing in DualModel

What we have. "No duals at all" is representable only as the all-nothing DualModel; nothing distinguishes "this solution kind never has duals" (flow-built) from "this solver produced no duals this time". Accessors and show branch at runtime per field.

• ✅ Zero migration. Already type-parametric per field (a DualModel{Nothing,…} is concretely typed).
• ❌ "Absent" is not a first-class state; every consumer re-derives it by testing 8 fields.
• ❌ dual(...) on an all-nothing model currently returns a functor wrapping nothing or fails downstream — the failure mode is late and obscure.

A2 — A trait on Solution (HasDuals / NoDuals type parameter)

Add a trait axis to Solution itself, dispatch accessors/printers on it.

• ❌ Contradicts the Handbook's criterion for trait parameters: "if an axis never propagates beyond one level of the hierarchy it is an implementation detail, not a structural axis." The duals axis lives entirely inside one component; state/control/costate accessors, plotting, export are all unaffected. Nothing in the Solution contract changes.
• ❌ The information would be duplicated: the DualModelType parameter already encodes it. Two sources of truth for one fact.
• ❌ Every constructor and downstream signature grows a parameter.

Rejected.

A3 — A separate type (ConstrainedSolution <: AbstractSolution)

Keep Solution dual-free; a subtype carries duals.

• ❌ Handbook, verbatim: "Two concrete types with identical fields that differ only by a flag should be one parametric type." The two types would share 9 of 10 fields.
• ❌ Combinatorial risk: the next orthogonal axis (multi-grid vs unified grid, GPU arrays, …) would demand 4 types, then 8.
• ❌ Every consumer (CTDirect, plots, serialization, export/import) must now accept both types or the abstract type — a large, ecosystem-wide migration for no layout benefit.
• ❌ It also draws the line in the wrong place: a solution from a constrained problem solved by an indirect method may legitimately have no stored duals. "Constrained" (problem property) ≠ "has duals" (solution content).

Rejected.

A4 — Sentinel component: EmptyDualModel <: AbstractDualModel ⭐ recommended

Encode absence where the data lives, exactly like EmptyTimeGridModel and EmptyDefinition already do:

struct EmptyDualModel <: AbstractDualModel end   # "this solution carries no duals"

has_duals(::EmptyDualModel) = false              # extractor on the component
has_duals(::DualModel)      = true
has_duals(sol::Solution)    = has_duals(sol.dual) # trait-style, compile-time

• Solution stays unique; the existing DualModelType parameter makes the empty/non-empty state dispatchable at compile time (Holy-trait pattern, no runtime branch).
• Consistent with the two sentinel precedents already in CTModels.Solutions — no new idiom to learn.
• build_solution constructs EmptyDualModel() when all dual kwargs are nothing (the current default), a DualModel otherwise. CTFlows' _build_ocp_solution then produces EmptyDualModel without changing a line on the CTFlows side.
• Partial presence (e.g. path duals but no box duals) stays inside DualModel with per-field nothing, as today. No per-field sentinel proliferation — one sentinel for "nothing at all", per-field nothing for "partially filled". This keeps the change small; per-field refinement remains possible later without breaking anything.

Migration cost: one new struct, one has_duals, a branch in build_solution, dispatch in show, and a guard in dual(...). CTDirect is untouched (it passes actual duals, so it keeps getting DualModel).

A.2 Accessor contract (roadmap question 2)

The v4 diagnosis is right: a global duals(sol) cannot be well-defined (several constraints, partial activity, ambiguous stacking). The supported query should remain the one that is already well-defined and already implemented:

dual(sol, model, label) — per-named-constraint lookup — is the sole public query.

Contract, precisely:

1. has_duals(sol) == false → dual(sol, model, label) throws a CTBase.Exceptions-style exception ("this solution carries no dual variables; it was produced by …"), regardless of the label. Rationale: the Handbook flags Union{T,Nothing} returns as an anti-pattern (conditional return types), and a silent nothing here propagates NaNs into plots and quadratures far from the cause. Loud and early beats quiet and late.
2. has_duals(sol) == true and the label exists → return what is returned today (functor of t, vector, or scalar per constraint kind).
3. Unknown label → IncorrectArgument, as today.
4. Callers that want to test first use the predicate has_duals(sol) (cheap, compile-time). A finer has_dual(sol, model, label) can be added if a real consumer needs it — not before (YAGNI).

On inactive constraints: CTModels stores what the producer gives it and does not interpret activity. A stored t -> μ(t) that is zero off the active set is the producer's convention (this is what NLP duals from CTDirect naturally look like). dual returns it as-is. The ambiguity the roadmap worries about ("which stored multiplier belongs to which constraint, and is it active?") is resolved by the label key for the first half and deliberately not resolved by CTModels for the second half — activity is not a Solution concept (next section).

A.3 Active-set bookkeeping (roadmap question 3)

Position: the active set is a solver-side concept and does not enter Solution or the flow.

• In indirect methods the arc structure (sequence of interior/boundary arcs) is an input fixed by the user before shooting; the flow just integrates each arc's dynamics. The flow never discovers activation — it is told.
• In direct methods activity is implicit in the NLP duals (complementarity); CTDirect already stores those functions without an activity flag, and no downstream consumer has asked for one.
• Storing an active set in Solution would therefore duplicate knowledge the producer has, in a form no current consumer reads. If a consumer appears (e.g. plot shading of boundary arcs), add a field or a metadata entry then, with that consumer's actual requirements on the table.

This answers the roadmap's question 3 with its own suspected answer: "maybe the flow only integrates the constrained dynamics and multipliers are never a flow output at all" — yes, and Part B makes that concrete.

A.4 Display (roadmap question 4)

Dispatch, don't branch:

_show_duals(io, ::EmptyDualModel, sol)  # prints nothing (or a one-line "no dual variables")
_show_duals(io, dm::DualModel, sol)     # current per-field behavior

The existing per-field nothing checks inside the DualModel method remain (partial presence). Graceful degradation becomes a property of the type, not of scattered ifs.

A.5 What CTModels needs to agree to (summary of Part A)

Change	Where	Size
EmptyDualModel <: AbstractDualModel + has_duals	Solutions/solution_types.jl, dual_model.jl	tiny
build_solution: all-nothing → EmptyDualModel()	Solutions/build_solution.jl	tiny
dual(...) throws on EmptyDualModel with context	Solutions/dual_model.jl	tiny
show dispatches on the dual component	Solutions/show.jl	small

No change to Solution's type or field list. No change to CTDirect. No change to CTFlows (it inherits EmptyDualModel through the build_solution defaults).

---

Part B — Constrained flows: constraints and multipliers

B.1 The mathematical setting (why this shapes the design)

For a path constraint g ≤ 0 treated by the indirect method, along a boundary arc the Hamiltonian gains a penalty term:

H(t, x, p, v) = p·f(t, x, u, v) + s p⁰ ℓ(t, x, u, v) + μ' g(t, x, u, v)

with, classically:

• the boundary control u determined by g = 0 (or d^q g/dt^q = 0 for a constraint of order q),
• the multiplier μ determined analytically from the stationarity along the arc,
• costate jumps at junction times τ: p(τ⁺) = p(τ⁻) − ν ∂g/∂x, with ν a shooting unknown,
• the arc structure fixed a priori; shooting solves for junction times and jump parameters.

Two structural consequences drive everything below:

1. If μ is a known function of (t, x, p, v), the constrained Hamiltonian is just another Hamiltonian. No state augmentation, no new integrator machinery, no DAE. The existing HamiltonianSystem / PseudoHamiltonianSystem pipeline applies verbatim to H + μ'g. The old makeH overloads are the proof: the entire "constrained" feature was one extra summand.
2. Junctions are already handled: costate jumps at switching times are exactly what MultiPhase concatenation supports today (jumps at switching_times). Tangency/entry conditions are shooting residuals, not flow code.

The hard problem the roadmap fears ("junction/activation handling … is the hard part") is hard only if the flow must discover activation itself. It shouldn't (A.3): in the indirect workflow, activation is structure, and structure is input.

B.2 The design axis: who provides μ, and is it bundled with g?

B1 — Separate carriers, user-supplied (g, μ), paired at the call site ⭐ recommended

The modernized old model. Two new nouns in CTBase.Data, next to ControlLaw (the symmetry the roadmap itself suggests):

# kind is an adjective → trait parameter, not three structs (Handbook rule)
struct PathConstraint{F, K<:ConstraintKind, TD, VD} <: AbstractPathConstraint{K, TD, VD}
    f::F                    # normalized call: g(t, x, u, v), like ControlLaw's (t,x,p,v)
end
# K ∈ StateKind | ControlKind | MixedKind  (drives which args are actually used)

struct Multiplier{F, TD, VD} <: AbstractMultiplier{TD, VD}
    f::F                    # normalized call: μ(t, x, p, v)
end

Same construction discipline as ControlLaw: parametric functor, trait extractors, typed constructors plus convenience constructors from plain functions (is_autonomous / is_variable kwargs). The pair is assembled at the Flow call:

fb = Flow(ocp, u_boundary; constraint = g, multiplier = μ)   # boundary arc
fi = Flow(ocp, u_interior)                                    # interior arc
f  = fi * (τ1, jump1, fb) * (τ2, jump2, fi)                   # existing MultiPhase syntax

• ✅ μ is an input exactly like the control law — same epistemic status (both derived analytically by the user from the arc structure), same wrapper philosophy.
• ✅ Each carrier is reusable alone: u_boundary without constraint on another arc, g fetched from the model (B.3), μ swapped during debugging.
• ✅ The flow stays a pure integrator; everything in A.3 holds.

B2 — Bundled arc descriptor (ConstrainedArc(u, g, μ))

One carrier for the triple, on the argument that on a boundary arc u and μ are both derived from g, so the triple is one mathematical object.

• ✅ Self-documenting; impossible to pair g with the wrong μ across arcs.
• ❌ Breaks the symmetry with ControlLaw, which is already a standalone citizen of CTBase.Data — the bundle would contain a ControlLaw, creating two ways to pass a control.
• ❌ Less composable, and the pairing is naturally (and visibly) expressed at the call site anyway.
• ❌ Control-free constrained problems (the old code supported Flow(ocp, g, μ) with a dummy control) would need a degenerate bundle.

Rejected as the primary API. Nothing prevents adding a convenience bundle later on top of separate carriers if call sites get noisy.

B3 — The flow computes μ automatically from the OCP's declared constraints

Activation detection, constraint order q, tangency conditions, boundary control from d^q g/dt^q = 0 — the constrained analogue of §1's implicit control, but substantially harder (event detection, structure discovery).

• ❌ Out of scope for v4; arguably solver territory forever (structure discovery is a homotopy/solver concern in practice, cf. how the team actually uses these flows in shooting).
• Note: §1's IFT machinery (∂²H̃/∂u², ∂²H̃/∂z∂u AD primitives) would be a prerequisite anyway — another reason §1 lands first, as the roadmap already orders it.

Deferred indefinitely. Design B1 does not preclude it (an automatic layer would produce (g, μ) pairs and feed the same API).

B.3 Linking to the OCP's declared constraints — by label

The roadmap asks for "a way for the flow to read the OCP's declared constraints (coordinate with CTModels)". The coordination already exists: constraints are labeled in ConstraintsModel, and constraint(model, label) returns a functor. So let constraint = accept either form:

Flow(ocp, u; constraint = :c1, multiplier = μ)          # by label: fetched from the model
Flow(ocp, u; constraint = g,   multiplier = μ)          # by function: ad-hoc, not in the model

The label form buys three things:

1. No re-declaration — dimension and kind come from the model; the carrier's K trait is inferred.
2. Validation — unknown label fails at flow construction with the model's own error machinery.
3. The door §3 keeps open — the flow now knows which model constraint it integrates. If one day a MultiPhase trajectory call should assemble path_constraints_dual (zeros off-arc, the user's μ on-arc), the label is exactly the key dual(sol, model, label) resolves. This is explicitly not implemented now (A.3: flow-built solutions carry EmptyDualModel), but the design costs nothing to keep it possible.

The function form is kept because the old code had it and it is genuinely useful (penalization experiments, constraints not worth declaring in the model).

B.4 API details and validations

• Keywords, not positional. constraint = …, multiplier = … rather than the old Flow(ocp, u, g, μ). Consistent with v4's keyword-first convention (variable, variable_costate, and §1's planned control = u0). Positional (g, μ) reads ambiguously once §1 adds more call shapes.
• Both or neither. constraint without multiplier (or vice-versa) → construction-time error, trait-style validation like variable for NonFixed.
• Control-free constrained flows stay supported: Flow(ocp; constraint = :c1, multiplier = μ) (pure state constraint, the old Flow(ocp, g, μ) case).
• Several simultaneous constraints: accept tuples — constraint = (:c1, :c2), multiplier = (μ1, μ2) — summing Σ μᵢ' gᵢ. Cheap to support in the H builder; matched lengths validated at construction. (Single pair remains the 95% case in indirect practice, since distinct constraints are rarely active on the same arc.)
• Both Hamiltonian paths get the term — with distinct AD semantics. The hamiltonian_type = :total | :partial option stays available on constrained flows:
  • :total — full composition: H(t,x,p,v) = p·f + sp0·ℓ + μ(t,x,p,v)'·g(t,x,u(t,x,p,v),v). AD differentiates through everything, including the (x,p)-dependence of both u and μ.
  • :partial — μ joins u as a frozen argument of the pseudo-Hamiltonian: H̃c(t,x,p,u,μ,v) = p·f + sp0·ℓ + μ'·g(t,x,u,v). Gradients are taken at fixed (u, μ): we differentiate through g (its explicit dependence on x, v, with u held fixed) but not through μ. This is mathematically consistent on a boundary arc: the neglected term (∂μ/∂(x,p))'·g vanishes there because g ≡ 0 along the arc — exactly the same mechanism that justifies freezing u (the neglected (∂u/∂(x,p))'·∂H̃/∂u vanishes when ∂H̃/∂u = 0). Implementation: the constrained pseudo-Hamiltonian functor gains a frozen μ slot alongside the frozen u slot; the PseudoHamiltonianSystem AD pipeline extends with one more pass-through argument.
• Keyword composition for NonFixed problems. constraint / multiplier are construction-time keywords — they define the integrated system: Flow(ocp, u; constraint = …, multiplier = …). variable / variable_costate remain call-time keywords, unchanged: f(t0, x0, p0, tf; variable = v, variable_costate = true). All four coexist on a NonFixed constrained problem; g and μ may depend on v (their VD trait). Under variable_costate, ṗᵥ = −∂H/∂v picks up the μ'g term with the same :total/:partial semantics: in :partial, ∂/∂v is taken at frozen μ, and the neglected (∂μ/∂v)'·g again vanishes on the arc. No interaction with the (nonexistent) dual storage.
• Sign/orientation convention (+μ'g, g ≤ 0, :min) must be documented once in the guide with the sp0 convention; the old code silently assumed it.

B.5 The §3→§2 gate, answered

The roadmap blocked §2 on §3 because "the storage question dictates the augmented-state layout". With B1:

Multipliers are inputs, not integrated state. The state stays (x, p) (+pv under variable_costate). There is no augmented-state slot for duals, hence nothing for the storage question to dictate.

The gate is lifted. §2 reduces to: two carrier types in CTBase.Data, one extra summand in the two Hamiltonian builders, keyword plumbing and validation, tests against an analytic reference (e.g. the constrained double integrator / Bryson–Denham problem, whose boundary-arc μ, ν and junction times are known in closed form).

---

C. Decisions summary

Decisions 1–9 were validated on July 8, 2026; decision 10 records the validated follow-up remark on the :total/:partial semantics.

#	Question	Decision	Alternative rejected because
1	Presence/absence of duals in Solution	Sentinel component EmptyDualModel + has_duals extractor; Solution type unchanged	Trait on Solution duplicates what DualModelType already encodes; subtype duplicates 9/10 fields and forks the ecosystem
2	Accessor contract	dual(sol, model, label) only; throws (with context) when has_duals is false; predicate has_duals(sol)	Returning nothing = conditional return type + silent NaN propagation
3	Active set	Solver-side concept; not in Solution, not in the flow	No consumer; duplicates producer knowledge; flow is told the structure, it doesn't discover it
4	Display	Dispatch on the dual component type	Runtime per-field branching for the "none at all" case
5	Constraints vs multipliers	Separate carriers (PathConstraint with kind-trait, Multiplier), in CTBase.Data next to ControlLaw; paired at the Flow call by keywords	Bundle hides the ControlLaw; auto-computed μ is a (much later, maybe never) layer on top
6	Who computes μ	The user (indirect-method practice), as a (t,x,p,v) function	Automatic computation = structure discovery = solver territory
7	Link to model constraints	constraint = :label fetches from ConstraintsModel (also accepts a plain function)	Re-declaration loses validation and the future per-label dual key
8	Flow output	Constrained flows return solutions with EmptyDualModel (v4 position: no duals from the flow); label-linked design keeps later assembly possible	—
9	Augmented state	None — (x, p) unchanged; junctions via existing MultiPhase jumps	DAE/augmentation only needed if μ were unknown, which it isn't under 6
10	:total vs :partial with constraints	Both supported. :partial freezes μ (like u) and differentiates through g; sound on-arc since the neglected (∂μ/∂z)'·g term vanishes where g ≡ 0	Differentiating through μ in :partial would contradict the frozen-argument semantics of the pseudo-Hamiltonian path and buys nothing on the arc

Impact per package

• CTModels — the four small changes of A.5. This note is the coordination artifact; agree before touching code.
• CTBase — Data.PathConstraint (kind as trait parameter), Data.Multiplier, traits + typed/convenience constructors, following the ControlLaw template.
• CTFlows — + μ'g in the :total and :partial H builders; constraint/multiplier keywords with both-or-neither validation and label resolution; tests (Bryson–Denham class reference).
• CTDirect — untouched (keeps filling DualModel through build_solution, whose non-default path is unchanged).

Suggested sequencing

1. CTModels PR (A.5) — small, self-contained, unblocks everything. Reviewable against this note.
2. CTBase PR — the two carriers.
3. CTFlows PR(s) — constrained :total path first (simplest: one summand in the composed Hamiltonian), then :partial (pseudo-Hamiltonian gains the frozen μ slot alongside the frozen u), then the control-free constrained case; analytic-reference tests per step, including NonFixed cases exercising variable / variable_costate alongside constraint / multiplier.
4. Guide chapter: constrained flows + the sign conventions + a full Bryson–Denham shooting example (arcs, jumps, junction residuals).

Open questions (deliberately left out of scope)

• Per-field sentinels inside DualModel (finer has_* granularity) — only if a consumer materializes.
• Box constraints as flow-side carriers (x − ub residuals) — expressible today via the function form; a dedicated label form for boxes can wait.
• Whether §1's implicit control and the boundary-arc control u_boundary from d^q g/dt^q = 0 should share machinery — revisit when §1's IFT primitives exist.
• The option-vs-strategy cross-cutting decision does not arise here: no new heavy dependency, no RHS branching beyond one summand → plain keywords, consistent with the roadmap's rule of thumb.

[ocots]

Plan — Constrained flows & duals (roadmap v4 §3→§2)

Context

The design note .reports/constraints-and-duals-design.md was validated on July 8, 2026:

• A4 — EmptyDualModel <: AbstractDualModel sentinel in CTModels; has_duals; dual(sol, model, label) throws with context when duals are absent; show handles the empty case gracefully.
• B1 — separate carriers PathConstraint{F,K,TD,VD} / Multiplier{F,TD,VD} in CTBase.Data (ControlLaw template); user-supplied μ; paired at the Flow call via construction-time keywords constraint = :label | g and multiplier = μ.
• :total differentiates through everything; :partial freezes μ alongside u but differentiates through g (sound on-arc: (∂μ/∂z)'·g vanishes where g ≡ 0).
• constraint/multiplier (construction-time) compose with variable/variable_costate (call-time) for NonFixed.
• No state augmentation; flow-built Solutions carry EmptyDualModel; active set stays solver-side.

User process decisions: 5 PRs (1 CTModels, 1 CTBase, 3 CTFlows), beta release after each merge so the chain pins registered betas instead of local dev wiring.

Process rules (every step)

• Branch from up-to-date main of the target repo, feat/... naming (repo convention, e.g. feat/hvf-getter-functors).
• Ask the user before any commit; never push or open a PR without explicit approval.
• One PR per step; PR body references the design note; plans-before-code and docstrings-last per CLAUDE.md.
• Local cross-repo testing before betas exist: Pkg.develop sibling repos + [compat] bump (Manifest gitignored, not committed). After each beta is registered, switch the downstream pin to the beta.
• Tests: full suite via ct-dev-mcp get_test_command (never subsets).
• Current versions: CTBase 0.26.3-beta, CTModels 0.14.1-beta (CTModels pins CTBase = "0.26").

---

PR 1 — CTModels: EmptyDualModel + has_duals + accessor/show contract

Branch feat/empty-dual-model from CTModels.jl main.

1. Sentinel + trait in src/Solutions/solution_types.jl (next to EmptyTimeGridModel, line ~175):

• struct EmptyDualModel <: AbstractDualModel end
• has_duals(::EmptyDualModel) = false, has_duals(::AbstractDualModel) = true (mirror the two-method is_empty pattern, lines 195/216)
• Solution-level has_duals(sol::Solution) = has_duals(dual_model(sol)) (dual_model at build_solution.jl:1453 already dispatches on the AbstractDualModel parameter — works unchanged).

2. Construction branch in src/Solutions/build_solution.jl: the single DualModel(...) construction site is lines 378–387 (legacy single-grid method forwards here; serialization reconstruction too — one site to patch). Branch: if the 8 assembled values (fpcd, boundary_constraints_dual, fscbd, fscud, fccbd, fccud, variable_constraints_lb_dual, variable_constraints_ub_dual) are all === nothing → EmptyDualModel(), else DualModel(...) as today.

3. Accessors on EmptyDualModel: the 8 DualModel-typed accessors in src/Solutions/dual_model.jl:135–344 would MethodError on the sentinel. Add one EmptyDualModel method per accessor returning nothing (keeps the Solution-level delegating accessors at build_solution.jl:1476–1548 working, and show's dim helpers too).

4. dual() guard: in dual(sol, model, label) (dual_model.jl:37–120), first line: has_duals(sol) || throw(CTBase.Exceptions.PreconditionError(...)) — PreconditionError (CTBase/src/Exceptions/types.jl:126) is the semantic fit ("arguments valid, call forbidden by state"); message says the solution carries no duals and names the producer context.

5. Show: src/Solutions/show.jl (single method, lines 18–159) branches on integer dims and reaches duals via the accessors — with step 3's nothing-returning methods, the existing nothing guards inside the dual blocks (lines 120–142, 145–158) must be checked: guard the two dual blocks with has_duals(sol) && so an EmptyDualModel solution prints cleanly without them.

6. Exports in src/Solutions/Solutions.jl export list (lines 84–108): add EmptyDualModel, has_duals next to DualModel, AbstractDualModel (line 85).

7. Tests (auto-discovered test/suite/solutions/test_*.jl, top-level test_<name>() shim convention):

• New test/suite/solutions/test_empty_dual_model.jl: sentinel construction, has_duals both ways, all 8 accessors return nothing, dual() throws PreconditionError, show of a dual-free Solution.
• Touch existing: test_dual_model.jl (all-nothing DualModel at line 119 — decide: still legal to construct directly, but build_solution now yields the sentinel; add a build-path assertion), test_solution_types.jl:74–76 (hand-built Solution with all-nothing DualModel — may switch to EmptyDualModel()).
• CTFlows side effect (no code change): _build_ocp_solution passes no dual kwargs → now gets EmptyDualModel automatically.

Release: after merge, bump CTModels version (next beta, e.g. 0.14.2-beta or per user's call) and register.

---

PR 2 — CTBase: Data.PathConstraint + Data.Multiplier + constraint-kind trait

Branch feat/path-constraint-carriers from CTBase main. No new deps.

1. New trait family src/Traits/constraint_kind.jl (pattern of feedback.jl): AbstractConstraintKind <: AbstractTrait, empty structs StateConstraintKind, ControlConstraintKind, MixedConstraintKind, accessor stub function constraint_kind end, predicates. Type-parameter-only contract (like AbstractFeedback — no has_ guard). Wire into Traits.jl includes (lines 78–88) and exports (94–116).

2. PathConstraint{F,K,TD,VD} — clone the ControlLaw template:

• src/Data/abstract_path_constraint.jl (from abstract_control_law.jl): AbstractPathConstraint{K,TD,VD}, trait accessors (Traits.constraint_kind, time_dependence, variable_dependence, predicates).
• src/Data/path_constraint.jl (from control_law.jl): struct + keyword ctor (is_autonomous/is_variable, inline flag→trait mapping as in control_law.jl:79–80) + typed ctor + per-kind convenience ctors (StateConstraint(g; ...), ControlConstraint(g; ...), MixedConstraint(g; ...) — analogous to OpenLoop/ClosedLoop/DynClosedLoop) + natural call operators per kind/traits + uniform call (t, x, u, v) + two-method Base.show.
• Natural arities: State g(x)/g(t,x)/g(x,v)/g(t,x,v); Control g(u)…; Mixed g(x,u)/g(t,x,u)/g(x,u,v)/g(t,x,u,v).

3. Multiplier{F,TD,VD} — no kind parameter; template is hamiltonian.jl (cleaner 3-param match, uniform call (t,x,p,v)): src/Data/abstract_multiplier.jl + src/Data/multiplier.jl.

4. Display helpers in src/Data/helpers.jl: _kind_label, _natural_sig_pc/_natural_args_pc/_uniform_sig_pc, _natural_sig_mult/_uniform_sig_mult (copy the ControlLaw block, lines 286–359).

5. Wiring: Data.jl includes after line 27 (abstract-then-concrete pairing) and export block (47–51 pattern): AbstractPathConstraint, PathConstraint, StateConstraint, ControlConstraint, MixedConstraint, AbstractMultiplier, Multiplier.

6. Tests (auto-discovered test/suite/data/test_*.jl, module-wrapped + trailing shim convention, cf. test_control_law.jl):

• test_path_constraint.jl, test_multiplier.jl (construction all trait combos, natural/uniform calls, extractors, show, subtyping).
• test_abstract_path_constraint.jl / test_abstract_multiplier.jl with top-level Fake* structs (repo convention: fakes at module top level).
• Update test_data_module.jl hardcoded symbol lists (EXPORTED_ABSTRACT_TYPES 36–42, EXPORTED_CONCRETE_TYPES 44–53, PRIVATE_SYMBOLS 55–77) or it fails.

Release: bump CTBase (next beta after 0.26.3-beta) and register. If it stays in the 0.26 line, CTModels' CTBase = "0.26" pin is untouched.

---

PR 3 — CTFlows: constrained :total path

Branch feat/constrained-flow-total from CTFlows.jl main. Pin the two betas in [compat] (dev wiring meanwhile). Wiring anchors: see "CTFlows wiring facts" below.

1. Keywords: declare :constraint and :multiplier as action OptionDefinitions in _flow_action_defs() (flow_routing.jl:107–117, default nothing) — accepted only where these action_defs are passed, auto-rejected elsewhere. Unwrap in _flow_from_ocp_control (flow.jl:456–476) with _unwrap_option.

2. Resolution helper (new, in flow.jl or optimal_control_flow.jl): constraint accepts a Symbol label → CTModels.Models.constraint(ocp, label), accept only kind === :path (reject :boundary/box kinds with IncorrectArgument), take the out-of-place g(t,x,u,v) functor; or a Data.PathConstraint (used as-is, uniform call); or a plain function (wrapped into Data.MixedConstraint with the OCP's TD/VD). multiplier accepts Data.Multiplier or plain function (wrapped, uniform call (t,x,p,v)).

3. Validation at construction: both-or-neither (IncorrectArgument if only one given); constrained + hamiltonian_type = :partial → IncorrectArgument "not yet supported" until PR 4 (keeps this PR's scope honest); single constraint only (tuples deferred to PR 5).

4. Constrained pseudo-Hamiltonian: new callable struct (e.g. ConstrainedPseudoHamiltonian{H,G,M} holding base h̃, g, μ) with call (t,x,p,u,v) = h̃(t,x,p,u,v) + dot(μ(t,x,p,v), g(t,x,u,v)), built next to _ocp_pseudo_hamiltonian (optimal_control_flow.jl:227–240). Selected by dispatch in _flow_from_ocp_control — the unconstrained hot path (OCPPseudoHamiltonianFunction, _ocp_pseudo_H at 182–195) is untouched. Fed into Data.ComposedHamiltonian exactly as today (flow.jl:290–294): correct for :total since the composed-H AD differentiates through everything inside h̃ — including μ(t,x,p,v) and the law, which is the :total semantics. Works for NonFixed too (∂H/∂v picks up both ∂μ/∂v·g and μ'∂g/∂v).

5. OptimalControlFlow carries g/μ (or nothing) for show; light show extension (abstract_flow.jl:466–486 and/or a dedicated method). Solution building unchanged (EmptyDualModel arrives via CTModels defaults).

6. Tests in test/suite/flows/test_ocp_control.jl (+ test_flow_routing.jl for the new options): fixtures at module top level per convention; constrained double integrator / Bryson–Denham-class boundary arc with closed-form μ, ν, junction times; keyword validation matrix (label vs function vs Data.PathConstraint, both-or-neither, non-:path label rejected, :partial guarded); MultiPhase concatenation with costate jump p(τ⁺) = p(τ⁻) − ν ∂g/∂x across the junction; NonFixed smoke test (variable at call time composes untouched, calling.jl:221–228).

PR 4 — CTFlows: constrained :partial path

Branch feat/constrained-flow-partial.

The correctness constraint (see wiring facts): reusing PR 3's baked-μ ConstrainedPseudoHamiltonian under the :partial gradients would be wrong — the closures in pseudo_hamiltonian_gradient (h̃_x(x_) = h̃(t, x_, p, u, v), CTBaseDifferentiationInterface.jl:141–142) differentiate through the internal μ(t, x_, p, v) evaluation. :partial semantics require μ frozen, g differentiated through.

1. Value-frozen wrapper (new, CTFlows-local, no CTBase change): callable struct FrozenMultiplierPseudoHamiltonian{H,G,T} holding base h̃, g, and the evaluated value μ_; call (t,x,p,u,v) = h̃(t,x,p,u,v) + dot(μ_, g(t,x,u,v)) — g is inside the AD closures (differentiated through), μ_ is a captured constant.

2. Constrained PseudoHamiltonianSystem variant (pseudo_hamiltonian_system.jl:36–46, builder building.jl:135–141) carrying (h̃, law, g, μ, backend); RHS functors extended by dispatch (unconstrained structs untouched): in the call (pseudo_hamiltonian_rhs_functors.jl:78–86, 114–121, 159–169), evaluate u_ = f.law(t, cx, cp, v) and μ_ = f.μ(t, cx, cp, v), wrap h̃c = FrozenMultiplierPseudoHamiltonian(f.h̃, f.g, μ_) (cheap immutable, per-call), then call the existing pseudo_hamiltonian_gradient(f.backend, h̃c, t, cx, cp, u_, v) / pseudo_variable_gradient unchanged — the augmented-variable path automatically gets ṗᵥ ∋ −μ_'∂g/∂v at frozen μ.

Fallback if the wrapper hurts AD/perf: CTBase's generic differentiate(backend, f, ::Val{slot}, active, consts...) (abstract_ad_backend.jl:368–379) already supports freeze-all-but-one-slot — a μ-slot primitive is expressible there, at the cost of a CTBase contract + extension change and a compat bump. Decide at implementation; wrapper first.

3. Lift the PR 3 guard: constrained + :partial now routes to the new system; _build_pseudo_flow(::Val{:partial}, …) (flow.jl:302–305) gains the constrained dispatch.

4. Tests: same analytic references as PR 3 under hamiltonian_type = :partial; cross-check :total ≡ :partial trajectories on-arc (where g ≡ 0, both at 1e-12 tolerances via _opts()); a deliberately state-dependent μ(t,x,p,v) off-stationarity case where :total and :partial differ, confirming μ is genuinely frozen; variable-costate quadrature check (_trapz) including the μ'g term.

PR 5 — CTFlows: control-free constrained case + NonFixed matrix + guide

Branch feat/constrained-flow-controlfree-docs.

• Flow(ocp; constraint=…, multiplier=…) (pure state constraint, no control law) routed on the control-free path.
• Multiple simultaneous constraints (tuples of constraint/multiplier, matched lengths — deferred from PR 3).
• NonFixed test matrix: variable / variable_costate × constraint / multiplier (ṗᵥ = −∂H/∂v picks up μ'g; :partial at frozen μ).
• Docs: guide chapter (constrained flows, sign conventions +μ'g, g ≤ 0, :min/sp0, full shooting example with arcs/jumps/junction residuals); docstrings last, api_reference sync.

---

CTFlows wiring facts (from exploration, July 8, 2026)

Call chain — src/Flows/flow.jl:

• Flow(ocp, law) (419–421) → _flow_from_ocp_control(Traits.feedback(law), …); DynClosedLoop branch 456–476: validates WithControl (462–469), routed = _route_flow_options(kwargs; action_defs=_flow_action_defs()) (470), ht = _unwrap_option(get(routed.action, :hamiltonian_type, nothing), :total) (472), h̃ = _ocp_pseudo_hamiltonian(ocp) (473), inner = _build_pseudo_flow(Val(ht), h̃, law, components) (474), OptimalControlFlow(inner, ocp, law) (475).
• Raw-function overload Flow(ocp, u::Function) 441–446 wraps into Data.DynClosedLoop with the OCP's TD/VD. Positional guard 521–531 (PreconditionError).
• _build_pseudo_flow(::Val{:total}, …) 290–294: H = Data.ComposedHamiltonian(h̃, law); Systems.build_system(H, backend) → HamiltonianSystem. ::Val{:partial} 302–305: build_system(h̃, law, backend) → PseudoHamiltonianSystem. Fallback 307–316 throws IncorrectArgument.
• Control-free Flow(ocp) 151–153 → _flow_from_ocp(ControlFree, …) 165–172; its Hamiltonian core _ocp_H is optimal_control_flow.jl:65–82 (PR 5 insertion point).

Options machinery — src/Flows/flow_routing.jl: _flow_action_defs() 107–117 holds the sole OptionDefinition (:hamiltonian_type, default :total — the default lives here, not in defaults.jl). New :constraint/:multiplier action options are declared here; they are then only accepted on the constructors that pass these action_defs (rejected elsewhere by route_all_options, which raises IncorrectArgument on unknown kwargs — _route_flow_options 84–95).

Pseudo-Hamiltonian & the + μ'g site — src/Flows/optimal_control_flow.jl: OCPPseudoHamiltonianFunction struct 165–170; core _ocp_pseudo_H 182–195 with the summation at 192–193 (val = sum(p .* r); … val += h.sp0 * h.lagrange(t, x, u, vv)); builder _ocp_pseudo_hamiltonian 227–240.

AD freeze mechanics (the critical fact) — ComposedHamiltonian core (CTBase/src/Data/composed_hamiltonian.jl:144–147) is u = law(t,x,p,v); h̃(t,x,p,u,v): everything inside h̃ is differentiated through. The :partial gradients (CTBase/ext/CTBaseDifferentiationInterface.jl:131–146) build closures h̃_x(x_) = h̃(t, x_, p, u, v) that capture u,v as constants — but an internal μ(t, x_, p, v) evaluation inside h̃ would be differentiated through. Consequence:

• Baking μ(t,x,p,v)'·g(t,x,u,v) into the pseudo-H functor is correct for :total only (differentiate through everything — desired).
• For :partial, μ must be evaluated in the RHS functor before the gradient call (like u_) and passed as a frozen value. RHS functors: src/Systems/pseudo_hamiltonian_rhs_functors.jl — PseudoHamIpRHS 68–86 (u_ = f.law(t, cx, cp, v) then pseudo_hamiltonian_gradient(f.backend, f.h̃, t, cx, cp, u_, v)), PseudoHamOoPRHS 104–121, PseudoHamIpAugRHS 148–169 (pseudo_variable_gradient). System struct: pseudo_hamiltonian_system.jl:36–46; builder Systems/building.jl:135–141.

Label resolution — CTModels.Models.constraint(model, label)::Tuple (CTModels src/Models/model.jl:1459–1549) returns (kind::Symbol, f, lb, ub) (not type-stable; fine at construction time). The :path branch yields an out-of-place g(t, x, u, v) (scalar if dim 1, else vector) — same arity as the dynamics, slots straight into _ocp_pseudo_H. Other kinds (:boundary, box :state/:control/:variable) return different functor shapes → accept only :path labels, reject the rest with IncorrectArgument.

Call-time validation pattern (for the variable/variable_costate composition) — src/Flows/calling.jl: _invoke_flow 221–228 dispatches on (variable_dependence(flow), typeof(variable)), PreconditionError for NonFixed-missing / Fixed-provided (303–409). Unchanged by this work — constrained flows compose for free.

Tests & show — primary test home test/suite/flows/test_ocp_control.jl (fixtures as top-level builders _build_lqr, _build_double_integrator, _build_variable_control…; _opts() = (; alg=Tsit5(), reltol=1e-12, abstol=1e-12); _trapz quadrature for variable-costate checks); option machinery in test_flow_routing.jl; system tests test/suite/systems/test_pseudo_hamiltonian_system.jl; auto-discovery suite/*/test_*. Flow show is generic (abstract_flow.jl:466–486, 501–511); PseudoHamiltonianSystem show at pseudo_hamiltonian_system.jl:191–213 — extend lightly to surface g/μ.

Deps — none new: CTBase.Differentiation primitives + CTModels accessors suffice. Project.toml: 0.10.0-beta, compat CTBase = "0.26.3", CTModels = "0.14".

---

Verification

• Per-repo: full test suite green (ct-dev-mcp get_test_command).
• End-to-end after PR 3: run the Bryson–Denham shooting script — interior/boundary/interior structure, costate jump at junctions via MultiPhase, compare (x, p, μ, ν, junction times) to closed-form; after PR 4 the same under :partial, :total ≡ :partial on-arc.
• After each beta registration: downstream repo resolves from the registry (no Pkg.develop) and its CI passes.