A two-phase computational scheme for solving bang-bang control problems

DOI 10.1007/s11081-006-0349-x

A two-phase computational scheme for solving

bang-bang control problems

Chih-Hung Huang· Ching-Huan Tseng

Received: 18 August 2003 / Revised: 6 July 2005 C

Springer Science+ Business Media, LLC 2006

Abstract This paper focuses on numerical methods for solving time-optimal con-trol problems using discrete-valued concon-trols. A numerical Two-Phase Scheme, which combines admissible optimal control problem formulation with enhanced branch-and-bound algorithms, is introduced to efficiently solve bang-bang control problems in the field of engineering. In Phase I, the discrete restrictions are relaxed, and the resulting continuous problem is solved by an existing optimal control solver. The information on switching times obtained in Phase I is then used in Phase II wherein the discrete-valued control problem is solved using the proposed algorithm. Two numerical examples, in-cluding a third-order system and the F-8 fighter aircraft control problem, are presented to demonstrate the use of this proposed scheme. Comparing to STC and CPET meth-ods proposed in the literature, the proposed scheme provides a novel method to find a different switching structure with a better minimum time for the F-8 fighter jet control problem.

Keywords Two-Phase Scheme . Time-optimal control problem . Bang-bang control . Enhanced branch-and-bound method . Admissible optimal control problem formulation

1 Introduction

Time-optimal control problems (TOCP) have attracted the interest of researchers in the area of optimal control because they often occur in practical applications. Thus, a

C.-H. Huang

Chief Technology Co., Ltd. No. 250, Sec. 4, Jhongsing Rd., Jhudong Township, Hsinchu County 310, Taiwan, R. O. C.

e-mail: tenor@ms9.url.com.tw; angushuang736@gmail.com C.-H. Tseng (_)

Department of Mechanical Engineering, National Chiao Tung University, Hsinchu 30056, Taiwan, R. O. C.

series of essential results has been published concerning applications of Pontryagin’s maximum principle to the time-optimal control of finite dimensional linear systems and low-order nonlinear systems. However, in the case of state- and/or control-constrained high-order nonlinear systems, solving the two-point boundary value problem that results from Pontryagin’s maximum principle is difficult. Moreover, analytic solutions are impractical if the dimension of a system exceeds three (Kirk, 1970). Therefore, in recent research, many numerical techniques have been developed and adopted to solve time-optimal control problems.

For a time-optimal control problem, one of the most common types of control function is the piecewise-constant function by which a sequence of constant inputs is used to control a given system with suitable switching times. Additionally, when the control is bounded, a very commonly encountered type of piecewise-constant control is bang-bang, which switches between the upper and lower bounds of the control input. When the controls are assumed to be of the bang-bang type, the time-optimal control problem becomes one of determining the switching times. Several methods for determining TOCP switching times have been extensively studied in the literature (see Kaya and Noakes, 1996; Bertrand and Epenoy, 2002; Simakov et al., 2002, for examples). However, in these methods, the number of switching times must be known before their algorithms can be applied. In most practical cases, however, the number of switching times is unknown before the control problems are solved. To overcome the numerical difficulties that arise during the process of finding the exact switching points, Lee et al. (1997) propose the Control Parameterization Enhancing Transform (CPET). CPET is also extended to deal with optimal discrete-valued control problems (Lee et al., 1999) and is applied to solving the sensor scheduling problem (Lee et al., 2001).

In like manner, this paper focuses on developing a numerical method to solve time-optimal control problems. This method consists of two computational phases: in the first, switching times are calculated using existing optimal control methods; and in the second, the resulting information is used to compute the discrete-valued control strategy. The proposed algorithm, which integrates the existing optimal control solver with an enhanced branch-and-bound method (Tseng et al., 1995), is implemented and applied to some example systems, including that of the F-8 fighter aircraft.

The remainder of the paper is organized as follows. Section 2 introduces the general formulation of the time-optimal control problem. Section 3 briefly discusses the theo-retical basis and the architectural framework of the admissible optimal control problem (AOCP) formulation. Section 4 presents an algorithm that integrates the AOCP with the enhanced branch-and-bound method. Section 5 presents the Two-Phase Scheme to overcome the difficulties associated with switching time computation methods. Section 6 then describes three nonlinear numerical examples taken from the literature and the numerical results obtained by applying the proposed algorithm. Section 7 concludes the paper.

2 Problem formulation

Consider a dynamical system described by the following nonlinear differential equa-tions on [0,tf]:

˙

x= f (t, b, x(t), u(t)) , t ∈0, tf 

with the initial condition

x(0)= x0, (2)

wheretf is the terminal time, b∈ Rkis the vector of design variables, u∈ Rm is a

vector of control functions and x∈ Rnis a vector of state variables. The function

f≡f1, f2,. . . , fn T

is assumed to be continuously differentiable with respect to all its arguments x0is a given vector inRn.

For a continuous control variable, any piecewise continuous function ucfrom [0,tf]

intoRmmay be taken as an admissible control function. For optimal discrete-valued control problems, a piecewise-constant function ud, ud :



0, tf 

→ Ud, may be taken

as an admissible control function, where Udis a finite set inRm. LetU be the class

of all such admissible control functions. Then a time-optimal control problem may be stated formally as follows: Given the dynamical system (1, 2), find u =[uc, ud]∈ U

such that the cost functional (performance index)

J0=

 tf 0

dt= tf (3)

is minimized subject to the constraint

Ji= i(b, x(tf), tf)+  tf 0 Li(b, u(t), x(t), t) dt  = 0; i = 1, . . . , Ne ≤ 0; i = Ne+ 1, . . . , NT (4) and the following continuous inequality constraint on the function of the state and control:

ψj(b, u(t), x(t), t) ≤ 0; j = 1, . . . , q, ∀t ∈ [0, tf]. (5)

wherei,Li, andψj are continuously differentiable with respct to their respective

arguments. Let this problem be referred to as Problem (PU). A control u∈ U is said to be a feasible control if it satisfies constraints (4) and (5). For a time-optimal control problem, the terminal time,tf, is not fixed and is treated as a design variable in b.

The differential equations that govern the system described by Eq. (1) are expressed in general first-order form. Equation (5) represents the mixed state and control con-straints, and the terminal conditions are treated as equality constraints in the first term of Eq. (4).

3 Admissible optimal control problem method

An admissible optimal control problem (AOCP) method was presented in earlier re-search (Sage and White, 1977) and used to solve continuous time-optimal control problems in Kaya and Noakes (1996); Lee et al. (1997). The core idea behind the AOCP method is to treat an optimal control problem as an initial value problem and

then adjust parameters and controls to satisfy optimality and constraints by applying iterative methods of nonlinear programming. SQP methods (Rao, 1996; Arora, 1989) are among the most widely used algorithms for solving general nonlinear program-ming problems and hence the SQP, which the BFGS method (Arora, 1989) applied to calculate the search direction, has been adopted in this paper. The approximate trajectory x(t ) is generated by solving the initial value problem defined in Eqs. (1) and

(2). Several effective and efficient numerical procedures can be used to integrate this initial value problem with the internal interpolation of the state variable x(t ). Hence, a

control parameterization technique is applied in which only control functions u(t ) are

discretized. The entire time interval [0,tf] is subdivided intoN unequal time intervals,

and the time grid points are defined in vector form as

T= [t1, t2, . . . , tN]T, (6)

where tN = tf. The control functions u(t ) can be discretized for each time

in-terval, and then the resulting control vector can be represented by the following equation:

U(D)=u(1)_{, u}(2)_{, . . . , u}(N )T

= [u1(t1), . . . , um(t1), u1(t2), . . . , um(t2), . . . , u1(tN), . . . , um(tN)]T (7)

where u(l) ∈ Rm, l = 1. . . N , is the vector of control variables for the l-th time

interval [tl, tl+1].

The time grid vector T and the discretized control vector U(D) are herein treated as design variables of the design variable vector P, and the terminal timetf is also

included as a design variable in T. To solve the system equation, the values of the control functions between two time grids are needed, which can be approximated by interpolating functions Iu, where Iu∈ Rm. The entire design variable vector can be

written as P = [bT_{, T}T_{, U}(D)T

]T_{. The admissible control functions are represented in}

the form u(P), and the state variable is expressed as x(P,t ) to emphasize that it is a

function of the design variable vector P.

In the AOCP method, the system equation, Eq. (1), and the initial condition, Eq. (2), constitute an initial value problem, and the corresponding state variables can be de-termined by using the values of the design variables obtained in each iteration of an SQP method, which the AOCP algorithm combines with the control discretization process. Figure 1 illustrates the architectural framework of the AOCP method (Huang and Tseng, 2003). This framework is based on an SQP method that cooperates with an optimal control problem solver responsible for solving the initial value problems and then using the solutions to calculate the values of cost functions (performance index) and constraints. An SQP algorithm uses these values to update the design variables and find the optimal solution.

Nonetheless, the AOCP method cannot deal with optimal discrete-valued con-trol problems; therefore, some mixed integer nonlinear programming techniques are needed. Hence, the algorithm developed in this paper is based on the AOCP but uses the enhanced branch-and-bound method (Tseng et al., 1995) to solve the optimal discrete-valued control problems. In each branching node of the branch-and-bound

Initialization H(0) = I, Initial guess P(0) Cost function J0(P(k)) Equality constraints h(P(k)₎ Inequality constraints g(P(k)) J0(P(k)) h(P(k)₎ g(P(k)₎ QP subproblem Descent direction d(k) Check convergence Line search

Search the step size (k)

Update Hessian matrix

H(k+1)

Update design variables

P(k+1) = P(k) + (k)d(k)

Show results

Discretization

Solve the initial value problems

Finite difference method

OCP Solver Design Variable P(k) Time Interval, N k = 0 k = k + 1 NO YES X(k)

Fig. 1 Conceptual flow chart of the AOCP algorithm

process, the AOCP is used as an optimal control problem solver to calculate continuous solutions.

4 Mixed-integer nonlinear programming algorithm for solving time-optimal control problems

Most discrete programming methods are based on the assumption that discontin-uous optimization problems are transformed into multiple contindiscontin-uous optimization sub-problems to take advantage of well-established continuous optimization algo-rithms. These continuous optimization problems are solved by imposing constraints on the discreteness of the design variables. The optimal discrete solution is taken

from among the continuous solutions obtained in the optimization sub-problems. However, the large number of discontinuous design variables greatly increases the number of the continuous optimization sub-problems. Tseng et al. (1995) presents an enhanced branch-and-bound method for reducing the number of executions of the continuous-optimization scheme by intelligently selecting the bounding route. Because such an enhanced branch-and-bound method dramatically reduces the to-tal number of continuous optimization runs executed and speeds up its convergence (Tseng et al., 1995), it is adopted herein and integrated with the AOCP to develop a mixed integer NLP algorithm for solving time-optimal control problems (TOCP). 4.1 Integrating the AOCP and enhanced branch-and-bound method

The algorithm developed in this paper consists of three major processes: branching, the AOCP, and bounding. Initially, all discrete-valued restrictions are relaxed and the resulting continuous NLP problem is solved using the AOCP. If the solution of con-tinuous optimum design problem occurs when all discrete-valued variable values are in the discrete set Ud, which is preset by the user to meet practical requirements,

then an optimal solution is determined and the procedure ends. Otherwise, the algo-rithm selects one of the discrete-valued variables whose value is not in the discrete set Ud—for example, the j -th design variable, Pj, with value ˆPj—and branches

on it.

4.1.1 Branching process

In the branching process, the original design domain is divided into three sub-domains by two allowable discrete values, ¯ui and ¯ui+1, that are nearest to the continuous

op-timum, as shown in Fig. 2. Of the three sub-domains, sub-domain II, included in the continuous solution but not in the feasible discontinuous solution, is dropped. In the other two sub-domains, called nodes, two new NLP problems are formed by adding simple bounds, ˆPj≤ ¯uiand ˆPj ≥ ¯ui+1, respectively, to the continuous NLP problems.

One of the two new NLP problems is selected and solved next. Many search methods based on tree searching—including depth-first search, breadth-first search and best-first search—can be applied to choose the next branching node. The branching process is repeated in each of the sub-domains until a feasible optimal solution is found in which all the discrete variables have allowable discrete values. Obviously, the number of sub-domains may grow exponentially so that a great deal of computing time is re-quired. Thus, in the enhanced branch-and-bound method (Tseng et al., 1995), multiple branching and unbalanced branching strategies have been developed to improve the method’s efficiency.

4.1.2 Bounding process:

In discrete optimization, the minimum cost is always greater than or equal to the cost of the original regular optimal design that was originally branched. This fact provides a guideline for when branching should be stopped. If the branching process yields a feasible discontinuous solution, then the corresponding cost value can be considered a bound. Any other sub-domain that imposes a continuous minimum cost larger than

Fig. 2 Conceptual layout of branching process

this bound does not need to be branched further. This bounding strategy can be used to select the branching route intelligently and avoid the need for a complete search over all the branches.

Algorithm: the AOCP and enhanced branch-and-bound

Initialization: Relax all discrete-valued restrictions and then place the resulting

con-tinuous NLP problem on the branching tree. Set the cost bound Jmax = ∞.

while (there are pending nodes in the branching tree) do 1. Select an unexplored node from the branching tree. 2. Control discretization.

3. Repeat (fork-th AOCP iteration )

(1) Solve the initial value problem for state variable x(k)_{of AOCP.}

(3) Solve the QP(k) _{problem by applying the BFGS method to obtain the}

descent direction d(k).

(4) if (QP(k)is feasible and convergent) then exit AOCP

(5) Find the step sizeα(k)of the SQP method by using the line search method. (6) Update the design variable vector: P(k+1)= P(k)+ α(k)d(k)_.

4. if (NLP is optimal) and (J0<Jmax) then

if (U(D)(k+1)is feasible ) then

Update the current best point by setting the cost boundJmax= J0. Add this node to the feasible node matrix.

else

Evaluate the values of criteria for selecting the branch node. Choose a discrete-valued variable U(_iD)(k+1) /∈ U and branch it. Add two new NLP problems into the branching tree.

Drop this node. endif

else

Stop branching on this node. endif

end while

5 Two-phase scheme for solving TOCP

The mixed integer NLP algorithm developed in this paper is one type of switch-ing time computation (STC) method. Most switchswitch-ing time computation methods [see, for example, (Kaya and Noakes, 1996; Lucas and Kaya, 2001; Simakov et al., 2002)] assume that the structure of the control is bang-bang and the number of switching times is known. Unfortunately, the information on the switchings of sev-eral practical time-optimal control problems is unknown and hard to compute using analytical methods. Hence, to overcome this difficulty, this paper proposes a Two-Phase Scheme that consists of the AOCP plus the mixed integer NLP method. In Phase I, the AOCP is used to calculate the information on switchings with rough time grids so that the information can be used in Phase II as the feasible ini-tial design of the mixed integer NLP method. This scheme is described briefly below.

Phase I

1. Follow the steps of the AOCP method proposed by Huang and Tseng (2003) to solve the time-optimal control problem using continuous controls.

2. Based on the numerical results, extract information about the switching times and terminal times,tf.

Phase II

3. Based on the information about switching times obtained in Phase I, the switchings are treated as design variables and added into the time grid vector T that is given in

Section 3. It should be noted that each interval between the upper and lower bounds on each of those design variables must include one switching.

4. Insert the terminal time,tf, into the design variable vector P.

5. Discretize each control variable into the number of switchings plus one. Then the discrete control vector, U(D), defined in Eq. (7) can be add to the design variable vector P and limit their corresponding upper and lower bounds by the original bounds of the controls.

6. Solve the problem by applying the mixed integer NLP method, and then find the optimal discrete-type control trajectories.

A third-order system shown in Section 6.1 is used to demonstrate those processes of this numerical scheme.

6 Numerical examples

The numerical results for the following examples are obtained on an Intel Celeron 1.2 GHz computer with 512 MB of RAM memory. The AOCP is coded in FORTRAN, and C language is used to implement the enhanced branch-and-bound method. The Visual C++ 5.0 and Visual FORTRAN 5.0 installed in a Windows 2000 operating system are adopted to compile the corresponding programs. The total CPU times for solving the F-8 fighter craft problem in Phase I and Phase II are 3.605 and 1.782 seconds, respectively.

6.1 Third-order system

The following system of differential equations is a model of the third-order system dynamics taken from Wu (1999).

˙ x1 = x2 (8) ˙ x2 = x3 (9) ˙ x3 = −10x3+ 10u (10)

The problem here is to find the control|u| ≤ 10 in order to bring the system from the initial state [−10, 0, 0]T to the final state [0, 0, 0]T in minimum time.

First, this problem is solved directly by the mixed integer NLP method. Assume four switching times (T1, T2, T3, T4) and five control arcs have values in the discrete set, Ud:

{−10, 10}. The terminal time, tf, is treated as a design variable, so the design variable

vector P can be expressed as [T1, T2, T3, T4,tf,Ud1, Ud2, Ud3, Ud4, Ud5]T. Notably,

the final conditions of the state variables are transferred to the equality constraints. The TOCP problem becomes one of determining the switching times. Figure 3(a) presents the continuous solution obtained by using the AOCP and the discrete solution determined by applying the mixed integer NLP method proposed herein. The results indicate that the control trajectory determined by the mixed integer NLP method is of the bang-bang type and the solution consistent with the results obtained by Wu (1999).

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 Time (sec.) -12 -10 -8 -6 -4 -2 0 2 4 6 8 10 12 14 Con tr ol (u )

Mixed Discrete NLP method (tf = 0.7813) AOCP (tf = 0.7826)

(a) AOCP vs. Mixed integer NLP method

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 Time (sec.) -12 -10 -8 -6 -4 -2 0 2 4 6 8 10 12 14 Con tr ol ( u ) Phase I (tf = 0.7864) Phase II (tf = 0.7813) (b) Phase I vs. Phase II Fig. 3 Control trajectories of the third-order system

As stated in Section 5, several assumptions must be made when the mixed integer NLP method is applied to solving TOCP directly. Unfortunately, these assumptions cannot be guaranteed to hold in practical cases. Consequently, the Two-Phase Scheme proposed in this paper is needed. For illustration, the third-order system is again solved using this Two-Phase Scheme. In Phase I, two switching times are found to be [0.330, 0.725]T and the terminal time tf is 0.7864. In the first phase, these switching data

need not be accurate because they are only used to help users decide on the number of switching times, the control arcs and their corresponding boundaries. Thus, in Phase II, the design variable vector P is re-formed as [T1, T2,tf,Ud1, Ud2, Ud3]T,

and the numerical result obtained by applying the mixed integer NLP method is as presented in Fig. 3(b). In Phase II, the switching times of discrete control input are [0.323, 0.713]T_{, and the terminal timet}

f is 0.7813 seconds. The control trajectory

also agrees with that obtained by Wu (1999). 6.2 F-8 fighter aircraft

The F-8 fighter aircraft has been considered in several pioneering studies (Kaya and Noakes, 1996; Banks and Mhana, 1992; Simakov et al., 2002), for example and has become a standard for testing various optimal control strategies. A nonlinear dynamic model of the F-8 fighter aircraft is considered below. The model is represented in state

0.0 1.0 2.0 3.0 4.0 5.0 6.0 Time(sec.) -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 0.00 0.01 0.02 0.03 0.04 0.05 0.06 Phase I (Tf = 5.74173) Phase II (Tf = 5.74216)

space by the following differential equations. ˙ x1= −0.877x1+ x3− 0.088x1x3+ 0.47x12− 0.019x 2 2− x 2 1x3+ 3.846x13 −0.215u + 0.28x2 1u− 0.47x1u2+ 0.63u3, (11)

Fig. 5 Trajectories of the states

and control input for the F-8 fighter aircraft

Table 1 Results of various methods for the F-8 fight aircraft problem

Accuracy of

Method tf Switching times terminal constraints

STC (Kaya and Noakes, 1996) 6.3867 0.0761, 5.4672, 5.8241, 6.3867 ≤ 10−5 CPET (Lee et al., 1999) 6.0350 2.188, 2.352, 5.233, 5.563 ≤ 10−10 Two-Phase Scheme 5.7422 0.098, 2.027, 2.199, 4.944, 5.265 ≤ 10−10 ˙ x2= x3, (12) ˙ x3= −4.208x1− 0.396x3− 0.47x12− 3.564x 3 1 − 20.967u +6.265x2 1u+ 46x1u2+ 61.4u3, (13)

wherex1is the angle of attack in radians,x2is the pitch angle,x3is the pitch rate and the control inputu represents the tail deflection angle. For convenience of comparison, the

standard settings (Kaya and Noakes, 1996; Lee et al., 1997) are used. The control|u| ≤ 0.05236 must be found that brings the system from its initial state26.7π/180, 0, 0T to the final state0, 0, 0T in minimum time.

When the Two-Phase Scheme is applied, as described in Section 5, the switch-ing times computed in Phase I are 0.115, 2.067, 2.239, 4.995, and 5.282, and the terminal time is tf = 5.7417. These switching data are used to set the

de-sign variables and their corresponding bounds, and then the problem is solved by the mixed integer NLP method. Finally, the switching times for the discrete con-trol input are 0.098, 2.027, 2.199, 4.944, and 5.265, and the terminal time tf is

5.74216. Figure 4 shows the comparison of the controls between Phase I and Phase II, while Fig. 5 shows the trajectories of the states and the control of Phase I and Phase II. This example is also solved by Kaya and Noakes (1996 using the switching time computation method and by Lee et al. (1997) using the Control Parameterization Enhancing Transform method. Table 1 shows the terminal time

tf, switching times and the accuracy of terminal constraints computed by

vari-ous methods for this problem. According to the numerical results, the Two-Phase Scheme provides a better solution, and the accuracy of the terminal constraints is acceptable.

7 Conclusion

This paper discusses a method for solving time-optimal control problems with discrete-type control inputs that include the bang-bang discrete-type most commonly encountered when the control is bounded. This two-phase computational scheme for finding a discrete optimal control for time-optimal control problems is novel because its discrete control can be more easily implemented than continuous control in practical applications. A simple example, a third-order system, is presented to demonstrate the usage of the proposed scheme. An F-8 fighter aircraft control problem is considered and solved by application of the proposed scheme. Numerical results are obtained efficiently and

accurately. The results reveal that the Two-Phase Scheme constitutes a viable method for solving time-optimal control problems with discrete-valued controls.

Acknowledgments The authors would like to thank the National Science Council of the Republic of China

for financially supporting this research under Contract No. NSC91-2212-1-009-022.

References

Arora JS (1989) Introduction to optimum design. McGraw-Hill, New York

Banks SP, Mhana KJ (1992) Optimal control and stabilization of nonlinear systems. IMA Journal of Mathematical Control and Information 9:179–196

Bertrand R, Epenoy R (2002) New smoothing techniques for solving bang-bang optimal control problems – numerical results and statistical interpretation. Optimal Control Applications and Methods 23:171–197 Huang CH, Tseng CH (2003) Computational algorithm for solving a class of optimal control problems. IASTED International Conference on Modelling, Identification, and Control (MIC2003), Innsbruck, Austria, pp 118–123

Kaya CY, Noakes JL (1996) Computations and time-optimal controls. Optimal Control Applications and Methods 17:171–185

Kirk DE (1970) Optimal control theory: an introduction. Prentice-Hall

Lee HWJ, Jennings LS, Teo KL, Rehbock V (1997) Control parameterization enhancing technique for time optimal control problems. Dynamic Systems and Applications 6:243–262

Lee HWJ, Teo KL, Rehbock V, Jennings LS (1999) Control parameterization enhancing technique for optimal discrete-valued control problems Automatica 35:1401–1407

Lee HWJ, Teo KL, Andrew EB Lim (2001) Sensor scheduling in continuous time. Automatica 37:2017– 2023

Lucas SK, Kaya CY (2001) Switching-time computation for bang-bang control laws. In: Proceedings of the 2001 American Control Conference, pp 176–181

Rao SS (1996) Engineering optimization: theory and practice. 3rd edn. Jonh Wiley & Sons Sage AP, White CC (1977) Optimum systems control, 2nd edn. Prentic Hall

Simakov ST, Kaya CY, Lucas SK (2002) Computations for time-optimal bang-bang control using a la-grangian formulation. 15th Triennial World Congress, Barcelona, Spain

Tseng CH, Wang LW, Ling SF (1995) Enhancing branch-and-bound method for structural optimization. Journal of Structural Engineering 121(5):831–837

Wu ST (1999) Time-optimal control and high-gain linear state feedback. International Journal of Control 72(9):764–772