BVProblem Interface Summary

[vyudu]

@parameters α=1.5 β=1.0 γ=3.0 δ=1.0 tf=1.0
@variables x(..) y(..)

eqs = [D(x(t)) ~ α * x(t) - β * x(t) * y(t),
           D(y(t)) ~ -γ * y(t) + δ * x(t) * y(t)]

u0map = []
tspan = (0.0, 1.0)
guesses = [x(t) => 4.0, y(t) => 2.]
constraints = [x(.6) ~ 3.5, x(.3) ~ 7., y(t) < 3]
cost = [norm(y(t) - x(t)) + tf, x(0.3)]
coalesce(u) = u[1] + 0.1*u[2]

@mtkbuild lksys = ODESystem(eqs, t; constraints, coalesce, cost)
BVProblem(lksys, [], (0.0, tf))

Constraints

Specified as a vector of Equation, gets codegened to the bc function. Constraints that have a numeric argument are interpreted to hold at that timestep, and constraints that have t as the argument are interpreted to hold over the whole interval.

Needs: minsol, maxsol on the BoundaryValueDiffEq side for the whole-interval constraints, and codegen should handle inequality constraints.

Cost Function

For optimal control, specified as a vector of expressions, gets codegened to a cost function. norm(x(t) - y(t)) is argument is interpreted as the integral norm, x(0.3) is just the value.

Needs: integralnorm, and codegen to cost function needs to be implemented in MTK.

Coalesce

A function that takes the vector returned by a cost function and returns a scalar (or a vector in the case of multi-objective). Passed directly to the ODESystem constructor.

To discuss?

Multiphase problems and minimum final time problems

[vyudu]

@ErikQQY For your reference, this is the planned interface for MTK's codegen to BVProblem, does this all seem reasonable?

[vyudu]

@ErikQQY What's the BVP interface for optimization-type problems going to look like? Are the inequality constraints passed as a separate argument, like BVProblem(f, bc, p, t; ineqs) or something?

Also what should the function for inequality constraints look like? If you have x(0.5) < 1, y(0.5) > 2, is it like

function ineqs(resid, sol, p, t)
    resid[1] = x(0.5) - 1
    resid[2] = 2 - y(0.5) 
end

and it ensures resid is less than zero?

[ErikQQY]

While there are no interfaces for optimization-type problems in BoundaryValueDiffEq.jl, that's something on the TODO list, I think we better follow the constraints function in Optimization.jl, e.g. in https://docs.sciml.ai/Optimization/stable/tutorials/constraints/, we can do like

cons(res, sol, p, t) = (res .= [x(0.5) - 1, 2 - y(0.5)])
fun = BVPFunction(f!, cons = cons)
prob = BVProblem(fun, bc, u0, p, lcons = [-Inf, 0.0], ucons = [0.0, Inf])