Nonlinear model predictive control (NMPC) can be applied to large nonlinear systems with complex objectives while enforcing constraints, but the cost of solving the underlying optimization problems limits real-time application for systems with fast dynamics or many states. Neural networks can reduce this online cost by approximating the control law, but may violate constraints and typically require many precomputed optimizer solutions for training. The learning-based iterative solver for constrained optimization (LISCO) instead learns to predict primal-dual steps toward satisfaction of the Karush-Kuhn-Tucker (KKT) conditions through self-supervised training without precomputed solutions, but its dense, monolithic architecture grows with the prediction horizon. We remove this horizon dependence by exploiting the temporal structure of the problem in the solver architecture, which makes the network size independent of the horizon while leaving the training loss and the solver iteration unchanged. A closed-loop resampling strategy further aligns the training distribution with closed-loop operation at no additional sampling cost. On an economic NMPC benchmark of a CSTR cascade, the proposed solver scales to larger problems than monolithic LISCO and achieves median speedups of up to a factor of four over IPOPT. In closed loop, it achieves comparable control performance to IPOPT with no observed constraint violations.
Figures & tables
Figure 1: Monolithic LISCO update map. From the parameters p and the current iterate z^k , the modified KKT residual ( 8 ) is evaluated and a single dense network maps it to the full primal-dual step Δz^k , which updates the iterate. Iterating this map drives the residual to the tolerance. The decomposed map proposed in this work is shown for comparison in Figure 2 .
Figure 2: LISCO update map with temporal decomposition. The global modified KKT residual ( 8 ) is split into overlapping stage windows. A single shared network is evaluated on all windows in parallel. The overlapping predictions are combined by the fixed averaging operator. The per-stage magnitude and the global step scaling γ are restored to form the primal-dual step.
Figure 3: Time decomposition of the horizon. The stage-partitioned primal-dual iterate z=(z0,…,zN−1) is covered by overlapping windows S0,…,SM−1 of m consecutive stages, shifted by the stride s , so consecutive windows overlap in m−s stages (red). A stage that lies in two windows, such as stage 2 , receives one prediction per window, and these are averaged into a single step by the merge operation ( 25 ). The number of windows is M=1+(N−m)/s .
Symbol
Value [Unit]
Description
nr
study dependent [–]
Number of reactors
Δt
1.0 [s]
Discretization step
VR
15 [L]
Total reactor volume
Vout
1 [L/s]
Volumetric flow rate
k1
2 [L/(mol s)]
Rate constant
k2
2 [L/(mol s)]
Rate constant
Table 1: Parameters of the CSTR cascade model. The first reaction A+B→R is denoted by subscript 1 and the second reaction 2A→S by subscript 2 .
Variable
Unit
Scaling range
Bounds
cA,n
mol/L
[0,0.1]
[0,∞)
cB,n
mol/L
[0,1]
unconstrained
cR,n
mol/L
[0,4]
unconstrained
cS,n
mol/L
[0,1]
unconstrained
Tr,n
∘ C
[20,80]
[20,80]
uA,n
mol/s
[0,1.5]
[0,1.5]
Table 2: Scaling ranges and decision variable bounds. The listed scaling ranges are mapped to [0,1] .
Group
Setting
Value
Network
hidden layers
2×1024 , GELU [ 19 ]
precision
bfloat16
Solver
step scaling γ
0.1
norm clamp ηmax
100
smoothing ϵ
10−1
penalization ρFB
0.8
Table 3: Shared configuration of all studies.
Figure 4: CPU solve time against the horizon length N at nr=5 for the decomposed solver and IPOPT with MA27, with the median (solid) and the maximum (dotted) over 256 instances per configuration. The monolithic solver is not shown, as it does not converge at nr=5 .
monolithic
decomposed
nr
N
npar [M]
ttrain [h]
conv.
npar [M]
ttrain [h]
conv.
2
10
1.84
0.11
1.000
1.47
0.14
1.000
2
20
2.62
0.23
0.945
1.47
0.27
1.000
2
50
4.95
0.54
0.000
1.47
0.77
1.000
2
100
–
–
–
1.47
2.05
1.000
3
10
2.23
0.12
0.996
1.67
0.15
1.000
Table 4: Training cost and convergence of the monolithic and the decomposed solver at equal hidden architecture and training budget for all reactor counts nr and horizon lengths N . Reported are the trainable parameters npar , the training times ttrain and the fraction of the 256 evaluation samples reaching ∥rϕ∥∞≤10−5 within 300 iterations, where 0.000 means that no sample reached the tolerance. Dashes mark configurations without a monolithic run.
tLISCO [ms]
tIPOPT [ms]
speedup [ × ]
nr
N
nz
med
p95
max
med
p95
max
min
med
max
2
10
380
0.89
1.29
1.66
2.65
2.92
3.39
1.64
2.99
4.14
2
20
760
1.69
1.99
6.21
3.47
3.77
4.42
0.52
2.03
2.83
2
50
1900
3.44
3.93
4.36
5.79
6.24
6.99
1.28
1.68
2.23
2
100
3800
6.24
9.91
10.66
10.97
11.63
12.00
1.00
1.72
2.02
3
10
570
1.08
1.46
2.05
3.50
3.81
4.43
1.65
3.23
4.86
Table 5: CPU solve times of the decomposed solver (LISCO) and IPOPT with MA27 over 256 instances per configuration, with the primal-dual dimension nz and the median, 95th percentile and maximum solve times. The speedup is the per-sample ratio tIPOPT/tLISCO on identical instances, reported by its minimum, median and maximum. Both solvers converge on every instance. The monolithic solver is omitted, as it does not converge at most configurations.
Figure 5: Fraction of the 256 samples reaching the residual tolerance ∥rϕ∥∞≤10−5 over the solver iterations for the smoothing parameters ϵ∈{10−1,…,10−6} at nr=2 and N=50 , with an iteration cap of 300.
Figure 6: Absolute relative objective gap ∣q−q∗∣/∣q∗∣ at the final iterate for the smoothing sweep at nr=2 and N=50 over 256 samples, where q∗ denotes the objective of the IPOPT reference solution at tolerance 10−8 . The line marks the median and the boxes span the 10th to 90th percentile.
Figure 7: 90th percentile of the iterations to the residual tolerance ∥rϕ∥∞≤10−5 for different window lengths m and strides s , with nr=5 reactors and N∈{49,50,51} . Each configuration is evaluated on 512 samples with a maximum number of 200 iterations. The grey cell indicates that fewer than 90% of the samples converge within the iteration limit. Dashes indicate combinations that violate the admissibility condition ( 18 ) and were not evaluated. The black outline marks the lowest observed 90th percentile.
metric
IPOPT
LISCO
approx. MPC
closed-loop cost (mean)
-358.4415
-358.8662
-358.4793
solve time median [ms]
18.61
5.20
0.08
solve time p90 [ms]
19.20
6.30
0.10
solve time max [ms]
21.41
11.62
0.43
input viol. (mean)
0.0000
0.0000
0.0000
input viol. (max)
0.0000
0.0000
0.0002
Table 6: Closed-loop comparison of IPOPT, the decomposed solver (LISCO) trained with closed-loop resampling and approximate MPC trained on closed-loop data, over 50 trajectories of 50 steps each, started from uniformly sampled initial states at nr=5 and N=50 . Reported are the mean closed-loop cost, the per-solve times and the mean and maximum constraint violations on the scaled variables for the input and state bounds, where violations below the tolerance 10−5 are not counted. The last two rows list the offline data sampling and training times.
step bin
closed-loop trained
uniform trained
steps 1–1
100.0%
100.0%
steps 2–5
100.0%
84.5%
steps 6–10
100.0%
0.4%
steps 11–25
100.0%
0.0%
steps 26–50
100.0%
0.0%
all 1–50
100.0%
8.8%
Table 7: Fraction of closed-loop solves reaching the residual tolerance ∥rϕ∥∞≤10−5 within 300 iterations, binned by the simulation step, for the decomposed solver trained with closed-loop resampling ( qreset=0.02 ) and with uniform resampling ( qreset=1 ), over the 50 trajectories of 50 steps of the closed-loop study.
Nonlinear Model Predictive Control requires solving a constrained nonlinear program (NLP) in real-time at every sampling instant, a computational bottleneck that limits deployment on resource-constrained hardware or at high sampling rates. We address this challenge for the broad class of input-affine nonlinear systems to show that the optimal control move can be approximated by a state-dependent quadratic program (QP) whose cost parameters depend on the current state and reference. We propose a single-network residual-corrector architecture: a state-dependent analytic baseline provides initial QP parameters, and the network learns only the corrections needed to match the full NLP solution; the QP is solved by a differentiable interior-point layer, guaranteeing constraint satisfaction for the first control action. The network is trained offline on data generated by an NLP solver using a hybrid loss that combines supervised imitation and KKT-residual penalties. We validate the approach on a three-link planar robotic arm with Cartesian end-effector tracking, demonstrating orders-of-magnitude speedup over the NLP solver while maintaining comparable tracking performance.
The practical deployment of nonlinear model predictive control (NMPC) is often limited by online computation: solving a nonlinear program at high control rates can be expensive on embedded hardware, especially when models are complex or horizons are long. Learning-based NMPC approximations shift this computation offline but typically demand large expert datasets and costly training. We propose Sequential-AMPC, a sequential neural policy that generates MPC candidate control sequences by sharing parameters across the prediction horizon. For deployment, we wrap the policy in a safety-augmented online evaluation and fallback mechanism, yielding Safe Sequential-AMPC. Compared to a naive feedforward policy baseline across several benchmarks, Sequential-AMPC requires substantially fewer expert MPC rollouts and yields candidate sequences with higher feasibility rates and improved closed-loop safety. On high-dimensional systems, it also exhibits better learning dynamics and performance in fewer epochs while maintaining stable validation improvement where the feedforward baseline can stagnate.
Institute of Computer Engineering, TU Wien, Vienna, Austria · AIT Austrian Institute of Technology, Vienna, Austria · Institute of Mechanics and Mechatronics, TU Wien, Vienna, Austria
Nonlinear Parametric Optimization Network (NLPOpt-Net) is an unsupervised learning architecture to solve constrained nonlinear programs (NLP). Given the structure of an NLP, it learns the parametric solution maps with guaranteed constraint satisfaction. The architecture consists of a backbone neural network (NN) followed by a multilayer (k-layered) projection. While the NN drives toward optimality through a loss function consisting of a modified Lagrangian augmented with a consistency loss, the projection ensures feasibility by projecting the NN predictions in the original constraint manifold. Instead of typical distance minimization, our projection exploits local quadratic approximations of the original NLP. Under certain conditions (such as convexity), the projection has a descent property, which improves the NN predictions further. NLPOpt-Net deploys an inversion-free, modified Chambolle-Pock algorithm to solve the constrained quadratic projections during the forward pass and uses the implicit function theorem for efficient backpropagation. The fixed structure of the projection further allows decoupling of the NN and the projection once the training is complete. NLPOpt-Net solves large-scale convex QP, QCQP, NLP, and nonconvex problems with near zero optimality gap and constraint violations reduced to machine precision. Additionally, it provides near accurate prediction of the active sets and corresponding dual variables, thereby enabling a scalable approach for multiparametric programming. Compiling the projection in C provides order of magnitude improvement in inference time compared to JAX. We provide the codes and NLPOpt-Net as a ready to use package that includes GPU support.