Skip to content

Repository files navigation

simple_casadi_mpc

build docs

Lightweight C++ utilities for building and solving MPC problems with CasADi. Includes runtime MPC, JIT-compiled MPC, and CMake-integrated compiled MPC solvers.

Dependencies

  • CasADi (install script)
  • IPOPT or FATROP (optional solver backends)
  • Eigen3
  • Python3 + NumPy + pybind11
  • matplotlibcpp17 (examples/benchmarks use it via pybind11)
  • doxygen (for docs; uses doc/doxygen-awesome-css submodule)

Build & install

mkdir build
cd build
cmake ..
make
sudo make install

CMake usage

find_package(simple_casadi_mpc REQUIRED)
target_link_libraries(my_target PRIVATE simple_casadi_mpc)

Solver overview

  • MPC: simplest runtime solver; easiest for quick validation.
  • JITMPC: JIT-compiles on the first solve for faster subsequent runs; expect a startup lag (cacheable with ccache).
  • CompiledMPC: builds solver code at CMake time; best steady-state speed with no runtime lag.

Limitation for CompiledMPC: the solver backend (IPOPT/FATROP/...) and its parameters are fixed at build time.

Problem format

The library solves the following discrete-time finite-horizon optimal control problem at every solve() call:

$$ \begin{aligned} \min_{\substack{x_0,\dots,x_N \cr u_0,\dots,u_{N-1}}} \quad & \sum_{k=0}^{N-1} L(x_k, u_k, k;,p) ;+; \Phi(x_N;,p) \cr \text{s.t.} \quad & x_0 = x_\text{init}, \cr & x_{k+1} = f_d(x_k, u_k;,p), && k = 0, \dots, N-1, \cr & g_\text{eq}^{(i)}(x_k, u_k;,p) = 0, && k = 0, \dots, N-1,; i = 1, \dots, n_\text{eq}, \cr & g_\text{ineq}^{(j)}(x_k, u_k;,p) \le 0, && k = 0, \dots, N-1,; j = 1, \dots, n_\text{ineq}, \cr & x_{\text{lb},k} \le x_k \le x_{\text{ub},k}, && k = 0, \dots, N, \cr & u_{\text{lb},k} \le u_k \le u_{\text{ub},k}, && k = 0, \dots, N-1. \end{aligned} $$

Symbol mapping:

Symbol Code
$x_k \in \mathbb{R}^{n_x}$ state at stage $k$ (Problem::nx())
$u_k \in \mathbb{R}^{n_u}$ control at stage $k$ (Problem::nu())
$N$ prediction horizon (Problem::horizon())
$L(x_k, u_k, k;,p)$ Problem::stage_cost(x, u, k)
$\Phi(x_N;,p)$ Problem::terminal_cost(x)
$f_d$ discretization of Problem::dynamics(x, u) per DynamicsType
$g_\text{eq}^{(i)},,g_\text{ineq}^{(j)}$ each call to Problem::add_constraint(Equality / Inequality, ...) adds one $g_\text{eq}^{(i)}$ or $g_\text{ineq}^{(j)}$
$x_{\text{lb},k},,x_{\text{ub},k}$ Problem::set_state_bound(...)
$u_{\text{lb},k},,u_{\text{ub},k}$ Problem::set_input_bound(...)
$p$ runtime parameters via Problem::parameter(...) / reference_trajectory(...)
$x_\text{init}$ first argument to MPC::solve(x, ...)

For DynamicsType::ContinuesForwardEuler, ContinuesModifiedEuler, and ContinuesRK4, $f_d$ is the corresponding one-step integrator applied to the user-supplied continuous dynamics $\frac{dx}{dt} = f(x, u)$ with step size $\Delta t$. For DynamicsType::Discretized the user supplies $f_d$ directly.

Soft constraints add slack $s \ge 0$ and a penalty term to the cost; see Soft path constraints below.

Defining a problem

Define your MPC problem by deriving from simple_casadi_mpc::Problem and overriding dynamics, stage_cost, and (optionally) terminal_cost:

#include "simple_casadi_mpc/simple_casadi_mpc.hpp"

class MyProblem : public simple_casadi_mpc::Problem {
public:
  MyProblem()
      // DynamicsType, nx, nu, horizon (N), dt
      : Problem(DynamicsType::ContinuesRK4, /*nx=*/2, /*nu=*/1, /*N=*/20, /*dt=*/0.05) {
    // Optional: per-stage input/state bounds
    set_input_bound(Eigen::VectorXd::Constant(1, -1.0),
                    Eigen::VectorXd::Constant(1,  1.0));
    // Optional: runtime-tunable parameters (updated each `solve` call)
    x_ref_ = reference_trajectory("x_ref"); // shape (nx, N)
  }

  // Continuous dynamics: return dx/dt for {ContinuesForwardEuler|ContinuesModifiedEuler|ContinuesRK4}
  // Discrete dynamics:   return x_{k+1} for {Discretized}
  casadi::MX dynamics(casadi::MX x, casadi::MX u) override {
    return casadi::MX::vertcat({x(1), u});
  }

  // Per-stage cost. `k` is the stage index in [0, horizon).
  casadi::MX stage_cost(casadi::MX x, casadi::MX u, size_t k) override {
    casadi::MX e = x - x_ref_(casadi::Slice(), k);
    return casadi::MX::mtimes(e.T(), e) + 0.1 * casadi::MX::mtimes(u.T(), u);
  }

  // Terminal cost on x_N (default returns 0 if not overridden).
  casadi::MX terminal_cost(casadi::MX x) override {
    return 10.0 * casadi::MX::mtimes(x.T(), x);
  }

private:
  casadi::MX x_ref_;
};

Solve loop:

auto prob = std::make_shared<MyProblem>();
simple_casadi_mpc::MPC mpc(prob); // or JITMPC / CompiledMPC

Eigen::VectorXd x = /* initial state */;
for (...) {
  // Update parameters declared via `parameter(...)` / `reference_trajectory(...)`
  casadi::DM x_ref_dm = /* (nx, N) DM */;
  Eigen::VectorXd u = mpc.solve(x, {{"x_ref", x_ref_dm}});
  x = prob->simulate(x, u, sim_dt);
}

Hard path constraints

Add equality / inequality path constraints applied at every stage:

// In your Problem subclass constructor:

// Inequality g(x, u) <= 0  (e.g. circular obstacle: r^2 - |xy - c|^2 <= 0)
add_constraint(ConstraintType::Inequality, [](casadi::MX x, casadi::MX u) {
  using namespace casadi;
  MX dx = x(0) - 1.0, dy = x(1) - 0.5;
  return MX::vertcat({0.5 * 0.5 - (dx * dx + dy * dy)});
});

// Equality g(x, u) = 0
add_constraint(ConstraintType::Equality, [](casadi::MX x, casadi::MX u) {
  return casadi::MX::vertcat({u(0) + u(1)}); // e.g. force a coupling
});

To apply a constraint only to a subset of stages, use add_constraint_at(type, func, start, end). Range semantics match set_input_bound:

// Only at stage 0 (single-stage form, omit `end`)
add_constraint_at(ConstraintType::Equality,
                  [&](casadi::MX x, casadi::MX) { return x - x_target; }, 0);

// Only on stages [N-3, N) (terminal-band)
add_constraint_at(ConstraintType::Inequality,
                  [&](casadi::MX x, casadi::MX) { return x - x_safe; }, N - 3, N);

Soft path constraints

add_soft_constraint(type, func, w1, w2) introduces non-negative per-stage slack variables $s \ge 0$ and adds the penalty $w_1 \mathbf{1}^\top s + \tfrac{1}{2} w_2 s^\top s$ to the cost. The original constraint is relaxed as:

  • Inequality $g(x, u) \le 0 ;\to; g - s \le 0$.
  • Equality $h(x, u) = 0 ;\to; |h| \le s$, i.e. $h - s \le 0$ and $-h - s \le 0$.
// Soft inequality with L1 weight 1e3 (default) and no L2 term.
add_soft_constraint(ConstraintType::Inequality, my_constraint);

// Mixed L1 + L2 penalty.
add_soft_constraint(ConstraintType::Inequality, my_constraint, /*w1=*/1e2, /*w2=*/1.0);

The default w2 = 0 gives a pure L1 (exact) penalty; w2 > 0 adds an L2 (smooth) term. Large w1 recovers hard-constraint behavior; small w1 lets the optimizer trade violation against the rest of the cost. Hard add_constraint and soft add_soft_constraint may be mixed on the same problem.

Stage-specific variants are also available — same start/end semantics as add_constraint_at:

add_soft_constraint_at(ConstraintType::Inequality, my_constraint,
                       /*w1=*/1e3, /*w2=*/0.0,
                       /*start=*/N - 3, /*end=*/N); // only the last 3 stages

Usage for CompiledMPC via CMake

find_package(simple_casadi_mpc REQUIRED)

# Generate a compiled solver (codegen step happens at build time)
add_simple_casadi_mpc_codegen(
  <solver_target_name>                  # e.g., my_problem
  <codegen_cpp>                         # e.g., my_problem_codegen.cpp (derives Problem)
  EXPORT_SOLVER_NAME <export_name>      # optional, default is <solver_target_name>_compiled_solver
  INCLUDE_DIRS ${CMAKE_CURRENT_SOURCE_DIR} # where your Problem header lives
  SOLVER_NAME <casadi_solver>           # optional; default is fatrop (e.g., ipopt/fatrop/...)
  # LINK_LIBS ...                       # optional; extra solver libs if needed
)

# Link your executable against the generated solver + simple_casadi_mpc
add_executable(<your_exe> main.cpp
                ${<solver_target_name>_COMPILED_SOLVER_CONFIG_SOURCE})
target_include_directories(<your_exe> PRIVATE
  ${CMAKE_CURRENT_SOURCE_DIR} ${<solver_target_name>_CODEGEN_DIR})
target_link_libraries(<your_exe> PRIVATE
  simple_casadi_mpc::simple_casadi_mpc
  ${<solver_target_name>_COMPILED_SOLVER})

Examples

double_integrator_mpc_example

Drives a frictionless point mass to the origin (position and velocity feedback).

From: example/double_integrator_mpc_example.cpp

cartpole_mpc_example

Cartpole swing-up and balance (problem setup from the linked gists).

From: example/cartpole_mpc_example.cpp

https://gist.github.com/mayataka/ef178130d52b5b06d4dd8bb2c8384c54 https://gist.github.com/mayataka/bc08faa63a94d8b48ceba77cc79c7ccc

inverted_pendulum_mpc_example

Rotary inverted pendulum swing-up with torque limits that force a multi-phase motion.

From: example/inverted_pendulum_mpc_example.cpp

diff_drive_mpc_example

Differential-drive robot from top-left to bottom-right while avoiding circular obstacles and respecting velocity limits.

From: example/diff_drive_mpc_example.cpp

diff_drive_soft_constraint_example

Same diff-drive setup with a single circular obstacle, comparing add_constraint (hard) and add_soft_constraint (soft) for the obstacle. With a large penalty weight the soft formulation matches the hard one; with a small weight the optimizer prefers cutting through the obstacle if the tracking gain dominates the violation cost.

From: example/diff_drive_soft_constraint_example.cpp

Tips

Choosing a solver backend

Backend Default config Best for
IPOPT default_ipopt_config() General-purpose NLP; default and easiest to debug
FATROP default_fatrop_config() OCP-structured problems; fastest for MPC
qpOASES default_qpoases_config() SQP method backed by qpOASES (linear-quadratic)

Pass a copy of one of these dicts as the third argument to MPC / JITMPC, mutating any keys you need:

auto cfg = simple_casadi_mpc::MPC::default_fatrop_config();
cfg["fatrop.max_iter"] = 100;
simple_casadi_mpc::MPC mpc(prob, "fatrop", cfg);

FATROP requires CasADi to be built with WITH_FATROP=ON and WITH_BLASFEO=ON; see install_casadi.sh.

Picking MPC vs JITMPC vs CompiledMPC

  • Reach for MPC while iterating on the problem itself; no compile step.
  • Switch to JITMPC once the model is stable. The first solve triggers a gcc -O3 -march=native compile; cache it with ccache (the default default_jit_options() already does so).
  • Use CompiledMPC (with the add_simple_casadi_mpc_codegen CMake helper) when you want zero runtime startup cost and the solver backend is fixed.

Customizing JIT compile options

Override the defaults from JITMPC::default_jit_options():

auto opts = simple_casadi_mpc::JITMPC::default_jit_options();
opts["compiler"] = "clang";
opts["flags"] = "-O2 -fno-fast-math";
opts["verbose"] = true;
simple_casadi_mpc::JITMPC mpc("my_prob", prob, "ipopt",
                               simple_casadi_mpc::MPC::default_ipopt_config(),
                               opts);

Runtime-tunable parameters

Problem::parameter(name, rows, cols) returns a symbolic MX; update it at each solve:

mpc.solve(x, {{"x_ref", x_ref_dm}, {"obstacle", obs_dm}});

reference_trajectory(name) is a shorthand for parameter(name, nx, horizon) whose column k is consumed at stage k.

Performance knobs (MPC / JITMPC config)

Two simple-casadi-mpc-specific options can be passed in the config dict (consumed before being forwarded to CasADi):

  • mapsum_stage_cost (default true): build the stage-cost sum via MapSum so first-order AD stays loop-shaped.
  • expand_inner_functions (default true): SX-expand per-stage F/L/G before mapping for faster JIT compilation.

Both default to true. If your stage_cost has stage-dependent branching beyond per-stage parameter slicing, the library auto-falls-back to a per-stage loop (warning emitted) and you can also set mapsum_stage_cost = false explicitly.

Variable time step

Problem(DynamicsType, nx, nu, horizon, dt) uses a uniform $\Delta t$. To use a per-stage step (e.g. coarse-to-fine schedules), pass a std::vector<double> of length horizon instead:

std::vector<double> dts(N);
for (size_t k = 0; k < N; ++k) dts[k] = 0.02 + 0.005 * k; // refine the prediction
MyProblem(N, dts);

The continuous integrator (f_d) is built with $\Delta t_k$ as a symbolic per-stage input, so the solver still uses one shared dynamics function. Problem::dt() returns the stage-0 step for backward compatibility, while Problem::dt(k) and Problem::dts() give per-stage access. Problem::has_uniform_dt() reports whether all stages share the same step.

Warm starting

MPC::solve caches the previous x, lam_x, lam_g internally and feeds them to the next solve, so closed-loop simulations naturally benefit. Solver-side warm start is also enabled in default_ipopt_config() (ipopt.warm_start_init_point = "yes").

Benchmarks

Runtime comparisons for cartpole MPC solver variants.

Documentation

  1. Fetch submodules:
git submodule update --init --recursive
  1. Generate docs:
cd doc
doxygen Doxyfile
  1. Open doc/build/html/index.html.

About

C++用CasADiでMPCする用ライブラリ

Resources

Stars

30 stars

Watchers

1 watching

Forks

Releases

Packages

Contributors

Languages