A Modular Approach for Diesel Engine Air Path Control Based on Nonlinear MPC

This article presents a systematic approach to realize highly dynamic control strategies for the air path of diesel engines. It is based on grey-box models of individual air path components, designed to be applied in nonlinear model predictive control (NMPC). Specifically, they are suited for algorithmic differentiation and gradient-based optimization. This modular approach allows to derive models for a variety of complex air path systems and can be identified with a low amount of measurement data. An NMPC structure, which is based on these models, enables the tracking of arbitrary air path reference signals and allows the introduction of further control objectives for overactuated systems. We demonstrate our approach’s general applicability and effectiveness on two different laboratory engines using rapid control prototyping hardware. For a turbocharged light-duty diesel engine with dual-loop exhaust gas recirculation (EGR), we apply the proposed control structure to track intake manifold gas conditions while simultaneously minimizing engine pumping losses. Experimental results show an excellent tracking performance and reduced pumping losses compared to other control strategies. For a turbocharged heavy-duty engine with high-pressure EGR, we experimentally demonstrate a superior tracking performance to that obtained with a reference controller. Due to the modular and systematic approach, the algorithm design is straightforward, and the experimental calibration effort is low.


I. INTRODUCTION
E NVIRONMENTAL concerns have pushed research to improve internal combustion engines (ICEs) in recent years and demand further development in the near future. With the aim to reduce fuel consumption and pollutant emissions, engine systems continuously increase in complexity. This specifically applies to the engine air path, as the combustion process of an ICE is highly sensitive to the gas conditions Manuscript  in the intake manifold. With the application of multiple turbochargers and exhaust gas recirculation (EGR) paths, modern applications allow these conditions to be set within a wide range and in a highly dynamic way. However, due to its complexity, the control of the air path is a challenging task. Multiple actuators are used to simultaneously control the conditions of the gas inside the intake manifold while respecting various component and actuator limits. In addition, the system is highly nonlinear and subject to cross-couplings. Therefore, the application of sophisticated methods to control the engine air path is a well-researched topic. Model-based control strategies, especially those based on model predictive control (MPC) principles, have been proposed to reduce the engine calibration effort and increase the control performance. Early publications proposed explicit MPC methods due to the limited computation power of embedded processors [1], [2], [3]. The control task investigated was the tracking of the intake manifold pressure and oxygen concentration of a diesel engine with a variable turbine geometry (VTG) turbocharger and a single EGR path. Later publications applied explicit MPC or linear MPC methods with a series of linear models to control the air paths of two-stage turbocharged diesel and gasoline engines [4], [5]. These air paths are overactuated systems because most intake manifold pressure set values can be achieved with various combinations of the two turbine bypass valve positions [6]. Another control task addressing an overactuated system is the tracking of the intake manifold pressure and oxygen concentration of a turbocharged diesel engine with dual-loop EGR. In [7] and [8], it is approached using a linear MPC algorithm with a series of models derived by linearizing a physics-based nonlinear model. While in the publications mentioned so far, the experimental validation was performed on processors with increased computational power, recent development efforts have made it possible to deploy linear MPC on production engine control units (ECUs). This was demonstrated in [9] for the air path of three gasoline engines with a turbocharger and variable camshaft timing.
Concerning linear MPC methods, the actuator signals derived are only optimal for the underlying linear model. Hence, any nonlinear effect present within the system is neglected, which can considerably affect the control performance. Therefore, nonlinear MPC (NMPC) has been proposed as a method with the potential to outperform linear MPC methods for air path control. The classic task to control the intake manifold gas pressure and oxygen concentration of a turbocharged diesel engine with EGR using NMPC is first addressed in a simulative study in [10]. An experimentally validated real-time feasible NMPC algorithm for the same control task is presented in [11]. An experimentally validated NMPC approach based on a physics-based model, which addresses the overactuated control task of tracking the intake manifold pressure of a gasoline engine with a two-stage turbocharger, is presented in [12]. An extended version of this algorithm, which also includes the control of the oxygen concentration using low-pressure EGR, is presented in [13].
For overactuated air path systems, experimentally validated approaches proposed in the scientific community either use the additional degree of freedom to reduce actuator movement [5], [12] or they introduce additional reference signals [4], [7], [8], [9], [13]. For the latter approach, the reference signals are derived offline under consideration of additional control objectives, such as the reduction of pumping losses. However, the determination of additional reference values in real-time or the direct integration of economic objectives (economic MPC) would offer extended tuning possibilities and is expected to improve the controller performance. Therefore, economic MPC formulations are a major topic of current engine control research [14]. An economic NMPC approach for a naturally aspirated gasoline engine with high-pressure EGR is experimentally demonstrated in [15]. Furthermore, simulative studies have proposed such controllers for overactuated air paths [16], [17], [18]. However, to the best of the authors' knowledge, experimentally validated control structures for overactuated air path systems have not been published as yet.
Concerning experimentally validated NMPC air path control, the literature presents the algorithms tailored to a specific control task, based on prediction models specifically designed for the algorithms proposed. Any reformulation of the control task or any changes to the engine air path configuration is thus associated with a labor-intensive redesign of the control algorithm and the underlying model.
Within this article, we focus on the control of complex air path systems of diesel engines. We present a library of grey-box models of common diesel air path components, which is specifically designed for algorithmic differentiation and gradient-based optimization. These components can be combined flexibly to represent a variety of air path configurations. Model parameters are identified for each component separately, yielding a series of small identification tasks that can be carried out with only a small set of measurement data. The resulting air path models achieve a high accuracy. Since the components are modeled individually, almost all measurable physical signals are included in the model and are available for the control design.
The resulting models are used within a control structure based on NMPC, which is designed to track reference signals under consideration of additional economic objectives and engine operation and component limits. A small set of algorithm tuning parameters allows to intuitively adjust the control behavior both in simulation and on the actual system in real time. Due to the high model accuracy and the possibility to tune these parameters in simulation, the calibration effort on the engine test bench is kept low. The control structure is  I   CHARACTERISTICS OF THE TWO TURBOCHARGED DIESEL RESEARCH  ENGINES AND LABORATORY MEASUREMENT HARDWARE efficiently implemented on a rapid control prototyping (RCP) unit using the open-source software package acados [19]. As a proof of concept and motivation for our approach, we investigate two engines that differ in size and air path configuration. For each of these engines, we design and evaluate two controllers based on the model library and control structure presented. A first controller is designed to achieve the best possible control performance given the computational power of the RCP. It is used as a benchmark controller and contains an increased number of tuning parameters, which can be adjusted in real time. We then derive a second controller of lower computational complexity and compare it with the benchmark controller in terms of computation time and control performance. This second controller can be seen as a first step toward developing a real-time capable controller for a production-type ECU with significantly lower computational power.
The remainder of this article is structured as follows.
The key characteristics of the experimental setup for both engines are described in Section II. Section III presents the air path component model library, while Section IV describes the workflow for the parameter identification. The control algorithm is defined in Section V, and its application to the two engines is presented in Sections VI and VII. Finally, Section VIII concludes this article and presents an outlook.

II. EXPERIMENTAL SETUP
The experimental controller evaluation is conducted on two laboratory diesel engines, further referred to as engines A and B. Table I lists their key characteristics and the measurement equipment used for this study. In addition, their air paths are schematically shown in Figs. 1 and 2.
The CO 2 concentration is measured in the intake and exhaust manifold and used to calculate the corresponding oxygen mass fractions. For engine A, an additional oxygen sensor allows to measure the oxygen mass fraction before the compressor.

III. MODEL LIBRARY
The core of the proposed control approach consists of mean-value models for common diesel air path components.  Depending on the arrangement of these components, various air path configurations can be represented.
In general, the components consist of grey-box models, which are based on physical principles whenever possible. Doing so reduces the calibration effort, facilitates an intuitive understanding, and allows model extrapolation outside of the identification data range to a certain extent. If physical principles are not available or do not lead to satisfactory model precision, analytic functions are implemented, which capture the basic component behavior. Our objective is to derive models with a precision high enough for economic optimization and a computational effort, which allows the implementation on an RCP unit at small sampling times. For the application in optimization-based control algorithms using algorithmic differentiation, the component models have to be twice differentiable.
Physical models often contain fractions or root functions that are continuously differentiable only for a specific domain. For components introduced in this section, pressures p, temperatures T , and rotational speeds ω need to be greater than zero. Except for the EGR, mass flowsṁ need to be greater than zero too. For simulation, various tests have shown that a model consisting of well-identified components and initialized with realistic values respects these limits. For optimization, such constraints could easily be included in the problem formulation. However, as constraints add computational complexity, these are only implemented if they increase the algorithm's robustness.
Many of the component models chosen are based on wellknown approaches. However, for their successful application in optimization-based control, we need to introduce a number of adaptations and respect certain details. If necessary, we saturate model signals or replace specific functions with alternatives that offer certain desired characteristics. In addition, we state boundaries on fitting parameters and model signals. These are considered within the parameter identification process and ensure that the individual air path component characteristics are correctly represented.

A. Common Analytic Functions
Some general functions are applied in various component models. Hence, we introduce these first.
If a model signal needs to be limited to values higher or lower than a predefined boundary, we apply the following smooth saturation functions [6], [20]: A smooth transition from the value y 1 to y 2 at a certain threshold u thresh is achieved by a sigmoid function of the following form: As an example, Fig. 3 shows saturation and sigmoid functions with various smoothing parameters. If physical principles are not available or unsuitable, a (multivariate) polynomial function with degree d might be used

B. Control Volumes and Gas Flow
Control volumes, such as the intake and exhaust manifold, have three state variables, i.e., the pressure p > 0, the temperature T > 0, and the oxygen mass fraction F ≥ 0 (or burnt-gas rate x bg ∈ [0, 1]). We choose a common adiabatic formulation [21] d As the number of state variables has a major influence on the computational complexity of optimization algorithms, individual states are replaced by their steady state (T = T in and F = F in ), a constant, or a measured value whenever possible.
The assumption of a homogeneous gas mixture inside a control volume may not apply equally to all state variables. Compared to changes in pressure, differences in oxygen concentration are compensated slowly. A change in F in measured with a sensor inside of the control volume might therefore better be described by a transport delay rather than by the first-order element behavior of (8). To account for this effect, we consider different volumes (V p , V F , and V T ) for each state, which represent the volumes in which they homogeneously mix. For the oxygen mass fraction, we further introduce a model for the transport delay τ of the gas. As schematically shown in Fig. 4, a gas mass flowṁ AB from V F,A of receiver A to V F,B of receiver B has to travel through the volumes V τ,A and V τ,B .
The time it requires to do so can be determined as follows: The application of such transport delays in optimization algorithms is not straightforward. If not decisive for the model precision, they are neglected. If decisive, a Padé approximation [21] allows to convert the delay into a number of states. However, this method could considerably increase the computational effort. In some cases, it is possible to reformulate the transport delay as a control input delay (as described in Section V), which allows to consider it without any increase of optimization variables.

C. Throttle Valves
For throttle valves, the literature commonly proposes the equations of an isenthalpic throttle with an adapted flow function thr to circumvent numerical problems and to increase the model precision for pressure ratios thr = ( p aThr / p bThr ) close to one [21], [22] We adapt this model approach by introducing own functions for the flow function approximation (14) and the relation between throttle actuation and mass flow (15) The gradient min = 0.05 guarantees a minimum sensitivity of the throttle position on the gas mass flow and therefore allows optimization-based controllers to derive unique solutions. Fig. 5 shows the identified normalized opening areasÃ for all the throttle valves of engines A and B. The complete throttle model with the two introduced functions consists of the five fitting parameters k , k thr , and k Ai .

D. Turbocharger
The turbocharger model can be divided into five submodels: the turbine mass flow and efficiency models, the compressor mass flow and efficiency models, and the dynamic turbocharger speed model. For systems with a variable geometry turbine, we further introduce a dynamic model for the guide vane position.
Both the compressor and turbine models are based on common approaches, which describe the component behavior for predefined and constant inlet gas conditions. Doing so allows compressor or turbine maps to be used for model identification. Equations (16) and (17) connect the standard and the operating conditions [23]

1) Compressor Mass Flow:
For the compressor mass flow, an approach first introduced in [24] is chosen. Although relevant for air path control, compressor surge detection and corresponding control actions are out of the scope of this article. Accordingly, neither a surge limit nor the effects of compressor surge are modeleḋ This approach states that, for a constant turbocharger speed ω tc,c,0 , the dimensionless compressor flow rate c,0 can be described by the function stated in (19), depending on the dimensionless head parameter c,0 , which is a monotonically increasing function of the pressure ratio c . As the shape of the function (19) itself depends on the turbocharger speed, the three parameters q i are modeled as a polynomial function of degree 2.
To be suitable for identification and optimization, we apply an upper limit c,0,max < q 2 to the dimensionless head parameter. This upper limit is derived based on a constant lower bound for the mass flowṁ c,0,min To ensure that the mass flow increases for decreasing pressure ratio or increasing turbocharger speed, the following conditions for the polynomials q i and the partial derivative of the mass flow are introduced: In total, this compressor mass flow model contains nine fitting parameters and the two tuning parametersṁ c,0,min and ε c . As an example, Fig. 6 shows its application to compressor map data of engine A. 2) Compressor Efficiency: As proposed in [24] and [25], the compressor efficiency is modeled as a multivariate polynomial function of order 2, depending on the dimensionless compressor flow rate c,0 and the turbocharger speed. In order to allow a model interpretation, the following notation is commonly used: This approach implies that at a certain dimensionless mass flow k cη,max and a turbocharger speed k cη,ωmax , the maximum compressor efficiency k cη,max is reached. Any deviation leads to a reduction in efficiency. Hence, the symmetric matrix cη ∈ R 2×2 needs to be positive definite. As the other fitting parameters k cη, represent physical quantities, these are constrained as well cη 0, k ηc,ηmax ∈ [0.1, 0.9], k cη, ≥ 0.
As the gradient of the turbocharger speed is a function of the inverse compressor efficiency, we apply a lower bound η c,min . Fig. 7 shows the efficiency model, which contains a total of six identification parameters, fitted to compressor map data of engine A.
3) Turbine Mass Flow: As reported in [25], the commonly used throttle model approach might not represent the turbine mass flow sufficiently well. We, therefore, introduce our own functionṁ The reduced turbine mass flowṁ t,red is described by the monotonically increasing product of a linear and an exponential function, which depends on the inverse turbine pressure ratio −1 t . The linear function has a constant gradient k t5 and the shape of the exponential function is defined by k t3 .  Both functions are zero for −1 t = k t4 . The effect of the turbocharger speed ω tc,t,0 on the turbine mass flow is accounted for via k t2 .
The turbine mass flow model consists of five fitting parameters k t,i . If, in addition, the position of the guide vanes s vtg can be adjusted, the first four parameters are modeled as a polynomial function of degree 2. In that case, the number of parameters increases to 13.
The constraints specified in (32) ensure that a decreasing inverse pressure ratio or an increasing turbocharger speed leads to a decrease of the reduced mass flow. Furthermore, we limit the partial derivative of the reduced mass flow with respect to s vtg to ensure that the mass flow strictly decreases with increasing VTG position. As an example, Fig. 8 shows the application of this model to turbine map data of engine A ∂ṁ t,red ∂s vtg < 0.

4) Turbine Efficiency:
Various publications report that the efficiency of a turbine can be modeled as a polynomial function depending on the blade speed ratio (BSR), which is a function of the pressure ratio and the turbocharger speed [21], [25]. We apply this approach and slightly adapt it with the introduction of a lower limit η t,min for the efficiency In order to increase the precision, we further extend this model by replacing its parameters by polynomial functions that depend on the turbocharger speed and, if existent, the position of the turbine guide vanes s vtg a tη = f poly ω tc,t,0 , 1 These additional functions increase the number of identification parameters from three to seven or, in the case of a VTG, to 13. To represent the basic characteristics of turbine efficiency, the following constraints are introduced:

5) Turbocharger Speed:
The turbocharger speed ω tc is modeled as a dynamic state using Newton's second law. The turbocharger inertia tc is used as a fitting parameter dω tc dt 6) VTG Guide Vane Dynamics: Depending on the turbine design, the connection between the guide vane actuator and the guide vanes themselves can be subject to considerable backlash. While investigating the turbine of engine A, we found that the backlash within the guide vane positioning mechanism can take values of up to ±10% abs of the measured electronic actuator position. We further found that, independent of the direction from which we approach it, the steady-state operation point is unique for every VTG actuator position. This observation indicates that, within their clearance, the guide vanes converge to a unique position, presumably due to vibrations and the force exerted by the airflow.
Given these observations, the steady state of the air path is not affected by backlash, but the dynamic behavior changes. Therefore, we introduce our own backlash model, which is based on ideas published in [26]. It consists of a mass-springdamper system with variable spring stiffness, described by two states and four fitting parameters Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.
Outside of the clearance (| s vtg | s vtg,lim ), the spring has a high stiffness k sp,max and the system is slightly underdamped with a damping ratio ζ ≤ 1, which leads to a fast convergence toward the clearance. Within the clearance (| s vtg | s vtg,lim ), the spring stiffness is significantly reduced to the value k sp,min , while the damping coefficient k dmp is kept constant. This leads to a strongly overdamped behavior and consequently to a slow convergence toward the actuator position u vtg at the center of the clearance. Two sigmoid functions provide a smooth transition between these two spring constants at the clearance limits (52)

E. Engine
As a component of the air path, the engine acts as a volumetric pump, which enforces a mass flow depending on the engine speed ω e , the pressure p be and the temperature T be upstream of the engine, and the pressure p ae downstream of the engine. The combustion of the injected fuelṁ fuel additionally increases the engine-out temperature T ae and decreases the oxygen mass fraction F ae . We adopt common modeling approaches [21] m ae =ṁ be +ṁ fuel (53) The equations for the mass flows and the engine out oxygen concentration contain only two fitting parameters for the volumetric efficiency (56). All other parameters are physical constants or geometric engine characteristics.
For the engine-out temperature, the heat added to the engine mass flow is described as a multivariate second-order polynomial function, which depends on the pressure and temperature of the gas upstream of the engine, and the energy released by the combustion of the fuel and the engine speed Q cyl = f poly E fuel,cyl , ω e , p be , T be , 2 E fuel,cyl =ṁ fuel,Ox 4π n cyl ω e H l,fuel .
Due to the data-based approach chosen, we limit the temperature to the range of the model identification data. The temperature model contains a total of 15 fitting parameters.

F. Heat Exchanger
Typically, heat exchangers are installed to cool the exhaust gas before its recirculation or the fresh air after its compression. Often, the temperature of the gas mass flow exiting the heat exchanger shows slow dynamics or is nearly constant. For these cases, the temperature can be either measured or set to a constant value. Otherwise, we use a second-order multivariate polynomial function that depends on the temperature of the in-flowing gas T he,in and the unidirectional gas mass floẇ m he . The resulting temperature T he,out is limited to temperature values between the cooling fluid T cl and the in-flowing gas IV. PARAMETER IDENTIFICATION The complete parameter identification is carried out in four successive steps.
1) Physical constants and other known quantities are set.
2) All parameters of the saturation functions, referred to as tuning parameters, are set manually according to the range of the signal to which they are applied. 3) A major share of the remaining parameters affect the model's steady state and is identified with steady-state measurement data as described in Section IV-A. 4) Finally, some parameters remain that only affect the model dynamics. These are a volume V for each pressure, temperature or oxygen mass fraction state, an inertia tc for each turbocharger, and four parameters to define the nonlinear spring stiffness k sp for each VTG. They are identified using dynamic measurement data as described in Section IV-B.

A. Steady-State Identification
We propose the following routine for successful steady-state identification of an entire air path model. 1) Every output signal of each individual component is fitted to steady-state measurement data, except for the turbine and compressor efficiency models. If map data is available for the turbine or the compressor mass flow, it can be combined with measurement data for the identification of individual components. Doing so increases the model precision for those turbocharger speed and pressure ratio combinations that are not included in engine steady-state measurement data. 2) Compressor and turbine efficiency are identified together and within the complete air path model using a constrained nonlinear least-squares problem, which simultaneously minimizes the model steady-state error of the turbocharger speed ω tc and the mass flowsṁ c andṁ t . With the separate treatment of each component, the number of identification parameters is kept small, which allows to work with small identification datasets and reduces the calibration effort. Furthermore, this strategy ensures that individual components do not compensate for any model errors of other components. Only the turbocharger efficiencies are identified together and with all air path components connected. Thus, high model accuracy and stability is achieved, and the identification can be performed without turbocharger map data. All steady-state parameter fitting tasks are carried out using the algorithmic differentiation and optimization software package CasADi [27] together with the nonlinear program (NLP) solver IPOPT [28]. Each individual task is formulated as a nonlinear least-squares optimization problem with constraints on fitting parameters and model signals as stated in Section III.
For the two engines used in the case studies, Table II provides an overview of the steady-state model precision for a selection of model states and outputs. The number of measurement data points is denoted by n D .

B. Dynamic Identification
After all steady-state fitting parameters have been identified, some parameters remain that only affect the model dynamics. These are identified by solving a nonlinear least-squares optimization problem with 20 min of dynamic measurement data using solvers provided by MATLAB. The measurement data are generated by applying step and chirp signals to a single actuator at a time while keeping the other actuators and the engine operation point constant.
The air path models for the two engines A and B contain a total of ten and seven dynamic parameters, respectively. Their identification turned out to be a straightforward and low-effort task. Besides the optimization-based approach chosen, manual tuning of these parameters in simulation has been shown to result in acceptable model accuracy too.

V. CONTROL ALGORITHM
The control algorithm presented here is designed to track predefined reference signals r by adapting a number of air path actuators, referred to as control inputs u. It allows the formulation of control tasks that combine reference tracking with the optimization of an economic cost and the consideration of system input and system output delays τ u and τ y . It is based on a common offset-free NMPC formulation consisting of an NLP derived by discretizing an optimal control problem (OCP). This OCP-NLP determines the dynamic (dyn) system actuation using a nonlinear least-squares cost for tracking the reference signals and, if required, a freely Fig. 9. Structure of the proposed control algorithm, used to create one-or two-layer NMPC formulations. definable economic cost term. An observer derives the initial states and model disturbance signals. As shown in Fig. 9, the required measured signals are the tracking signals y meas and known disturbances υ.
Depending on the control task and the available computational capacity of the real-time hardware, the proposed one-layer formulation can be extended to a two-layer NMPC structure [6], [29]. An additional NLP, called target selector NLP (TS-NLP), is used to determine the optimal feasible system steady state (ss) with respect to reference values and possibly an economic objective. The OCP-NLP introduced in the one-layer structure then determines the dynamic system actuation to track this optimal steady state. This two-layer structure introduces additional tuning parameters and can increase the control performance.
The TS-NLP and the OCP-NLP, as well as the observer (obs), are based on an air path model, which describes the evolution of the state x by a system of nonlinear ordinary differential equations (ODEs) that are denoted as follows: In addition, the model provides the analytic functions to calculate arbitrary model output signals y. For the control algorithm, we introduce output functions for the reference signals r ( f y,ss ), the dynamic reference r dyn ( f y,dyn ), and for signals that are constrained ( f y,g ) If there exists a measurement for a model output signal y ,i , it is denoted as y meas,,i ∈ y meas .

A. Observer: Offset-Free NMPC
To obtain an offset-free NMPC formulation, we apply one of the most widely used methods in industry [30]. As this method was commonly applied in dynamic matrix control (DMC), it is referred to as the DMC-like correction scheme [30], [31], [32]. It is based on the assumption that there exist additive output disturbances only. Hence, only the model output functions f y, (·) are corrected with an offset d , while the model dynamics remain unchanged. The states x obs can thus be derived using a numerical simulation of the air path model (open-loop observer). An integrator f int (·), which applies the implicit midpoint method n obs times, is used to propagate the states over the period [t, t + t s ]. The control input u and the known disturbances υ are kept constant during this period Actuators and measured signals of the air path system can be subject to time delays. Therefore, we introduce vectors τ u ∈ R n u ≥0 and τ y, ∈ R n y, ≥0 , which contain the times by which the components of u and y are delayed. The observer accounts for those delays by shifting the control input applied to the simulation and by shifting the output signals y ,obs before they are used to calculate the model errord ,i . Finally, the disturbances d ,i are obtained by filtering the model errors. For each disturbance, a dynamic filter function f filt,,i (·) is introduced, which allows to set a tradeoff between disturbance rejection and noise reduction . . , n y, , and ∈ {g, ss}.
The model error is calculated for the reference signals (d ss ) and the constrained signals (d g ). If there does not exist a measurement for a specific signal, the error is set to zero.

B. Target Selector NLP Formulation
Engine air path control algorithms typically use tracking signals that are subject to strong cross-couplings, such as the intake manifold pressure and the burnt-gas ratio. For such control tasks, finding feasible reference values while considering disturbances can be challenging. If needed, the TS-NLP can be used to compute a feasible steady state where some reference signals can be prioritized over others min w ss f y,ss (w ss , υ) + d ss − r 2 Q ss +J eco,ss s.t. f x (w ss , υ)= 0 f y,g (w ss , υ) + d g ∈ [ y g,min , y g,max ].
The diagonal weighting matrix Q ss is a tuning parameter that determines the prioritization of the reference signals. If the control task is overactuated, the additional degrees of freedom can be used to further optimize the system behavior by adding an economic cost term J eco,ss . Thus, (71) computes an economically optimal and feasible steady state w ss . The OCP-NLP can then track this steady state by choosing at least as many dynamic tracking references r dyn as there are control inputs [29]. If neither the prioritization of tracking signals for nonreachable reference values nor the derivation of an economically optimal steady state is essential for the control task, the TS-NLP can be omitted. In this case, the dynamic reference r dyn is set to the reference r itself, corrected by the model error d ss dyn (w ss , υ), two-layer structure r − d ss , one-layer structure. (72)

C. Optimal-Control-Structured NLP Formulation
The OCP-NLP (73) is based on a standard multiple-shooting discretization [33] of an optimal control problem with a nonlinear least-squares cost formulation. For the sake of readability, the subscripts dyn, which indicate that the states and inputs belong to the OCP-NLP, are omitted In the context of classic optimal control formulations, w andu can be considered as the state and the control input, respectively. The optimal control problem is discretized on a possibly nonuniform time grid with time steps t 0 , . . . , t N −1 . The state propagation between shooting nodes is derived with the same integrator function f int (·) as applied for the observer introduced in Section V-A. With a nonuniform time grid, the number of intermediate integration steps n int,i is defined for each time step individually as well. If the control task contains an economic objective that has to be considered dynamically, the nonlinear least-squares cost formulation can be extended by a term J eco . As the air path actuators are often speedlimited, their signal derivative is considered within the tracking cost as well as in the constraints.
The prediction horizon N, the nonuniform time grid, the grid-dependent number of intermediate integration steps, and the diagonal weighting matrices Q and R are the main tuning possibilities of the OCP-NLP.
1) Consideration of Input Delays: The actuator delays described in Section V-A also need to be considered in the OCP-NLP formulation. For this purpose, we constrain the control inputs within the prediction horizon up to the times τ u to the values that were previously applied to the actual system. Thus, we recreate the actual system behavior with delayed inputs within the OCP-NLP. The solution to the OCP-NLP thus includes the equality-constrained control inputs that have already been applied to the system. Hence, we apply the first unconstrained actuator values of the solution to the NLP as the next control inputs.

D. NLP Solver and Implementation
For computation-time-critical control tasks, sequential quadratic programming (SQP) with a Gauss-Newton Hessian approximation is a well-suited method [34]. Instead of solving the complete NLP at once, a sequence of quadratic programs (QPs) that approximate the NLP at the current solution are solved in each time step. For the OCP-NLP stated in (73), we solve a single QP per sample time, which results in the well-known real-time iteration (RTI) scheme [35]. As the TS-NLP (71) provides the reference for the OCP-NLP, we apply multiple SQP iterations to increase the precision of its solution. In addition, the TS-NLP has to be robust for a large range of initial optimization parameters and reference values. Therefore, we use an SQP approach with a globalization line search strategy based on merit functions as suggested in [36,Sec. 3.5.1], [37].
With an efficient solver method defined, our objective is to efficiently implement the control structure on our real-time hardware: a dSPACE MicroLabBox RCP system. For this purpose, we use the modular open-source software package acados [19], which allows the creation of fast embedded nonlinear solvers based on a number of further open-source software tools. The QPs that approximate the NLP in each SQP step are derived by the application of the algorithmic differentiation software CasADi [27]. In order to solve these QPs, the interior-point QP solver HPIPM, which is based on the linear algebra library BLASFEO, is used [38], [39]. For the QPs derived from the OCP-NLP, we use HPIPM with its efficient partial condensing implementation to fully exploit the optimal control structure within the QP [40]. The QPs derived from the TS-NLP are solved using the dense solver implementation of HPIPM.
To realize a control structure as it is proposed here, the acados interface was adapted and extended. Specifically, the definition of TS-NLPs and the deployment of the resulting NLPs on the RCP system were subject of further development. The new functionalities thus derived are now publicly available as part of the acados software package.
The main solver tuning parameters for both NLPs are the solver tolerances and the Hessian regularization with a Levenberg-Marquardt term. For the TS-NLP, the line search parameters and the number of SQP steps constitute further tuning possibilities.

VI. CASE STUDY A: PUMPING-LOSS-OPTIMAL DUAL-LOOP EGR CONTROL
In this first case study, we design an air path controller, which tracks the pressure and the burnt gas ratio in the intake manifold of engine A using the VTG and the throttles of the high-and low-pressure EGR paths. The resulting control problem is overactuated. The desired burnt-gas ratio can be tracked with either the low-pressure or the high-pressure EGR path. Table III lists two steady-state measurement data points of engine A, one of which uses high-pressure EGR, while the other one uses low-pressure EGR. Although both data points are measured at the same engine speed and injected fuel mass and nearly identical intake manifold conditions, the engine efficiency increases by 1% abs if we switch from high-pressure EGR to low-pressure EGR. A comparison of the indicated mean effective pressure (IMEP) for the gas exchange and the combustion cycle individually indicates that the efficiency reduction originates from the pumping losses. These are considerably higher with high-pressure EGR for the selected intake manifold pressure and burnt-gas ratio.
Given these observations, we extend our control task to include an economic objective. Specifically, we want to track  III   RESULTS OF TWO MEASUREMENT DATA POINTS GATHERED AT EQUAL  ENGINE SPEED AND FUEL MASS FLOW AND NEARLY IDENTICAL  INTAKE MANIFOLD GAS CONDITIONS. ONCE THE HIGH-PRESSURE  EGR PATH AND ONCE THE LOW-PRESSURE EGR PATH IS USED  TO REACH THE DESIRED INTAKE MANIFOLD BURNT-GAS RATIO Fig. 10. Dynamic measurement data of a low-pressure EGR throttle step (blue) and a high-pressure EGR throttle step (green) at t = 1 s. intake manifold gas conditions with the lowest possible pumping losses p em − p im . The dynamic behavior of the two EGR paths is plotted in Fig. 10. It shows the pressure and the burnt-gas ratio in the intake manifold, once for a step with the high-pressure EGR throttle and once for a step with the low-pressure EGR throttle. The high-pressure path has a considerably shorter response time but reduces the mass flow through the turbine, resulting in a reduction of the intake manifold pressure. The low-pressure EGR path has no cross-couplings with the intake-manifold pressure but has a considerable time delay, which must be considered in a control strategy that aims for high tracking performance. Fig. 1 shows the air path model with all of its states, control inputs, and measured disturbances. It contains the following simplifications and assumptions.

A. Air Path Model
1) The low-pressure EGR path also contains a throttle in the air path u lp,air to increase the maximum EGR mass flow. We eliminate this actuation as a model and control input by setting it to u lp,air = 1 − 0.7 u lp,egr . This is a common approach, which ensures that the fresh-air throttle only closes with the opening of the low-pressure EGR throttle and is never closed more than 70%. Therefore, it prevents a pressure drop in front of the compressor, which would increase pumping losses. 2) Given the introduced strategy for the actuation of the throttle u lp,air , the effect of the low-pressure EGR path on the pressure before the compressor is negligible. Also, the cooled low-pressure EGR does not significantly change the temperature before the compressor. Hence, except for the oxygen mass fraction, the gas conditions can be modeled as a disturbance. 3) We simplify the model of the low-pressure EGR path by assuming that the pressure and the temperature of the cooled low-pressure EGR are equal to the conditions of fresh air and that the flow functions of the two corresponding throttles have the same exponent k . Doing so allows us to model the oxygen mass fraction before the compressor (9) as a function of F at and u lp,egr only. 4) As the dynamics of the heat exchanger after the compressor (intercooler) are slow, the temperature of the out-flowing gas is not modeled but measured and treated as a disturbance. 5) As the low-pressure EGR actuation only affects the oxygen mass fraction before the compressor, we can model the transport delay of the gas from before the compressor to the intake manifold as an actuator delay and include it in the dynamic NLP, as described in Section V-C1.

B. Definition of Two Control Algorithms
Based on the air path model of engine A, two control algorithms are designed, tuned, and experimentally validated. They are derived in sequence and serve different purposes.
First, we design a controller that achieves the best possible control performance given the full computational power of the RCP system. It is used to analyze the control behavior in simulation and as an experimentally applicable benchmark for other control algorithms. It offers extended tuning possibilities during experiments, as it is robust for a wide range of cost weights and other tuning parameters. We denote this controller as C bench .
While real time capable on our RCP system, the C bench controller would not be real-time feasible on today's ECUs. These offer significantly less computing power and must perform other real-time tasks in addition to air path control. However, the derived benchmark controller can be used as a starting point for designing a computationally less complex controller. The execution time can be significantly reduced by eliminating online tuning capabilities while making minimal sacrifices in control performance. Therefore, starting with C bench , we design a second controller C rt that is tuned for both control performance and computational effort. This second controller is used to evaluate the real-time suitability of our approach for production ECUs.
Table IV defines the two controllers C bench and C rt by stating all the signals needed for the control structure presented in Section V.
1) Benchmark Controller: We first focus on the C bench controller, which consists of a TS-NLP and an OCP-NLP. The purpose of the TS-NLP is to calculate the pumping-loss optimal steady state for the set values r. Hence, besides the tracking cost for the intake manifold pressure and burntgas ratio, it further contains an economic cost term J eco,ss . This term is used to penalize the exhaust pressure, causing the TS-NLP to find a steady state that achieves the reference values with the lowest possible exhaust pressure. Since reducing the exhaust pressure is equivalent to reducing the pumping losses for a given intake pressure, this formulation achieves the desired control behavior. The constant offset value p em,0 < p em , as well as the scalar weighting factor Z ss , is the tuning parameters of the TS-NLP.
By applying constraints ( f y,g ), we limit all control inputs to their physical range u ∈ [0, 1]. Furthermore, we also restrict p em and T em to physically reachable values, which leads to increased robustness of the TS-NLP for an extended range of initial conditions. The solver uses eight SQP steps per sample time to accurately determine the optimal steady state.
Given the pumping-loss-optimal steady state, the OCP-NLP can be formulated as a pure tracking task ( J eco = 0). Its objective is to control the pressure and the burnt-gas ratio in the intake manifold and the burnt-gas ratio before the compressor to the values obtained by the TS-NLP. With this definition of the reference signals r dyn , the amount of lowand high-pressure EGR is uniquely defined.
A prediction horizon of 3.3 s was found to provide good control performance for a wide range of cost weights and further tuning parameters. For the OCP discretization, we use a nonuniform time grid, consisting of 30 shooting intervals with time steps t 0 , . . . , t 11 = 0.05 s and t 12 , . . . , t 29 = 0.15 s. As we consider the transport delay of the oxygen mass fraction upstream of the compressor using an input delay, the smaller time steps for the first 0.6 s lead to increased precision of the discretized transport delay.
Regarding numerical integration, two intermediate steps for t 0 , . . . , t 11 and three steps for t 12 , . . . , t 29 provided a good tradeoff between robustness, precision, and computational effort. As for the TS-NLP, constraints are used for all control inputs and the two states p em and T em .
The resulting NLP is solved using an RTI scheme with full SQP steps, and together with the TS-NLP, it is solved every t s = 0.05 s. This sampling time was determined in simulations and is in good agreement with various studies on diesel air path control [1], [2], [7], [10], [41]. Although recent publications propose up to five times faster sample times [14], [42], a decrease of t s did not considerably increase the control performance for the engine and control task investigated in this case study.
2) Reduction of the Computational Effort: To reduce the computational complexity, the C rt controller omits the TS-NLP. The OCP-NLP directly tracks the intake manifold pressure and burnt-gas ratio without knowledge of the pumping-loss optimal steady state. The economic objective is directly integrated into the OCP-NLP, with the addition of the cost term J eco , which slightly penalizes the exhaust pressure over the prediction horizon.
Besides omitting the TS-NLP, a further reduction of the computational complexity is achieved by considerably reducing the prediction horizon from 3.3 to 1.4 s. It is discretized with a nonuniform time grid with the time steps t 0 = 0.05 s, t 1 = 0.1 s, t 2,...,4 = 0.15 s, and t 5,6 = 0.4 s. With this discretization, the gas transport delay of the low-pressure EGR path can still be represented, albeit with lower accuracy.
The resulting controller consists of a single NLP with 109 optimization variables, 29 linear and 56 nonlinear equality constraints, and 70 linear inequality constraints.

C. Simulative Steady-State Control Analysis
The TS-NLP derived for the C bench controller is used in simulation to compare the steady-state pumping losses of the dual-loop EGR strategy with the pumping losses of a strategy using single-loop EGR only. The results are shown in Fig. 11.
In general, pumping losses increase with higher intake manifold pressure. Using high-pressure EGR only, an increase of the burnt-gas ratio can lead to a decrease or an increase in pumping losses, depending on the intake manifold pressure. On the other hand, the application of low-pressure EGR only increases the intake manifold burnt-gas ratio without altering the pumping losses. The pumping-loss optimal combination of these two EGR paths drastically increases the range of burnt-gas ratio references at which the lowest possible pumping losses are reached.

D. Reduction of Pumping Losses
To dynamically evaluate the C bench controller in terms of pumping losses, we compare it to a high-pressure and a lowpressure EGR strategy, both of which are obtained by adapting the actuator constraints of C bench . Given a reference trajectory that contains a variation of the burnt-gas ratio at constant intake manifold pressure, Fig. 12 shows the effectiveness of our benchmark controller in simulation, and Fig. 13 shows the corresponding experimental results obtained on engine A. The behavior of the algorithm in simulation translates well to the actual system. Also, the mean performance values of the algorithm stated in Table V show the same tendencies in simulation as in the experimental results. Compared to both  single-loop EGR strategies, the C bench controller successfully manages to reduce the pumping losses. Compared to the simulation, the tracking performance on the actual system is slightly reduced due to model mismatches, disturbances, and measurement noise.
Experimental results are also obtained for the C rt controller and illustrated together with the results of C bench in Fig. 14 and Table V.   Compared to C bench , the computationally less demanding algorithm has a slightly increased tracking error at similar pumping losses.

E. Tracking Performance
The high model precision required for the economic optimization of the pumping losses also enables excellent  To follow the negative reference step of the burnt-gas ratio, the dual-loop strategy uses the high-pressure EGR path to temporarily follow the reference until the delayed reaction of the low-pressure EGR path starts to affect the intake manifold conditions. Hence, the consideration of the low-pressure EGR transport delay clearly improves the tracking performance. Another example is given in Fig. 15, which shows the tracking of a step in the intake manifold pressure at a constant burnt-gas ratio with the C bench controller, once in simulation and once measured. The control algorithm exploits the maximum actuator speed and range of the VTG and the maximum speed of the high-pressure EGR throttle to reduce the tracking error as fast as possible. Table VI lists the computation times required by the RCP unit to calculate the results shown in Fig. 14. With only minor sacrifices in tracking accuracy, the C rt controller manages to reduce the computation time by a factor of 5 compared to the benchmark controller.

F. Computational Effort
Production ECUs provide considerably less computational power and memory. If equipped with a 256-MHz processor [42], their clock speed is eight times lower compared to our RCP. However, with a sample time t s = 50 ms, the C rt controller utilizes only a fraction of the computing power provided by our RCP system. If we were to assume that clock speed is proportional to computation time in first approximation [42], C rt would run on a production ECU but would use up to 90% (mean: 70%) of the processor capacity. However, this linear assumption is questionable, and with a size of 1.1 MB, the controller would most likely exceed the memory of production ECUs. We conclude that further development toward lower computational effort and memory usage is required for the algorithm to run on production hardware.

VII. CASE STUDY B: TRANSFER AND TUNING
OF THE ALGORITHMS Given an air path component library able to represent various air path systems and the NMPC formulation stated in Section V, the derivation of new controllers for different engines is straightforward. We demonstrate this by transferring the controllers derived in Section VI to engine B. Compared to engine A, it has no low-pressure EGR path and no variable geometry turbocharger, but an exhaust throttle after the turbine. Without any adaption of the model functions, we reidentify the air path component models with measurement data of engine B and connect them to represent the new air path configuration, as shown in Fig. 2.

A. Definition of the Algorithms
Analogous to engine A, the objective is to track the intake manifold pressure and the burnt-gas ratio. Given the new air path structure, our control inputs are the high-pressure EGR throttle and the exhaust throttle. Since we track as many signals as we have control inputs, this control task is a pure tracking problem and we can neglect the economic cost term. Table VII summarizes the characteristics of the two controllers for engine B.
Since a complete closure of the exhaust throttle is not desired, we use the constraints u eb ∈ [0.2, 1] in all NLPs. As for engine A, we also constrain u hp,egr ∈ [0, 1] and restrict the two states p em and T em to physically reachable values.
The prediction horizon of the benchmark controller OCP-NLP is set to 3.75 s and is discretized by a uniform time grid consisting of t 0 , . . . , t 24 = 0.15 s. Numeric integration has proven to be robust and precise with three intermediate steps (n int,i = 3). For C rt , it was found that a short horizon of 0.85 s, discretized with a nonuniform time grid with three shooting intervals, is sufficient.
The sample time t s and the number of SQP iterations of the TS-NLP are taken from the control structure of engine A without any adaption. The weighting matrices Q, R, and Q ss are chosen using simulations.

B. Experimental Results: Tracking Performance
To put the control performance of the algorithm derived into perspective, we compare it to a linear reference controller, which is designed and tuned for one specific engine operation point. It is based on a linear quadratic Gaussian control approach with integral action (LQGI) and uses the intake manifold pressure and burnt-gas ratio as control feedback. On the RCP unit, it is implemented with a sample time of 1 ms. Fig. 16 shows the tracking results of the two nonlinear controllers C bench and C rt and the linear controller. Although the linear controller is reevaluated 50 times faster and is specifically tuned for the selected engine operation point, it is outperformed by both nonlinear controllers. Crosscoupling effects are less pronounced, and the tracking of reference steps is generally faster with a lower overshoot. These results show that, even at a constant engine operation point, it is advantageous to consider system nonlinearities within the control approach.  The two nonlinear controllers C bench and C rt have an almost identical control behavior. Overshoot and cross-coupling effects are slightly increased with the C rt controller but still reduced compared to the linear control approach.

C. Experimental Results: Controller Tuning
For the C bench controller, the weighting matrix Q ss of the TS-NLP can be set to prioritize tracking signals for nonreachable reference values. For changes in engine speed and load, Fig. 17 contrasts the reference tracking results of C bench for two different choices of the matrix entries of Q ss , once these values are set to prioritize the burnt-gas ratio and once to prioritize the intake manifold pressure.
For reachable set values, the two controllers react similarly, and they focus on their prioritized signal for infeasible reference values, which leads to a completely different control behavior.

D. Computational Effort
Table VIII lists the computation times of the nonlinear controllers for the results shown in Fig. 16. With only minor sacrifices in tracking accuracy, the C rt controller manages to reduce the computation time by a factor of 10 compared to the benchmark controller.
As for engine A, we assume that an ECU with 256 MHz would need about eight times longer computation times. With a sample time t s = 50 ms, the algorithm should therefore be feasible to run on a production ECU and would take up to 30% (mean: 20%) of the processor's capacity.

VIII. CONCLUSION AND OUTLOOK
We present a systematic and modular approach for the control of diesel engine air paths. It consists of a library of grey-box models for common air path components and a two-layer NMPC structure. We show a tracking performance superior to that obtained with linear control strategies, and we demonstrate the applicability of the algorithm to various air path configurations and control tasks with a low calibration effort. We further demonstrate a model precision high enough to derive control strategies with an economic objective and experimentally validate the resulting algorithms with RCP hardware.
The experimental evaluation performed within the two case studies motivates the approach presented and acts as a proof of concept. Further research should focus on the extended experimental evaluation of the framework and the reduction of the computational complexity. Specifically, the following extensions are worth further investigation: 1) the evaluation of the controller performance and robustness for complete regulatory drive cycles under varying ambient conditions; 2) a comparison of the control performance with the control strategies of a fully calibrated series production engine; 3) the derivation and evaluation of linear time-varying MPC based on the model library presented as an extension of the framework; 4) the evaluation of the potential of more sophisticated methods for offset-free NMPC; 5) further development toward lower computational complexity combined with the algorithm implementation on a production ECU. These additional developments would result in a software framework enabling the systematic generation of controllers for complex diesel engine air paths, both for experimental benchmarking and for use in production-type control units.