# Closed-loop training of static output feedback neural network controllers
for large systems: A distillation case study

Evren Mert Turan  Johannes Jäschke

###### Abstract

The online implementation of model predictive control for constrained
multivariate systems has two main disadvantages: it requires an estimate of
the entire model state and an optimisation problem must be solved online.
These issues have typically been treated separately. This work proposes an
integrated approach for the offline training of an output feedback neural
network controller in closed loop. Online this neural network controller
computers the plant inputs cheaply using noisy measurements. In addition, the
controller can be trained to only make use of certain predefined measurements.
Further, a heuristic approach is proposed to perform the automatic selection
of important measurements. The proposed method is demonstrated by extensive
simulations using a non-linear distillation column model of 50 states.

###### keywords:

Optimal control , Neural networks , Closed-loop optimal controllers , output
feedback controller, measurement noise, Explicit MPC, End-to-end learning

[inst1]organization=Department of Chemical Engineering, Norwegian University
of Science and Technology (NTNU),addressline=Sem Sælandsvei 4, city=Trondheim,
postcode=7491, country=Norway

##  1 Introduction

Due to the difficulty of explicitly defining a control law for multivariate
interacting systems with constraints, a popular approach is to use model
predictive control (MPC) to implicitly define a feedback control law. Given an
estimate of a system’s state, MPC uses a dynamic model to forecast the
behaviour of the system subject to a sequence of control actions, optimises to
find this sequence to minimise an objective and satisfy system constraints.
This sequence of control actions results from an open-loop optimisation, i.e.
the third control action does not take into account the sequence of
disturbances that may have occurred. It is well established that use of open
loop optimisation can have dramatic effects due to stochasticity or
innacuracies in the system forecast. MPC incorporates feedback by only
implementing the first control action, before repeating the optimisation with
a new estimate of the system state.

Thus, MPC incorporates feedback due to measurements that are used to estimate
the system state. However, each time a control action is determined, MPC
requires (a) the current state of the system to be estimated from
measurements, and (b) the solution of an optimisation problem. For many
systems one can use a linear(-ised) model which allows for computing the state
estimate and optimal solution reliably and efficiently. Some systems are
inherently non-linear and the use of linearised models can result in poor
controller performance. However, non-linear models generally incur a
significantly larger computational cost, and as such there may be a non-
negligible delay between measurement and control action. This is especially
true for large systems, or if uncertainty is considered. This delay reduces
performance and can (if significant enough) destabilise the system [1, 2].

Various approaches have been proposed to reduce the computational cost
associated with MPC. This can be done by general approaches, such as improving
the optimisation algorithm, or more tailored approaches such as using the time
between measurements to compute the control sequence in advance, and then
adjust this based on the new measurement [3, 4, 5]. Many authors have taken
motivation from recent advances in machine learning to reduce the
computational cost of MPC. The main idea is to approximate some control
relevant mapping that is typically unavailable or expensive to compute.
Examples of these mapping include: terms in the objective function [6, 7],
steady state inverses of processes [8], and offline control policies [9, 10,
11, 12, 13, 14]. This work focuses on the latter, where the aim is to learn an
explicit output feedback control policy offline instead of solving an optimal
control problem online.

For a “standard” linear MPC problem the optimal control law is piecewise
affine on polyhedrals, and can be computed offline by solving a parametric
program [15]. This has become commonly known as ”explicit MPC”. Unfortunately,
the online computational requirements of this method grows exponentially in
the number of states and length of control horizon, and hence this approach is
not suitable for large systems. Thus, a more generally applicable approach is
to train an approximator to describe the control policy. Typically neural
networks are used as they are flexible and cheap to evaluate – e.g. one can
find a neural network that is exactly equivalent to the explicit linear MPC
solution, while avoiding the exponential scaling of the online evaluation
[10].

When training a neural network control policy one may either (a) separate the
training from the control problem by first collecting a database of states and
control actions, and then train the controller on this database (optimise-
then-learn, or imitation learning) [9, 10, 11], (b) combine the two problems
by embedding the neural network in the control problem and train the
controller on the control objective (optimise-and-learn) [14, 13], or (c) use
a reinforcement learning approach. As we wish to find a static, offline policy
we do not consider reinforcement learning approaches in this work.

Although the imitation learning approach is conceptually and
implementationally simpler, it is restrictive as it requires data for the
controller to imitate. In addition, as the control problem is solved
separately to the training, a controller that does well on the
training/validation data may have poor closed-loop performance [14]. In the
optimise-and-learn approach the controller is directly trained based on the
controller performance, and hence avoids problems occured by separating the
control and training problems [14, 16]. In addition, if the controller is
directly embedded into the system dynamics, e.g. Turan and Jäschke [17, 14],
then the training occurs in closed-loop. Closed-loop training allows for
flexible handling of uncertainty and also allows the use of a smaller
controller compared to the open-loop alternative, e.g. [13]. In general, the
optimise-and-learn approach is a more complicated training problem, however it
allows for a more flexible specification of the resulting controller.

In this work we find an explicit output-feedback neural network control
policy, that can be cheaply evaluated, using noisy measurements, to give
control actions. This policy is trained in a closed-loop optimise-and-learn
formulation [14], which allows considering measurement noise, and pragmatic
usage of a subset of only some measurements. We show that one can perform
measurement selection as part of the training of the control policy, and
highlight some potential issues of this joint approach. The control policy is
demonstrated on simulations of a non-linear distillation column. This is a
considerably larger problem (50 states) than all other examples in the
optimise-and-learn literature.

Prior literature has predominantly considered relatively small scale examples,
with the notable exception of Kumar et al. [11] which considered an imitation
learning approach applied to large linear MPC problem for control of
distillation column.

To the authors knowledge the only work considering offline optimisation of an
output feedback neural network policy is [13], which considers a small single-
input single-output system. Unlike this prior work, we incorporate knowledge
of the measurement noise in the closed loop training of the policy and
consider control policies that do not use all available measurements. In
addition, we consider a much larger nonlinear model of a distillation column
consisting of 50 states. A key challenge of the training is the size of the
state space. Any practical implementation requires restriction of attention to
a smaller region of the state space. We sample the state space to find the
typical operationally relevant region [11, 12], and train the policy to start
from that region. Unlike an imitation which would requires sampling along
trajectories starting in this region [12], the control policy is optimised in
closed loop and hence these samples are “implicitly generated” during training
[14].

This work is structured as follows: we briefly summarise the background and
relevant literature in Section 2, and proceed to formulate the output-feedback
training problem in Section 3. In Section 4 we introduce the distillation
problem. In Section 5 we construct the control policies, and present extensive
numerical results of their closed-loop performance in Section 6. Lastly in
Section 7 we discuss important aspects of the results and potential future
work.

##  2 Background

###  2.1 Model predictive control

In a continuous time formulation of model predictive control (MPC) we wish to
solve the time dynamic optimisation problem111Note that in practice, most MPC
problem are given in a discrete time formulation, but we focus on the
continuous time formulation as we are interested in obtaining a continuous
time feedback control law.:

$\displaystyle u^{*}_{MPC}(t,z_{0})=\arg\min_{u(t)}\>\>$ $\displaystyle
J(z_{0},u)$ (1a) $\displaystyle
J(z_{0},u)=\int_{t_{0}}^{t_{f}}l(z,u,t)\>dt+V_{f}\left(z(t_{f}),u(t_{f})\right)$
(1b) $\displaystyle 0=f\left(\dot{z}(t),z(t),u(t),p\right)$ (1c)
$\displaystyle z(t)\in\mathcal{Z}$ (1d) $\displaystyle u(t)\in\mathcal{U}$
(1e) $\displaystyle z(t_{0})=z_{0}$ (1f)

where $t$ is time, $t_{0}$ is the initial time, $t_{f}$ the final time,
$z\in\mathcal{Z}\in\mathbb{R}^{n_{z}}$ is the system state, with time
derivative $\dot{z}$, $u\in\mathcal{U}\in\mathbb{R}^{n_{u}}$ is the control
input, $u^{*}_{MPC}$ is the optimal control input, $f$ is an implicit
differential equation, $z_{0}$ is the state at initial time, $l$ is the stage
cost, $V_{f}$ is the terminal cost, and $\mathcal{Z}$ and $\mathcal{U}$ are
constraint sets. In process control problems these constraints sets are
typically defined by upper and lower bounds constraints on $z$ and $u$.

To incorporate feedback MPC is implemented in a receding-horizon approach.
Given a state estimate, $z_{0}$, (1) is solved, and the control
$u^{*}_{MPC}(t_{0},z_{0})$ is implemented. Then at $t_{0}+\Delta t$, given a
new estimate the problem is resolved to find a new control action. Thus, MPC
implicitly defines the control policy:

$\kappa_{MPC}(z_{0})=u^{*}_{MPC}(t_{0},z_{0})$ (2)

where $\kappa_{MPC}$ returns the first control action of the open-loop
solution. We emphasise that to solve (1) the full initial state of the system,
$z_{0}$, needs to be specified. In a typical process many of the process
variables are unmeasured, and hence have to be estimated from the available
noisy measurements.

###  2.2 Neural network and control policies

When using MPC, due to the need to perform state estimation, and then solve an
optimisation problem (1) there will be a delay between receiving the system
measurement and sending the optimal input to the plant. If large enough this
delay can reduce the controller performance and can even destabilise the
system. Although there are efficient optimisation algorithms for MPC, managing
the computational delay can be challenging for large systems, especially when
non-linear dynamics and uncertainty is considered. We consider the use of
neural network control policies to eliminate this computational delay and
provide fast online evaluations for the control actions.

Authors dating back to 1995 [9] have proposed the use of neural networks to
learn control policies that otherwise require an expensive evaluation such as
MPC, as they are universal approximators, i.e. sufficiently large feed-forward
neural networks are able to approximate bounded, continuous functions defined
on a compact subset of $\mathbb{R}^{n_{z}}$ to an arbitrarily low tolerance
[18, 19]. More recently authors have shown that neural networks with a
specific architecture can exactly represent the control policies of linear MPC
[10], and under some assumptions can approximate nonlinear MPC policies to
arbitrary precision [14].

Consider the feed-forward neural network, $\kappa_{NN}$:

$\displaystyle\kappa_{NN}(z,\theta)=\zeta^{(N_{L})},\qquad\zeta^{0}=z$ (3a)
$\displaystyle\zeta^{(i+1)}=\alpha^{(i)}(W^{(i)}\zeta^{(i)}+b^{(i)}),\qquad
i=0,\dots,N_{L}-1$ (3b) $\displaystyle\theta=\texttt{vec}(W^{(0)},\ b^{(0)},\
\dots,\ W^{(N_{L}-1)},\ b^{(N_{L}-1)})$ (3c)

where $N_{L}$ is the number of layers, $\zeta^{(i)}\in\mathbb{R}^{w^{(i)}}$ is
the latent state of layer $i$, $\alpha^{(i)}$ is an activation function, and
$W^{(i)}$ and $b^{(i)}$ are weights and biases that are collected in the
vector $\theta$. In this work $\kappa_{NN}$ is trained in closed-loop to yield
a feedback control policy.

####  2.2.1 Optimise-and-learn formulation

In this paper, we use an optimise-and-learn approach to train the control
policy in closed loop [14, 17, 13, 16]. In this approach the neural network is
embedded into the system dynamics to form a single dynamic optimisation
problem. This problem is not computationally feasible to solve online, however
one can find an offline policy to be used for controlling the system starting
in some region $\mathcal{Z}_{0}$ by solving:

$\displaystyle\min_{\theta}\>\>$
$\displaystyle\mathop{\mathbb{E}}_{\pi_{z_{0}}}\left[J(z(t_{0}),u)\right]$
(4a) $\displaystyle 1$
$\displaystyle=\mathop{\mathbb{P}}_{\pi_{z_{0}}}\left[z(t)\in\mathcal{Z}\right]$
(4b) $\displaystyle 0$ $\displaystyle=f(\dot{z}(t),z(t),u(t),p)$ (4c)
$\displaystyle u(t)$ $\displaystyle=\kappa_{NN}(z(t),\theta),\qquad
u(t)\in\mathcal{U}$ (4d) $\displaystyle z(t_{0})$
$\displaystyle\sim\pi_{{z}_{0}}$ (4e)

where $\pi_{z_{0}}$ is some non-zero probability distribution defined on
$\mathcal{Z}_{0}$, $\mathop{\mathbb{E}}_{\pi_{z_{0}}}$ is the expectation with
respect to $\pi_{z_{0}}$, and $\mathop{\mathbb{P}}_{\pi_{z_{0}}}$ is the
probability with respect to $\pi_{z_{0}}$. If the control policy defined by
the related MPC problem is continuous in $z$ then under mild conditions,
minimisation of (4) yields the a policy equivalent in performance to the MPC
policy [14]. Typically $u$ is constrained between upper and lower bounds, and
thus (4d) can be enforced through the activation function of the final layer.
As this is the standard case we assume that the network architecture is chosen
to satisfy the constraint.

The expectation of the objective and probability constraint need to be
approximated in some way to yield a tractable problem. At each iteration we
consider a stochastic approximation of (4) by evaluating:

$\displaystyle\phi$
$\displaystyle=\sum_{s=1}^{n_{s}}\omega_{s}J(z_{s}(t_{0}),u_{s}(t))+\rho(z_{s}(t_{0}),u_{s}(t))$
(5a) $\displaystyle 0$ $\displaystyle=f(\dot{z_{s}}(t),z_{s}(t),u_{s}(t),p)$
(5b) $\displaystyle u_{s}(t)$
$\displaystyle=\kappa_{NN}(z_{s}(t),\theta),\qquad u_{s}(t)\in\mathcal{U}$
(5c) $\displaystyle z_{s}(t_{0})$
$\displaystyle=\texttt{sample}(\pi_{{z}_{0}}),\qquad s=1,\dots,n_{s}$ (5d)

where $n_{s}$ is the number of samples, $z_{s}(t)$ is the trajectory of sample
$s$ taken from $\pi_{{z}_{0}}$, $\omega_{s}>0$ weights the different
trajectories, and $\rho$ is some penalty function used to enforce the
constraint (1d) [14]. Note that the $s^{th}$ contribution of
$\nabla_{\theta}\phi$ can be calculated by a standard single shooting
approach, i.e. $J(z_{s}(t_{0}),u_{s}(t))$ and $\rho(z_{s}(t_{0}),u_{s}(t))$
can be evaluated by solving a differential equation with the current $\theta$,
and the gradient with respect to $\theta$ can then be found by algorithmic
differentiation. As each contribution to $\nabla_{\theta}\phi$ is independent,
this computation can be done in parallel, greatly reducing the computational
cost of taking multiple samples.

The key aspect that makes (5) appealing is that if random samples can be taken
such that $\phi$ is an unbiased estimator of the penalised objective of (4),
then $\nabla_{\theta}\phi$ is a stochastic approximation of the gradient of
(4). Thus, $\nabla_{\theta}\phi$ can be used in a stochastic optimization
algorithm, e.g. stochastic gradient descent, to solve (4). There are various
options on how to select the samples, weights and number of samples. At one
extreme a single sample can be randomly taken at each iteration, and at the
other extreme the number of samples, weight, and sample locations can be
adapted to maintain a certain approximation tolerance of the expectation and
probabilistic constraint (a form of adaptive mini-batching). Although the
latter gives a better quality estimation of the gradient at each iteration,
each iteration of the former is likely to be much more computationally
efficient.

Note that in some cases the constraint (1d) may be enforceable implicitly by
the system dynamics, network architecture or similar. If this is not the case
then a penalty approach must be used and adequate satisfaction of the
constraint should checked after training.

####  2.2.2 Selection of initial conditions

An important choice in the optimisation-based design of off-line control
policies is the region $\mathcal{Z}_{0}$ on which the policy is developed for.
For small systems one may consider selecting $\mathcal{Z}_{0}=\mathcal{Z}$ ,
however this becomes computationally impractical for systems with more than a
few states. The majority of chemical processes operate in a relatively small
region of the total feasible state space [11, 12]. If a model is available
then this region can be found by simulating closed-loop trajectories of the
controlled system subject to assumed disturbances, e.g. changes in set-points
and disturbances in the feed. Then from the simulations one can directly
estimate $\mathcal{Z}_{0}$ and $\pi_{z_{0}}$.

##  3 Optimisation of an output feedback policy

The previous literature [9, 17, 14, 11, 10, 12] has predominantly focused on
neural network control policies based on the assumption of noiseless full-
state feedback. In general the current state of the system is not available,
instead we have measurements $\mathscr{y}\in\mathbb{R}^{n_{\mathscr{y}}}$
which are related the state by:

$\mathscr{y}(t)=h(z(t),u(t),p)+\eta(t)$ (6)

where $h$ is a (potentially nonlinear) measurement equation, and
$\eta\in\mathbb{R}^{n_{\mathscr{y}}}$ is a vector of noise and/or biases that
influence the measurements.

In this work, we leverage the optimise-and-learn framework to directly
optimise a control policy that:

  1. (a)

Incorporates knowledge of the expected measurement noise into the optimisation
of the feedback policy.

  2. (b)

Only uses a subset of the potentially available measurements to control the
system.

The advantages of (a) is that the control policy will be less sensitive to
small perturbations because of noise, which should yield a more robust
controller. The advantage of (b) is more system dependent, but resolves around
controlling the desired behaviour of the policy. These two points are
developed in the following sections:

###  3.1 Noise and uncertainty

Consider measurements of some of the state subject to additive noise $\eta$ as
in (6) with $\eta$ distributed by $\pi_{\eta}$, some compactly defined
multivariate probability distribution. If the noise is significant it should
be included in the training of the policy. This can be done by adapting the
approach for parametric uncertainty from Turan and Jäschke [14]. Consider the
augmented state vector $\bar{z}=[z,\ \eta]$, where $\eta$ has zero dynamics.
Let $\eta(t_{0})$ be distributed by $\pi_{\eta}$. During training let the
control actions be given by:

$u=\kappa_{NN}(\mathscr{y}(t),\theta)$ (7)

As $\eta$ is measurement noise not a stochastic input, it should only
influence the controller, i.e. in training the dynamics, objective, etc.
should still use $z$. The objective to be minimised is:

$\mathop{\mathbb{E}}_{\pi_{z_{0}},\pi_{\eta}}\left[J(\bar{z}(t_{0}),u)\right]=\mathop{\mathbb{E}}_{\pi_{z_{0}}}\left[\mathop{\mathbb{E}}_{\pi_{\eta}}\left[J(\bar{z}(t_{0}),u)\right]\right]$
(8)

where for clarity we separate the expectations on the right. Along each
trajectory in the training $\eta$ has zero dynamics, i.e.
$\eta(t)=\eta(t_{0})$, and hence acts as a constant measurement bias. However,
the controller cannot learn to “see” this constant bias as the controller has
no memory, and for any state in $\mathcal{Z}_{0}$ any $\eta$ can be sampled.

###  3.2 Measurement selection

In general the entire state of system is rarely measured, and instead
typically only a subset of the states are measured. Thus, it is a natural
desire to have a controller that only uses this subset of states. In addition,
we may also select a subset of these measurements for use. One may interested
selecting only some of the available/possible measurements when:

  * •

Reliable measurements are expensive.

  * •

Certain states can only be estimated with high uncertainty.

  * •

We wish to influence the behaviour of the control policy, e.g. only specific
measurements should influence the output.

  * •

We wish to identify which measurements are more important for control.

  * •

The system under consideration consists of interacting sections, or is
distributed. In this case we may wish that specific control actions are only
dependent on a sparse selection of the total measurements.

Consider the control policy:

$u=\kappa_{NN}(H\mathscr{y},\theta)$ (9)

where $H\in\mathbb{R}^{n_{d}\times n_{\mathscr{y}}}$ is a matrix containing
only zero entries and $n_{d}$ ones on the main diagonal. By specifying $H$ one
specifies which measurements the control policy should use, which can be done
using engineering judgement. Another other approach is to include the
selection of measurements (i.e. selecting the non-zero entries of $H$) in the
training problem.

One cannot directly optimise for which $n_{d}$ entries to select as this is a
computationally intractable problem. Instead we use elastic net regularisation
as an heuristic. For simplicity of presentation we assume that
$n_{d}=n_{\mathscr{y}}$, i.e. all measurements are available candidates to be
used for feedback in the neural network controller. Let $H$ be a diagonal
matrix, and let $\bar{\theta}=\texttt{vec}(H,\ \theta)$. Consider the use of
the regularised objective:

$\min_{\bar{\theta}}\>\>J(\bar{\theta})+\lambda_{1}(\lambda_{2}\|\bar{\theta}\|_{1}+0.5(1-\lambda_{2})\|\bar{\theta}\|_{2}^{2})$
(10)

where $\lambda_{1}$ controls the regularisation strength, and
$0\leq\lambda_{2}\leq 1$ controls the type of regularisation. $\lambda_{2}=1$
and $\lambda_{2}=0$ correspond to $l_{1}$ (lasso) and $l_{2}$ (ridge)
regularisation respectively. This regularisation is called elastic net
regularisation, and it promotes small values of $\bar{\theta}$. The $l_{1}$
term values promotes sparsity in $\bar{\theta}$, while the $l_{2}$ term
penalises large values. Importantly as $H$ is included in $\bar{\theta}$ this
regularised objective can be used for input selection as a measurement is
unimportant it’s corresponding entry in $H$ will be zero or close to zero.
Note that the regularisation has to be applied to all of $\bar{\theta}$ as
otherwise the network can adjust to use large values of $\theta$ to compensate
for the small values in $H$. After training one can select the important
inputs by selecting the entries of $H$ larger than some magnitude. Lastly,
note that although this heuristic approach can be used to find a policy that
uses $n$ measurements, there are no guarantees that the best $n$ measurements
will be selected. However, due to the intractability of the exact measurement
selection problem this is an attractive heuristic that has been used in
similar context [20].

##  4 Distillation problem formulation

In this work we consider the separation of an ideal binary mixture by a
distillation column of 25 theoretical stages (including reboiler and total
condenser, shown in fig. 1). The desired product composition (for the low
boiling component) is 0.99 mole fraction in the top and 0.01 in the bottom. We
use the model of Skogestad [21], with parameters shown in Table 1. Although
the model is relatively simple, it captures the typical dynamic behaviour of a
distillation column. For completeness the model is briefly described in
Section 4.1. The major assumptions of the model are: constant molar flow
rates, constant pressure, constant relative volatility, vapour liquid
equilibrium, negligible vapour holdup, and linearised liquid dynamics. In
section 4.2 aspects of the control problem formulation are described.

Importantly apart from the constraints on the control usage, when integrating
the system the inequality and equality constraints of the system are naturally
satisfied by the physics of the model. This removes the issue of ensuring
constraint satisfaction of (1d).

Figure 1: Sketch of a distillation column with LV-configuration and external
flows shown.

###  4.1 Distillation model

The model primarily consists of total and component material balances for each
stage of the column. The stages are indexed by $i$, with the bottom stage
assigned $i=1$, and the top $i=N_{T}$. The balance equations are described
first for the trays, then the condenser and last the reboiler. This is
followed by the consititutive equations of the column.

####  4.1.1 Material balances equations

Excluding the feed stage, each tray in the column is described by the total
mole balance:

$\frac{dM_{i}}{dt}=L_{i+1}-L_{i}+V_{i-1}-V_{i}$ (11)

where $M_{i}$ is the liquid holdup on stage $i$, $L_{i}$ is the liquid flow
from stage $i$, and $V_{i}$ is the liquid flow from stage $i$. Additionally,
the material balance for the lighter component is:

$\displaystyle\frac{dM_{i}x_{i}}{dt}$
$\displaystyle=L_{i+1}x_{i+1}+V_{i-1}y_{i-1}-L_{i}x_{i}-V_{i}y_{i}$ (12a)
which combined with (11) yields: $\displaystyle\frac{dx_{i}}{dt}$
$\displaystyle=\frac{L_{i+1}x_{i+1}+V_{i-1}y_{i-1}-L_{i}x_{i}-V_{i}y_{i}-x_{i}(L_{i+1}-L_{i}+V_{i-1}-V_{i})}{M_{i}}$
(12b)

where $x_{i}$ and $y_{i}$ are the mole fraction of the lighter component in
the liquid and vapour phases on stage $i$.

Assuming that the feed is mixed directly into feed stage, the total and
component mass balances of the feed stage are:

$\displaystyle\frac{dM_{N_{F}}}{dt}=$ $\displaystyle
L_{{N_{F}}+1}-L_{N_{F}}+V_{{N_{F}}-1}-V_{N_{F}}+F$ (13a)
$\displaystyle\frac{dx_{N_{F}}}{dt}=$
$\displaystyle\frac{L_{{N_{F}}+1}x_{{N_{F}}+1}+V_{{N_{F}}-1}y_{{N_{F}}-1}-L_{N_{F}}x_{N_{F}}-V_{N_{F}}y_{N_{F}}+Fz_{F}}{M_{N_{F}}}$
$\displaystyle-
x_{N_{F}}\frac{L_{{N_{F}}+1}-L_{N_{F}}+V_{{N_{F}}-1}-V_{N_{F}}+F}{M_{N_{F}}}$
(13b)

where $F$ is the feed, $N_{F}$ is the feed stage, and $z_{F}$ is the mole
fraction of the lighter component in the feed. Similarly, the total and
component balances of the condenser are:

$\displaystyle\frac{dM_{N_{T}}}{dt}$
$\displaystyle=-L_{N_{T}}+V_{{N_{T}}-1}-D$ (14a)
$\displaystyle\frac{dx_{N_{T}}}{dt}$
$\displaystyle=\frac{V_{N_{T}}y_{{N_{T}}-1}-L_{N_{T}}x_{N_{T}}-Dx_{N_{T}}}{M_{N_{T}}}-x_{N_{T}}\frac{-L_{N_{T}}+V_{{N_{T}}-1}-D}{M_{N_{T}}}$
(14b)

where $D$ is the distillate flow rate. And lastly, the total and component
balances of the reboiler are:

$\displaystyle\frac{dM_{1}}{dt}$ $\displaystyle=L_{2}-V_{1}-B$ (15a)
$\displaystyle\frac{dx_{1}}{dt}$
$\displaystyle=\frac{L_{2}x_{2}-V_{1}y_{1}-Bx_{1}}{M_{1}}-x_{1}\frac{L_{2}-V_{1}-B}{M_{1}}$
(15b)

where $B$ is the bottoms flow rate.

####  4.1.2 Constitutive equations

By assumption, the liquid flow rates are governed by linearised dynamics:

$\displaystyle L_{i}$
$\displaystyle=L^{0}_{i}+\frac{M_{i}-M^{0}_{i}}{\tau_{l}}+\lambda(V_{i-1}-V^{0}_{i-1}),\qquad
i=1,\dots,N_{T-1}$ (16a) $\displaystyle L_{NT}$ $\displaystyle=L_{T}$ (16b)

where $L_{i}^{0}$ and $M_{i}^{0}$ are nominal values of the liquid flow and
holdup on stage $i$, $\tau_{l}$ is the time constant for the liquid flow
dynamics and $\lambda$ describes how the vapour flow rate influences the
liquid flow rate (the K2-effect [21]). The liquid flow rate from the top tray
is simply given by the reflux flow rate, $L_{T}$.

By the assumption of constant molar flows and negligible vapour hold-up the
vapour flow rates are given by:

$\displaystyle V_{i}$ $\displaystyle=V_{i-1}\qquad
i=2,\dots,N_{T-1},\>i\not=N_{F}$ (17a) $\displaystyle V_{N_{F}}$
$\displaystyle=V_{N_{F}-1}+(1-q_{F})F$ (17b) $\displaystyle V_{1}$
$\displaystyle=V_{B}$ (17c)

where $V_{B}$ is the boilup flow rate, and $q_{F}$ is the liquid fraction of
the feed.

The vapour liquid equilibrium is assumed to given by:

$y_{i}=\frac{\alpha x_{i}}{1+(\alpha-1)x_{i}}$ (18)

where $\alpha$ is the relative volatility.

To calculate tray temperatures a linear relationship based on the pure
component boiling points is assumed:

$T_{i}=x_{i}T_{b,L}+(1-x_{i})T_{b,H}$ (19)

where $T_{b,L}$ and $T_{b,H}$ are the boiling points of the light and heavy
components.

Table 1: Summary of nominal column parameters. Parameter | Description | Value
---|---|---
$F$ | Feed rate [kmol/min] | 1.0
$z_{F}$ | Feed composition | 0.5
$q_{F}$ | Feed liquid fraction | 1.0
$\alpha$ | Relative volatility | 1.75
$\tau_{L}$ | Time constant for liquid flow dynamics [min] | 0.063
$\lambda$ | Constant describing the K2-effect | 0.0
$M^{0}_{i}$ | Nominal liquid holdup[kmol] | 0.5
$L_{i}^{0}$ | Nominal liquid flow rate[kmol/min] | $\begin{cases}3.564i\leq N_{F}\\\ 2.564i<N_{F}\end{cases}$
$V_{i}^{0}$ | Nominal vapour flow rate [kmol/min] | 3.065
$T_{b,L}$ | Light boiling point [K] | 341.9
$T_{b,H}$ | Heavy boiling point [K] | 357.4
$K_{D}$ | P-controller tuning for distillate | 10.0
$K_{B}$ | P-controller tuning for distillate | 10.0
$D^{0}$ | Nominal distillate flow [kmol/min] | 0.5
$B^{0}$ | Nominal boilup flow [kmol/min] | 0.5

Lastly, we specify that the column is operated in LV-configuration which means
that the liquid levels in the reboiler and condenser ($M_{1}$, $M_{N_{T}}$)
are controlled by the product flows, $B$ and $D$:

$\displaystyle D$ $\displaystyle=D_{s}+K_{D}(M_{N_{T}}-M^{0}_{N_{T}})$ (20a)
$\displaystyle B$ $\displaystyle=B_{s}+K_{B}(M_{1}-M^{0}_{1})$ (20b)

The reflux and boilup, $L_{T}$ and $V_{B}$, remain as control variables hence
the name LV-configuration. This is a reasonable assumption as it is the
“conventional” choice for distillation columns [22]. Note that although this
controller stabilises the liquid levels, the column itself remains unstable.

###  4.2 Control problem formulation

In the following we describe details of the control problem – the objective,
constraints, assumed noise of the measurements and MPC parameters.

####  4.2.1 Objective and constraints

As the objective we consider regulating the product and distillate
compositions to their set-points with a small penalty on moving the inputs
from the nominal values:

$\displaystyle J=\int_{t_{0}}^{t_{f}}$
$\displaystyle(x_{1}(t)-0.01)^{2}+(x_{N_{T}}(t)-0.99)^{2}$
$\displaystyle+0.001\left(\left(V_{B}(t)-V_{B}^{0}\right)^{2}+\left(L_{T}(t)-L_{T}^{0}\right)^{2}\right)\>dt$
(21)

In addition we impose the following inequality constraints:

$\displaystyle 0$ $\displaystyle\leq x(t)\leq 1$ (22a) $\displaystyle 0$
$\displaystyle\leq y(t)\leq 1$ (22b) $\displaystyle 0$ $\displaystyle\leq
M(t)$ (22c) $\displaystyle 0$ $\displaystyle\leq V_{B}(t)\leq 3.25$ (22d)
$\displaystyle 0$ $\displaystyle\leq L_{T}(t)\leq 2.75$ (22e)

Note that, apart from apart from eqs. 22d and 22e, these constraints are
implicitly enforced by the system dynamics, i.e. adaptive time-stepping of a
differential equation solver can ensure that the constraints are satisfied.
This is similarly the case if the MPC problem is solved with single shooting
[23].

####  4.2.2 Disturbances, measurement noise and operating region

We assume that potential disturbances, and the ranges, are the column feed
flow-rate $[0.8\ 1.2]$, feed composition $[0.4\ 0.6]$ and feed liquid fraction
$[0.8\ 1.0]$. Using these disturbance ranges we generate a multi-level pseudo-
random sequence of disturbances, and by simulating the distillation column
controlled by MPC (see section 5.1) find the effective operating region of the
column. 100 disturbances are used with the time between disturbances randomly
selected as one of 10 linearly spaced levels between $0.5$ and $10$ minutes,
yielding a total simulation time of 485 minutes. At each disturbance time-
point one of the disturbances are randomly selected and changed to one of 15
levels equally spaced between the upper and lower bounds of the respective
disturbance. The column is initialised with the nominal feed composition and
nominal liquid holdups on each stage and after 15 minutes the disturbance
sequence starts. The same approach is used to generate a disturbance sequence
to test different controllers against each other, this sequence is shown in
Figure 2.

Figure 2: Disturbance profile used in the comparison of the control policies.
Figure 3: Violin plot showing the temperature range of the distillation column
subject to disturbances, to find the relevant operating region of the column.

The resulting temperature operating range of the column is shown in Figure 3.
Due to the controller the ends of the column are not significantly influenced
by the disturbances while in the middle of the column there is considerable
variation. Importantly, from this data (and the corresponding hold-up data)
one can see that a very small region of the feasible operating space is
visited in standard operation. This reduction of the feasible state space is
key aspect of computational feasibility when considering high dimensional
systems.

We assume that there are 30 candidate column measurements available: the
temperatures on each stage ($T_{1:N_{T}}$), feed flow rate ($F$), feed
temperature ($T_{F}$), feed liquid fraction ($q_{F}$), and the liquid holdups
in the reboiler and condenser ($M_{1}$ and $M_{N_{T}}$). Because composition
measurements are unreliable and typically subject to delays, temperature
measurements are normally used instead. All measurements are assumed to be
subject to measurement noise as summarised in Table 2 where
$\mathcal{N}(\mu,\sigma)$ denotes the normal distribution with mean and
standard deviation $\mu$ and $\sigma$, and the distributions are truncated at
the specified min and max values.

Table 2: Summary of assumed measurement noise Measurement | Distribution | (min, max)
---|---|---
Temperature | $\mathcal{N}(0,\ 0.015)*(T_{bH}-T_{bL})$ | $(-0.775,0.775)$
Flow rate | $\mathcal{N}(0,\ 0.03)$ | $(-0.1,0.1)$
Liquid fraction | $\mathcal{N}(0,\ 0.03)$ | $(-0.1,0.1)$
Holdup | $\mathcal{N}(0,\ 0.03)$ | $(-0.1,0.1)$

##  5 Controller optimisation and training

The proposed method is implemented in Julia 1.7, with major use of the
following packages: Flux.jl [24], Zygote.jl [25], DifferentialEquations.jl
[26], and JuMP.jl [27]. For optimisation of the MPC and neural networks Ipopt
[28] and RMSProp [29] are used.

The control policies we consider are summarised in Table 3. We construct an
MPC policy $\kappa_{mpc}$ and four neural network policies. The first two
policy are trained using all of the available measurements (note this excludes
the tray holdups). We consider training this controller with and without
measurement noise included in the training, yielding $\kappa_{all}$ and
$\kappa_{all}^{nonoise}$ respectively. The second policy, $\kappa_{reg}$, also
has all the available measurements, but is trained with the elastic net
regularisation term (10). The training will thus select a reduced set of
inputs that the policy will use. The third policy, $\kappa_{sel}$, is trained
with only four user specified measurements as inputs using the unregularised
objective.

The inputs of $\kappa_{sel}$ are chosen as
$\zeta=[\bar{T}_{5},\bar{T}_{10},\bar{T}_{16},\bar{T}_{21}]$. These
temperatures are chosen in emulation of classic control configurations for
distillation columns. In such a control scheme instead of controlling the
product purities directly, the deviation of a temperature in the bottom and
top section of the column from some setpoint is controlled. This is because
the temperatures away from the top and bottom are more sensitive to changes in
the inputs than the product temperatures [21, 22]. Although we don’t directly
specify a set-point for these trays, we also do not give $\kappa_{sel}$ the
product temperatures. Any control policy that $\kappa_{sel}$ defines must
implicitly transform and combine the temperatures and regulate this
combination to some set-point.

After training of the neural network control policies, all the control
policies are tested on their performance subject to the disturbance sequence
shown in Figure 4.2.2. This comparison is detailed in section 6.

Table 3: Summary of controllers compared in the distillation case study Controller | Objective | Controller inputs | Training noise
---|---|---|---
$\kappa_{mpc}$ | Discretised form of (4.2.1) | 55 (all states) | N/A
$\kappa_{all}^{no\ noise}$ | (4.2.1) | 30 | No
$\kappa_{all}$ | (4.2.1) | 30 | Yes
$\kappa_{reg}$ | (4.2.1) with elastic net penalty | 30 | Yes
$\kappa_{sel}$ | (4.2.1) | 4 temperatures | Yes

###  5.1 The benchmark MPC policy

We construct a MPC policy to be used a base-line, “best-case” control policy.
As such we conservatively assume that the MPC has perfect full state feedback.
In a more realistic comparison a state estimator would be used with the noisy
measurements that are provided to the other control policies. In the MPC
problem we use implicit Euler to discretise the dynamic system. Note that in
the discretised problem the objective and constraints are only evaluated at
the discrete time points. A discretisation time of 30 seconds, and control
horizon of 20 minutes are chosen, thus there are $N_{H}=41$ points in the time
discretisation. The temperature profile of the distillation column when
controlled by MPC is shown in Figure 4. Note that although the column
temperatures are in continuous time, the control output of the MPC is
piecewise-constant on 30 second intervals.

Figure 4: Temperature profiles of distillation column controlled by MPC with
test disturbances. The dashed black lines indicate the pure component boiling
points, and the red dash-dotted line indicates the temperature in the reboiler
and condenser.

###  5.2 Training control policies

####  5.2.1 Specification of neural policies

We consider optimisation of a two layer feed-forward neural network control
policy of the form:

$\displaystyle[L_{T},\ V_{B}]$ $\displaystyle=\frac{1}{2}[2.75,\
3.25]\circ(1+\kappa_{NN}(H\zeta,\theta))$ (23a)
$\displaystyle\alpha^{(i)}(\zeta)$ $\displaystyle=(1+e^{-\zeta})^{-1},\qquad
i=0,1$ (23b) $\displaystyle\bar{\theta}$ $\displaystyle=\texttt{vec}(H,\
\theta)$ (23c)

where $\circ$ is the element-wise product and $H\in\mathbb{R}^{n_{\zeta}\times
n_{\zeta}}$ is a diagonal matrix. As discussed in section 3.2, the use of $H$
with a regularised objective promotes a sparse selection of the inputs during
training. For consistency, when comparing between the different policies we
still use $H$ when training with the un-regularised objective. The sigmoid
function is used as the activation function for all layers of the network
including the neural network output. This choice of activation function means
that the controller inequality constraints (22) are satisfied for all
$\bar{\theta}$.

For $\kappa_{all}$ and $\kappa_{reg}$ the policy takes the input vector:

$\displaystyle\zeta$ $\displaystyle=[\bar{T}_{1},\
\dots,\bar{T}_{N_{T}},F,\bar{T}_{F},q_{F},M_{1},M_{N_{T}}]$ (24a)
$\displaystyle\bar{T}_{i}$ $\displaystyle=\frac{T_{i}-T_{bL}}{T_{bH}-T_{bL}}$
(24b)

where $\bar{T}$ is the normalised temperature. A normalised temperature is
used as then all elements of the input vector are the same magnitude. For both
$\kappa_{all}$ and $\kappa_{reg}$ we use a network of width 30, which means
that there are $\sim 1000$ network parameters.

The inputs of $\kappa_{sel}$ are chosen as
$\zeta=[\bar{T}_{5},\bar{T}_{10},\bar{T}_{16},\bar{T}_{21}]$, in emulation of
classic control configurations for distillation columns. We use the same
structure of the network but with only 4 nodes in the input layer. To keep the
network a similar size ($\sim 1000$ parameters) we chose a width of 150.

We note that any control policy found by $\kappa_{sel}$ is within the scope of
$\kappa_{reg}$ which is itself within the scope of $\kappa_{all}$. Thus one
may expect that $\kappa_{all}$ should perform the best. However, this is not
entirely to be expected. Although inclusion of extra information can
theoretically only be used to improve the control, it is well known that
despite neural networks being universal approximators they can have worse
performance without judicial selection and manipulation of inputs (“feature
engineering”).

####  5.2.2 Embedded training of the control policy

The control policies are embedded in the dynamical system and trained in the
optimise-and-learn formulation (5) using the control objective (4.2.1) (and
regularisation term for $\kappa_{reg}$). RMSProp [29] is used for training,
and was chosen by benchmarking several choices from Flux.jl [24] on the
problem.

Due to the relatively high dimensionality of the system (50 states) and as the
training was performed on a laptop (16 GB RAM, i5-101310U processor) to
perform the training efficiently we followed the following strategy:

  1. 1.

$\pi_{z_{0}}$ is estimated from the closed-loop data.
As described in section 4.2.2 we simulate the column under control of the MPC
policy to find the typical operating region of the column (see Figure 3).
Using this data we fit a multi-variate normal distribution to the temperatures
of the column and a separate multi-variate normal distribution to the hold-ups
of the column222Separation of the temperatures and hold-ups is conservative,
and is to ensure that the controller doesn’t somehow learn the initial hold-up
from the initial temperature profile.. As a compact probability distribution
is required these are truncated at $\pm 3$ standard deviations from the mean.

  2. 2.

$\pi_{z_{0}}$ is discretised to yield $\pi^{d}_{z_{0}}$ by taking 1000 quasi-
random samples.
We wish to use quasi-random samples and not random samples as they cover the
domain more evenly at low number of samples. Generating quasi-random samples
with correlations is difficult to do directly. We first use a Sobol sequence
to generate points from the marginalised cumulative probability distribution
of each state. Correlations between these points are then induced by the Iman-
Connover method [30].

  3. 3.

Noise is similarly sampled from the assumed distributions to give
$\pi^{d}_{\eta}$
Note that as the noise is uncorrelated, correlations are not induced by the
Iman-Connover method.

  4. 4.

RMSProp is used to train the controller. At each iteration $n_{s}$ points are
independently sampled from $\pi^{d}_{z_{0}}$ and $\pi^{d}_{\eta}$.
For the first 2000 iterations we use $n_{s}=1$, and then perform 750
iterations with $n_{s}=2$ (with $w_{s}=0.5$). The number of iterations is
chosen by picking a point on a plot of the objective function versus the
number of iterations at which the objective stopped improving. The last 750
iterations serve to reduce the variability of the estimate of
$\bar{\theta}^{*}$. The objective gradient is calculated by reverse-mode
automatic-differentiation on the differential equation solver [25].

Lastly if a regularised objective is used then after 4. the entries of $D$
that are smaller than 0.001 are set to zero. Then the training is repeated
without regularisation with $D$ frozen to its new values, and $\theta$ set
back to the initial guess. This is because the elastic net regularisation also
influences the expressivity of the network, which is an undesired behaviour.

##  6 Results

We now compare the policies trained in the previous section, on controlling
the column subject to the disturbance sequence shown in Figure 4.2.2. Section
6.1 compares the nominal performance of the policies with and without
measurement noise. In Section 6.2 we compare how the policies handle model
mismatch, first with respect to the specified noise distribution, and then
with respect to the implemented control to the column. Although the control
policies are not directly interpretable, by examining their outputs we aim to
provide some interpretations of their behaviour.

###  6.1 Comparison of policies with no model-mismatch

In this section the policies are compared based on their performance of
controlling the column subject to a test sequence of disturbances, prepared in
the same way as when determining the operating region of the column. Note that
as before the column starts from initially equimolar composition on each
stage. This “start-up” period is outside of the operating region used in the
policy optimisation (see Figure 3). Thus, the control policies $\kappa_{all}$,
$\kappa_{reg}$, and $\kappa_{sel}$ will perform poorly in this region. The
performance of the control policies are summarised in table 4, with further
details given in the following.

Table 4: Cumulative closed loop objective value from using the trained policies on the test sequence of disturbances excluding the start-up portion. Noise refers to constant measurement noise along a trajectory. Unlike the other policies, $\kappa_{mpc}$ uses full state feedback and hence does not have noise entries. | Cumulative objective
---|---
| no noise | with noise | averaged with noise
$\kappa_{mpc}$ | 0.0076 | - | -
$\kappa_{all}$ | 0.0092 | 0.0079 | 0.0096
$\kappa_{sel}$ | 0.0087 | 0.0103 | 0.0094
$\kappa_{reg}$ | 0.0142 | 0.0175 | 0.0157

####  6.1.1 $\kappa_{all}$ and $\kappa_{mpc}$

Figure 5: Temperature profile using $\kappa_{all}$. Black dashed lines
indicated $T_{bL}$ and $T_{bH}$, red dash-doted lines indicates $T_{1}$ and
$T_{N_{T}}$, and the green lines indicate column temperatures.

We first consider the case where $\kappa_{mpc}$ and $\kappa_{all}$ are used to
control the column with no model-mismatch or noise in the measurements, with
the resulting temperature profiles shown in Figures 4 and 5. By inspection the
behaviour of the system is qualitatively similar, although there are
differences. For example, between 30-60 minutes $\kappa_{mpc}$ is able to
stabilise the system very quickly, while $\kappa_{all}$ has a much larger
transient response.

This is to be expected as with no measurement noise or model error
$\kappa_{mpc}$ is able to exactly compensate for (some) incoming disturbances.
On the other-hand, $\kappa_{all}$ does not use the hold-ups of the internal
stages of the column (which the MPC does), and thus worse performance should
be expected. However the actual difference in the control objective is minor –
the cumulative objective value (excluding the start-up portion) is 0.0076 for
$\kappa_{mpc}$ and 0.0092 for $\kappa_{all}$ (table 4).

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 6: Closed loop profile of $\kappa_{all}$ with normally distributed
measurement noise. In fig. 6(a) a constant realisation of measurement noise is
used across the horizon. In fig. 6(b) the solid line indicates the nominal
control output, while shaded region shows the range of possible control
outputs due to normally distributed noise of the measurements.

To gain further understanding of the behaviour of $\kappa_{all}$ we can look
at the actual control output and how this is influenced by measurement noise,
as shown in Figure 6. An interesting thing to note in that in the control
output of $\kappa_{all}$ (fig. 6(b)) there are sharp spikes at the points of
some disturbances – this means that the controller is directly making use of
the feed (disturbance) measurements for feed-forward control. Although this
comes with improved performance it does mean that the controller can be more
influenced by model error. Despite this, fig. 6(a) shows that in the presence
of a constant measurement bias $\kappa_{all}$ is able to maintain reasonable
control of the system.

To demonstrate the role of including noise in the training process, we train a
controller set-up equivalently to $\kappa_{all}$ but without noise in the
training, $\kappa^{no\ noise}_{all}$. Figure 7 shows the closed loop
performance of using $\kappa^{no\ noise}_{all}$. This controller output is
clearly much more sensitive to noise, and in addition, the spikes in fig. 7(b)
are much more pronounced. This difference is expected – $\kappa^{no\
noise}_{all}$ is more aggressive in its feedforward behaviour at it “expects”
perfect measurements of the disturbance and column. However, we would also
like to note that despite the higher sensitivity of $\kappa^{no\
noise}_{all}$, the controllers have very similar performance when examining
the temperature profiles.

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 7: Closed loop profiles using $\kappa^{no\ noise}_{all}$ with normally
distributed measurement noise. Green lines indicate column temperatures used
by the controller, and grey lines indicate unused column temperatures.

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 8: Closed loop profile of $\kappa_{reg}$ with normally distributed
measurement noise.

####  6.1.2 Regularised input selection: $\kappa_{reg}$

$\kappa_{reg}$ is optimised with the elastic net parameters
$\lambda_{1}=0.01$, $\lambda_{2}=0.99$. After the first stage of optimisation,
measurement selection is performed by checking for entries of $D$ larger in
magnitude than 0.001. These hyper-parameters were chosen as being
demonstrative of the kind of measurement selection that can occur, and were
not selected based on the closed loop performance. The training yielded nine
selected measurements: $T_{1}$, $T_{3}$, $T_{5}$, $T_{6}$, $T_{7}$, $T_{10}$,
$T_{13}$, $F$, and $q_{F}$. Upon first glance this selection may seem strange
– the seven temperature measurements are from the bottom section of the
column, and the other two measurements relate to the feed.

The repeated bottom temperature measurements may allow the controller to
combine measurements for less sensitivity to noise, but this doesn’t explain
why only bottom measurements are used. Note that as the feed is mainly liquid
(even with the disturbances to $q_{F}$), and upon entering the column most of
the feed immediately proceeds into the bottom section. Combined with the
liquid flow dynamics this means that multiple temperature measurements in the
bottom section of the column provides a short term “record” of previous
disturbances, and are thus very informative of the overall column state
(assuming that there is no model mismatch). By this reasoning one should
expect $\kappa_{reg}$ behaviour substantially in a feedforward manner as the
disturbances are measured exactly, and indirectly the column state and
disturbance history are measured. Considering that the training problem is to
reject disturbances while penalising the use of many measurements, this
behaviour is not surprising as feedforward control will yield tighter control
than pure feedback control in the absence of model error.

The closed loop simulations of $\kappa_{reg}$ is shown in Figure 8. The
feedforward nature of $\kappa_{reg}$ can be seen by examining the nominal
control profile (Figure 8(b)) which varies very slowly between disturbances,
while having sharp changes whenever there is a large disturbance in the feed.
In general, $\kappa_{reg}$ shows much less variation than $\kappa_{all}$ both
in the nominal control response, and under the influence of noise.

The temperature profile has relatively sharp and rapid changes compared to
compared to Figure 6. Visually in some places $\kappa_{reg}$ appears to
perform better than $\kappa_{all}$, however the cumulative control objective
over the test horizon using $\kappa_{reg}$ is worse than using $\kappa_{all}$
(see table 4).

####  6.1.3 Control policy with manual input selection: $\kappa_{sel}$

Unlike $\kappa_{reg}$ where we regularise for input selection, for
$\kappa_{sel}$ the inputs are chosen as
$\zeta=[\bar{T}_{5},\bar{T}_{10},\bar{T}_{16},\bar{T}_{21}]$. As discussed
earlier, the motivation behind this highly reduced set of temperature
measurements is that these are similar temperatures that may be used in
classic distillation control strategies. In addition note that this removes
the possibility of the controller performing feedforward control on the
disturbances. The resulting temperature and control profiles using
$\kappa_{sel}$ are shown in Figure 9.

In general there is more variation in the control output than $\kappa_{reg}$ ,
but less than $\kappa_{all}$. However, the temperature profiles are much less
sensitivity to noise than $\kappa_{reg}$. Interestingly, the temperature and
control profile of $\kappa_{sel}$ and $\kappa_{all}$ are very similar. In
addition, the cost of using these controllers is very similar (table 4). In
particular, $\kappa_{sel}$ seems to be giving a smoother version of
$\kappa_{all}$, as it lacks the initial kick from the feed measurements
provided to $\kappa_{all}$. What this suggests is that $\kappa_{sel}$ contains
the key elements of the “feedback part” of $\kappa_{all}$ while neglecting the
feedforward part on the disturbances. Due to the time scale, the slightly
slower response of $\kappa_{sel}$ compared to $\kappa_{all}$ does not make a
significant difference to the cumulative objective.

(a) Temperature profile (constant bias)

(b) Controller output (varying noise).

Figure 9: Closed loop profile of $\kappa_{sel}$ with normally distributed
measurement noise.

###  6.2 Controller sensitivity to mismatch

It is well known that optimal control policies can potentially be sensitive to
mismatch between the model used in their solution and the true system. In this
section we compare the controllers when: (1) the measurement noise is
mispecified and (2) when there is a multiplicative output disturbance. We find
that of the trained controllers $\kappa_{sel}$ is much more robust to these
disturbances and is better able to regulate the system. However, as all
control schemes do not have an integrator, offset free control is not
achieved.

####  6.2.1 Mispecification of noise distribution

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 10: Closed loop profile of $\kappa_{all}$ with extreme noise
realisations.

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 11: Closed loop profile of $\kappa_{reg}$ with extreme noise
realisations.

(a) Temperature profile (constant bias)

(b) Controller output (varying noise)

Figure 12: Closed loop profile of $\kappa_{sel}$ with extreme noise
realisations.

figs. 10, 12 and 11 show the closed loop response of the system when instead
of normally distributed noise, the noise realisations are exclusively at the
extremes of the truncated probability distributions. This is an incredibly
unlikely possibility with respect the normally distributed noise and serves as
an example of severe mispecification of the noise characteristics of the
system.

$\kappa_{all}$, fig. 10, has the most variable controller output due to noise
and also has the worst temperature regulation. A potential explanation is that
as $\kappa_{all}$ has access to many measurements $\kappa_{all}$ may
compensate for the noise by combining measurements together. When the
assumption behind the noise is false then this compensation does not work well
and results in a large variance in the controller output, as the controller
did not learn to be less sensitive to variations. Additional support for this
is noting that $\kappa_{reg}$ and $\kappa_{sel}$ are much less sensitive to
noise when examining both the controller output and temperature profiles.

Comparing $\kappa_{reg}$ and $\kappa_{sel}$ it is interesting to note that
despite $\kappa_{reg}$ showing more variation with noise than $\kappa_{sel}$
(figs. 9(b) and 11(b)), in terms of regulation $\kappa_{sel}$ performs better.
As noted in the previous section this is likely because $\kappa_{reg}$
primarily makes changes based on the feed measurements, and only small changes
based on the temperatures.

####  6.2.2 Multiplicative model error

(a) $\kappa_{mpc}$

(b) $\kappa_{all}$

(c) $\kappa_{reg}$

(d) $\kappa_{sel}$

Figure 13: Closed loop temperature profiles of the system, with a constant
multiplicative output disturbance.

We consider the same model, but assume that due to model-mismatch the reflux
and boil up are respectively 10% more and less than what the policy specifies,
i.e.

$[L_{T},\ V_{B}]_{actual}=[1.1L_{T},\ 0.9V_{B}]$ (25)

This is a considerably large error, and as the control policies considered in
this work are static and not integrating we can reasonably expect significant
steady state offset. However, offset free control is not the primary goal,
instead we can see how robust the policies are at managing this large,
unmodelled disturbance.

The closed loop response of the system, without any measurement noise, and
using policies $\kappa_{mpc},\ \kappa_{all},\ \kappa_{reg},$ and
$\kappa_{sel}$ is shown in fig. 13. Immediately it is clear that the best
performing controller is $\kappa_{mpc}$ followed by $\kappa_{sel}$. The good
performance of $\kappa_{sel}$ is unsurprising given the previous results, as
the controller is forced to perform feedback control and robustness is
exacerbated by feedforward control.

##  7 Conclusion

We have proposed and demonstrated the closed loop training of (static)
measurement based control policies on a distillation case study. This case
study is significantly larger compared to prior literature. The key features
of the proposed approach are: (1) the controller is trained in closed loop
using an optimise-and-learn approach, (2) the controller is trained in an
operationally relevant region of the state-space to reduce the computational
demands of the training, and (3) using this approach the controller can be
trained to use a selection of measurements. Three controllers are trained
using different selections of measurements, and their performance is compared
for the nominal system (with and without measurement noise) and also for a
system with model mismatch. A controller using only 4 measurements along the
column is able to perform well, achieving close to MPC performance on the
nominal system while being more robust than the other controllers in the
example with significant model mismatch.

####  7.0.1 Further work

There are key aspects of the approach that need to be further developed are
the selection of measurements for the controller and the incorporation of
integration in the control policy.

We have demonstrated that one can set up a training problem to automatically
select important measurements for use by a control policy. Through numerical
simulation we show that such a control policy ($\kappa_{reg}$) achieves good
nominal performance. However, this example also shows that doing so can lead
to fragile policies that perform poorly beyond the nominal system (fig. 13).
This is because regularisation for measurement selection removes measurements
that are not necessary for control of the (assumed perfect) model, which can
result in a dramatically worse performance when there is model mismatch. A
potential approach to fix this problem is to incorporate uncertainty into the
training. However, doing so requires in a meaningful manner requires an
accurate description of the uncertainty, and can also dramatically increase
the computational complexity of the training.

Another choice is to use engineering judgement to select measurements that can
be used to find good control actions ($\kappa_{sel}$). Despite using
considerably less measurements, this controller was significantly more robust
to the others in the numerical examples. Compared to classical control
techniques that pair single measurement to control variables, $\kappa_{sel}$
can be trained with multi-variable pairings without specifying the structure
of this relationship. However, the central issue of selecting which
measurements remains [31].

Without integration, a static controller cannot achieve offset free control by
itself. A potential option is to directly include the potential for
integration in the controller by allowing it to use an additional state as
“memory”. The challenge with this approach is for the controller to not over-
fit during training, and to generalise beyond the model-mismatch, measurement
bias, and similar used in the training. Alternatively, integration may be
performed separately, e.g. by a disturbance model. In this approach the
controller would be trained with the augmented model, and online the
disturbance parameters would be provided to the controller to achieve offset
free control.

## References

  * Santos et al. [2001] L. O. Santos, P. A. Afonso, J. A. Castro, N. M. Oliveira, L. T. Biegler,  On-line implementation of nonlinear mpc: an experimental case study,  Control Engineering Practice 9 (2001) 847–857.
  * Findeisen and Allgöwer [2004] R. Findeisen, F. Allgöwer,  Computational delay in nonlinear model predictive control,  IFAC Proceedings Volumes 37 (2004) 427–432.
  * Diehl et al. [2005] M. Diehl, H. G. Bock, J. P. Schlöder,  A real-time iteration scheme for nonlinear optimization in optimal feedback control,  SIAM Journal on control and optimization 43 (2005) 1714–1736.
  * Zavala and Biegler [2009] V. M. Zavala, L. T. Biegler,  The advanced-step nmpc controller: Optimality, stability and robustness,  Automatica 45 (2009) 86–93.
  * Jäschke et al. [2014] J. Jäschke, X. Yang, L. T. Biegler,  Fast economic model predictive control based on nlp-sensitivities,  Journal of Process Control 24 (2014) 1260–1272.
  * Seel et al. [2022] K. Seel, A. B. Kordabad, S. Gros, J. T. Gravdahl,  Convex neural network-based cost modifications for learning model predictive control,  IEEE Open Journal of Control Systems 1 (2022) 366–379.
  * Turan et al. [2023] E. M. Turan, Z. Mdoe, J. Jäschke,  Learning a convex cost-to-go for single step model predictive control (2023). URL: http://arxiv.org/abs/2312.02650. arXiv:2312.02650.
  * Ramchandran and Rhinehart [1995] S. Ramchandran, R. R. Rhinehart,  A very simple structure for neural network control of distillation,  Journal of Process Control 5 (1995) 115–128.
  * Parisini and Zoppoli [1995] T. Parisini, R. Zoppoli,  A receding-horizon regulator for nonlinear systems and a neural approximation,  Automatica 31 (1995) 1443–1451. URL: https://linkinghub.elsevier.com/retrieve/pii/000510989500044W. doi:10.1016/0005-1098(95)00044-W.
  * Karg and Lucia [2020] B. Karg, S. Lucia,  Efficient representation and approximation of model predictive control laws via deep learning,  IEEE Transactions on Cybernetics 50 (2020) 3866–3878. doi:10.1109/TCYB.2020.2999556. arXiv:1806.10644.
  * Kumar et al. [2021] P. Kumar, J. B. Rawlings, S. J. Wright,  Industrial, large-scale model predictive control with structured neural networks,  Computers and Chemical Engineering 150 (2021). doi:10.1016/j.compchemeng.2021.107291.
  * Bonzanini et al. [2021] A. D. Bonzanini, J. A. Paulson, G. Makrygiorgos, A. Mesbah,  Fast approximate learning-based multistage nonlinear model predictive control using gaussian processes and deep neural networks,  Computers & Chemical Engineering 145 (2021) 107174.
  * Drgoňa et al. [2022] J. Drgoňa, K. Kiš, A. Tuor, D. Vrabie, M. Klaučo,  Differentiable predictive control: Deep learning alternative to explicit model predictive control for unknown nonlinear systems,  Journal of Process Control 116 (2022) 80–92. URL: https://doi.org/10.1016/j.jprocont.2022.06.001. doi:10.1016/j.jprocont.2022.06.001.
  * Turan and Jäschke [2024] E. M. Turan, J. Jäschke,  Closed-loop optimisation of neural networks for the design of feedback policies under uncertainty,  Journal of Process Control 133 (2024) 103144. doi:10.1016/j.jprocont.2023.103144.
  * Bemporad et al. [2002] A. Bemporad, M. Morari, V. Dua, E. N. Pistikopoulos,  The explicit linear quadratic regulator for constrained systems,  Automatica 38 (2002) 3–20. doi:10.1016/S0005-1098(01)00174-1.
  * Li et al. [2022] Y. Li, K. Hua, Y. Cao,  Using stochastic programming to train neural network approximation of nonlinear mpc laws,  Automatica 146 (2022) 110665.
  * Turan and Jäschke [2022] E. M. Turan, J. Jäschke,  Designing neural network control policies under parametric uncertainty: A koopman operator approach,  IFAC-PapersOnLine 55 (2022) 392–399.
  * Hornik et al. [1989] K. Hornik, M. Stinchcombe, H. White,  Multilayer feedforward networks are universal approximators,  Neural networks 2 (1989) 359–366.
  * Cybenko [1989] G. Cybenko,  Approximation by superpositions of a sigmoidal function,  Mathematics of control, signals and systems 2 (1989) 303–314.
  * Feng and Simon [2017] J. Feng, N. Simon,  Sparse-input neural networks for high-dimensional nonparametric regression and classification,  arXiv preprint arXiv:1711.07592 (2017).
  * Skogestad [1997] S. Skogestad,  Dynamics and control of distillation columns: A tutorial introduction,  Chemical Engineering Research and Design 75 (1997) 539–562.
  * Skogestad [2007] S. Skogestad,  The dos and don’ts of distillation column control,  Chemical Engineering Research and Design 85 (2007) 13–23.
  * Biegler [2010] L. T. Biegler, Nonlinear programming: concepts, algorithms, and applications to chemical processes, Society for Industrial and Applied Mathematics, 2010\. URL: http://epubs.siam.org/doi/book/10.1137/1.9780898719383. doi:10.1137/1.9780898719383.
  * Innes et al. [2018] M. Innes, E. Saba, K. Fischer, D. Gandhi, M. C. Rudilosso, N. M. Joy, T. Karmali, A. Pal, V. Shah,  Fashionable modelling with flux,  CoRR abs/1811.01457 (2018). URL: https://arxiv.org/abs/1811.01457. arXiv:1811.01457.
  * Innes [2018] M. Innes,  Don’t unroll adjoint: Differentiating ssa-form programs,  arXiv preprint arXiv:1810.07951 (2018).
  * Rackauckas and Nie [2017] C. Rackauckas, Q. Nie,  DifferentialEquations.jl–a performant and feature-rich ecosystem for solving differential equations in Julia,  Journal of Open Research Software 5 (2017).
  * Lubin et al. [2023] M. Lubin, O. Dowson, J. D. Garcia, J. Huchette, B. Legat, J. P. Vielma,  JuMP 1.0: recent improvements to a modeling language for mathematical optimization,  Mathematical Programming Computation (2023). doi:10.1007/s12532-023-00239-3. arXiv:2206.03866.
  * Wächter and Biegler [2006] A. Wächter, L. T. Biegler,  On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming,  Mathematical Programming 106 (2006) 25–57. URL: http://link.springer.com/10.1007/s10107-004-0559-y. doi:10.1007/s10107-004-0559-y.
  * Hinton et al. [2012] G. Hinton, N. Srivastava, K. Swersky, Overview of mini-batch gradient descent, http://www.cs.toronto.edu/~tijmen/csc321/slides/lecture_slides_lec6.pdf, 2012\. Accessed: 2023–12-02.
  * Iman and Conover [1982] R. L. Iman, W.-J. Conover,  A distribution-free approach to inducing rank correlation among input variables,  Communications in Statistics-Simulation and Computation 11 (1982) 311–334.
  * Foss [1973] A. Foss,  Critique of chemical process control theory,  IEEE Transactions on Automatic Control 18 (1973) 642–652.