Learning Instance-Solution Operator For Optimal Control
Authors: Mingquan Feng
Zhijie Chen
Junchi Yan
Venue: Submitted to ICLR 2023 (2023)

Abstract: Optimal control problems (OCPs) aim at finding a control function for a dynamical system such that a cost functional is optimized. These problems are central to physical system research in both academia and industry. In this paper, we propose a novel instance-solution operator learning perspective, which solves OCPs in a one-shot manner with no dependence on the explicit expression of dynamics or iterative optimization processes. The design is in principle endowed with substantial speedup in running time, and the model reusability is guaranteed by high-quality in- and out-of-distribution generalization. We theoretically validate the perspective by presenting the approximation bounds for the instance-solution operator learning.  Extensive experiments on 6 physical systems verify the effectiveness and efficiency of our approach. The source code will be made publicly available.

Under review as a conference paper at ICLR 2023
LEARNING INSTANCE-SOLUTION OPERATOR FOR OP- TIMAL CONTROL
Anonymous authors Paper under double-blind review
ABSTRACT
Optimal control problems (OCPs) involves finding a control function for dynam- ical system such that cost functional is optimized, which are central to physical system research both academia and industry. this paper, we propose a novel instance-solution operator learning perspective, which solves OCPs in a one-shot manner with no dependence on the explicit expression of dynamics or iterative op- timization processes. The design is in principle endowed with substantial speedup in running time, and the model reusability is guaranteed by high-quality in- and distiiibibtttoionn generalization. We theoretically validate the perspective by presenting the approximation bounds for the instance-solution operator learning. Experiments on synthetic environments and real-world dataset verify the ef- fectiveness and efficiency of our approach. The source code will be made publicly available.
INTRODUCTION
The explosion of data for embedding the physical world is reshaping the ways we understand, model, and control dynamical systems. Though control theory has been classically rooted in model-based design and solving paradigm, the demands of model reusability, and the opacity of complex dynam- ical systems call for rapprochement of modern control theory, machine learning, and optimization. Recent years have witnessed emerging trends of control theories with successful applications to en- gineering and scientific research, such as robotics (Krimsky & Collins, 2020), aerospace technology (He et al., 2019), and economics and management (Lapin et al., 2019) etc.
We consider the establstall formulation of optimal control (Kirk, 2004) in finite time horizon T = [40, 'r] Denote X and U as two vector-valued function sets, representing state functions and control functions respectively. Functions in X (resp. U) are defined over 7 and have their outputs in dd (resp. dd. State functions æ E X and control functions u E U governed by differential equation. The optimal control problem (OCP) is targeted at finding a control function that minimizes the cost functional f (Lions, 1992; Kirk, 2004; Vinter & Vinter, 2010; Lewis et al., 2012):
f(x, u) / p(x(t), u(t)) dt + h(*(t)))
min WEU
(1a)
s.t.
+(t) d(=((),u(()))
(Ib)
xtto) = Zo.
(1c)
where d is the dynamics of differential equations; p evaluates the cost alongside the dynamics and h evaluates the cost at the termination state <=ty); and zu is the initial state. We restrict our discussion optimal control problems, leaving the control problems in stochas-
to differential tic networks (Dai & Gluzman, 2022), inventory management (Abdolazimi et al., 2021), etc. out of the scope of this paper. The analytic solution of Eq. is usually unavailable, especially for complex dynamical systems. Thus, there has been a wealth of research towards accurate, efficient, and scal- able numerical OCP solvers (Rao, 2009) and neural network based solvers (Kiumarsi et al., 2017) in recent years. However, both classic and modern numerical OCP solvers are facing challenges, especially emerging in the big data era, which we briefly discuss follows.
1) Opacity of Dynamical Systems. Existing works (Böhme Frank, 2017a; Effati & Pakdaman, 2013; Jin et al., 2020) assume the dynamical systems priori and exploit their explicit forms to ease

---

Under review as a conference paper at ICLR 2023
Table I: Comparison of modern optimal control approaches. The proposed OptCtrlOp naturally
covers all the merits in the sense of performing a single-phase mippcng paradigm that does not rely on known system dynamics, and supports arbitrary input-domain queries.
Methods
Phase
Continuity
Dynamics
Reusability
Paradigm
Direct (Böhme Frank, 2017a)
Single
Discrete
Required
No
Iterative
Two-Phase (Hwang et al., 2022)
Two
Discrete
Dispensable!
Partial
Iterative
PDP (Jin al., 2020)
Single
Discrete
Required
No
Tterative
DP (Tassa et al., 2014)
Single
Discrete
Required
No
Iterative
OptCtriOP (ours)
Single
Continuous
Dispensable
Yes
One-pass
if phase-1 uses PINN loss (Wang et al., 2021a): required.
only phase- is reusable.
the optimization. However, the real-world dynamical systems can be unknown and hard to model. It raises serious challenges in data collection and system inference (Schmidt al., 2021; Ghosh al., 2021), where special care is required for error reduction.
2) Model Reusability. Model reusability conceptually measured by the capability of utilizing historical data when facing an unprecedented problem instance. Since solving an individual instance of Eq. from scratch expensive, a reusable model that can be well adapted to systems simi- lar forms is welcomed for practical usage. This point traces to the sensitivity analysis of optimal control (Oniki, 1973) yet absent in recent works.
3) Running Paradigm. As of typical paradigms adopted in optimization, numerical optimal control solvers traditionally use iterative methods before picking the control solution, and thus introducing multiplicative term regarding the iteration in the running time complexity. This sheds light on the high computational cost of solving a single optimal control problem.
4) Control Solution Continuity. Control functions defined continuous domain (typically time), although being intractable for numerical solvers. Hence resorting to discretization in the input domain gives rise to the trade-off in the precision of the control solution and the computational cost, as well as the truncation errors. While the discrete solution can give point-wise queries, learning a control function for arbitrary time queries is much more valued.
5) Running Phase. A two-phase model (Chen al., 2018; Wang al., 2021a; Hwang al., 2022) can (partially) overcome the above issues the cost of introducing an auxiliary dynamics inference phase. This thread of works first ap- proximates the state dynamics by a differentiable surrogate model and then, in its second phase, solves an optimization problem for control vari- ables (more explanation in Appx. B). However, the two-phase paradigm increases computational cost and manifests inconsistency between the two phases. A motivating example in Fig. 1 shows the two-phase paradigm leads crucial failures. When the domain of phase-2 optimization goes outside the training distribution in phase-1, this method might collapse.
10
2x
True Cost
100
10
10000
10500
Epechs
Phase-2 cost curves of two failed instances of two-phase control on Pendulum system. The control function gradually moves outside the training data distribution of phase- 1. As result, the control function converges w.r.t. the cost predicted by the surrogate model (blue), but diverges w.r.t. true cost (red).
Figure 1:
Table compares the methods regarding the above aspects. We propose an instance-solution oper- ator perspective for learning to solve OCPs, thereby tackling the issues above. The contributions are:
1. We propose the operator perspective and solve OCPs by learning direct mappings from OCPs to their solutions. The design holds the following merits. The system dynamics is implicitly learned during the training, which relies on neither any explicit form of the system nor the optimization process at test time. As such the operator can be reused and generalized to himilarly-formed OCPs without retraining, and such generalization ability even missing for learning- free solvers. The single-phase direct mapping paradigm avoids iterative processes with substantial speedup.
2

---

Under review as a conference paper at ICLR 2023
2. We theoretically validate the instance-solution mapping perspective by leveraging Pontryagin's Maximum Principle, thereby converting Eq. 1 to a boundary value problem. We instantiate a neural solver: OptCtrlOP (Optimal Control OPerator), and derive bounds on its approximation enror.
3. Experiments on both synthetic and real systems show that OptCtrlOp is versatile for various forms of OCPs. It achieves about 100x speedup against MLP baseline, and 10Kx speedup (on synthetic environments) over classical direct method solvers. It also generalizes well on both in- and out-of-distribution OCP instances.
Related Works. Most OCPs can not be solved analytically, and numerical methods are developed, which can be mainly categorized into three groups: direct methods, indirect methods, and dynamic programming. Our work bears resemblance to indirect methods. Due to page limitation, we leave the detailed discussion on related works to Appendix B.
2 METHODOLOGY
In this section, we present the instance-solution operator perspective for solving OCPs. Inspired by Pontryagin's 's Maximum Principle (PMP) (Pontryagin, 1987), OCPs (Eq. 1) can be converted to boundary value problem (BVP), the solution of which is typical function operator. The operator can further be approximated accurately by neural networks with theoretical guarantees. Hence, we propose OptCtrlOP, an end-to-end OCP neural solver that learns the underlying infiniteddimensional operator G governed by BVP. The input of operator G is the cost functional f. and the output the optimal control function u", i.e. G F + U. OptCtrlOP learns such operator from paired data by minimizing the error between predicted and optimal control, without explicit knowledge of the dynamics. Our theoretical analysis guarantees that for arbitrary small approximation error tolerance, there exists an instance of OptCtrlOP that satisfies the tolerance with bounded size and depth.
2.1 ASSANCE-SOLUTTR OPERATOR PERSPECTIVE OF OCP
This section presents the pre-condition of our model, novel perspective that converts OCPs into an operator that maps problem instances (defined by cost functionals and dynamics) into their solutions. The conversion has two steps: 1) convert OCP into BVP by PMP: 2) define the operator based on BVP, as explained below. Here we assume the dynamics are deterministic DE, and the derivation of stochastic DE (SDE) (Carius et al., 2022; Kishida & Ogura, 2022) is similar but based on Hamilton- Cecobi-Bellman (HJB) equation (Yong & Zhou, 1999), which is left for future work.
We elaborate how PMP converts OCP into BVP, following the derivation of the first-order necessary optimality conditions in the calculus of variations. To begin with, define the Hamiltonian H. of Eq. 1:
H(x, u, A) == p(a, u) + A d( (x, u).
(2)
Then suppose u* U is the optimal control function, and -X is the optimal state trajectory.
The PMP asserts that there exists a co-state (adjoint) function A : + Rd such that the following BVP satisfied (Pontryagin, 1987):
dynamics:
4² H(r"" u² X²),
(3a)
0x
δ
co-state:
X H(x u² ('))
(3b)
dx
optimal control:
0
H(w"", u" Xt),
(3c)
Du
boundary conditions:
X² dz h(a (1z)), x: (to) = Cenit
(3d)
A BVP is composed of set of DEs along with boundary conditions. The BVP defines an implicit mapping from cost functional and dynamics, both of which can be represented by functions, to the optimal control function, with state and co-state being intermediate variables.
Based on BVP, one can define an operator according to concrete problem settings or protocols. In the theoretic analysis of this paper, we assume that: the initial state Dinit and the dynamics X x U + R are unknown but constant (Hwang et al., 2022). Consequently, the state æ is uniquely determined by control u, thus the cost functional X x U RR can be rewritten
3

---

Under review as a conference paper at ICLR 2023
interchangeably fa U J mapping from cost functional to optimal control.
R. Under the protocol, the operator is defined G
F U,
In this way, have converted the OCPs into operator, which directly maps OCP instances to their solutions. Theoretically, the operator able to solve all OCPs of the same type (e.g. governed by the same dynamics d), thus highly reusable and efficient. Since in practice such non-linear operator can hardly have a closed-form expression, we propose the following neural model.
2.2 ARCHITECTURE OF THE PROPOSED NEURAL MODEL
From the operator perspective discussed above, constructing an OCP solver is equivalent to an op- erator learning problem. Inspired by DeepONet (Lu al., 2021), we propose Operator Learning based Control Network (OpiCtrIOP) to solve the OCP by learning the non-linear operator ç di- rectly, and theoretically estimate the error bounds. Specifically, we decompose OptCurOP into three components, as suggested in (Lanthaler et al., 2022):
1. Encoder. We define the following mapping as the encoder: & F IR" which converts the cost functional f into finite-dimensional vector e that can be processed by neural networks. For demonstration, we now let the encoder be an ad-hoc function that returns the target state Egaat- Such a simplified encoder is still effective, since the form of the cost functional assumed fixed, and the only parameter is the target state (Eq. 12, following Jin al. (2020)). We term such encoder as physical priors based. Otherwise, the encoder can be implemented by functional approximations, c.g. point-wise evaluation, basis expansion, or neural network. For ex- ample, in the Heating system (Appendix E.6) the encoder output is the coefficients of trigonometric polynomial.
Approximator. The encoded vector e mapped to another finite dimensional vector a by the approximator mapping defined as: A IRm + IR".
The approximator is implemented by neural network termed branch net B, which is ReLU activated multi-layer perceptron (MLP) (Goodfellow et al., 2016). The output vector a = ((e) will serve the coefficient vector of basis functions produced by the following reconstructor.
3. Reconstructor. It reconstructs the control function u by the mapping: R RR + U.
The reconstructor firstly constructs p + basis functions via another MLP named trunk net, + RP+1, a parametric mapping from time to basis function values. Note that the origi- nal DeepONet trunk net output is p-dimensional. we add an additional bias function dimension as suggested by (Lanthaler et al., 2022). Then, the control function u = R(a) is reconstructed by an affine combination of basis functions, where the coefficient is produced by the approximator:
R(a)= == To + a Tl:p-
(4)
Combining three components together, OptCtrlOp maps the cost functional to control function, : - U, and can be viewed as a composition of the above three mappings:
: : R a AoE.
(5)
The architecture of OptCtrlOP is summarized in Fig. 2 and the pseudo code of forward step is described Alg. 1 in Appendix A.
2.3 APPROXIMATION ERROR ESTIMATION
We give the estimation for the approximation error of the proposed model. The theoretic result guar- antees that there exists neural network instance of OptCurlPP architecture (Eq. 5) approximating the operator G arbitrary error tolerance. Furthermore, the size and depth of such a network are upper-bounded. The technical line of our analysis is partly inspired by (Lanthaler et al., 2022) which provides error estimation for DeepONets (Lu al., 2021).
The error of interest in OCP is the distance between the cost of the solution and the optimal cost. Let μ be the probability measure of cost functional, with E P2(F), where P2(F) is the set of probability measures with finite second moments. We assume the control dimension d, 1, and fa > 0, w.l.o.g. Then the approximation error OptCtrlOP can be defined as:
fa N()) fa 9(f) dn(s)-
=
(6)
fa G()
4

---

Under review as a conference paper at ICLR 2023
Encoder
Approximator A Branch Net
Physical Priors
Functional
R
R
fEF
Point-wise Evaluation
Neural Network
Reconstructor R
Trunk Net
Optimal Control
Time Index tET
RRP1
Eq.4
a (t) E Rd
The architecture of OptCtrlOP. The network takes 2 inputs: cost functional f and time index The input is processed by encoder, approximator and reconstructor, as introduced Sec- tion 2.2. And finally it outputs itt), the estimation of optimal control at time
Figure 2:
Following the decomposition described Sec. 2.2, the approximation error of OptCtrlOp can also be decomposed into three parts: 1) encoder error, 2) approximator error, and 3) reconstructor error.
The error of our simplified encoder zero, since there exists one-to-one mapping between the target state (encoder output) and the cost functional (encoder input). In other words, for given OptCtrlOP encoder E, there exists an inverse mapping E-1, such that 8-1 a Id. the encoder is implemented by functional approximations, e.g. point-wise evaluation (Lu et al., 2021), I), then the encoder error should be considered (Appx D).
For reconstructor R. its error &R is estimated by the mismatch between R and its approximate inverse mapping, projector P, weighted by push-forward measure Gyu(u) == u(g-1(u)).
&R =
P(u) - uulled (G+#) (u) (1RO
P
argmin CR, s.t. Po R Id.
Intuitively, such reconstructor error quantifies the information loss induced by R. An ideal with- out any information loss should be invertible, i.e. its optimal inverse P is exactly R-1, thus we have &R O.
Given the encoder and reconstructor, and denote the encoder output is R" the error & induced by approximator A is defined as the distance between the approximator output and the optimal coefficient vector, weighted push-forward measure ( ((E-
E = ( [IA(e)) - - (e))l (E (e)
(8)
With the definitions above, we can estimate the approximation error of our OpICtrlOP by the error of each of its components, as stated in the following theorem (see detailed proof in Appendix C.I):
Theorem (Decomposition of OptCtrIOP Approximation Error). Suppose the cost functional fa is Lipschitz continuous, with Lipschitz constant Lip(fa). Define constant C' sup JaoG(T
Lipifa)
The approximation error & (Eq. 6) of OptCtrlOp N R A E is upper-bounded by & < &R).
(9)
The estimation of reconstructor error & is analyzed in the previous work (Lanthaler et al., 2022), which gives detailed discussion of the error estimation of DeepONet, and the reconstructor error component is the same as that of OptCtrlOP. We cite the result to establish the following theorem:

---

Under review as a conference paper at ICLR 2023
Theorem 2 (Reconstructor Error (Lanthaler et al., 2022)). If 9 defines Lipschitz mapping G F + (T) for some 8 >> 0, > O, with Se [l]((f)) , S M < oc, then there exists a constans C = Css, M) > O, such that for any p E N, there exists a trunk net T T RR+1 (with bias term TO and the associated reconstruction R. _P- +P satisfies:
size(T) S Cp (1 (1+1og())))
depth(() S C (opp((t)=<((++og(())) (1+1og())))
8R. S Cp
(10)
Furthermore, the reconstruction R. and the optimal projection P satisfy Lip(R),, Lip(P) S 2
Note that size(-) ) is defined as the number of trainable parameters of a neural network, and depth(.) denotes the number of hidden layers. H is Sobolev space with degrees of regularity and I norm. The proof (omitted here) is based an observation that reconstructor minimized when trunk net outputs {Ti, .2...., } are linearly independent. The observation is consistent with intuition, since trunk net outputs are regarded as basis functions (Eq. 4).
....
Next, the error bound of the approximator will be presented. Recall that the approximator learns a mapping between vector space by MLP, whose error estimation well studied in deep learning theory. One of the existing works (Gühring al., 2020) derives the estimation based on the Sobolev regularity of the mapping. We extend their result to form the following theorem (proof given in Appendix C.2):
Theorem 3 (Approximator Error). Given operator 9 F U, encoder F + Rm and recon- structor R E" + U. let P denote the corresponding projector. for some E Nz2. M > 0, the following bound is satisfied: (Pog. o - ((())))))) S M < co, then there exists constant
C' = =C(mm, s, M) > O, such that for any E (0, 1)) there exists an approximator A: Rm + R':
Aize(A) S Cppc- log 1
depth(A) S Colog (c-
E S VTc.
(11)
In summary, we have proved that the approximation error is bounded by the sum of the reconstructor and approximator errors. Those errors can hold under arbitrary small tolerance, with bounded size and depth. The theorems hold as long as some integratable and continuous conditions are satisfied, which is trivial in real-world OCPs.
3
EXPERIMENTS
3.1 SYNTHETIC CONTROL SYSTEMS
3.1.1 CONTROL SYSTEMS AND DATA GENERATION
We evaluate OptCtrlOP on six representative optimal control systems by following the same protocol of (Jin et al., 2020; Hwang et al., 2022), as summarized in Table 4. We postpone the rest systems to Appendix E, and only describe the details of the Quadrotor control here:
1 q 1= )lWWaa
p m
0 mgg
R (g)
0
Jō = Tu- w x Jw.
1 2
This system describes the dynamies of helicopter with four rotors. The state x 1p TT w T Iº consists of parts: position p. velocity v, and angular velocity W. The control E i is the thrusts the four rotating propellers of the quadrotor. q E IR4 is the unit quaternion (Jia, 2019) representing the attitude(spacial rotation) of quadrotor w.r.t. the inertial frame. J is the moment of inertia in quadrotor's s frame, and Tu is the torque applied to the quadrotor. Our setting is similar to (Jin et al., 2020, Appx. E.I,. but we exclude the quaternion from the state. We set the initial state Zinit = [-8, - 6, 9, 0 T the initial quaternion Jinst = O. The matrices D(w), R(q), T are coefficient matrices, see definition in Appx. E.4. The cost functional is defined as following, with
coefficients ee = 1, C. 0.1.
min / ['c((ett) cT(((()) xgoat)² 2 + c.. cull(x(t)]]² dt
(12)

---

Under review as a conference paper at ICLR 2023
33
1077 Reddkkmhhokknn Cartask Quadstor Roclet
Pandakurs Grrihile Quadroter leddet
Pandhfamm RobotAre Carthile Guadertrr Rocist
(a) Running Time (s/problem)
(h) ID MAPE
(c) OOD MAPE
Figure 3: The rumning time and mean absolute percentage (MAPE) on in-distribution (ID) and out-of-distribution (OOD) benchmark problem sets. Compared with baselines, OptCtrlOP (red bars) achieves higher comparable accuracy, with over 100x speedup.
With the settings above given, the solution of OCP only depends on the target state good. There- fore, we generate datasets (for model training/validation) and benchmarks (for model testing) by sampling target states from pre-defined distribution. To fully evaluate the generalization ability, we define both in-distribution (ID) and distrf-distriuution (OOD) (Shen et al., 2021). Specifically, we design two random variables, Xgoot := base + Cins and Xooa == a base + Couts where a: achase baseline goal state, and out are different noise applied ID and OÔD. In Quadrotor problems for example, we set goal 0, and uniform noise Ein ~ ((.11, 1.1), and Eout -U(-0.1,0.1)
gool is
The training data are sampled from ID, while validation data and benchmark sets are both sam- pled from ID and OOD separately. The data generation process is shown in Alg. in Appendix. For given distribution,
sample target state gooat, and construct the corresponding cost func- tional f and OCP. Then define 100 time indices uniformly spaced in time horizon = [0, tf], If - (11 1.01). The length of T is slightly perturbed that time indices fall in the whole horizon, instead of fixed points. Then we solve the resulting OCP by the DM solver, and get the optimal con- trol u² at those time indices. After that we sample 10 time indices {ti]isigio, creating 10 triplets [(f,t, u² (t:)))ssggio and adding them to the dataset. Repeat the process, until the size meets the requirement. The benchmark set is generated in the same way, but we store (f, Jopt) pair (Jopt is the optimal cost) for each problem instead.
3.1.2 IMPLEMENTATION AND BASELINES
For all systems and all neural models, the learning rate starts from 0.01, decaying every 1,000 epochs at a rate of 0.9, The batch size is 10,000, and the optimizer is Adam (Kingma & Ba, 2014). The loss is the mean squared error defined below, with a dataset of N samples: L= NV Eijen [t)) wz [4 For comparison, choose the following baselines. Other details of implementation and baselines are recorded in Appendix F.
1) Direct Method (DM): classical direct OCP solver, with Interior Point OPTimizer (IPOPT) (Biegler& & Zavala, 2009) backend NL.P solver.
2) Pontryagin Differentiable Programming (PDP) (Jin et al., 2020): an adjoint-based indirect method, differentiating through PMP, and optimized by gradient descent.
3) Multi-layer Perceptron (MLP): a fully connected counterpart (Alg. 4) OptCtrlOP.
4) Fourier Neural Operator (FNO) (Li et al., 2020): A neural operator consists of consecutive Fourier transform layers.
5) Graph Element Network (GEN) (Alet et al., 2019): A graph neural operator with graph convo- lution backbone.
3.1.3 RESULTS AND DISCUSSION
We present the numerical results the six systems to evaluate the efficiency and accuracy of OptC- trIOP. The metrics of interest are 1) the running time of solving problems; 2) the quality of solution, measured by mean absolute percentage emrrr (MAPE) between the true optimal cost and the pre- dicted cost, which defined as the mean of [JJppt - Juat)/ Jopels where Jopt is the optimal cost generated by DM (regarded as ground truth), and Juot is the cost of the solution produced by the
7

---

Under review as a conference paper at ICLR 2023
Table 2: Performance of Quadrotor environment.
Model
Time (sec./instance)
ID MAPE OOD MAPE
ID-OOD Gap
DM (Böhme & Frank, 2017a)
9.23 10-2
OptCtrlOp
7.79 x 10-7
1.74 x 10
2.13 x 10
1.96 x 10-1
MLP (Alg. 4)
3.04 x 10-6
4.12 x 10-4
6.09 x 10-4
1.97 x 10-1
GEN (Alet al., 2019)
6.15 10
1.40
10
7.31
10-4
5.91 10-4
FNO (Li et al., 2020)
6.38 x 10-4
9.87
10
1.28 X 10-3
1.27 X 10-3
PDP (Jin et al., 2020)
7.25 x 10a
1.24 X 10-1
1.16 X 10
8.74 X 10-6
101
10a
10
GEN
20-1
PNO
10
10"4
100
200
300
400
500
100
200
500
600
500
100
200
300
400
500
100
200
300
400
500
(a) Loss
(b) MAPE
Figure 4: (a) Quadrotor Loss curve. (b) Quadrotor cost MAPE (mean absolute percentage error) curve. For ID (train and ID_tes (Dttes)) OptCtrlOP (red curves) performs competitively, or better than others. It also outperforms all neural baselines on OOD.
model. The MAPE is calculated on ID/OOD benchmarks respectively, and the running time is aver- aged for 2,000 random problems. The results of consttrnid OCP are visualized in Fig. 3.
First, the comparison of running time shown Fig. 3a. which shows that the neural operator solver is much faster than the classic solver. For example, OptCtrIOP achieves over 10 times speedup against the DM solver. The acceleration can be reasoned in two aspects: 1) the neural operator solvers produce the output by single forward propagation, while the classic methods need to iterate between forward and backward pass multiple times; 2) the neural solver calculation is highly paralleled. Moreover, OptCtrlOP is 100 times faster than MLP, although both of them are neural operator models. The difference rooted in OptCtrIOP architecture, where two inputs (i.e. cost functional and time index) are processed by branch net and trunk net separately. Consequently, OptCtrOPP runs O (M + .N) forward passes, if the benchmark set consists of N problems, and the time horizon is discretized into M indices. In comparison, MLP entangles requiring O(MN) forward passes. The difference is also given in Alg. 3-4 in Appendix. The FNO follows similar diagram as MLP, but 10+ times slower than MLP, since involves computationally intensive Fourier transformations. The GEN processes two inputs in a disentangled manner, similar that of OptCtrlOP. But GEN is 100+ times slower than OptCtr1OP, probably because of the intrinsic complexity of graph construction and graph convolution.
inputs together,
Next, the accuracy on in- and out-of-distribuutin benchmark sets compared in Fig. 3b-3c. Com- pared with other neural models, OptCtrlOP achieves better or comparable accuracy in general. In addition, OptCtrlOp outperforms classical PDP on more than half of the benchmarks. crete example, we investigate the performance of Quadrotor environment. From Tab. 2 (and Fig. 4 for learning trajectory of neural models), one can observe that MAPE of OptCtrlOP on ID the second lowest among all models, with slight disadvantage compared to FNO. On OOD, however, OptCtrlOP outperforms FNO by a clear margin, and the performance is close to that of the classical method PDP. Among all neural models, OptCtrlOP achieves the lowest OOD MAPE well as the smallest ID-OOD gap (defined as the absolute distance between ID and OOD MAPE). We conjec- ture that the OOD generalization ability of OptCtrlOp results from its architecture, where branch net and trunk net (coefficients and basis functions) are explicitly disentangled. Such a structure may inherit the inductive bias from numerical basis expansion methods (c.g. Kafash al. (2014)), thus being more robust to distribution shifts. Due to space limitations, the numerical results of other systems are given in Tab. 6-10 Appendix.
a con-
3.2 REAL-WORLD CONTROL DATASET OF PLANAR PUSHING
8

---

Under review as a conference paper at ICLR 2023
Table 3: Results of Pushing environment.
Model
Time(sec./instacce)
ID MAPE
OOD MAPE
OptCtrlOP
6.07 x 10-6
0.109
0.136
MLP
8.85 x 10-1
0.128
0.149
GEN
6.59 x 10-5
0.201
0.207
FNO
1.34 x 10-3
0.135
0.144
This section presents how to learn optimal control of robot arm for pushing objects of varying shapes var- ious planar surfaces. We use the Pushing dataset Yu et al. (2016), challenging dataset consisting of noisy, real-world data produced by an ABB IRB 120 industrial robotic arm (Fig.5, right part). The robot arm is controlled to push objects starting from various contact positions and 9 angles (initial state), along different trajectories (target state functions), with 11 object shapes and 4 surface ma- terials (dynamics). The control function is represented by the force exerted on to object, measured by the force sensor. Left of Fig.5 gives a compact overview of input variables, and more details are given in Appx. E.8.
trajectory
robot arm
s
friction
PT sensor
pusher
shape
object
Figure 5: Pushing environment, com- posed of images from Yu et al. (2016).
In our experiment, we apply OptCtrlOp to learn map- ping from a pushing OCP instance (represented by vari- ables above) to the optimal control function. The input now is no longer cost functional f only, but fa with abuse of notation, where subscript d denotes dynamics and initial conditions. And the encoder is realized by different techniques for different inputs, such as Savitzk--Golay smoothing (Savitzky & Golay, 1964) and down- wwwplng for trajectories, mean and standard value extraction for friction map, and Convolution (etWork(CNN) (LeCun et al., 1989) for shape images. The en- coder error now is not negligible, but the analysis framework of approximation error can be extended to include this error (Appx. D). The estimation of encoder error itself depends on the specific choice of the encoder, which is beyond the scope of this paper and left future work.
We extract training data from ID, validation and test data from both ID/OOD, The ID/OOD is dis- tinguished by different initial contact positions. The accuracy metric MAPE is now defined as || - [1/l|=""] |]. We compare OpiCtrIOP performance only with neural baselines, since the explicit expression of pushing OCP is unavailable. All neural models share the same encoder structure, while the parameters of CNN are trained end-to-end individually for each model. The results are displayed Tab. 3, from which one can observe OptCtrlOP outperforms all baselines in both running time and ID/OOD accuracy. Notice that the performance of FNO and GEN degrades compared with that of synthetic data. The reason might be that their architecture is not suitable for complex OCP tasks like pushing. For example, the essence of FNO is to learn parametric transformation in Fourier space. For the pushing dataset, however, the input OCP instance a composition of several functions and environment parameters, of which the Fourier transform does not have clear physical meaning.
4 CONCLUSION
Future Works. This paper studies an effective way to solve general OCPs with data-driven ap- proaches. We do not specify the forms of problem instances or investigate sophisticated models for specific problems, which calls for careful designs and exploitation of the problem structures. We also leave rigorous analysis OOD performance and encoder error future works.
Conclusion. We have proposed a novel instance-solution operator perspective of OCPs, where the operator directly maps cost functionals to optimal control functions. Based on this perspective, present a neural operator OptCtrIOP, with a theoretic guarantee on its approximation capability. Extensive experiments
various OCP system benchmarks show the outstanding generalization ability and efficiency of OptCurOPP, on both ID and OOD settings. We envision the proposed model will be beneficial in solving numerous high-dimensional problems in the learning and control fields.

---

Under review as a conference paper at ICLR 2023
REFERENCES
Omid Abdolazimi, Davood Shishebori, Fariha Goodarzian, Peiman Ghasemi, and Andrea Appol- loni. Designing new mathematical model based on abc analysis for inventory control problem: A real case study. RAIRO-operations research, 55(4):2309-2335, 2021.
Asma Al-Tamimi, Frank L Lewis, and Murad Abu-Khalaf. Discrete-time nonlinear hjb solution using approximate dynamic programming: Convergence proof. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics), 38(4):943-949, 2008.
Ferran Alet, Adarsh Keshav Jeewajee, Maria Bauza Villalonga, Alberto Rodriguez, Tomas Lozano- Perez, and Leslie Kaelbling. Graph element networks: adaptive, structured computation and memory. In International Conference on Machine Learning, pp. 212-222. PMLR, 2019.
Anima Anandkumar, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Nikola Kovachki, Zongyi Li, Burigede Liu, and Andrew Stuart. Neural operator: Graph kernel network for partial differ- ential equations. In ICLR 2020 Workshop on Integration of Deep Neural Models and Differential Equations, 2020.
Joel A E Andersson, Joris Gillis, Greg Horn, James B Rawlings, and Moritz Diehl. CasADi A software framework for nonlinear optimization and optimal control. Mathematical Programming Computation, 11(1):1-36, 2019. doi: 10. 10.1007/>12332-018-4
Mokhtar Bazaraa, Hanif Sherali, and Chitharanjan M Shetty. Nonlinear programming: theory and algorithms. John Wiley & Sons, 2013.
Richard Bellman and Robert Kalaba. Dynamic programming and adaptive processes: mathematical foundation. IRE transactions on automatic control, (1):5-10, 1960.
Lorenz T Biegler and Victor M Zavala. Large-scale nonlinear programming using ipopt: An integrat- ing framework for enterprise-wide dynamic optimization. Computers & Chemical Engineering, 33(3):575-588 2009.
Hans Georg Bock and Karl-Josef Plitt. A multiple shooting algorithm for direct solution of optimal control problems. IFAC Proceedings Volumes, 17(2):1003-1088 1984.
Paul T Boggs and Jon W Tolle. Sequential quadratic programming. Acta numerica, 4:1-51, 1995.
Thomas J. Böhme and Benjamin Frank. Direct Methods for Optimal Control, pp. 233-273. Springer International Publishing, Cham, 2017a. ISBN 978-3-319-5131711. doi: 10.
URL https: //doi. org/10 1107/998--33199-5555
Thomas J. Böhme and Benjamin Frank. Indirect Methods for Optimal Control, pp. 215-231. Springer International Publishing, Cham, 2017b. ISBN 978-3-319-51317-1
doi: 10.1007/
978-3-319-51317-1 URL
1007 11w.11107977---11---11
Greg Brockman, Vicki Cheung. Ludwig Pettersson, Jonas Schneider, John Schulman, Jie Tang, and Wojciech Zaremba. Openai gym, 2016.
Jan Carius, René Ranftl, Farbod Farshidian, and Marco Hutter. Constrained stochastic optimal control with learned importance sampling: A path integral approach. The International Journal of Robotics Research, (1(2):18922)) 2022.
Tianping Chen and Hong Chen. Universal approximation to nonlinear operators by neural networks with arbitrary activation functions and its application to dynamical systems. IEEE Transactions on Newral Networks, 6(4):911-917,, 1995.
Yize Chen, Yuanyuan Shi, and Baosen Zhang. Optimal control via neural networks: A convex approach. arXiv preprint arXiv: 11835, 2018.
Lin Cheng. Zhenbo Wang, Yu Song, and Fanghua Jiang. Real-time optimal control for irregular asteroid landings using deep neural networks. Acta Astronautica, 170:66-79, 2020.
Jim G Dai and Mark Gluzman, Queueing network controls via deep reinforcement learning. Stochas- fie Systems, 12(1) 2(11):0-67, 2022.
10

---

Under review as a conference paper at ICLR 2023
Andrea D' a'mbbooio,, Enrico Schiassi, Fabio Curti, and Roberto Furfaro. Pontryagin neural networks with functional interpolation for optimal intercept problems. Mathematics, 99))999, 2021.
Sohrab Effati and Morteza Pakdaman. Optimal control problem via neural networks. Neural Com- puting and Applications, 23(7):2093-2100, 2013.
Wendell Fleming and Raymond Rishel. Deterministic and stochastic optimal control, vol- ume 1. Springer Science & Business Media, 2012.
Sanmitra Ghosh, Paul Birrell, and Daniela De Angelis. Variational inference for nonlinear ordinary differential equations. In International Conference on Artificial Intelligence and Statistics, pp. 2719-2727. PMLR, 2021.
Ian Goodfellow, Yoshua Bengio, and Aaron Courville. Deep learning. MIT press, 2016.
Ingo Gühring, Gitta Kutyniok, and Philipp Petersen. Error bounds for approximations with deep relu neural networks in p norms. Analysis and Applications, 18(05):803-859, 2020.
Xiaoxiao Guo, Wei Li, and Francesco lorio. Convolutional neural networks for steady flow ap- proximation. In Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining, pp. 481 1-49, 2016.
Xiaolong He, Jichao Li, Charles A Mader, Anil Yildirim, and Joaquim RRA Martins, Robust aero- dynamic shape optimization- from a circle airfoil. Aerospace Science and Technology, 87: 48 61, 2019.
Lukas Hewing, Kim P Wabersich, Marcel Menner, and Melanie N Zeilinger. Learning-based model predictive control: Toward safe learning in control. Annual Review of Control, Robotics, and Autonomous Systems, 3:269 296, 2020.
Rakhoon Hwang, Jae Yong Lee, Jin Young Shin, and Hyung Ju Hwang. Solving pde-constrained control problems using operator learning. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 36, pp. 4504 4512, 2022.
Yan-Bin Jia. Quaternions. Com S, 477:577, 2019.
Wanxin Jin, Zhaoran Wang, Zhuoran Yang, and Shaoshuai Mou. Pontryagin differentiable pro- gramming: An end-to-end learning and control framework. Advances Neural Information Processing Systems, 33:7979-7992 2020.
B Kafash, A Delavarkhalafi, and SM Karbassi. A numerical approach for solving optimal control problems using the boubaker polynomials expansion scheme. Interpolat. Approx. Sci. Comput, 3:1-18, 2014.
Yuehaw Khoo, Jianfeng Lu, and Lexing Ying. Solving parametric pde problems with artificial neural networks. European Journal of Applied Mathematics, 32(3):421-435, 2021.
Diederik P Kingma and Jimmy Ba. Adam: A method for stochastic optimization, arXiv preprint arXiv:1412.6988. 2014.
Thomas N Kipf and Max Welling, Semi-supervised classification with graph convolutional net- works. arXiv preprint arXiv: www.li.09909900 2016.
Donald E Kirk. Optimal control theory: an introduction. Courier Corporation, 2004.
Masako Kishida and Masaki Ogura. Temporal deep unfolding for constrained nonlinear stochastic optimal control. IET Control Theory & Applications, 16(2):139-150, 2022.
Bahare Kiumarsi, Kyriakos G Vamvoudakis, Hamidreza Modares, and Frank L Lewis. Optimal and autonomous control using reinforcement learning: A survey. IEEE transactions neural networks and learning systems, 29(6):042-2208 2017.
Nikola Kovachki, Samuel Lanthaler, and Siddhartha Mishra. On universal approximation and error bounds for fourier neural operators. Journal of Machine Learning Research, Art-Nt, 2021.

---

Under review as a conference paper at ICLR 2023
Erez Krimsky and Steven H Collins. Optimal control of eeergy-recyclln actuator for mo- bile robotics applications. In 2020 IEEE International Conference on Robotics and Automation (ICRA), pp. 3559-3565 IEEE, 2020.
Samuel Lanthaler, Siddhartha Mishra, and George E. Karniadakis. Error estimates for deeponets: A deep learning framework in infinite dimensions. Transactions of Mathematics and Its Applica- tions, 6 6(1)tmac001, 2022 2022.
Alexander Lapin, Shuhua Zhang, and Sergey Lapin. Numerical solution of parabolic optimal control problem arising in economics and management. Applied Mathematics and Computation, 361:715-729, 2019.
A Lasota. A discrete boundary value problem. In Annales Polonici Mathematici, volume 20, pp- 183-190. Institute of Mathematics Polish Academy of Sciences, 1968.
Yann LeCun, Bernhard Boser, John Denker, Donnie Henderson, Richard Howard, Wayne Hubbard, and Lawrence Jackel. Handwritten digit recognition with back- pack-propagation network. Advances in neural information processing systems, 2, 1989.
Frank Lewis, Draguna Vrabie, and Vassilis Syrmos. Optimal control. John Wiley & Sons, 2012.
Weiwei Li and Emanuel Todorov. Iterative linear quadratic regulator design for nonlinear biological movement systems. In ICINCO (I), pp. 222-229. Citeseer, 2004.
Zongyi Li. Nikola Borislavov Kovachki. Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, Anima Anandkumar, et al. Fourier neural operator for parametric partial differential equa- tions. In International Conference Learning Representations, 2020.
Marten Lienen and Stephan Günnemann. Learning the dynamics of physical systems from sparse observations with finite element networks. In International Conference on Learning Representa- fions, 2021.
PL Lions.
Optimal control. In ICIAM 91: Proceedings of the Second International Conference Industrial and Applied Mathematics, volume 61, pp. 182. SIAM, 1992.
Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis. Learning nonlinear operators via deeponet based the universal approximation theorem of operators. Nature Machine Intelligence, 3(3):218-22) 2021.
Frederik Baymler Mathiesen, Bin Yang, and Jilin Hu. Hyperverlet: symplectic hypersolver for hamiltonian systems. In Proceedings of the AAAI Conference Artificial Intelligence, vol- ume 36, pp. 4575-4582, 2022.
Sanjay Mehrotra. On the implementation of a primal-dual interior point method. SIAM Journal on optimization, 2(4):575-601, 1992.
Julian Nubert, Johannes Köhler, Vincent Berenz, Frank Allgöwer, and Sebastian Trimpe. Safe and fast tracking on a robot manipulator: Robust mpc and neural network control. IEEE Robotics and Automation Letters, 5(2):3050-3057, 2020.
Hajime Oniki. Comparative dynamics (sensitivity analysis) in optimal control theory. Journal of Economic Theory, 6(3):265-283, 1973.
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, high- performance deep learning library. In Advances in neural information processing systems, 2019.
Michael Poli, Stefano Massaroli, Federico Berto, Jinkyoo Park, Tri Dao, Christopher Re, and Ste- fano Ermon. Transform once: Efficient operator learning in frequency domain. In ICML 2022 2nd Al for Science Workshop, 2022.
Lev Semenovich Pontryagin. Mathematical theory of optimal processes. CRC press, 1987.
Antonin Raffin.
RI basclines3
zoo.
https://github. com DDR-RDLRRRRR RM/
rl-basel - 2020.
12

---

Under review as a conference paper at ICLR 2023
Antonin Raffin, Ashley Hill, Adam Gleave, Anssi Kanervisto, Maximilian Ernestus, and Noah Dormann. Stable-baselines Reliable reinforcement learning implementations. Journal of Machine Learning Research, 22(268):1-88 2021. URL http::/jmlz... org/ wttp://.....papps///222. papers/v22/ 20- 1364. html.
Maziar Raissi, Paris Perdikaris, and George Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal of Computational physics, 378:686-707, 2019.
Maziar Raissi, Alireza Yazdani, and George Em Karniadakis. Hidden fluid mechanics: Learning velocity and pressure fields from flow visualizations. Science, (67(118)) 1030, 2020.
Anil Rao. A survey of numerical methods for optimal control. Advances the Astronautical Sciences, 135(1):497-528, 2009.
Abraham Savitzky and Marcel JE Golay, Smoothing and differentiation of data by simplified least squares procedures. Analytical chemistry, 36(8):1627-1639, 1964.
Jonathan Schmidt, Nicholas Krämer, and Philipp Hennig. A probabilistic state space model for joint inference from differential equations and data. Advances in Neural Information Processing Systems, 34: 2214-1 1385,, 2021.
John Schulman, Filip Wolski, Prafulla Dhariwal, Alec Radford, and Oleg Klimov. Proximal policy optimization algorithms. arXiv preprint arXiv:1707.06347, 2017.
Zheyan Shen, Jiashuo Liu, Yue He, Xingxuan Zhang, Renzhe Xu, Han Yu, and Peng Cui. Towards out-of cuttrf--dittibution generalization: A survey. arXiv preprint arXiv:2 10..2/08..300000 3324,, 2021.
Justin Sirignano and Konstantinos Spiliopoulos. Dgm: A deep learning algorithm for solving partial differential equations. Journal of computational physics, 375:1339-1364 2018.
Mark W Spong, Seth Hutchinson, and Mathukumalli Vidyasagar. Robot modeling and control. John Wiley & Sons, 2020.
Yuval Tassa, Nicolas Mansard, and Emo Todorov. Control-limited differential dynamic program- ming. In 2014 IEEE International Conference on Robolics and Automation (ICRA), pp. 1168- 1175. IEEE, 2014.
Russ Tedrake. Underactuated Robotics. 2022. URL URLLtttp:/underacctuttet
edu.
Fredi Tröltzsch. Optimal control of partial differential equations: theory, methods, and applications, volume 112. American Mathematical Soc., 2010.
Clifford A Truesdell. A First Course in Rational Continuum Mechanics VI. Academic Press, 1992.
George E Uhlenbeck and Leonard S Ornstein. On the theory of the brownian motion. Physical review, 36(5):823, 1930.
Richard Vinter and RB Vinter. Optimal control. Springer, 2010.
Sifan Wang, Mohamed Aziz Bhouri, and Paris Perdikaris. Fast pde-constrained optimization via self-supervised operator learning. arXiv preprint 10. 13297, 2021a.
Sifan Wang, Yujun Teng, and Paris Perdikaris. Understanding and mitigating gradient flow patholo- gies in physics-informed neural networks. SIAM Journal on Scientific Computing, 43(55: A3081, 2021b.
Shaobo Xie, Xiaosong Hu, Shanwei Qi, and Kun Lang. An artificial neural network-enhanced energy management strategy for plug-in hybrid electric vehicles. Energy, 163:837-848, 2018.
Dongbin Xiu and Jan S Hesthaven. High-order collocation methods for differential equations with random inputs. SIAM Journal on Scientific Computing. 1189)1118- 2005.
Jiongmin Yong and Xun Yu Zhou. Stochastic controls: Hamiltonian systems and HJB equations, volume 43. Springer Science & Business Media, 1999.
13

---

Under review as a conference paper at ICLR 2023
Bing Yu et al. The deep ritz method: deep learning-basse numerical algorithm for solving varia- tional problems. Communications in Mathematics and Statistics, 6(1):1-12, 2018.
KT Yu, M Bauza, N Fazeli, and A Rodriguez. More than a million ways to be pushed. A High- Fidelity Experimental Data Set of Planar Pushing In: IEEE/RSJ IROS, 2016.
Yinhao Zhu and Nicholas Zabaras. Bayesian deep convolutional encoder-decoder networks for surrogate modeling and uncertainty quantification. Journal of Computational Physics, 366:415- 447, 2018.
Yinhao Zhu, Nicholas Zabaras, Phacdon-Stelion Koutsourelakis, and Paris Perdikaris. Physics- constrained deep learning for high-dimensional surrogate modeling and uncertainty quantification without labeled data. Journal of Computational Physics, 394:56 81, 2019.
14

---

Under review as a conference paper at ICLR 2023
A
ALGORITHMS
Algorithm 1: OptCtrlOP forward step (single sample version)
Input: cost functional and time index
Output: estimated control u(4) = (N((4)
< /// Encoder
< A(e) :// Approximator
r < T(4): // trunk net (Reconstructor step 1)
a(1.) +- -ro- Ttppp/ affine combination (Reconstructor step 2)
return ú(t;)
Algorithm 2: Data generation
Input:
Distribution of target state Xgoal
Output: Dataset
N < Number of samples per trajectory;
Dataset < empty set;
while not Dataset is full do
sample a goal from Xgool; construct cost functional f (e.g. Eq. 12) with Egoat;
construct OCP (Eq. 1) with cost functional f;
u" OC So..Soove(OOP)) // Any solver is applicable. We choose DM.
sample [4itisss from time horizon 7 of OCP;
add triplets {(f.t, tg u (t.))ssigN to Dataset;
return Dataser
Algorithm 3: OptCtrlOP inference on benchmark set
Input: Benchmark set of cost functionals F (fylisn, time indices t = (t, Jig
Output: estimated controls Ú = (()(4s))
E < E(F); // Encoder Approximator, O(N)
A < A(E):
R r(t); // trunk net (Reconstructor step 1), O(M)
for j t to N do
L (t) + Ro + A] Ripi // affine combination (Reconstructor step 2)
return Ú
B BACKGROUND AND RELATED WORK B. OCP SOLVERS
Traditional numerical solvers well developed over the decades, which learning-free and often involve tedious optimization iterations to find an optimal solution.
Direct methods (Böhme & Frank, 2017a) reformulate OCP inite-dimensional nonlinear pro- gramming (NLP) (Bazaraa et al., 2013), and solve the problem by NLP algorithms, e.g. sequential quadratic programming (Boggs & Tolle, 1995) and interior-point method (Mehrotra, 1992). The reformulation essentially constructs surrogate models, where the state and control function (infinite dimension) is replaced by polynomial or piece-wise constant functions. The dynamics constraint is discretized into equality constraints. The direct methods optimize the surrogate models, thus the so- lution is not guaranteed to be optimal for the origin problem. Likewise, typical direct neural solvers (Chen al., 2018; Wang al., 2021a: Hwang al., 2022), termed as Two-Phase models, consist of two phases: 1) approximating the dynamics by a neural network (surrogate model); 2) solving the NI J. via gradient descent, by differentiating through the network. The advantage of Two-Phase
15

---

Under review as a conference paper at ICLR 2023
Algorithm 4: MLP inference on benchmark set
Input: Benchmark set of cost functionals F (ffiSN, time indices t = {ts (ts
Output: estimated controls U = = (t)))ss
for 4- 1 to .M do
for j - 1 to N do
ej < 8()):// Encoder
(4) - MLP(e,;4): // repeat O(MN) times
return
models against traditional direct methods is computational efficiency, especially in high-dimensional cases. However, the two-phase method is sensitive to distribution shift (see Fig. 1).
Indirect methods (Böhme & Frank, 2017b) are based on Pontryagin's Maximum Principle (PMP) (Pontryagin, 1987). By PMP, indirect methods convert OCP (Eq. 1) into a boundary-value problem (BVP) (Lasota, 1968), which is then solved by numerical methods such as shooting method (Bock & Plitt, 1984), collocation method (Xiu & Hesthaven, 2005), adjoint-basee gradient descend (Effati & Pakdaman, 2013; Jin et al., 2020). These numerical methods are sensitive to the ini- tial guesses of the solution. Some indirect methods based neural solvers approximate the finite- dimensional mapping from state x (t) E Rd= to control u (t) E Rd (Cheng et al., 2020), or to co-state X (t) E R (Xie et al., 2018). The full trajectory of the control function is obtained by repeatedly applying the mapping and getting feedback from the system, and such sequential nature is the efficiency bottleneck. Another work (D' a'mbroooi et al., 2021) proposes to solve the BVP via a PINN, thus its trained network works only for one specific OCP instance. Distinctive from all operator that maps cost
these models, OptCtrlOP solves BVP by learning an functional f E to control function u E U. The trained model is available for parallel queries on different OCP instances, with high efficiency and accuracy.
Dynamic programming (DP) is an alternative, based on Bellman's principle of optimality (Bell- man & Kalaba, 1960). It offers rule divide high-dimensional optimization problem with long time horizon into smaller, easier-to-solve auxiliary optimization problems. Typical methods are: Saailtoon-Jaoob--ellllaa (HJB) equation (Al-Tamimi et al., 2008), differential dynamical pro- gramming (DDP) (Tassa et al., 2014), which assumes quadratic dynamics and value function, and iterative linear quadratic regulator (iLQR) (Li & Todorov, 2004), which assumes linear dynamics and quadratic value function. Similar to dynamic programming, model predictive control (MPC) synthesizes the approximate control function via the repeated solution of finite overlapping hori- the number
zons (Hewing et al., 2020). The main drawback of DP is the curse of dimensionality and complexity of the auxiliary problem. MPC alleviates this problem at the expense optimality. Yet fast implementation of MPC is still under exploration and remains open (Nubert et al., 2020).
B.2 DIFFERENTIAL EQUATION NEURAL SOLVERS
A variety of networks have been developed to solve DE, including Physics-informed neural networks (PINNs) (Raissi et al., 2019), neural operators (Lu et al., 2021), hybrid models (Mathiesen et al., 2022; Lienen & Günnemann, 2021), and frequency domain models (Li et al., 2020; Poli et al., 2022). We will briefly introduce the first two models for their close relevance to our work.
PINNs parameterize the DE solution as neural network, and learn the parameters by minimizing the residual loss and boundary condition loss (Yu et al., 2018; Raissi et al., 2019; Sirignano & Spiliopoulos, 2018). PINNs are similar to those numerical methods e.g. the finite element method in that they replace the linear span of a finite set of local basis functions with neural networks. PINNs usually have simple architectures (e.g. MLP), although they have produced remarkable results across a wide range of areas in computational science and engineering (Raissi et al., 2020; Zhu et al., 2019). However, these models are limited to learning the solution one specific DE instance, but not the operator. In other words, if the coefficients of the DE instance slightly change, then new neural network needs to be trained accordingly, which is time-consuming Another major drawback of PINNs, as pointed out by (Wang al., 2021b). is that the magnitude of two loss terms (i.e.residual loss and boundary condition loss) is inherently imbalanced, leading to heavily biased predictions even for simple linear equations.
16

---

Under review as a conference paper at ICLR 2023
Neural operators regards DE operator that maps the imput to the solution. Learning operators using neural networks was introduced in the seminal work (Chen & Chen, 1995). It proposes the universal approximation theorem for operator learning, i.c. network with single hidden layer can approximate any nonlinear continuous operator. Lu et al. 12021) follows this theorem by designing a deep architecture named DeepONet, which consists of two networks: branch net for input functions and a trunk net for the querying locations in the output space. We choose this architecture as our OCP operator learner, and our analysis of optimal control error bound is partly inspired by Lanthaler et al. (2022), providing error estimation of DeepONet. Note that DeepONet is designed to solve DE operators, and not ready handle optimization tasks e.g. OC.
In addition, there exist many neural operators with other architectures. For example, another type of neural operator to parameterize the operator as a convolutional neural network (CNN) between the fnite-dimensional data meshes (Guo et al., 2016; Zhu & Zabaras, 2018; Khoo et al., 2021). The major weakness these models is that it impossible to query solutions at off-grid points. More- over, graph neural networks (GNNs) (Kipf & Welling, 2016) are also applied in operator learning (Alet et al., 2019; Anandkumar et al., 2020). The key idea is to construct the spacial mesh as graph, where the nodes are discretized spatial locations, and the edges link neighboring locations. Compared with CNN-based models, the graph operator model is less sensitive to resolution, and is capable of inferencing off-grid points by adding tional cost is still high, growing quadratically with the number nodes. Another category of neural operators is Fourier transform based (Li al., 2020; Kovachki et al., 2021). The models learn para- metric linear functiins in the frequency domain, along with nonlinear functions in the time domain. The conversion between those two domains realized by discrete Fourier transformation.
nodes to the graph. However, its computa-
C PROOFS
C.I PROOF OF THEOREM 1
Firstly, extract the constant, and decompose the error by triangle inequality (subscript of norm omitted):
Proof.
8 f NN
Saog
Lip(fa)
S sup
Ja G
SEF
|| -NII = [ - R &
[| Ro Po
S
|| - R. 9ll PoooPoooRooo
The first term is exactly the reconstructor error &R, by definition of push-forward:
|| R - P 9ll lld R Pl 1.(()u)) ER
And the second term is related to the approximator error EA:
|| P o Pog- Ro RoAoEl =
30 R o AoE oo AoE (P(m)
< Lip(R)|PoGoE - LA(H)
Pop 3F - AoE (7())))
Lip(R)
C.2 PROOF OF THEOREM 3
Proof. Our estimation of approximator error is based on the approximation rates of deep ReLU neural networks derived in Gühring et al. (2020) (notation modified for consistency):
17

---

Under review as a conference paper at ICLR 2023
Theorem 4. Let Letme N 00, C = C(m, s, 9, M.M, and onediimensional output in subsets of the Sobolev space W*a;
O, and 0 S S 1, then there exists a constant M,n), n). with the following properties: For any function f with d-dimensional input
I/'l S M,
and any E the exists ReLU MLP N such that:
I - fllwai S,
and:
size(N)
Ce Ce-m/((on) log (-0/(o-n))
depth(A)
C log s/((-n))
Such error bounds can not be directly applied to the approximator A in our framework, since A is implemented by single branch net |3 with p-dimensional output. It is different with stacking p independent me-dimensioonll output networks (Wihssspp and concatenating the outputs. The key difference lies in the parameter sharing of hidden layers. To fill the gap, design special structure of the branch net 33 without parameter sharing, as explained below.
Given p independent one-dimensional output networks (Ng : I R hisss denote the weight matrix of the i-th layer the N₁ as Wij. The weight matrix of i-th layer of the branch net, w can be constructed as:
W1.1
W
0
W W
0
W 122
W
Wi.
W
The weight of first layer W is vertical concatenation of (W1,)) (Wi.hhsssp And the weight w of any remaining layer 2 2 block diagonal matrix, with the main-diagonal blocks being
(W,, It is easy to verify that such approximator is computationally equivalent stacking of p independent one-dimensional output networks stack((A) [ hisss))
Let q 0. and = Go 8-1 then by Theorem 4 the upproximator error bounded by:
ε | P PoGo
1/2
Σ IING - PoG0
S (pe 112 Vpc.
And the depth and size of B3 can be calculated by comparing with any Ns: (epth(AA) = depth((,)) < C'log (s-1),
sizo(A) sizc(W) + size (W:1)
psize(Wllj) p² size (W.,)
p² (Wi))
i=1
p² size(N)))
< - m/s $log (E (1-
18

---

Under review as a conference paper at ICLR 2023
Table 4: Typical control systems used in literature and our experiments and their dimensions.
System
Description
Dynamics
#
#
Pendulum
single pendulum
ODE
2
RobotArm
two-link robotic arm
ODE
4
CartPole
one pendulum with cart helicopter with 4 rotors 6-DoF rocket
ODE
4
Quadrotor
ODE
9
Rocket
ODE
9
Heating
source beating on 2D plane
PDE
1
Pushing
pushing objects on various surfaces single pendulum with stochastic dynamics
ODE
3
StoPendulum
SDE
2
D EXTENSION TO NON-ZERO ENCODER ERROR
The approximation error estimation can be naturally extended to model non-zero encoder error, as derived in Lanthaler al. (2022). For encoder &, its error &g is estimated by the distance to its optimal approximate inverse mapping, decoder D, weighted by measure
& (D - -Illdu(()
(13)
D argmin &c. EOD ld.
D
Similar to reconstructor error, this emor quantifies the information loss during encoding. An ideal encoder should be invertible, iee the decoder D E-1, thus we have 6. 0.
When the encoder error is non-zero, the definition of the approximator error & (Eq. 8) should be modified accordingly, by replacing 8-1 to D:
$
8 1 - Pag ((e))ll (E (e)
(14)
Also, the approximation error bound (Eq. 10) is changed to (proof similar to C.I, omitted here):
8 & (Lip(( o P)δ + Lip(R) &xxx)
(15)
E EXPERIMENT ENVIRONMENTS
E.1 PENDULUM
S"" /f ca (a(()) - Dgoaa) 2 + (1)
min
(16)
x2(() s.t. a(t) = [i((t) g Isin x. (t) b (1)/T
x(0) = 2 Enit
where state = [1 x2] denoting the angle and angular velocity of the pendulum respectively, and control is the external torque. The initial condition Vinit j0. 0] i.e. the pendulum starts from the lowest position with zero velocity. The cost functional consists of two parts: state mismatching penalty and control function regularization, and C [10, 1] C, 0.1 are balancing coefficients. Other constants are: 1, 10, 1, 1/3.
E.2 ROBOT ARM
M(x)= = mpr? + 1 + m2 + + 24, cos(æ2)) m2 (r + 1172 c05(x2)) +12
ma(r + h cos(12)) + In męr + da
19

---

Under review as a conference paper at ICLR 2023
-mahr2sin(z2)( + xa)) maling 222
m22lry (22())
((x)
0
mirig cos(xc)) + m2 cos(2) + x2) + I cos(x
((x)
m2gr2 cos(21 + x2)
E3
4/3
E3 2r4
M
(x)
+ C(x)
+ g(x)..
The RobotArm (also named Acrobot) is planar two-link robotic arm in the vertical plane, with an actuator at the elbow. The state is æ [m] I2, TyT1a where 21 is the shoulder joint angle, and ø2 is the elbow (relative) joint angle, 23, 24 denotes their angular velocity respectively. The control the torque at the elbow. Note that the last equation the manipulator equation, where M is the inertia matrix, C' captures Coriolis forces, and g is the gravity vector. The details of the derivation can be found in (Spong et al., 2020, Sec. 6.4).
The initial condition Finit [x/4, π 0, T]T, and the target state baseline Tigoal = /m/2, 0, O The cost functional coefficients of two links m 1,2 = 1. gravitational acceleration
c2 [0.1, 0.1, 0.1, 0.1] Cu. 0.1. Other constants are: mass 0. links length 4.2 = 1, distance from joint to
the center of mass 2 2 0.5, moment of inertia / 1/3.
E.3 CART-POLE
41 T3
4ig T4+
433 me + mp sin 2
[u + my sin (lx + cos(x2))]
1
41
cos(((()) mplai (22 sin(x2) (m, + mp) g sin()
(m, + mp sin (x2))
In the CartPole system, an unactuated joint connects a pole(pendulum) to a cart that moves along frictionless track. The pendulum is initially positioned upright on the cart, and the goal is to balance the pendulum by applying horizontal forces to the cart. The state x [ar 22, T3344] where æ1 is the horizontal position of the cart, Z2 is the counter-clockwise angle of the pendulum, and angular velocity of cart and pendulum respectively. We refer the reader to of (Tedrake, 2022, Sec. 3.2) for derivation of above equations.
velocity
The initial condition Tinit = [0, 0, and the target state baseline Byoot = [0, n, 0, 0J The cost functional coefficients are Ca [0.1, 0.6, 0.1] Cu 0.3. Other constants are: mass of cart and pole me,p 0, gravitational acceleration g 10, pole length / 1.
E.4 QUADROTOR
0
0
+ R T
(q)
0
[mg]
1 u
4 ((())) Jw Tu w x Jw.
This system describes the dynamics of helicopter with four rotors. The state x= ip w T Rº9 consists of three parts: position p, velocity v, and angular velocity The control E R4 is the thrusts of the four rotating propellers of the quadrotor. q E R the unit quaternion (Jia, 2019) representing the Aititude(spaciaa rotation) of quadrotor w.n.t. the inertial frame. J is the moment of inertia in quadrotor's 's frame, and Tu the torque applied to the quadrotor. Our setting is similar to (Jin et al., 2020, Appx. E. 1), but exclude the quaternion from the state.
20

---

Under review as a conference paper at ICLR 2023
The derivation is straightforward. The first two equations Newton's laws of motion, the third equation is time-derivative of quaternion (Jia, 2019, Appx. B), and the last equation is Euler's rotation equation (Truesdell, 1992, Sec. L10). And the coefficient matrices and operators used in the equations are defined as follows:
0
W1
W3
(ww)
W1
0
W3
W/
W2
W3
0
W1
W3
W2
W1
(13 + 1
(1293 - q4(1)
[g244 + (2(1)
R(g))
2 (9293 + 9441)
- 2 (02 + g
(g344 9291
2 (9291
9391)
(9=91 - 4291)
1 (g3 +
0
112
1/27
T
- 1/2 c
0
1/2
0
C
We the initial state Dinit [[-8, 6, 6. 0,, 00,0],, the initial quaternion qinit O, and the target
state baseline gool = 0. Cost functional coefficients ce = 1, Cu = 0.1. Other constants are
configured as: mass wing length 0.4, moment of inertia 1. z-axis torque constant
= 0.01.
E.5 ROCKET
ÿ v,
[mg
0
+ R (g)u.
1 1iw)g
2
Jw Tu x Jw.
The rocket system models a 6-DoF rocket in 3D space. The formulation is very close to Quadrotor mentioned above. The state æ Ip u l.w'l " E R2 is same as that of Quadrotor, but the control IR3 is slightly different. Here a denotes the total thrust dimensions. Accordingly, the torque Tu is changed to:
2
Tu
0
X u
We set the initial state Vinit [10, 1,,,-855--0 - the initial quatenion qinit
(0.7(0.55), 0, 0, (in77.)5) and the target state baseline goot = 0. The cost functional coeffi- cients C2 1, C. 0.4. Other constants are configured as: mass m rocket length | the moment of inertia J = diag([0.5, 1])
E.6 HEATING
Ar u in
T 0 an
The heating system mimics 2D plane Ω whose temperature is controlled by heating source u, such plane heated by electromagnetic induction or microwaves (Tröltzsch, 2010, Chap. 1.2.1). The dynamics is Poisson equation with zero boundary, where V 22 is the Laplace operator. The objective is thus double integral over Ω:
min [ (x(() Fyout (s))² c. (.) ds.
The target state here a function, not a vector. To simplify the problem, we let (s) = a sin (TS)) sin(xT2)) + b bin((Ts)) + csin((s2) with [a. b, c) being the parameters. Let the baseline 1. and the in and out of distribution problem
parameters to be (a, b, c| noise ein and Eout to the baseline. The cost functional coefficients are e, 1, C2 = 10-4 During testing, the area Ω [0, 1|2 is evenly divided into 100 areas, with 121 grid points.
generated by adding
21

---

Under review as a conference paper at ICLR 2023
Shape
rectl, rect2, rect 3, hex, el allipl, ellip2, el1 ellip3 butter tri tri2, tri3
Surface
abs, der in, polywoed, pu
Speed (mm/s)
10, 20, 50, 75, 100, 150, 200, 300, 400, 500
Acceleration (ms ")
0, 0.1, 0.2, 0.5, 0.75, 1, 1,5, 2, 2.5
Initial contact
33 points for tri 1-333333 and hex, 40 el ellipp and butter, and 44 for rect1-3
Initial push direction
0", 20", 40", 60%. 80", -20", 40" -60",
Figure List of variables explored Pushing dataset, credited to Yu al. (2016).
E.7 STOCHASTIC PENDULUM
dr (t) = x2 (t df,
dx² (1) |ur [u({(()) mgl sin on (1) + odB(1),
sin # (t)
(t) cos x (t)
(1)
This system similar to the simple Pendulum system introduced in the main text, but with Brownian motion B() (Uhlenbeck & Ornstein, 1930) involved in the dynamics, resulting in stochastic OCP (Fleming Rishel, 2012). The state function is now a random process, thus the cost functional is defined as expectation:
2s E f (x(()) xgoal) + cu (o(t)) dt.
min
The state is defined by angle 21 and angular velocity T2, and control u the torque applied the pendulum. We follow the convention of Gym environment Brockman et al. (2016) by adding the observation state æ. And further constrain the scale of state and control as |u] S 2, [222| < 8. The initial state is Winit = 0, and the target state baseline is goot [πr, 0] Cost functional coefficients
Ca 1, cu = 0.001. Other constants are: mass m = 1, length = 1, gravitational accelcration 0.01.
g 10, scale of noise
We solve the problem in a closed-loop optimal control scheme, where the model takes the current state from the environment input at each time index, then outputs the control the environment. This setting is the same as RL, we define the reward RL as the negative cost of OCP.
E.8 PUSHING
In this dataset, the robot executes open-loop straight push along straight line of 5 cm, with different shapes of objects, materials of surface, velocity and accelerations, and contact positions and angles, see Fig. for details.
The control and state trajectories are recorded at 250 Hz. The length of the recorded time horizon varies among samples, due the difference in velocity and accelcration. We select acceleration = 0.5ms 2 with initial velocity 0. Then define time horizon T = 0.44s, and extract 110 time indices per instance.
The input to the encoder is 4167-dim, including 768-dim gray-scale image of the shape, 3280- dim friction map matrix, 110-dim trajectory, and 9-dim other parameters (e.g. mass and moment of inertia). The input to the neural network (including CNN) is 801-dim, and the encoded vector is 44-dim.
F DETAILS ON IMPLEMENTATIONS
For all the systems, have 2,000 samples in each ID/OOD validation set, and 100 problems in each ID/OOD benchmark set. The size of the training set varies among systems, and is roughly proportional the number of trainable parameters, displayed in Table 5.
22

---

Under review as a conference paper at ICLR 2023
Table 5: Hyper-paramelee settings of the proposed OptCtrlOp for different systems.
System
Depth
#Params
#Train Data
#Epochs
Pendulum
7681
10
10
RobotArm
3601
X 10a
X 10a
CartPole
6881
X 10a
X 103
Quadrotor
17284
10a
5 10a
Rocket
13883
X 10a
5 10a
Heating
3601
X 104
X 104
Pushing
12202
10*
103
StoPendulum
6881
X 10
1 1 10*
PDP, DM, and synthetic control systems are implemented CasADi (Andersson et al., 2019), which adapted from the code repository'. For classical methods, the dynamics discretized by Euler method. To limit the running time, we set the maximum number iterations PDP to 2,500.
OptCtrlOP and other neural models are implemented in PyTorch (Paszke et al., 2019).
For MLP, hyper-parametern are the same or as close as to that of OptCtrIOP. We adjust the width number of parameters. For FNO2, set the number of Fourier layers to as suggested in the open-source codes, and tune the network width such that the number of parameters is in the same order as that of OptCtrIOP. Notice that the original FNO outputs function values at fixed time indices, which is inconsistent with our experiment setting. Thus we slightly modify it by adding time indices its input. For GEN set graph nodes uniformly spaced in time (or space) horizon, and perform 3 graph convolution steps on them. The input function initializes the node features at time index
of MLP layer reach almost the
(multiplied by weights). For any time index ,, the GEN output is defined the weighted average of all node features. Both input/output weights are softmax negative distances between / and node positions.
For stochastic OCP (i.e. Stochastic Pendulum), we choose reinforcement learning algorithm named Proximal Policy Optimization(PPO) (Schulman et al., 2017) as the ground truth closed-loop OCP solver, by defining the reward negative cost. The PPO implementation is credited to open-source package Stable Baselines (Raffin et al., 2021), with hyper-parameter settings cited from Raffin (2020).
For fairness, all training/testing cases are executed on an Intel 19-10920X CPU, without GPU.
G MORE EXPERIMENT RESULTS
Table Results Pendulum environment.
Time (sec./instance)
ID MAPE OOD MAPE
ID-OOD Gap
DM
3.80 x 10-2
Ours
3.33 x 10-7
3.05 x 10-4
4.09
10 10-3
3.79
10 10-3
MLP
3.40 x 10-5
4.96 x 10-4
1.67 x 10-2
1.62
10
GEN
6.19 x 10-5
5.61 x 10-3
4.46
10-2
3.90
10-2
FNO
3.47 x 10-4
3.62 x 10-4
3.28
10-3
2.91
10-3
PDP
5.29 x 10
2.79 x 10-1
1.04 10
1.75
10 10-1
http::/www......aaa.aaaaccaaaaaa 'https: / / http://th...eeeeeraaaaagggggggga
23

---

Under review as a conference paper at ICLR 2023
Table 7: Results of RobotArm environment.
Time (sec//instance)
ID MAPE OOD MAPE
ID-OOD Gap
DM
4.58 x 10-2
Ours
1.88 x 10-?
1.24 x 10-5
5.95
10 10-5
4.70
10-5
10
MLP
1.56 x 10
1.16 x 10
1.43
10
1.42
x
10
GEN
6.24 x 10-5
1.53 x 10-5
1.69
10-4
1.53
10-4
FNO
3.48 x 10-4
8.86 x 10-6
1.29
10-4
1.20
10-4
PDP
5.62 x 10¹
4.73 x 10
1.43
10
9.55
10
Table 8: Results of CartPole environment.
Time (sec//instance)
ID MAPE OOD MAPE
ID-OOD Gap
DM
4.63 x 10-2
Ours
1.99 x 10-7
1.60 x 10
4.79
10
3.19
10 5
MLP
2.02 x 10-5
3.79 x 10-5
8.97
10-5
5.18
10-5
GEN
6.41 x 10-5
2.35 x 10
9.09
10 10
8.85
10 10-4
FNO
3.49 x 10-4
1.82 x 10
4.98
10-5
3.16
10
PDP
5.61 x 10
5.96 x 10-4
4.52
10-5
5.50
10-4
Table 9: Results Rocket environment.
Time (sec./instance)
ID MAPE OOD MAPE
ID-OOD Gap
DM
9.71 x 10-2
Ours
3.37 X 10-7
4.42 10
5.24
10
4.80
10
MLP
2.87 X 10-5
5.71 x 10
5.83
10
5.26
10
GEN
6.16 x 10-5
3.57 x 10-5
6.00
10
5.65
10 10-4
FNO
6.33 X 10-4
7.03 10-5
1.81
10
1.11
10-4
PDP
7.01 x 10'
1.80 x 10
1.52
10
2.79
10
Table 10: Results of Heating environment.
Time(s/problem)
ID MAPE OOD MAPE
ID-OOD Gap
DM
5.87 10-2
Ours
2.14 x 10-7
1.07 X 10
1.14 X 10
7.00 x 10
MLP
2.05 x 10-5
1.93 x 10-2
1.67
x 10-2
2.60 x 10-3
GEN
6.70 10-5
1.67 x 10-1
2.20
10
5.30
10
FNO
4.37 x 10-4
2.43 x 10
1.00
10
1.43 x 10
Results of closed-loop control on Stonhusticeendulum environment. result is generated by PPO.
Table 11:
The ground truth
Time (sec./instance)
In Dist. MAPE
Out of Dist. MAPE
PPO
3.81 x (11.6 x 10-1)
MI.P
2.10 x 10-1(+1.87 x 10-6)
9.13 10 2 10-2)
4.69 x 10-2 x 10-2)
OptCtrlOP
3.33 x 10-4 x 10-6)
9.19 10 10 10-2)
3.04 x 10 10 x 10-2)
24

---

Under review as a conference paper at ICLR 2023
Ours train
Ours __st__ttrrg valid MLP_out valid
Ours in_ valid
MLP_ train
MLP in valid
MSE
10-1
10-2
2500
5000 Epochs
7500
10000
0
2500
5000 Epochs
7500
10000
(a) Pendulum
Ours train
Ours _atss_tttttr valid MLP_ out valid
10
Ours in valid
10
MLP train
MLP_in valid
MSE
10-2
10-3
2500
5000
7500
10000
0
2500
5000 Epochs
7500
10000
Epochs
(b) RobotArm
Ours train
Ours out valid
10
Ours in valid
MLP out valid
MLP train
2
MLP in_ valid
MSE
10-3
10
500
1000 Epochs
1500
2000
500
1000
1500
2000
Epochs
(c) CartPole
Ours train
Ours Ours_out _ut__out____ valid
10
Ours in valid
MLP MLP_out valid
MLP_ train
MLP_train
10-1
MLP in in_iinvvvo valid
MSE
10-2
10-3
100
200
300
400
500
200
300
400
500
Epochs
Epochs
(d) Rocket
Figure 7: The loss eurves of systems the training set, in- and out-of-distribution testing sets. All curves are visualized after exponential moving average with weight=0.55
25

---

Under review as a conference paper at ICLR 2023
10a
Ours out test MLP out test
10
10-2
Ours in test
10-3
MLP in test
10-4
2500
5000 Epochs
7500
10000
2500
5000 Epochs
7500
10000
(a) Pendulum
Ours in test
Ours Ours_out test
MLP in test
MLP MiP__out_e test
MAPE
10-2
10-4
2500
5000 Epochs
7500
10000
2500
5000 Epochs
7500
10000
(b) RobotArm
Ours in_iiinttt_e test
Ours out test MLP_out_ test
10-2
MLP MLP_in_tee test
10
10-4
10
500
1000
2000
Epochs
Epochs
(c) CartPole
Ours iest
Ours_ purs_out__ test MLP _uttottt_hce test
10-2
MLP in_ test
MAPE
10-3
10
10-5
100
200
300
400
500
0
100
200
300
400
500
Epochs
Epochs
(d) Rocket
The cost MAPE (mean absolute percentage error) curves of 4 systems on the in- and out-of-distribution testing sets. All curves are visualized after exponential moving average with weight=0.5
Figure 8:
26

---

Under review as a conference paper at ICLR 2023
DM
5.0
5.0
5
5.0
2.5
2.5
MLP
0.0
PDP
1.0
0.5
2
t
X1
(a) Pendulum
DM
1.5
1.7
1.4
Ours
1.4
1.6
1.0
1.0
1.0
1.0
0.80 0.82
0.84
0.0
0.74
0.76
0.78
0.0
0.85
0.90
0.0
0.85
0.0
0.5
(b) RobotArm
DM
0.4
0.50
0.4
0.4
0.6
0.2
0.2
0.
PDP
0.0
0.00
0.0
0.0
10
-151 -13110..00 5o.50.
0.5
0.5
0.5
-1
-1
t
X2
(c) CantPole
-5.6
-5.8
5.8
55.
X
Ours
5.8
5.9
-5.9
5.9
MLP
1.0
1.0
10
-8.0
7.8
0.0
-8.0 8
67.67.
7.4
0.0
-8.0
0.0
-8.0
7.8
0.0
(d) Quadrotor
DM
7.0
7.0
6.5
7.0
7.0
7.5
775
PDP
8.0
8.0
8.0
8.0
1.0
1.0
1.0
1.0
0.5
0.5
0.5
0.5
10
0.0
8
10
8
10
X;
8
10
0.0
t
(e) Rocket
Figure 9: The state trajectories of five systems on four randomly sampled OCPs. The left two columns are distribution problems, and the right two columns are out distribution problems. The three dimensions consist of time and the first two dimensions of state Z2,
27

---

Under review as a conference paper at ICLR 2023
10A
MLP_ MLP_in_ted Jn_test
Ours_in In test
10-1
MAPE
_oout test
Ours_ out test
10
(x10%)
(x10P)
(a) Pendulum
MLP in test
10
Ours in test iestes
10-3
10
MLP_oue duP__utt_re test
Ours_out, test
10-4
6
10
10
#Train Data (x10')
#Train Data (x10')
(b) RobotArm
MLP instes
Ours in
10-4
MLP
MLP_out out test
Ours out test
10-1
18
24
#Train Data (x10?)
#Train Data (x10%)
(c) CartPole
10-2
MLP In test
Ours, in test
10-3
10-1
MLP out test
Ours out test
10-1
80
100
60
80
100
#Train (x107)
#Train Data (x10P)
(d) Quadrotor
10
Ours test
10-4
MLP out test
10
Ours out test
20
60
80
100
40
60
80
100
#Train Data (x10')
#Train Data (x10")
(e) Rocket
Figure 10: The MAPE (mean absolute percentage error) w.r.t. the number of training samples curves of systems on the in- and
testing sets.
28


Reviews:
[{"review_id": "cLdZHhQ45Wk", "reviewer": "Reviewer_p9mD", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "In this paper, the authors propose a new application of neural operators. They target a well-studied problem of optimal control where there is a dynamical system and the goal is to find a control function that minimizes a certain cost specified by the control system. \n\nTo do this, they first transform the OCP (optimal control problem) into a BVP by utilizing a known technique in optimal control known as the Pontryagins Maximum Principle. Then they apply a DeepONet-based operator on the resulting BVP. The authors then use previous work on DeepONet approximation bounds to provide a new theorem on the approximation error of their approach. They then provide simulations\nof multiple systems as their set of experiments (mainly focused on a pendulum).", "summary_of_the_review": "This paper proposes a new application of neural operators in learning the optimal deterministic control operator derived from the Pontryagins Maximum Principle. The problem of the study is interesting. The main limitation/weakness of the proposed paper is as follows: \n\n1) the theory analysis in Section 3 does not match the actual method. The authors should include the encoder approximation error to the approximation error bound; \n\n2) Lack of necessary comparison to baselines to justify the choice of the model architecture. The authors need to compare against FNO, GNO, and PINN in terms of computational efficiency as well as approximation accuracy (in both in-distribution data and OOD samples) to benchmark the proposed approach and argue why such a design is desired; \n\n3) Except for the pendulum example, other experiments are not discussed in necessary detail - the authors may want to re-assess the distribution of space to include at least one more real-world like control example and comparison in the main text beside the pendulum example. ", "strength_and_weaknesses": "Notes on Each section:\n\nRelated Work and Baseline:\n\nWeakness:  When talking about neural operators, only the DeepONet is mentioned but another related line of work Fourier Neural Operator (https://arxiv.org/pdf/2010.08895.pdf) is not mentioned, nor compared in the experiment section. In addition, the authors mentioned about another neural operator work Graph neural operator https://arxiv.org/pdf/2003.03485.pdf, still, no comparison of the GNO to DeepONet is conducted when learning the optimal control operators. The authors mentioned GNO may has high computational complexity this slower computational time, but a comparison between DeepONet, FNO, and GNO in learning the optimal control operator is needed to justify the choice of network structure.\n\nMethodology: \n\nStrength: The methodology section is well explained and easy to follow. \n\nWeakness: \n\n* One of the major weaknesses of this paper is that the approximation error bound in Section 3 does not match the experiment setting in Section 4. The author mentioned that “if the encoder is implemented by functional approximations e.g. point-wise evaluation, then the error should be considered”. In fact, their experiments do use a discretization grid/point-wise evaluation, so the approximation bound in Section 3 should match the actual setting and include the error from the encoder.\n\n* V needs to be defined in Eqn. 8 \n\nExperiments: \n\nWeakness: \n\n* My biggest concern is regarding the comparison to necessary baselines in order to benchmark the proposed approach. The authors only compared the DeepONet against MLP. They need to compare against other neural operators including both FNO and Graph-based Neural Operator, in terms of both accuracy and computation speed. \n\n* Related to the first weakness point, another useful baseline to compare against is PINN. Though the authors mentioned in related work that PINNs need to solve instance by instance, while DeepONet learns the entire operator. However, how about learning accuracy? Similarly, as my first comment, such a comparison between the proposed method and PINN in terms of both accuracy and computational efficiency is necessary to justify the choice of network structure. Especially on data that is unseen during the training stage or completely Out of Distribution test data, PINN might have better accuracy than neural operators since it solves the equation rather than suffering from generalization error. I would love to see such a comparison on test data that are sufficiently different from the training case and see the comparison between the proposed approach and PINN.\n\n* How can we tell which method is better from the state trajectories plot (e.g., Fig. 4 and Fig. 7)? Why is one trajectory better than another?\n\n* I don’t understand how the data for the stochastic pendulum was created using RL. “For stochastic OCP (i.e. Stochastic Pendulum), we choose a reinforcement learning algorithm named Proximal Policy Optimization(PPO) (Schulman et al., 2017) as the ground truth OCP solver. ” Page 21. The authors need to provide sufficient descriptions and details if they want to include this part of the results.\n\n* Another thing was I think the authors ran out of space. There is no extensive discussion about the real-world examples provided.\n\nSmall questions/weaknesses:\n* Why vary the time horizon? What did this achieve in demonstrating?\n* Why are the noise distributions for the pendulum different? One is positive and the other is negative (epsilon_in and epsilon_out)\n* Why is the batch size 10k for a data set of 1k (Pendulum example) - which seems like a typo\n\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "FJCyrqIpBNg", "reviewer": "Reviewer_Yqn8", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposes a neural network architecture to estimate a mapping from cost functionals to optimal inputs in order to solve optimal control problems. The authors extend the work involving DeepONet (Lu et al., 2021 and Lanthaler et al., 2022) to optimal control problems and also extend the analysis to theoretically bound the error of their method, OptCtrlOP. The performance of OptCtrlOP was shown through simulation on six different systems ranging in states from one to nine dimensions and input dimensionality ranging from one to four. The provided results show an improvement in running time and accuracy against two baselines.", "summary_of_the_review": "The concept of using a neural network as an operator mapping from cost functionals to control functionals is an important area of research. The paper is well written and the inclusion of both theoretical analysis and empirical results is a big strength. My primary concern with the paper as stands is the lack of support throughout for the very strong claims. Additionally, based on my current understanding of the approach, I do not see the novelty and how this approach would apply to control problems (if limited to a single initial condition). If this approach is the initial step towards work to provide a mechanism for a feedback/replanning process, then I recommend stating that more explicitly.", "strength_and_weaknesses": "The paper does a good job of providing background information and providing a theoretical analysis of their approach. I also appreciate how the paper is organized and includes both theoretical and empirical results. However, I think the largest weakness is the lack of evidence to support the strong claims of the paper. \n\n- Claim: \"[OptCtrlOP] solves OCPs in a one-shot manner with no dependence on the explicit expression of dynamics or iterative optimization processes.\"\n\n  - While I understand how this claim can be reached based on the wording, I believe it is a bit misleading. As presented, the method requires quite a bit of data involving the optimal control of the system. While the execution of the operator is not depending on the expression of the dynamics or an iterative nature, developing the mapping and training certainly do.\n\n- Claim: \"The design is in principle endowed with substantial speedup in running time...\"\n\n  - Results from the simulations support this claim, but I do not believe the baselines support such a claim. The speed shown for inference of OptCtrlOP is impressive and I believe would still be the better performer. However, the performance shown for the direct method appears to be slow compared to existing methods in that category. I understand the use of IPOPT as a generic non-linear solver; however, most implementations of a DM involve more specific implementations. I recommend a comparison of the systems to SCP-based implementations closer aligned with current state-of-the-art implementations. Additionally, do you expect this performance disparity to scale as the problem complexity increases? In Figure 2, there it appears there is more than a 2x increase in run time when going from problems with an input size of one to four (pendulum/robot arm/cart pole to quadrotor). This leads me to wonder how the inference time scales as the network size needs to grow.\n\n- Claim: \"... the model reusability is guaranteed by high-quality in- and out-of-distribution generalization.\"\n  - The discussion about model reusability was confusing for me to follow and see the connection to this claim. The experiments did provide a look at what happens with the goal location was outside of the training distribution. However, what about other model parameters (e.g. length of the pendulum, weight, other physical parameters)? Would you expect the same performance with modifications to other parameters? Providing a discussion on the connection of what is demonstrated or experiments expanding on this would add support for this claim.\n\n- Claim: \"Extensive experiments on 6 physical systems verify the effectiveness and efficiency of our approach.\"\n  - While I do appreciate the experiments conducted on different systems. I feel like the use of \"extensive\" and \"physical\" are misleading with simulation-only approaches involving only two modest baselines. I would recommend reevaluating this claim or supporting the experiments with more baselines and hardware.\n\nMy understanding of operators in the control space is that they are global operators mapping cost functionals to another functional space. The authors do a great job of discussing that concept in section 3.1. However, I believe there is a disconnect when generating data and training the operator with a single initial condition. This process appears to produce an approximation to the operator that is only valid in a local region (starting from the initial condition or near). Does this process work from arbitrary initial conditions or is a separate trained network needed? My understanding is that it does not work from arbitrary initial conditions. If it does, I recommend adding a discussion supporting that and/or experiments that also show this attribute. If it does not work from arbitrary initial conditions, I am struggling to see the applicability to a control approach. For example, when using direct methods in an MPC fashion, we do solve for an open loop control sequence from an initial condition. However, we replan at each step. So while a lot of MPC-based methods assume deterministic dynamics when solving an open-loop problem, a lot of uncertainty is handled by the replanning at the new state at each time step.\n\nDoes the frequency of training data influence the performance at different control frequencies? For example, if we are wanting to control something at a rate much slower than its dynamics, would this approach still be valid to estimate the operator? Would the rate of data gathered for training influence this performance?\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}, {"review_id": "7WIkNOgt0Y", "reviewer": "Reviewer_vSM1", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "In this paper, the authors propose a new, data-driven method to solve optimal control problems that is significantly faster than existing methods at execution time and produces similar suboptimality of solution when compared to other approaches. The method is explained through the lens of learning an operator that maps an optimal control problem to a control function. Specifically, the authors learn a neural network that maps the goal state and the time index to a control input for the dynamical system. The network architecture is inspired by the DeepONet paper.", "summary_of_the_review": "To summarize, I think that this is exciting work but the current presentation does not do it justice. Were this an option I would give the recommendation of: accept conditioned on major revision. ", "strength_and_weaknesses": "The main strength of this paper is that the method is incredibly fast – the method is many orders of magnitude faster than more traditional methods. Speed is one of the main bottlenecks of optimal control approaches, so methods that obtain similar performance with much less computational cost are valuable. \n\nThe main weakness of this paper is that the writing presents the content hastily; with little exposition, explanation, or context, making it difficult to follow the logical argumentation. Despite being someone who does research in the area of control theory, I had to read several other papers just to understand this paper. The audience will never have exactly the same background and expertise as the authors and the paper should be written with that expectation. This is especially relevant because this is a controls paper submitted to a machine learning conference – the vast majority of registrants to ICLR have little to no background in control theory but it should be accessible to them, as they are the intended audience. \n\nSpecific points to address:\n1. The second paragraph of the introduction jumps immediately into \"Lebesgue measurable function space\" with no definitions. I would ask – is this reference to the foundations of optimal control theory necessary? Or can it be omitted and instead a plain-language description of the optimal control problem added, after which the paper transitions directly to the formulation of the optimal control problem (Eq. 1)? This is the type of hasty writing without exposition, explanation, or introduction that appears throughout the paper and makes it difficult to read and understand.\n\n2. The concept of model reusability is not defined but should be defined in the paragraph where it is mentioned on page 2. \n\n3. The concept of two-phase optimal control solvers is not defined but should be defined in the paragraph where it is mentioned on page 2.\n\n4. In the discussion of PINNs in section 2.1, the paper mentions that \"the magnitude of the two loss terms is inherently imbalanced\" but this has little meaning for the reader without the appropriate context. Which two loss terms? Why should the reader care about this imbalance? As a note for here and throughout the paper: the paper should summarize other works with enough context so that the reader can appreciate the points the paper is making without having to re-read every reference. \n\n5.  It says in section 2.1 that for Direct methods of solving optimal control problems, \"The\nreformulation essentially constructs surrogate models, where the state and control function (infinite dimension) is replaced by polynomial or piece-wise constant functions.\" I do not understand what reformulation the paper is referring to. To my knowledge, using a direct method to solve an optimal control problem does not necessitate using a surrogate model. What am I missing here?\n\n6. At the end of section 3.1 the paper claims that the approach is highly reusable, but as I mentioned above, there is no definition of reusability that the paper has given so it is difficult to evaluate this claim. On my first read through this paper, I understood this to mean that the cost function could change even if the dynamics had to remain the same. However, later in the methodology the paper defines the encoder as simply returning the target state and effectively learn with respect to a fixed cost function. It seems like a major simplification. Does this claim of reusability just refer to the fact that the goal point can change? (because it's an input to the network)? I think this point needs major clarification. If possible, it would be nice to see an experiment where the encode of the cost function is more sophisticated. \n\n7. The main algorithm of the paper (Alg. 1) appears in the appendix. It would be nice to see it appear in the text, although I know space is limited. In contrast, the definition of a neural network in section 3.2 list item 2 (Approximator) is not necessary for the readers of ICLR papers who are all familiar with machine learning. It should be moved to the appendix. \n\n8. I had difficulty following the theoretical proof in 3.3. Much more intuition, definition, and explanation could be added. The terms in eq. 6 are not well defined. What is a push forward measure? Even using Lip() to denote Lipschitz constant (I presume) should be defined. \n\n9. For the pendulum example appearing in the text – I would encourage you to switch this out for something with more interesting dynamics that you have already run experiments on. It seems that every controls paper I read uses an inverted pendulum and it does not demonstrate the capability of the method. One would only use an involved learning based control algorithm like the one in this paper for a more complex system. Also, I believe there may be an error in the way in-distribution and out-of-distribution noise was defined for the pendulum – they appear to be switched, I think? It would make more sense for the in distribution noise to be centered around the goal point. \n\n10.  As a high-level point of feedback – I believe that a more concrete description of the method is in order. E.g. perhaps a diagram of network architecture similar to the DeepONet paper as well as a straightforward description like this line from the DeepONet paper: \"Let G be an operator taking an input function u, with G(u) being the corresponding output function. For any point y in the domain of G(u), the output G(u)(y) is a real number. Hence, the network takes inputs composed of two parts: u and y, and outputs G(u)(y) (Fig. 1a).\"\n\n11. Because the scale of the plots is logarithmic, the reader cannot easily verify a claim like the method is \"50% faster\" by e.g. looking for a bar that is 50% larger than another. \n\n12. The paper discusses the performance gap between ID and OOD performance but does not give numbers and also plots these on separate graphs, so again, it is hard for the reader to validate these claims. \n\n13. The paper gives a *hypothesis* as to why your architecture is better at handling OOD shift, but it is stated more that you *know.* I would clarify that this is an explanatory hypothesis.\n\n14. In the abstract, six experiments are mentioned, but the plots in Fig. 2 omit the heat problem. Why?\n\n15. It is difficult to extract information from fig. 4. Perhaps you could use 2D plot views instead? Thinner lines? Euclidean distance from the optimal path?\n\n16. The paper does not mentioned limitations of the work in the conclusion, but should. \n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works."}, "extra_text": {}}]

Author rebuttal:

Decision: Reject