Toward Scalable Multirobot Control: Fast Policy Learning in Distributed MPC

In plain words

A group of robots has to keep a formation, such as a line or a rectangle, while moving, avoiding each other and avoiding obstacles. The usual method is distributed model predictive control. At every time step, each robot solves a small optimisation problem with a numerical solver. That gets expensive when the robots are many and their motion is nonlinear, and small onboard computers struggle to keep up.

This paper replaces the solver with a learned rule. Each robot has a small neural network, the actor, that turns the position errors of itself and its neighbours into a steering command. A second network, the critic, estimates how the future cost depends on each robot's state. Both networks are trained a little at each step of the planning window, moving forward in time. The weights from one planning window are the starting point for the next. Once trained, the rule is a simple formula that a robot can run without a solver. Because each robot only looks at its neighbours, a rule learned on two robots can be reused on many more robots.

The paper also gives a second version of the rule for obstacles and limits. It adds barrier terms that push the robot away from the limits, in the spirit of interior-point methods for constrained optimisation.

TL;DR

The paper replaces the numerical solve in distributed MPC with explicit closed-loop policies. Each robot's actor and critic are trained online, incrementally, in every prediction interval §IV-B. The unconstrained method is stability-checked through a terminal penalty and a runtime cost-bound condition §IV-A, §IV-C. A barrier-based variant handles state and control limits §V; its steps and proofs are only in an auxiliary file that is not part of this source package. The repository's main Python script does not reproduce several of the paper's stated settings (see Paper vs code).

Problem

• Distributed MPC solves a finite-horizon optimisation at each step for each robot, usually with a numerical solver. For nonlinear robots this is computationally heavy and hard to run on small onboard hardware §I.
• Only the first control action is applied; the optimisation is then repeated at the next step §III-B.
• Existing RL approaches for multirobot control often lack online learning with stability guarantees §II. Cost shaping for safety can cause weight divergence in actor-critic learning §I, §V.
• Assumptions:
  • the robots follow the discrete local error model of eq. (3), with f_i(0) = 0 §III-A;
  • a local stabilising feedback exists for the whole formation (Assumption 1) §III-B;
  • the communication network is time-invariant and delay-free (Assumption 2) §III-B.
• Scope: the stability argument is for quadratic costs only; time-varying networks and communication delays are future work §VII.

Method

Model and cost (unconstrained case, Section IV)

• Robot i has state q_i = (p_{x,i}, p_{y,i}, \theta_i, v_i) and input (a_i, \omega_i) eq. 1. Its formation error e_i is built from its neighbours and the leader eq. 2.
• The local error model is input-affine, e_i(k+1) = f_i(e_{\mathcal N_i}(k)) + g_i(e_i(k))\,u_i(k) eq. 3. Here \mathcal N_i is the neighbour set including i, and e_{\mathcal N_i} stacks the neighbours' errors.
• The cooperative cost is J(e(k)) = \sum_i J_i(e_{\mathcal N_i}(k)), with J_i = \sum_j r_i + \|e_i(k+N)\|^2_{P_i} and stage cost r_i = \|e_{\mathcal N_i}\|^2_{Q_i} + \|u_i\|^2_{R_i} eq. 6, eq. 7. The paper requires Q_i, R_i, P_i \succ 0 §III-B.

Distributed policy learning

• The optimal control of robot i depends only on the one-step-ahead optimal costs of the robots that have i as a neighbour, \bar{\mathcal N}_i eq. 8, §IV-A.
• Within each interval, a value update and a policy update are iterated in parallel across robots eq. 9a, eq. 9b:
  • value: J_i^{t+1}(\tau) = r_i(\tau) + J_i^t(\tau+1) eq. 9a;
  • policy: u_i^{t+1} = \arg\min \{ r_i + \sum_{j \in \bar{\mathcal N}_i} J_j^{t+1}(\tau+1) \} eq. 9b.
• Terminal penalty: P_i solves a Lyapunov-type equation with F_i = A_{\mathcal N_i} + B_i K_{\mathcal N_i} and a tuning factor \beta_i that absorbs the linearisation error eq. 10, §IV-A.

F_i^\top P_i F_i - \bar P_i = -\beta_i\big(Q_i + K_{\mathcal N_i}^\top R_i K_{\mathcal N_i}\big) + \Gamma_{\mathcal N_i} \qquad [\text{eq. 10}]

Actor-critic implementation

• Critic: the costate is approximated as \hat\lambda_i = W_{c,i}^\top \sigma_{c,i}(e_{\mathcal N_i}, \tau) eq. 11. The basis \sigma may be polynomial, RBF, sigmoid, tanh, or others §IV-B.
• Desired costate, from the value update eq. 12:

\lambda_i^d(\tau) = 2 Q_i e_{\mathcal N_i}(\tau) + \sum_{j \in \bar{\mathcal N}_i} \Big(\frac{\partial f_j}{\partial e_{\mathcal N_i}}\Big)^\top \hat\lambda_i^{[j]}(\tau+1) \qquad [\text{eq. 12}]

• Critic update, with local learning rate \gamma_{c,i} eq. 13.
• Actor: \hat u_i = W_{a,i}^\top \sigma_{a,i}(e_{\mathcal N_i}, \tau) eq. 14. Let u_{o,i} = 2R_i \hat u_i and define the target u_{o,i}^d(\tau) = -\sum_{j \in \bar{\mathcal N}_i} g_i^\top(e_i(\tau))\,\hat\lambda_j^{[i]}(\tau+1) eq. 15. The actor update uses rate \gamma_{a,i} eq. 16.
• Policy used on the robot: the actor output at the current state, applied for one step. The network is then updated again in the next interval, so the policy keeps improving while it runs §IV-B.

Algorithm 1 (online, unconstrained) alg. 1

Per robot i, the paper's Algorithm 1 is an outer loop over time steps k with an inner loop over sweeps t until the change in weights, Err, is below \epsilon or t passes t_{\max}. Each sweep runs forward over the interval \tau = k, \dots, k+N-1. The pseudocode below follows the paper's steps. The paper's loop condition reads "Err \ge \epsilon \vee t \le t_{\max}" with an "or" (see Pitfalls).

# One forward sweep over the prediction interval  [alg. 1, l.6-11]
# errors: neighbourhood errors of all robots at time tau

def sweep(errors, W_c, W_a, horizon, gamma_c, gamma_a):
    for tau in range(horizon):
        # Actor output u_hat_i(tau) from the neighbourhood errors  [eq. 14]
        actions = actor(errors, W_a)
        # One-step-ahead errors e_i(tau+1) from the model  [alg. 1, l.7; eq. 3]
        next_errors = model_step(errors, actions)
        # Costate estimates at tau and tau+1  [alg. 1, l.8; eq. 11]
        costate = critic(errors, W_c)
        costate_next = critic(next_errors, W_c)
        # Desired costate from the value update  [alg. 1, l.9; eq. 12]
        costate_des = desired_costate(errors, costate_next)
        # Desired action target from the policy update  [alg. 1, l.9; eq. 15]
        action_des = desired_action(costate_next)
        # Gradient steps on the squared critic and actor errors  [alg. 1, l.10]
        W_c = W_c - gamma_c * grad_critic(costate_des - costate)  # eq. 13
        W_a = W_a - gamma_a * grad_actor(action_des - 2 * R @ actions)  # eq. 16
        errors = next_errors
    return W_c, W_a

# Outer loop at time step k  [alg. 1, l.3-5, l.12-20]
# The paper's l.5 reads "Err >= eps or t <= t_max"; "and" is used here.

def control_step(errors, W_c, W_a, horizon, tol, t_max):
    # Start from the weights of the previous interval  [alg. 1, l.4]
    # Repeat sweeps until the weight change is below tol  [alg. 1, l.5, l.13]
    err = tol + 1.0
    t = 0
    while err >= tol and t <= t_max:
        W_c_new, W_a_new = sweep(errors, W_c, W_a, horizon)
        err = weight_change(W_c_new, W_c, W_a_new, W_a)  # Err
        W_c, W_a = W_c_new, W_a_new
        t = t + 1
    # Apply the first action of the learned actor  [alg. 1, l.20; eq. 14]
    return actor(errors, W_a), W_c, W_a

# Actor: linear in basis functions of the neighbourhood error  [eq. 14]

def actor(errors, W_a):
    basis = basis_a(errors)     # sigma_a, n_u = 4  [tab. I]
    return W_a.T @ basis        # one control per robot

Practical stability check §IV-C

• The learned policy may not decrease the overall cost monotonically, so the usual MPC argument does not apply directly §IV-C.
• The paper therefore requires a baseline stabilising policy with a decreasing cost J^b eq. 17, and checks at every step that J(e(k)) \le J^b(e^b(k)) eq. 18. If this fails, the weights are re-initialised and the interval is repeated alg. 1, l.16-18.
• Two ways to build J^b are given [Remark 4]: an optimal first sequence with a terminal feedback, or any monotonically decreasing function with J^b(\infty) = 0.
• If C_{c,i} = \gamma_{c,i}\Gamma_{c,i}(\sigma_{c,i}) < 1 and C_{a,i} = 2\lambda_{\max}(R_i)\,\gamma_{a,i}\|\sigma_{a,i}\|^2 < 1, the weight-estimation errors are uniformly ultimately bounded app. A. The constant \Gamma_{c,i} depends on the model bounds \bar f from the appendix's model assumption, which the paper does not give numerically.

Safe policy learning (Section V)

• Barrier functions: for a set \mathcal Z_i, the log-barrier is \mathcal B^o_{z,i}(z_i) = -\sum_j \log(-\Xi^j_{z,i}(z_i)) on the interior of the set. It is recentred at z_{c,i} and relaxed for points close to the boundary [def. 1, eq. 19].
• Cost with barrier terms: \bar r_i = r_i + \mu(\mathcal B_{e,i} + \mathcal B_{u,i}) and a terminal barrier term eq. 20.
• Force-field policy: the control is a nominal part (objective force) plus gradient terms of the barriers (constraint forces) eq. 21:

\bar u_i(e_{\mathcal N_i}) = \nu_i(e_{\mathcal N_i}) + L_i\cdot\big(\nabla\mathcal B_{e,i}(e_{\mathcal N_i}),\ \nabla\mathcal B_{\nu,i}(\nu_i(e_{\mathcal N_i}))\big) \qquad [\text{eq. 21}]

• Terminal penalty for the barrier cost: H_{z,i} bounds the Hessian of the relaxed barrier eq. 22, and P_i solves eq. (23) §V-A.
• Actor and critic: the critic is built from the basis and the barrier gradient eq. 24, with target eq. (25) and update eq. (26). The actor has three weight blocks eq. 27, eq. 28, eq. 29 §V-B.
• The steps of the safe algorithm, its learning-rate conditions and its stability proof are "omitted due to space limitations" and left to an auxiliary file §V-B, §VI.

Robustness (Appendix A)

• With an additive disturbance bounded by \varepsilon_w and a Lipschitz constant L_p, the deviation between the real and nominal state grows at most as \frac{L_p^j - 1}{L_p - 1}\varepsilon_w app. A.
• The closed-loop state then converges to a bounded set \mathcal D^\infty_{\varepsilon_w} defined in the appendix app. A.

Hyperparameters

Two-robot actor weights from the straight-line run, to check an implementation against §VI:

w_1 = \begin{bmatrix} 0.65 & 0.4 & 0.3 & 1.58 \\ -0.26 & 0.47 & 1.2 & -0.16 \end{bmatrix}, \quad w_2 = \begin{bmatrix} -0.01 & -0.2 & 0.05 & -0.23 \\ 0.05 & 0.15 & 0.12 & -0.19 \end{bmatrix}

Here W_{a,1} = [w_1\ w_2]^\top for the first robot, and W_{a,i} = [w_1\ w_2\ w_2]^\top for the others when the policy is reused on 4, 200 or 1,000 robots §VI.

Collision-avoidance setup (16 robots): obstacles of diameter 0.4 m at centres (0,4), (0,2), (30,0), (30,-2); desired distances 1 m and 2 m; d_i = 0.2 §VI.

Results

Author-reported. Setup: the simulations used MATLAB on a laptop §VI-A. The Python timings in fig. 8 ran on a laptop and a Raspberry Pi 5. The Pi was about four times slower, and both grew linearly with the number of robots fig. 8. The baseline DMPC used IPOPT (via ALADIN-\alpha and CasADi) for [45], and fmincon for [14] §VI. The cost V is \frac{1}{M}\sum_j r_i over N_{\mathrm{sim}} = 180 steps tab. III, footnote 1. The computational time of the learning column is measured on one centralised machine tab. III, footnote 3.

Table III: comparison with DMPC tab. III

"Inf" means infeasible tab. III, footnote 4. Timings are per time step for all robots.

Other headline results

• Learning at scale. The formation forms from a disordered start, with a transient of about 3 s even at 10,000 robots §VI, fig. 5A.
• Transfer. The 2-robot straight-line policy was deployed directly to 4, 200 and 1,000 robots in straight-line and circular formations §VI, fig. 5B. The paper gives no error numbers for these runs.
• Cost comparison, two metrics disagree. Table III's per-step cost V is lower for the proposed method without obstacles. The cumulative cost J_c over N_{\mathrm{sim}} = 180 steps at M=2 is 8.9 for DMPC and 16.9 for the proposed method §VI. The paper argues that closed-loop performance matters more than per-interval optimality §VI.
• Collision avoidance versus cost-shaping RL (100 tests, 2 robots) §VI, fig. 7. Read from the plot: the proposed method succeeds in 100% of runs at 1 and 1.5 m/s, then about 91% at 2 m/s and 77% at 2.5 m/s. Cost-shaping RL succeeds about 12% of the time at 1 m/s and 0% at 2 m/s.
• Learned formation with obstacles, 16 robots §VI, tab. II:

• Multirotor drones in Gazebo (XTDrone, PX4) with M = 6, 18, 40. Horizontal control only; height from a baseline controller; roll and pitch treated as disturbances §VI-B, fig. 9. For six drones, state errors E_x and E_y were lower than the authors' baseline formation controller at leader speeds from 0.5 to 3 m/s fig. 10. The bars are read from the plot.
• Wheeled robots, real world (Δt = 0.1 s) §VI-C, tab. IV:

Maximum absolute errors, as reported. The paper states per-robot deployment cost of about 0.02 ms per step §VI-C.

Paper vs code

Repository. The paper's multimedia section links xinglongzhangnudt/policy-learning-for-distributed-mpc §Multimedia. The owner is the first author's GitHub account, and CODEOWNERS lists the same handle. Treat it as official. We cloned the default branch and pinned code_commit to its HEAD, 2e8accb.

Licence. The LICENSE file is headed "GPLv4 License", but its body is the GNU GPL, version 3 or later. We record GPL-3.0-or-later. This is copyleft, so the card cites code by path and line only and quotes none of it. Other folders (DSLC_training/, DSLC_xtdrone6/, DSLC_xtdrone18/, DSLC_solving_one_robot_control/) were only listed, not read. Barrier-function MATLAB files exist (for example DSLC_training/Barriera.m). I did not audit whether they implement the safe variant of Section V. The MATLAB twin of the scale script is dlpc_online_train_scales_to_10000.m.

Compared files. Main script dlpc_online_train_scales_to_10000.py and its MATLAB twin, both unconstrained.

Discrepancies between paper and code (paper first, code second):

• Discount factor. The paper's targets, eq. (9a) to (15), have no discount. The code uses \Gamma = 0.95 in the costate and actor targets dlpc_online_train_scales_to_10000.py.
• Prediction horizon. Paper: N = 20 §VI-A. Code: PreHor_len = 10 dlpc_online_train_scales_to_10000.py. The MATLAB twin also uses 10 dlpc_online_train_scales_to_10000.m.
• Stage cost scale. Paper: Q_i = I, R_i = 0.5 I_2 tab. I. Code: a scale of 0.01 multiplies Q, R and P dlpc_online_train_scales_to_10000.py.
• Barrier weight. Paper: \mu = 0.02 tab. I. Code: \mu = 0.001 dlpc_online_train_scales_to_10000.py. In the unconstrained script \mu is only used in a commented-out line dlpc_online_train_scales_to_10000.py.
• Weight initialisation. Paper: [0, 0.1] for the 10,000-robot run §VI. Code: U(-0.5, 0.5) for both actor and critic dlpc_online_train_scales_to_10000.py.
• Network form. Paper: one hidden layer with tanh §VI-A. Code: output is W\tanh(x) + b, with HiddenUnitNum = 5 declared but never used dlpc_online_train_scales_to_10000.py. Inputs are divided by 2 before the network dlpc_online_train_scales_to_10000.py. Neither the scaling nor the tanh-only form is in the paper.
• Learning loop. Paper: sweeps repeat until the weight change Err is below \epsilon or t > t_{\max} alg. 1, l.5, l.13. Code: one forward sweep per interval. If robot 1's error exceeds its bound, the sweep is retried and networks are reset after 30 failed attempts dlpc_online_train_scales_to_10000.py. Err is never computed; Err1 and Err2 are set but never used dlpc_online_train_scales_to_10000.py.
• Stability check. Paper: J(e(k)) \le J^b(e^b(k)) is checked each step, eq. (18). Code: no baseline J^b. The per-interval cost is accumulated only when timesteps == 2, is never used, and its output line is commented out dlpc_online_train_scales_to_10000.py.
• Terminal penalty. Paper: P_i solves eq. (10) with \beta_i = 1.1 and the gains K_{\mathcal N_i} §IV-A, tab. I. Code: six hard-coded matrices, scaled by 0.01. Robots i \ge 6 reuse P_6 dlpc_online_train_scales_to_10000.py. The code does not solve eq. (10).
• Default robot count. Paper: up to 10,000 robots. Code: num = 7500 dlpc_online_train_scales_to_10000.py. MATLAB: num = 20 dlpc_online_train_scales_to_10000.m.
• Simulation length. Paper: N_{\mathrm{sim}} = 180 for Table III tab. III. Code: 20 steps dlpc_online_train_scales_to_10000.py. MATLAB: 200 dlpc_online_train_scales_to_10000.m. The comment next to each default says 180.
• Constraint check. The loop stops when robot 1's error exceeds a bound dlpc_online_train_scales_to_10000.py. The variable state_bound is 5 at start dlpc_online_train_scales_to_10000.py. draw_fig2Scale declares it global and then sets it to 2 dlpc_online_train_scales_to_10000.py, so the bound changes after the first call. MATLAB has the same check dlpc_online_train_scales_to_10000.m.
• Parallelism. The "distributed" update is a serial loop over robots dlpc_online_train_scales_to_10000.py. Time grows linearly with M by construction. The timings in Table III are from one machine tab. III, footnote 3.
• Initial positions. Paper: "disordered initialisation" §VI. Code: two rows of robots, offset by a uniform random term in [0, 1) dlpc_online_train_scales_to_10000.py.

Consistent with the paper: \Delta t = 0.05 s dlpc_online_train_scales_to_10000.py and the learning rates \gamma_c = 0.4, \gamma_a = 0.2 dlpc_online_train_scales_to_10000.py. Desired formation offsets and the leader speed of 1 m/s are in the code dlpc_online_train_scales_to_10000.py.

Pitfalls

• Algorithm 1 loop condition. Read literally, "Err \ge \epsilon \vee t \le t_{\max}" keeps looping while t \le t_{\max} even when Err is small. The sensible reading is "and". The paper does not say which was meant alg. 1, l.5.
• Re-initialisation. Algorithm 1 says to "repeat steps 3-15" after a failed check alg. 1, l.17-18. It does not say whether this means the whole interval or only the sweep. The released code restarts the sweep and resets networks after a fixed count dlpc_online_train_scales_to_10000.py.
• Learning rates were tuned per experiment. Table I's note says the parameters were "dynamically tuned in the simulation studies" tab. I. No general recipe is given, and the convergence conditions are not checked in code.
• Convergence conditions are hard to check. C_{a,i} < 1 with R_i = 0.5I_2 and \gamma_{a,i} = 0.2 requires \|\sigma_{a,i}\|^2 < 5 (from 2 \cdot 0.5 \cdot 0.2 \cdot \|\sigma\|^2 < 1). Tanh features with four outputs satisfy this. The critic condition depends on the model bound \bar f and is not given numerically app. A.
• Safety degrades with speed. In collision avoidance, success falls above 1.5 m/s. The authors attribute this to a violated learning-rate condition, as the barrier gradient grows quickly with speed §VI, fig. 7. Higher speeds need smaller rates, which the paper does not quantify.
• Collision-avoidance cost. At M=2 the proposed method's cost V (13.5) is higher than IPOPT's (5.0). The DMPC baseline is infeasible at M=4 tab. III. Do not read the collision-avoidance rows as a win in cost.
• Cost claims. Table III favours the proposed method on V without obstacles. The text's cumulative J_c at M=2 favours DMPC (8.9 vs 16.9) §VI. The two metrics are different and the paper does not test a third.
• Horizon mismatch in the comparison. The baselines use N = 10, the proposed method N = 20 §VI. The table compares unequal horizons.
• Auxiliary material missing. The safe-learning algorithm, its learning-rate conditions, the stability proof (including the theorem cited as Theorem 8) and the 20 repeated tests are in "auxiliary-results.pdf" in "auxiliary-material.zip". That file is not in the source package, so the safe variant cannot be implemented from this card alone §V-B, §VI, §VI-A.
• Sim-to-real sampling. Simulation uses \Delta t = 0.05 s and real robots use 0.1 s tab. I, §VI-C. The paper does not say whether the weights were retrained for the new step.
• Nominal model only. The policy is learned from the nominal kinematic model. Model mismatch is handled only by a bounded-disturbance argument and a Lipschitz assumption §VII, app. A. Quadratic cost only §VII.
• Timing claims are not end-to-end. The deployment time of about 0.02 ms per robot per step §VI-C is for evaluating the learned policy, not for learning. The learning times in Table III are for all robots on one machine.
• Equation numbers. The equation and algorithm numbers above were read from the compiled PDF, not counted in the LaTeX. The appendix uses its own numbering, cited here as app. A.

Lineage

• Builds on: the cooperative DMPC formulation of [45], which this paper extends with a nonlinear terminal penalty §III-B, §IV-A. RL-based MPC with receding horizon in actor-critic form [16, 40] §II. None of these has a card in this wiki.
• The paper extends the authors' earlier conference version §I; that version is not identified here.
• Compared against: numerical DMPC of [45] (IPOPT via ALADIN-\alpha and CasADi) and [14] (fmincon, platoon control) §VI, tab. III; cost-shaping RL of [32] §VI, fig. 7; the authors' baseline formation controller for drones §VI-B, fig. 10.
• Superseded by: not known from this wiki's content.