Solution for state constrained optimal control problems applied to power split control for hybrid vehicles

Automatica 50 (2014) 187–192

Contents lists available at ScienceDirect

Automatica journal homepage: www.elsevier.com/locate/automatica

Brief paper

Solution for state constrained optimal control problems applied to power split control for hybrid vehicles✩ Thijs van Keulen 1 , Jan Gillot, Bram de Jager, Maarten Steinbuch Department of Mechanical Engineering, Eindhoven University of Technology, P.O. Box 513, 5600 MB Eindhoven, The Netherlands

article

info

Article history: Received 23 July 2011 Received in revised form 24 July 2013 Accepted 20 September 2013 Available online 28 October 2013 Keywords: Optimal control Maximum principle Trajectories Hybrid vehicles

abstract This paper presents a numerical solution for scalar state constrained optimal control problems. The algorithm rewrites the constrained optimal control problem as a sequence of unconstrained optimal control problems which can be solved recursively as a two point boundary value problem. The solution is obtained without quantization of the state and control space. The approach is applied to the power split control for hybrid vehicles for a predefined power and velocity trajectory and is compared with a Dynamic Programming solution. The computational time is at least one order of magnitude less than that for the Dynamic Programming algorithm for a superior accuracy. © 2013 Elsevier Ltd. All rights reserved.

1. Introduction This paper concerns numerical solutions for convex scalar state constrained optimal control problems. In the literature, several approaches are known to solve problems of this type. Firstly, the problem can be solved ‘‘directly’’ with, e.g., multiple shooting (Sager, 2005), or other approximate methods that require a quantization, see Gerdts (2008), Gerdts and Kunkel (2008), Wang, Gui, Teo, Loxton, and Yang (2009) and Yu, Li, Loxton, and Teo (2012). Dynamic Programming (DP) (Bertsekas, 2000) is another approach often applied to optimal control problems of low dimension. Drawback of the quantization of control and state variables is that it only leads to an approximation of the original problem. This results in the classical trade-off between accuracy and computational demand induced by the quantization grid size. Secondly, indirect methods derive the optimal solution using a two step procedure; as a first step the Pontryagin Minimum Principle (PMP) (Pontryagin, Boltyanskii, Gamkrelidze, & Mischenko, 1962) is applied to derive the necessary conditions for optimality in the form of a differential equation on the costate variable and the

✩ The material in this paper was not presented at any conference. This paper was recommended for publication in revised form by Associate Editor Kok Lay Teo under the direction of Editor Ian R. Petersen. E-mail addresses: [email protected] (T. van Keulen), [email protected] (J. Gillot), [email protected] (B. de Jager), [email protected] (M. Steinbuch). 1 Tel.: +31 402478390; fax: +31 402461418.

0005-1098/$ – see front matter © 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.automatica.2013.09.039

static optimization of the Hamiltonian function. These analytical results allow one to write the unconstrained optimal control problem into a boundary value problem. As a second step, the remaining boundary value problem can then be solved, without quantization of the state and control space, but generally using the discretization of the time space, using information about the state and costate at the boundaries. Extensions to the state constrained optimal control case can be found (Fabien, 1996; Jacobson & Lele, 1967). These solutions involve, e.g., a penalty function which comes with the cost of an increased state dimension, however. This contribution presents a novel numerical approach for convex scalar optimal control problems with ‘‘pure’’ state constraints which has superior results in terms of computational demand and accuracy compared to other known numerical techniques for problems of this type. It takes advantage of the PMP, it does not require quantization of the state and control space, it includes state constraints, state dependent losses and non-smooth cost function descriptions, all without the introduction of a penalty function. A proof for optimality of the solution is included. The novel approach is applied to the power split control problem in hybrid vehicles and is bench marked, for computation time and accuracy, with a DP algorithm. The influence of battery voltage increase as a function of the battery state-of-energy is also evaluated. The paper is organized as follows. In Section 2, the novel numerical algorithm based on the PMP is given. Section 3 introduces the power split control problem for hybrid vehicles. A comparison of the novel algorithm with the DP algorithm is presented in Section 4. Conclusions can be found in Section 5.

188

T. van Keulen et al. / Automatica 50 (2014) 187–192

2. Problem description

• the Hamiltonian H has a global minimum with respect to control u:

Consider an optimal control problem with a scalar state and control:

u∗ (t ) = arg min H (t , x∗ (t ), u(t ), p∗ (t )),

 t1    F (t , x(t ), u(t ))dt min  u∈U   t0  subject to: x˙ = f (t , x(t ), u(t )), (P0 )  1, x(t0 ) = x0 ,  x(t1 ) = x     x ( t ) 0 x ( t ) −  h(t , x(t )) = ≤ , x(t ) − x(t ) 0

where x∗ (t ) is the optimal state trajectory, u∗ (t ) the optimal control trajectory, p∗ (t ) the corresponding adjoint multiplier function. 2.2. Numerical solution for the unconstrained problem

here, u is the control variable, U is a closed convex set of admissible controls for every t and x, t0 the initial time, t1 the end time, F the time dependent cost function, x the state variable, f the state dynamics equation, x0 and x1 the boundary conditions, h the inequality constraint on the state, x the upper state constraint, and x the lower state constraint. It is assumed that the functions F , f and h are continuous in all their arguments, and continuously differentiable in x. 2.1. Necessary conditions of optimality In this section, the necessary conditions for optimality for problem P0 are given (Hartl, Sethi, & Vickson, 1995; Maurer, 1977; Seierstad & Sydsæter, 1987). The Hamiltonian is defined with: H (t , x(t ), u(t ), p(t )) = F (t , x(t ), u(t )) − p(t )f (t , x(t ), u(t )).

(1)

The state inequality constraints can be adjoined to the Hamiltonian to form the following Lagrangian: L(t , x(t ), u(t ), p(t ), λ(t ))

= H (t , x(t ), u(t ), p(t )) + λ⊤ (t )h(t , x(t )).

(2)

Applying the PMP as in Theorem 9.3.1 of Vinter (2000, p. 339), it follows that if the control is optimal, then there exists a nontrivial piecewise continuous multiplier function p(t ) ̸≡ 0 such that the following conditions are satisfied:

• the differential equation on the adjoint multiplier function: p˙ (t ) =

∂L ∂x

• the complementary slackness condition:   λ1 (t ) = 0 for t ∈ v : x∗ (v) < x ,   λ2 (t ) = 0 for t ∈ v : x∗ (v) > x ,

(3)

(4) (5)

• the condition on the adjoint multiplier, see also Hartl et al. (1995, Theorem 4, p. 186), for ta < tb in [t0 , t1 ]:  tb  ∂h p(tb+ ) − p(ta+ ) = p˙ (t )dt + dξ1 (t ) ta (ta ,tb ] ∂ x  ∂h − dξ2 (t ), (6) (ta ,tb ] ∂ x where ξ1 and ξ2 are of bounded variation, non-increasing, constant on intervals where x < x < x, right continuous and have left-sided limits everywhere. The multiplier trajectory p has a discontinuity given by the following jump condition: p(τ + ) = p(τ − ) + µ1 (τ ) − µ2 (τ ),

(8)

u∈U

(7)

with µ1 ≥ 0 and µ2 ≥ 0. Under the assumption that ξ1 and ξ2 have a piecewise continuous derivative, it is possible to set λ1 (t ) = ξ˙1 (t ), λ2 (t ) = ξ˙2 (t ), for every t for which ξ1 and ξ2 exist and µ1 (τ ) = ξ1 (τ − ) − ξ1 (τ + ), µ2 (τ ) = ξ2 (τ − ) − ξ2 (τ + ), for all τ ∈ [t0 , t1 ] where ξ1 and ξ2 are not differentiable,

In this section, a numerical solution for the unconstrained problem is discussed. Here ‘‘unconstrained’’ refers to a problem of type P0 without the inequality constraint. To derive a well defined two point boundary value problem from conditions (3) and (8), the following property is required. Lemma 1. Let H = F − pf with p ≥ 0, x a scalar variable, F a convex function in u, and f a strictly concave and strictly monotonic decreasing function in u, additionally, assume the differential equation p˙ = ∂∂Hx , with ∂∂Hx locally Lipschitz in p on a domain defined by U, then the solution u∗ of (8) is a monotonic decreasing function of p and there is a monotonic increasing relation between the initial value of the multiplier p(t0 ) and the final state x(t1 ) and between the inverse relation of the final state x(t1 ) and the initial value of the multiplier p(t0 ). Proof. The proof uses the convexity properties and the existence and uniqueness of the solution of the differential equation (3). If f is strictly concave in u, and p > 0, then the function −pf is strictly convex in u. The sum of two convex functions is also convex. So, if F is convex and −pf strictly convex, then H is strictly convex in u. If H is strictly convex in u, then it has a unique minimum defined ∂f by ∂∂ Fu − p ∂ u = 0. If f is strictly concave and strictly monotonic decreasing, then ∂f < 0 and also strictly monotonic decreasing. In the optimum it ∂u ∂f

∂f

holds that p = ∂∂ Fu / ∂ u with ∂ u < 0 and p ≥ 0, then ∂∂ Fu ≤ 0, and, because F is convex, ∂∂ Fu is also monotonic increasing. Again using ∂f

∂f

p = ∂∂ Fu / ∂ u with ∂∂ Fu ≤ 0 and monotonic increasing and ∂ u < 0 and strictly monotonic decreasing, it follows that the minimum u∗ of H is a monotonic decreasing function of p. If ∂∂Hx is Lipschitz continuous in p, it follows that p is a unique solution of (3), see Khalil (2002, Theorem 3.1), hence pa (t0 ) > pb (t0 ) implies p0a (t ) > p0b (t ), where p0 (t ) denotes the trajectory resulting from p(t0 ). Using that the minimum u∗ of H is a monotonic decreasing function of p, it follows that u0a (t ) ≤ u0b (t ) if p0a (t ) > p0b (t ), where u0 denotes the control trajectory resulting from p(t0 ). Finally, given the state dynamics x˙ = f with f monotonic decreasing in u, a monotonic increasing relation is found between p(t0 ) and x(t1 ) and likewise between x(t1 ) and p(t0 ).  Given the necessary conditions of optimality, the following boundary value problem is obtained:

 ∗ ˙  x(t ) = f (t , x(t ), u (t∗)), ∂ H (t , x(t ), u (t ), p(t )) (PBVP ) p˙ (t ) = ,  ∂x  x(t0 ) = x0 , x(t1 ) = x1 , in which u∗ (t ) is the solution of (8). Using single shooting, an initial value problem associated with this boundary value problem can be derived. Generally, we have to sample the time space and apply numerical integration methods to solve this initial value problem, e.g., the Euler scheme, where the discrete time sample moments are indicated by variable k = [1, . . . , n] with length n ∈ N defined by t1 − t0 with equidistant step size ∆t:

∆t =

t1 − t0 n−1

.

(9)

T. van Keulen et al. / Automatica 50 (2014) 187–192

189

The discretized initial value problem becomes:

  x(k + 1) = x(k) + f (k, x(k), u(k)), ∂ H (k, x(k), u(k), p(k)) (PIVP ) p(k + 1) = p(k) + ,  ∂x  x(1) = x0 , p(1) = p0 . The solution of PIVP solves also problem PBVP if the initial condition p0 is chosen such that the boundary condition x(n) = x1 is met. The boundary value problem can thus be solved if a root for the following nonlinear function is found:

ψ(p0 ) ≡ x(n, p0 ) − x1 = 0.

(10)

Here x(n, p0 ) is obtained by solving PIVP with p0 . A root finding algorithm, e.g., bisection, can be applied to find the solution p0 . 2.3. Numerical solution for the constrained problem Next, a method is proposed to find each boundary interval or contact time with the state constraint h(t , x(t )) based on the times where the constraint is exceeded the most in the unconstrained optimal trajectory such that the problem can be split in two subproblems that can again be solved with the approach outlined in Section 2.2. This procedure is repeated with a recursive scheme until all state constraints are met. Algorithm 1. The optimal multiplier p∗ (t ) and state x∗ (t ) trajectories for the state constrained optimal control problem with scalar state are found by the following sequence:

• compute the unconstrained optimal solution defined by the initial value of the costate pτ 0 , i.e., solve the two point boundary value problem as in Section 2.2. If a state constraint is violated: repeat • find the sample time moment τi where the state boundary is exceeded the most,

τi = arg max (x(k) − x(k, pτ i ), x(k, pτ i ) − x(k)),

(11)

k∈[kai kbi ]

where kai and kbi are the initial and final sample times of the sub-trajectory in which the state-boundary is exceeded the most, pτ i the initial costate of the sub-trajectory at kai , and i is the ith iteration of this recursive scheme, • split the initial (sub-)trajectory in two sub-trajectories: [kai τi ] and [τi kbi ]. – in the case the upper state constraint is exceeded, use x(τi ) := x(τi ) as the end point constraint for the interval [kai τi ] and as the initial condition for the interval [τi kbi ], then solve the two unconstrained boundary value problems, i.e., solve (10) + resulting in p− τ i and pτ i at the sub-trajectories, – in the case the lower state constraint is exceeded, use x(τi ) := x(τi ) as the end point constraint for the interval [kai τi ] and as the initial condition for the interval [τi kbi ], then solve the two unconstrained boundary value problems, i.e., solve (10) + resulting in p− τ i and pτ i at the sub-trajectories, • until maxk∈[kai kbi ] (x(k) − x(k, pτ i ), x(k, pτ i ) − x(k)) ≤ 0 for all i sub-trajectories. Proof. The proof is based on the following two observations:

• from the jump conditions of the PMP, the following properties of the optimal solution are obtained: p(τi+ ) > p(τi− ) if the upper constraint is reached, p(τi+ ) < p(τi− ) if the lower constraint is reached,

• using the invertible monotonic increasing relation between p(1) and x(n) described in Lemma 1 it follows that the time at which the boundary is exceeded the most is a contact point or part of the boundary interval, and, therefore, part of the optimal constrained solution.

Fig. 1. Sketch of solution. The upper sub-figure displays the state of energy solution, and the lower sub-figure the optimal multiplier values, in which the solid line is the unconstrained solution, the dashed line is the solution of the subtrajectories after the first iteration, and the dash-dotted line the solution after the third iteration.

All solutions on (sub-)trajectories can be reduced to three situations. The first situation is obtained when the unconstrained solution on a (sub-)trajectory exceeds a constraint, here it is immaterial whether the lower of upper constraint is exceeded. The second and third situation are needed to check if the result of the first situation is not invalidated by successive iterations. The second and third situations are obtained when one of the sub-trajectories resulting from the first situation, with a boundary condition on the constraint that was exceeded the most in the first situation, exceeds either the upper or lower constraint. Other situations can be mapped on these three situations. A proof for optimality for these three situations is given below. Without loss of generality, it is assumed that the maximum violation of the state boundaries of the unconstrained solution is at the upper bound, see solid lines in Fig. 1. Let us now denote by x∗ (k, p∗τ0 ) the unconstrained optimal state trajectory and p∗τ0 the multiplier trajectory of the unconstrained optimal solution. The time instance of the maximum violation is given by:

τ1 = max (x(k) − x(k, p∗τ0 ), x(k, p∗τ0 ) − x(k)).

(12)

k∈[1 n]

Let us compute the two sub-trajectories x(k, p− τ1 ) on the interval − [1 τ1 ] and x(k, p+ τ1 ) on the interval [τ1 n] with x(τ1 , pτ1 ) − x(τ1 ) = + x(τ1 , pτ1 ) − x(τ1 ) = 0. At the first sub-trajectory, the final state has to decrease, so, by Lemma 1 it follows that p has to decrease as well. At the second sub-trajectory, the initial state has a lower value, hence, by Lemma 1, it follows that p has to increase. Note that, pτ1 is discontinuous at τ1 and the jump direction meets the necessary conditions for optimality. Next, one of the two sub-trajectories is considered. The subtrajectory x(k, p− τ1 ) does not necessarily have a global maximum at contact time τ1 and also the lower bound can be exceeded. These two cases are now discussed one by one. First, it can also happen that the maximum violation of the boundary of the sub-trajectory is at the lower bound, for instance at τ2 < τ1 . Then, another iteration is necessary, splitting the original trajectory into three parts. By Lemma 1 it can be checked that − + + − p− τ2 > pτ1 > pτ2 and, again the jump direction pτ2 < pτ2 meets the necessary conditions for optimality. Second, it can happen that in a sub-trajectory the upper bound is still exceeded, so x(τ3 ) < x(τ3 , p− τ1 ), see the dashed line in Fig. 1. Here τ3 is the time where the unconstrained solution exceeds the upper boundary the most in the sub-trajectory. If we set x(τ3 ) = x(τ3 ) as the initial condition and x(τ1 ) = x(τ1 ) as the final condition for the sub-trajectory of another iteration, the state trajectory on the interval [τ3 τ1 ] has to decrease since x(τ3 ) < x(τ3 , p∗τ0 ) ≤ x(τ1 , p∗τ0 ). Lemma 1 is used once more to prove that the multiplier value on the interval [τ3 τ1 ] has to decrease compared to p∗τ0 . The − + ∗ + following inequality holds p− τ3 < pτ1 < pτ3 ≤ pτ0 < pτ1 .

190

T. van Keulen et al. / Automatica 50 (2014) 187–192

This proves that not only the necessary conditions for optimality are met, but that τ1 is also a contact point of the constrained optimal solution, i.e., by using the sample time moment at which the boundary is exceeded the most, results in a contact point of part of the boundary interval. The above reasoning holds also on all the sub-trajectories and a recursive strategy will converge to a solution that fulfills all necessary conditions of optimality. Given the discrete time space, the algorithm will finish in a finite number of steps.  3. The power split problem for hybrid vehicles Hybrid Electric Vehicles (HEVs) employ an electric machine, in combination with a battery, to recover and to store energy and also to optimize the operating condition of the internal combustion engine. Fuel consumption can be optimized by controlling the power split between the engine and the electric machine, using a power request from the driver as the input. In Tate and Boyd (2000), the power split problem is approximated and rewritten as a linear programming problem. Dynamic Programming is also applied to optimize the power split for prescribed trajectories; the interested reader may refer e.g., de Jager and van Keulen (2013) for details on DP in this context. In Delprat, Guerra, and Rimaux (2002), the necessary conditions of optimality are used to set up a boundary value problem, approximating the multiplier function as a constant scalar optimization variable, calculating the unconstrained optimal stateof-energy trajectory with a root finding algorithm. The algorithm of Delprat et al. (2002) is extended here to include state constraints, a non-smooth cost function description, as well as the state dependent battery losses. Below a model description for a hybrid drive train is introduced which fulfills the requirements of Lemma 1. The mechanical output power of the internal combustion engine Pp and electric machine Pm matches the power request Pr : Pr = Pp + Pm .

(13)

The fuel cost of the engine, at rotational velocity ω, can be approximated with a piecewise affine relation. The conversion characteristics of the electric machine are approximated with two piecewise quadratic functions with a non-smooth convex union which are monotonically increasing in the domain of interest. The battery open circuit voltage is approximated with an affine function of the state-of-energy and the battery loss power is approximated with an internal resistance model. More information about the models of the individual components is provided in de Jager, van Keulen, and Kessels (2013, Ch. 2). The objective for power split control is to minimize the fuel consumption Ef for a known power request Pr (t ) and velocity trajectory ω(t ) with length t1 − t0 :



t1

Ef =

Pf (Pp , ω)dt .

(14)

t0

Here, Pm can be used as control variable by substituting Pp with Pr − Pm using (13). The optimization (14) is subject to the state dynamics: E˙ s (t ) = −Ps (Es , Pm , ω),

(15)

and to the power throughput limitations on the drive train components and the state constraints E s ≤ Es (t ) ≤ E s ,

Fig. 2. Power request and velocity input trajectories.

(16)

where Es is the energy stored in the battery and Ps is the corresponding power flow. The boundary states of the battery are also fixed; Es (t0 ) = Es0 and Es (t1 ) = Es1 .

Fig. 3. Optimal state-of-energy and multiplier trajectories. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Eqs. (13)–(16) describe an optimal control problem of the form

P0 . An explicit expression for the multiplier function (3) can be derived which is also locally Lipschitz continuous in p and continuous in t and the minimization of the Hamiltonian (8) provides an analytical expression for the optimal control value, see again de Jager, van Keulen, and Kessels (2013, pp. 63–64). 4. Evaluation of the novel algorithm In this section, the algorithm presented in Section 2.3 is evaluated on accuracy and computational demand. Besides, the influence of the battery state dependent losses is evaluated. The DP algorithm presented in de Jager and van Keulen (2013) is used as benchmark. The power and velocity input trajectories are displayed in Fig. 2. The input data is sampled with a constant sampling time of 1 s. The computed state and multiplier results are depicted in Fig. 3. Both the state dependent (blue) and state independent (red) results of the novel algorithm are compared with the state dependent DP solution (green-dashed) hereby using a grid size of ∆Ps = 1 kW. The constrained optimal solution is calculated in 3 iterations. The DP results overlap the results of the PMP algorithm (greendashed on top of blue line). The state dependent model has a state trajectory which is slightly higher than the state trajectory resulting from the state independent model. The battery open

T. van Keulen et al. / Automatica 50 (2014) 187–192

191

circuit voltage increases with the state-of-energy Es , thus, a power request Pb at a slightly higher state-of-energy requires a lower current than at low state-of-energy and therefore has lower losses. ∂P

Given the optimal control by (8) it follows that ∂∂PH = ∂ P f + m m −∂ P p ∂∂PPs = 0 which leads to the condition p = ∂ P f . From the DP m s solution (on a grid) it is possible to provide bounds on the multiplier value by computing p = −∂ Pf /∂ Ps ≈ −δ Pf /δ Ps . By using both positive and negative values for δ Ps we get an upper and lower bound for p due to the convexity of Pf (Ps ). Pf (k, Ps∗ (k)) − Pf (k, Ps∗ (k) + ∆Ps )

∆Ps

≤

≤ p(k)

Pf (k, Ps∗ (k)) − Pf (k, Ps∗ (k) − ∆Ps )

−∆Ps

.

(17)

Note that, (17) is not valid in the case control constraints are active. It can be seen in Fig. 3 that the solution of the novel algorithm lies between the upper and lower DP bounds if no constraints are active. Here, the pink squares indicate the upper bound and the green circles the lower bound. The accuracy in fuel consumption is depicted in Fig. 4. The DP solution is obtained for different grid sizes. Typically, the DP provides an upper bound in the calculated fuel, by increasing the number of grid points the fuel consumption lowers due to the increasing number of possible state trajectories. Note that, the solution obtained with the indirect approach lies at the point which is expected to be reached if the DP grid size approaches 0. So, it is concluded that a superior accuracy in fuel consumption and multiplier estimation is obtained with the novel algorithm. The computation times for DP, for different grid sizes, and the novel algorithm, both for the state dependent and independent situation, are also depicted in Fig. 4. All computations where performed on a standard R laptop with a 2.00 GHz Intel dual core chip using Matlab⃝ 2011a. The computation time of the proposed ‘‘indirect’’ algorithm depends on the location and the number of times the constraints are reached. In the worst case the computation time of the novel algorithm increases linearly with the number of times a constraint is reached. The length of each sub-trajectory reduces after each bound that is reached, such that the computation time increases less than linearly. The computation time of the novel algorithm compared to the DP algorithm, for a comparable case with a grid size of 1 [kW] (length Ps = 96 and length Es = 2280), is in the same order of magnitude in the case the state dynamic losses are accounted for and, at least, one order of magnitude less for the algorithm ignoring the state dependent losses. 5. Conclusions In this paper, a numerical solution for scalar convex optimal control problems with ‘‘pure’’ state constraints is developed. Using this solution, the optimal control of power split in hybrid vehicles based on a known power and velocity trajectory has been addressed. The novel algorithm calculates the global optimal power split trajectory without quantization of the state and control variables, handles state constraints and non-smooth characteristics, all with a computational effort, at least, one order of magnitude less than for Dynamic Programming. The computation accuracy is superior to that of Dynamic Programming. The novel algorithm can also include the state dependent battery losses. Acknowledgments The research leading to these results has received funding from the ENIAC Joint Undertaking and from Senter-Novem in the

Fig. 4. Accuracy and computation time vs DP grid size.

Netherlands under Grant Agreement number 120009, and is part of a more extensive project in the development of advanced energy management control for urban distribution trucks which has been made possible by TNO Business Unit Automotive in cooperation with DAF Trucks NV. References Bertsekas, D. P. (2000). Dynamic programming and optimal control. Belmont, USA: Athena Scientific. de Jager, B., & van Keulen, T. Receding horizon real-time energy management strategy for hybrid vehicles. Accepted for the ASME 2013 International Mechanical Engineering Congress & Exposition, San Diego, CA, USA, 2013. de Jager, B., van Keulen, T., & Kessels, J. (2013). Optimal control of hybrid vehicles. In Advances in industrial control. London Heidelberg New York Dordrecht: Springer. Delprat, S., Guerra, T.M., & Rimaux, J. Optimal control of a parallel powertrain: from global optimization to real time control strategy. In Proc. of the IEEE Vehicular Technology Conference (pp. 2082–2088), Birmingham, Al., USA, 2002. Fabien, B. C. (1996). Numerical solution of constrained optimal control problems with parameters. Applied Mathematics and Computation, 80, 43–62. Gerdts, M. (2008). A nonsmooth Newton’s method for control-state constrained optimal control problems. Mathematics and Computers in Simulation, 79, 925–936. Gerdts, M., & Kunkel, M. (2008). A non-smooth Newton’s method for discretised optimal control problems with state and control constraints. Journal of Industrial and Management Optimization, 4, 247–270. Hartl, R. F., Sethi, S. P., & Vickson, R. G. (1995). A survey of the maximum principles for optimal control problems with state constraints. SIAM Review, 37, 181–218. Jacobson, D. H., & Lele, M. M. (1967). A transformation technique for optimal control problems with a state variable inequality constraint. IEEE Transactions on Automatic Control, 14, 457–464. Khalil, H. K. (2002). Nonlinear systems (3rd ed.). New Jersey, USA: Prentice-Hall, Inc. Maurer, H. (1977). On optimal control problems with bounded state variables and control appearing linearly. SIAM Journal on Control and Optimization, 15, 345–362. Pontryagin, L. S., Boltyanskii, V. G., Gamkrelidze, R. V., & Mischenko, E. F. (1962). In K. N. Tririgoff, & L. W. Neustadt (Eds.), The mathematical theory of optimal processes. New York, USA: Wiley, Transl. Sager, S. Numerical methods for mixed integer optimal control problems. Ph.D. Thesis, Universität Heidelberg, Heidelberg, Germany, 2005. Seierstad, A., & Sydsæter, K. (1987). Advanced textbooks in economics: vol. 24. Optimal control theory with economic applications. Amsterdam, the Netherlands: Elsevier. Tate, E.D., & Boyd, S.P. (2000). Finding ultimate limits of performance for hybrid electric vehicles. In Proc. of society of automotive engineers 2000 future transportation technology conference (p. 12), paper 2000-01-3099, 00FTT-50. Vinter, R. (2000). Systems & control: foundations & applications. In Optimal control. New York, USA: Springer. Wang, L. Y., Gui, W. H., Teo, K. L., Loxton, R. C., & Yang, C. H. (2009). Timedelayed optimal control problems with multiple characteristic time points: computation and industrial applications. Journal of industrial and Management Optimization, 5, 705–718. Yu, C., Li, B., Loxton, R., & Teo, K. L. (2012). Optimal discrete-valued control computation. Journal of Global Optimization, 1–16.

192

T. van Keulen et al. / Automatica 50 (2014) 187–192

Thijs van Keulen received the M.Sc. degree in Mechanical Engineering from Delft University of Technology, Delft, The Netherlands, and the Ph.D. degree from Eindhoven University of Technology, Eindhoven, The Netherlands, in 2007 and 2011, respectively. Since 2011 he is with the Engine Test Center at DAF Trucks N.V., Eindhoven, The Netherlands. Since 2013, he also holds a position as a part time Assistant Professor in the Control Systems Technology group at Eindhoven University of Technology, Eindhoven, The Netherlands. His research interests include the optimal control theory, application of optimal control in hybrid vehicles and heavy-duty powertrains, and identification and control of heavy-duty diesel engines. Jan Gillot received the B.Eng. degree in Automotive Engineering from HAN University of Professional Education, Arnhem, The Netherlands, and the M.Sc. degree in Mechanical Engineering from Eindhoven University of Technology, Eindhoven, The Netherlands, in 2005 and 2010, respectively. Since 2011 he is with BRACE Automotive B.V., Eindhoven, The Netherlands. His research interests include On Board Diagnostics and optimization of Energy Management Strategies in hybrid vehicles.

Bram de Jager received the M.Sc. degree in mechanical engineering from Delft University of Technology, Delft, The Netherlands, and the Ph.D. degree from the Technische Universiteit Eindhoven, Eindhoven, The Netherlands. Currently, he is with the Technische Universiteit Eindhoven. He was with the Delft University of Technology and with Stork Boilers BV, Hengelo, The Netherlands. His research interests include robust control of (nonlinear) mechanical systems, integrated control and structural design, control of fluidic systems, control structure design, and applications of (nonlinear) optimal control.

Maarten Steinbuch received the Ph.D. degree from Delft University of Technology, Delft, The Netherlands, in 1989. Since 1999, he has been a Full Professor with the Control Systems Technology Group, Department of Mechanical Engineering, Eindhoven University of Technology, Eindhoven, The Netherlands. He is the Editor-in-Chief of IFAC Mechatronics and the Scientific Director of the Centre of Competence High Tech Systems of the Federation of Dutch Technical Universities. His research interests are in the modeling, design, and control of motion systems and automotive powertrains.