So I have a system I’m trying to control which has a valve (can be open or closed, not in between) .(should be extended with multiple valves, but for now just simulating with 1.)
I force it to be 0 or 1 through constraints (as below). The issue is it just follows the initial conditions and seems to get “stuck” in that condition, i.e. the optimizer doesn’t change it.
I’ve tried to put very high cost terms related to the input (starting at 1 and then trying to penalize being 1), but even then it seems to just be following initial conditions. The solver also errors when initializing to 0.5.
any known approaches how to handle this?
# Constraint to ensure valve is binary: 0 or 1
ocp.model.con_h_expr = model.u[valve_idx] * (model.u[valve_idx] - 1)
ocp.constraints.lh = np.array([0])
ocp.constraints.uh = np.array([0])
I was trying to use the linear cost functions, hoping to keep it light. Not sure if it would help going more complex.
You are modelling a binary variable with a nonconvex expression which creates a complementarity constraint. This type of constraints violates standard constraint qualification, so the solver might not converge to a local optimal solution.
Some robust solver like ipopt, can converge to a solution but in general it is heavily biased by your initial guess.
I would suggest you to formulate your problem with binary variables as a mixed-integer program (MIP).
Unfortunately, acados does not support MIP formulations, but you can use Casadi if you are familiar with it.
You need to use Casadi qpsol, set the solver option discrete as a list that take 1 if the corresponding variable is discrete or 0 otherwise. Lastly, choose a solver that can handle discrete variables, for instance the open-source solver HiGHS. If you have access to commercial solver you can use Gurobi, CPLEX, Mosek, …
thanks for the suggestion! I would not have found this in the casadi documentation :).
I played a bit around but seems like the problem is harder to solve than I thought … also the highs solver errors out using qpsol in casadi.
If I interpret correctly the highs solver does its mixed integer programming for linear programs but not for qp.
Running HiGHS 1.7.2 (git hash: 5ce7a275): Copyright (c) 2024 HiGHS under MIT licence terms
WARNING: Ignored 6 entries of Hessian in opposite triangle
Coefficient ranges:
Matrix [2e-08, 1e+00]
Cost [3e-02, 4e+02]
Bound [1e+00, 4e+01]
RHS [3e-04, 2e+01]
ERROR: Cannot solve MIQP problems with HiGHS
I’ll probably need to revise how to formulate my problem and what simplifications to do (or just close my eyes and let ipopt do whatever it wants to do … )
Yes, HiGHS does not solve MIQP yet.
If you want to stay within casadi I think you need to use a commercial solver like gurobi.
Otherwise on the open-source side you can use SCIP but it is not interfaced with casadi (afaik), so you should change you modelling framework to something else like AMPL or Pyomo.
Thanks for the suggestions. In casadi I’ve tried “bonmin” as a mixed integer solver.
But this was actually performing significantly worse compared to ipopt (where integers were handled by constraints). I’ll see if I can get c-code out of it through casadi using ipopt solver, before jumping yet again to another framework … .