Skip to content

[WIP] Improvements to the dynamic optimization interface#4394

Draft
SebastianM-C wants to merge 14 commits into
SciML:masterfrom
SebastianM-C:smc/dyn_opt
Draft

[WIP] Improvements to the dynamic optimization interface#4394
SebastianM-C wants to merge 14 commits into
SciML:masterfrom
SebastianM-C:smc/dyn_opt

Conversation

@SebastianM-C

@SebastianM-C SebastianM-C commented Mar 16, 2026

Copy link
Copy Markdown
Member

Checklist

  • Appropriate tests were added
  • Any code changes were done in a way that does not break public API
  • All documentation related to code changes were updated
  • The new code follows the
    contributor guidelines, in particular the SciML Style Guide and
    COLPRAC.
  • Any new documentation only uses public API

Additional context

I will be splitting this in several PRs

@SebastianM-C SebastianM-C changed the title Improvements to the dynamic optimization interface [WIP] Improvements to the dynamic optimization interface Mar 16, 2026
@SebastianM-C
SebastianM-C marked this pull request as draft March 16, 2026 12:26
SebastianM-C and others added 9 commits July 2, 2026 17:23
FunctionWrapper parameters (e.g., interpolators) fail during fixpoint_sub
with fold=Val(true) because they can't accept JuMP types. Register them
as JuMP nonlinear operators instead, with automatic derivative detection
via Symbolics' derivative_rule system.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Allow users to provide constraint scaling, variable bounds, and
function-valued start trajectories at problem construction time.
All keys are normalized via default_toterm for consistent lookup.
Thread kwargs through all backend extensions (InfiniteOpt, CasADi, Pyomo).

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Replace inequality-constraint-based bounds with JuMP set_lower_bound/
set_upper_bound via dispatched set_variable_bounds!. Extract bounds
logic into extract_variable_bounds with support for user-provided
bounds that override metadata. This gives Ipopt better barrier method
structure for optimal control problems.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Allow users to provide function-valued initial trajectories for state
variables via the initial_trajectory kwarg. Functions are applied via
dispatched set_initial_trajectory! after variable creation.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Formulate dynamics as (∂x - tₛ*f(x)) / scale == 0 instead of scaling
each side independently. This prevents degenerate tₛ → 0 solutions and
improves Ipopt convergence for problems with multi-scale state variables.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Substitute observed variables before parameters so that observed
equations referencing parameters get fully resolved during fixpoint_sub.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Set print_level after set_optimizer to prevent the optimizer reset
from discarding the attribute.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Extend the optimal control bounds system to support observed variables:
- `extract_variable_bounds` now returns `observed_bounds` from variable
  metadata and user `bounds` dict (user bounds take priority)
- InfiniteOpt backend lifts observed bounds into auxiliary bounded
  decision variables with equality constraints, which interior-point
  solvers handle much more efficiently than nonlinear inequalities
- Auxiliary variables get proper start values computed from the initial
  state, and equality constraints are scaled via `scales` dict or
  variable metadata
- Change `add_initial_constraints!` from `fix` to `@constraint` to
  preserve variable bounds at the initial time point

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
SebastianM-C and others added 5 commits July 2, 2026 18:15
process_DynamicOptProblem relied on the ODEInputFunction default
inputs = unbound_inputs(sys), which is empty for compiled systems whose
inputs are bound through connector equations (e.g. flattened component
hierarchies). generate_control_function then substituted the controls at
their operating-point values, producing control-free dynamics: none of
the input variables appeared in the collocation stage equations, so the
NLP could not actuate the system and Ipopt converged to pseudo-optimal
constant-input trajectories in a handful of iterations.

Pass ctrls = inputs(sys) explicitly, matching how the rest of the
dynamic optimization interface identifies the controls.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
unbound_inputs classifies an input as bound if it appears in an equation
together with variables from another namespace. That connection-structure
heuristic is only meaningful before flattening: after mtkcompile,
connection and output equations are eliminated, so any effective input
co-occurs with variables from other components and unbound_inputs is
empty for every compiled hierarchical system.

Four codegen entry points used unbound_inputs(sys) as their default
input set and silently degraded on compiled hierarchical systems:
- ODEInputFunction generated dynamics with the inputs bound to their
  operating-point values (control-free dynamics in dynamic optimization)
- BVProblem computed n_controls = 0 and mis-shaped f_prototype
- calculate_control_jacobian returned an empty jacobian
- generate_control_function warned and bound the inputs when handed an
  already-scheduled system

Add default_codegen_inputs(sys), which returns the inputs declared to
mtkcompile for scheduled systems (stored in the inputs field) and falls
back to the unbound_inputs heuristic for unscheduled hierarchies, and
use it at all four sites.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
DiffEqDevTools v3 removed the construct* tableau functions (moved to
OrdinaryDiffEqExplicitTableaus / OrdinaryDiffEqImplicitTableaus when it
joined the OrdinaryDiffEq monorepo), and DiffEqDevTools v2's DiffEqBase
cap conflicts with ModelingToolkitBase's DiffEqBase >= 7.2 requirement,
making the whole Optimization test group unresolvable.

Depend on the two tableau packages instead and use the new names
(constructRK4() -> ExplicitTableaus.RK4(), etc.). Return types are
unchanged (DiffEqBase.ExplicitRKTableau), so no source changes.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
The dynamic optimization problem structs collected kwargs as positional
varargs, storing a Tuple in the kwargs field. Current SciMLBase remake
asserts values(prob.kwargs)::NamedTuple and reconstructs the problem
with T(f, u0, tspan, p, wrapped_model; kwargs...), and ConstructionBase
setproperties (used by remake internals) reconstructs with all six
fields positionally — both of which the varargs form breaks. This
surfaces as a TypeError from solve() -> remake(prob, p = new_p) for any
problem with tunable parameters, e.g. the parameter estimation tests.

Give each problem struct an exact six-argument positional constructor
(for setproperties) and a keyword-varargs constructor (for remake), and
pass kwargs through process_DynamicOptProblem as keywords.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
`InfiniteOptDynamicOptProblem` threw `Symbolic * GeneralVariableRef` from
`add_equational_constraints!` whenever a callable parameter was stored
unwrapped (e.g. a bare `DataInterpolations` object) rather than as a
`FunctionWrapper`. Two independent gaps caused it:

1. Operator registration was gated on `val isa FunctionWrapper`, so a bare
   callable was never turned into a JuMP nonlinear operator. Its substitution
   stayed `p => raw_callable`, producing `callable(symbolic_t)` which
   `unwrap_const` cannot collapse, and the subsequent `residual * tₛ` mixed a
   `BasicSymbolic` with a `GeneralVariableRef`.
2. The bare independent variable was never lowered onto the backend time
   variable, so even a registered `p(t)` kept a symbolic inner `t`.

Fix:
- Collect callable-parameter arities from the equations/observed/constraints/
  costs (`callable_parameter_arities`) and register every called parameter,
  reading its arity from the call site and falling back to a `FunctionWrapper`'s
  own argument-tuple type. Parameters that are never called are untouched.
- Add a `lowered_time_variable(model)` hook (default `nothing`,
  `model.model[:t]` for InfiniteOpt) and a `t => lowered_time_variable(model)`
  rule in `get_model_vars_substitution_rules!`. Because Symbolics substitution
  replaces the largest matching subtree first, `x(t) => U[i]` is applied before
  the walk reaches its inner `t`, so this only fires on genuinely-bare `t`.
- Factor the InfiniteOpt `register_operator!` body into a shared
  `_register_callable_operator!(m, dim, underlying, val, name)`; the bare and
  `FunctionWrapper` methods differ only in `underlying` (`val` vs `val.obj[]`).
  The analytic gradient is still resolved through Symbolics'
  `_derivative_rule_proxy` registry, so a registered derivative (e.g. from a
  DataInterpolations→Symbolics extension) is used with no direct dependency.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant