Pseudospectral optimal train control

Delft University of Technology

Pseudospectral optimal train control

Goverde, Rob M.P.; Scheepmaker, Gerben M.; Wang, Pengling

DOI

10.1016/j.ejor.2020.10.018

Publication date

2021

Document Version

Final published version

Published in

European Journal of Operational Research

Citation (APA)

Goverde, R. M. P., Scheepmaker, G. M., & Wang, P. (2021). Pseudospectral optimal train control.

European Journal of Operational Research, 292(1), 353-375. https://doi.org/10.1016/j.ejor.2020.10.018

Important note

To cite this publication, please use the final published version (if applicable).

Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons. Takedown policy

Please contact us and provide details if you believe this document breaches copyrights. We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

ContentslistsavailableatScienceDirect

European

Journal

of

Operational

Research

journalhomepage:www.elsevier.com/locate/ejor

Interfaces

with

Other

Disciplines

Pseudospectral

optimal

train

control

Rob

M.

P.

Goverde

a,∗

_,

_Gerben

_M.

_Scheepmaker

a,b

_,

_Pengling

_Wang

a

a Department of Transport and Planning, Delft University of Technology, Delft, the Netherlands

b Netherlands Railways, Department of Performance management and Innovation, Utrecht, the Netherlands

a

r

t

i

c

l

e

i

n

f

o

Article history: Received 8 August 2019 Accepted 13 October 2020 Available online 30 October 2020

Keywords:

Optimal train control Train trajectory optimization Pontryagin’s Maximum Principle Pseudospectral method Singular solution

a

b

s

t

r

a

c

t

Inthelastdecade,pseudospectralmethodshavebecomepopularforsolvingoptimalcontrolproblems. Pseudospectralmethodsdonotneedpriorknowledgeabouttheoptimalcontrolstructureandarethus veryflexibleforproblemswithcomplexpathconstraints,whicharecommoninoptimaltraincontrol,or traintrajectoryoptimization.Practicaloptimaltraincontrolproblemsarenonsmoothwithdiscontinuities inthedynamicequationsandpathconstraintscorrespondingtogradientsandspeedlimitsvaryingalong thetrack. Moreover,optimaltraincontrol problemstypicallyincludesingularsolutions witha vanish-ingHessianoftheassociatedHamiltonian.Thesecharacteristicsmaketheseproblemshardtosolveand alsoleadtoconvergenceissuesinpseudospectralmethods.Weproposeacomputationalframeworkthat connectspseudospectralmethodswithPontryagin’sMaximumPrincipleallowingflexiblecomputations, verificationandvalidationofthenumericalapproximations,andimprovementsofthecontinuous solu-tionaccuracy.Weapplytheframeworktotwobasicproblemsinoptimaltraincontrol:minimum-time traincontrolandenergy-efficienttraincontrol,andconsidercaseswithshort-distanceregionaltrainsand long-distanceintercitytrainsforvariousscenariosincludingvaryinggradients,speedlimits,and sched-uledrunning timesupplements. The frameworkconfirms theflexibility ofthepseudospectral method with regardstostate, control and mixed algebraicinequalitypathconstraints, and is ableto identify conditionsthatleadtoinconsistenciesbetweenthenecessaryoptimalityconditionsand thenumerical approximationsofthestates,costates,andcontrols.Anewapproachisproposedtocorrectthediscrete approximationsbyincorporatingimplicitequationsfromtheoptimalityconditions.Inparticular,the is-sueofoscillationsinthesingularsolutionforenergy-efficientdrivingascomputedbythepseudospectral methodhasbeensolved.

© 2020TheAuthor(s).PublishedbyElsevierB.V. ThisisanopenaccessarticleundertheCCBYlicense(http://creativecommons.org/licenses/by/4.0/ )

1. Introduction

Optimal control theory is widely applied indifferent fieldsto find the controls that minimize a cost functional subject to dy-namic constraints,pathconstraintsandboundaryconditions.One oftheapplicationsofoptimalcontroltheoryisminimum-timetrain control (MTTC) with the aim to compute the minimum time be-tween two train stops. Another important application is energy-efficient train control (EETC) or also referred to asenergy-efficient traintrajectoryoptimization,withtheaimtominimizetotaltraction energy. The train controls consist of traction and braking. These optimaltraincontrolproblemsarehighlycomplexinpracticewith discontinuousdynamicequations,purestateandcontrolalgebraic inequality constraints, andmixed state-control algebraic

inequal-∗_{Corresponding author.}

E-mail addresses: r.m.p.goverde@tudelft.nl (R.M.P. Goverde),

g.m.scheepmaker@tudelft.nl (G.M. Scheepmaker), p.l.wang@tudelft.nl (P. Wang).

ityconstraints. Moreover, the solution to the EETC problem gen-erallycontains singulararcs.Hence, solving optimal traincontrol problemsisquitechallenging,andspecificallyEETCproblems.Most research on optimal train control has focused on the application of Pontryagin’s Maximum Principle (PMP) (Pontryagin, Boltyanskii, Gamkrelidze,& Mishchenko,1962).The applicationofthistheory leadsto the optimal driving regimes consistingof maximum ac-celeration,cruising,coastingandmaximumbraking.Thechallenge is to determine the optimal switching structure with the exact switching points aswell as the optimal sequence of the driving regimes.Scheepmaker, Goverde,andKroon (2017) provided a re-cent extensiveliterature review on EETC. One ofthe conclusions was that the pseudospectral method was a promising approach to solve generic EETC problemswith varying trackgradients and speed limits. In this paper, we further explore the Radau pseu-dospectral method to optimal train control based on structured numericalexperiments thatare assessed onconsistency withthe optimalityconditionsobtainedbythePMP.Weidentifynumerical

https://doi.org/10.1016/j.ejor.2020.10.018

inaccuraciesduetoconvergenceissues,discretizationand oscillat-ing singularcontrol,andpropose anewapproachbased on com-biningunderdeterminedimplicitequationsfromthePMPwiththe pseudospectral approximate solutions to get accurate continuous optimalcontrolsolutions.

Optimal control problemscanbe solvedby indirectanddirect methods (Betts, 2010;Rao, 2014). Anindirectmethod derivesthe necessary optimality conditions based on PMP, which leads to a two-point boundary value problem (BVP) in the states and ad-jointcostatesthatneedstobesolvedbynumericalmethods(Ross, 2005). These methods thus indirectly solve the original optimal controlproblembysolvinganassociatedBVP.Ingeneral,however, theBVPishardtosolve,asitisverysensitivetotheinitialguesses fortheunknown(costate)boundaryconditionsandapriori knowl-edge oftheswitchingstructure oftheinequality pathconstraints is required.Therefore,mostindirectmethods use(heuristic) con-structivemethodstofind(sub)optimaldrivingstrategiesbasedon the PMPoptimality conditions,with (linearity) simplifications or assumptions on drivingregimes (e.g.nocoasting,no cruising be-low the speed limit). Mostpapers on optimal train control used indirect methods based on PMP optimality conditions (Albrecht, Howlett, Pudney, Vu, & Zhou, 2016a; 2016b; Howlett & Pudney, 1995;Khmelnitsky,2000;Liu&Golovitcher,2003).Forthe energy-efficient train control problem, algebraicformulaefor thecostate along sections withconstant gradient can be derived, which led to efficient real-time algorithms (Albrecht et al., 2016a; 2016b; Howlett,Pudney,&Vu,2009;Liu&Golovitcher,2003).Inaddition, many heuristicmethods havebeen applied usingimplicit knowl-edgeoftheoptimalcontrolstructure,forexampleChevrier, Pelle-grini,andRodriguez(2013),Sicre,Cucala,andFernández-Cardador (2014)andHaahr, Pisinger, andSabbaghian(2017).Formoredetails andreferences,seeScheepmakeretal.(2017).

On the other hand, direct solution methods transcribe the infinite-dimensional optimal control problem into a finite-dimensionalnonlinear programming(NLP) problemusing colloca-tion of thedifferential equations andthe integral objective func-tional. The resulting NLP problem is then solved using efficient nonlinear optimization algorithms (Betts, 2010). This direct ap-proachcanbeusedforhighlycomplexproblems,sincethereisno needtoderivethefirst-ordernecessaryoptimalityconditions,they arelesssensitivetoinitialsolutions,andnoaprioryknowledgeis neededoftheactiveandinactiveinequalitypathconstraints.Inthe optimaltraincontrolliteraturedirectmethodshavebeenusedonly recently. Most approaches are based on pseudospectral methods (Scheepmaker & Goverde, 2016;Wang & Goverde, 2016a; 2016b; 2017; 2019; Wang, De secondchutter, Van den Boom, & Ning, 2013; Ye& Liu,2016; 2017). Wangetal. (2013); Wang, DeSchut-ter,Van denBoom,andNing(2014) alsoapplied amixed-integer linear programming(MILP)approach usingpiecewise affine func-tions. At this stage the computationtimes forthe direct method are not competitive withthe computation times forthe indirect methodsthatarecurrentlyusedforon-boardcalculations.

Over the past two decades, pseudospectral methods have be-come popularforsolvingoptimalcontrol problemswithin partic-ular the aerospace domain (Garg et al., 2010; Ross & Karpenko, 2012). A pseudospectral method discretizes thestate andcontrol and then transcribes the continuous-time optimal control prob-lem to a finite-dimensional NLPthat is solved using established NLPsolvers.Inapseudospectralmethodthestateandcontrolare approximated by global polynomials at collocation points using a basis of Lagrange (or Chebyshev)polynomials. Typical colloca-tionpointsareLegendre-Gauss(LG),Legendre-Gauss-Lobatto(LGL), andLegendre-Gauss-Radau(LGR)points(Gargetal.,2010;Ross& Karpenko,2012).Thesepointsareobtainedfromtherootsofa Leg-endrepolynomialand/oritsderivatives.Theyarealldefinedonthe domain [−1,1] butdiffer inthattheLG pointsincludeneitherof

theendpoints, theLGRpointsincludeone endpoint,andtheLGL points include both endpoints. Different pseudospectral methods have beendeveloped based on the collocation points: theGauss (LG)pseudospectralmethods (Benson,Huntington,Thorvaldsen,& Rao, 2006; Rao et al., 2010), the Legendre (LGL) pseudospectral methods(Elnagar,Kazemi,&Razzaghi,1995;Fahroo&Ross,2001; Gong,Fahroo,&Ross,2008a;Ross&Fahroo,2004),andtheRadau (LGR)speudospectral methods (Garg etal., 2011). Thedifferences betweenthemethods lie inthethedegree oftheLagrange poly-nomial, the boundary conditionand the typical problemhorizon (Gargetal.,2009;Ross&Karpenko,2012).Theuseofglobal poly-nomialstogetherwithGaussianquadraturecollocationpoints pro-videsaccurateapproximationsthatconvergeexponentiallyfastfor smoothproblems(Gargetal.,2010).Foroptimalcontrolproblems with discontinuities in the constraints the problem can be par-titioned into multiple phases, where the phase boundaries(also calledmesh pointsorknots) canbe chosen at thepoints of dis-continuities (Betts, 2010). Inthesemultiple-phase optimalcontrol problems,each phase issolved witha separate setofcollocation points,whileadditionallinkingconditionsgluethevariablesofthe adjacentphasestogether(Darby,Garg,&Rao,2011a;Darby,Hager, &Rao,2011b;Patterson,Hager,&Rao,2015;Patterson&Rao,2014; Raoetal.,2010;Ross&Fahroo,2004).Hence,withmultiplephases thestateisapproximatedbymultiplepolynomialsthatarelinked atthephaseboundaries.

The pseudospectral methods can also compute estimations of thecontinuouscostatesbasedontheassociatedLagrange multipli-erscomputedforthediscretizedNLPproblem(Darbyetal.,2011a; Fahroo & Ross, 2001;Garg etal., 2011; 2010). FortheGauss and Radaupseudospectralmethodsthediscreteapproximationsare ob-taineddirectlybyatransformationoftheLagrangemutipliers us-ingthequadratureweightsofthecollocationpoints.Forthe Legen-drepseudospectralmethodsthediscretecostatesarenotuniquely determined andneeda closureconditionrelatedto the transver-sality conditions of the necessary optimality conditions for the continuous optimal control problem (Fahroo & Ross, 2001; Garg et al., 2010; Gong, Ross, Kang, & Fahroo, 2008b). The covector mapping principle forthe Legendre Pseudospectral method then statesthat these multipliersofthe discretized NLPproblem con-vergeto thecostates ofthe discretizednecessaryoptimality con-ditions(Gongetal., 2008a;Gong etal.,2008b; Ross&Karpenko, 2012).ComparativestudiesshowedthattheRadaupseudospectral methodprovidesthemostaccurateresults,includingthecostates, witha fastconvergencerate(Garg,2011;Gargetal.,2009; Hunt-ington,2007).Moreover,theLGRcollocation pointsaremost suit-able for multiple-phase optimal control problems since they in-cludeone endpointper phase, andthusall phase boundariesare collocation points corresponding to discontinuities in the input data.Together,theLGRcollocationpointsthuscoverallphasesand theirboundariesexcepttheterminalpointwhichishowever esti-matedwellbythestateapproximations.TheRadaupseudospectral method is implemented in the MATLAB toolbox GPOPS (General PseudospectralOPtimalControlSoftware)(Raoetal.,2010),which weusedinthispaperfortheexperiments.

The pseudospectral method works well forsmooth problems, buttheoptimaltraincontrolproblemsmaybenonsmoothin prac-tice due to infrastructural constraints. In particular, the varying slopesoftherailwaytracksareusuallymodelledaspiecewise con-stantgradientsresultinginadiscontinuousdynamicequation,and varying speedlimits alongthe trackresultindiscontinuousstate inequalitypathconstraintsthatmayleadtodiscontinuouscostate functions (Bryson & Ho, 1975). Moreover, the optimal train con-trolproblemshaveaHamiltionianthatislinearinthecontroland thus has a zero Hessian.This may causesingular arcsfor which noexplicit analytical expressions canbe derived fromthe neces-sary optimality conditions and for which the pseudospectral

so-lutions showoscillatingcontrol andstate behaviour.Furthermore, stop distances maybe very long withpossibly manyfluctuations in gradients andspeed limits leadingto manyphases ofthe re-sultingmultiple-phaseoptimalcontrolproblemwithasmany inte-riorstate equalityconstraints.Theliteratureonoptimaltrain con-trol using pseudospectral methods indeed showsirregularities in the control and/or state profile plots that need tobe understood before these methods can be used inpractice (Wang & Goverde, 2016a;2019;Ye&Liu,2016;2017).YeandLiu(2017)comparedthe pseudospectralmethodwithaheuristicsolutionmethodconsisting ofan NLPproblemformulationbasedonclosed-formexpressions for each driving regime obtainedwithsome simplifying assump-tions.Foronecasethepseudospectralmethodfoundthesame en-ergyconsumption,despitesomefluctuationsinthesingular cruis-ing regime. In asecond case studytheheuristic methodfound a solutionwithlessenergyconsumptionandtheyconcludedthatit isquitedifficulttochooseappropriateparametersofGPOPSto en-sureconvergencewithin acceptabletolerances. Wangetal.(2013, 2014) compared a pseudospectral method withan MILP method. The pseudospectral methodprovided solutionswiththe least en-ergy but wasmuch slower than the MILP method.However, the MILPmodelshowedaroughapproximationwhilethe pseudospec-tral method gave smooth results,which raises the questionhow accurate thepseudospectral mustbe:withlesscollocation points the pseudospectral method is faster and may still perform bet-ter, or theother wayaround, enforcing a higher accuracy ofthe MILPmodelwillincreaseitscomputationtime.WangandGoverde (2016a)andYeandLiu(2016,2017)observedstrongoscillating be-haviourinboththecontrolandstateforthesingularsolution cor-respondingtoanonsmoothapproximationofthecruisingregime. Zhong,Lin,Loxton,andTeo(2019)modelledtheoptimaltrain con-trolproblemasanoptimalswitchingcontrolproblem.Theydivide thetrackintoafinitenumberofsegmentsofconstantgradientand speed limitsimilartoamulti-phaseoptimalcontrol problem,and thenusecontrolparametrizationandtime-scalingtoobtainan ap-proximatefinite-dimensionalnonlinearprogrammingproblem.The traction and braking control are approximated aspiecewise con-stantfunctionswheretheswitchingtimesandjumpsarethe deci-sionvariables.Theycomparedtheir algorithmwitha pseudospec-tral methodusingGPOPS.Thecomputationtimesare comparable andtheir solutiondoesnot show singularcontroloscillations, al-thoughthetractionandbrakingcontrolarenowapproximatestep functions.Instead,thepseudospectralmethodshowsaccurate con-tinuoustractioncontrolexceptforthesingularcontroloscillations. Chen and Biegler (2016) proposed a nested direct transcription optimization method for solving singular optimal control prob-lems. The nonlinearprogrammingproblemresultingfromthe di-recttranscriptionisdecomposedintoaninnerandouterproblem. Intheinnerproblemmovingfiniteelementsareusedtofind accu-rateswitchingtimesandadditionalmonotoniccontrolconstraints oneachfiniteelementguaranteealow-ordercontrol. ‘Pseudomul-tipliers are introduced that reconstruct the necessary optimality conditionsforthesingularoptimalcontrolintheouterproblem.

In this paper, we propose a computational evaluation frame-work where the PMP is applied to verify, validate and improve pseudospectral solutions to optimal train control problems. It is based on the costate approximations that the pseudospectral methods provide nextto the controlsandstates.Thiscan be ex-ploited to analyse thepseudospectral results using thenecessary optimality conditions of the original continuous optimal control problem. It is shown that there are inconsistencies due to the discretization ofthe continuousproblem.The solutions can how-ever be corrected using analytical results from the PMP analysis ina postprocessingstepresultinginfeasiblecontinuoussolutions satisfying the optimality conditions. We consider both the MTTC and EETC problems as examples, and exclude regenerative

brak-ing. However, the proposed approach can be applied to any op-timaltrain control problem formulation.The paperthereby gives thefollowingcontributionstotheliterature:

1. A computational framework connects the pseudospectral methodwithPontryagin’sMaximumPrinciple.

2. Theframeworkisappliedtocompute,validateandimprove so-lutionstooptimaltraincontrolproblems.

3. Structuredexperimentsillustrate theimpact ofstop distances, varying speed limitsand gradients, andrunning time supple-ments.

4. Convergence issuesare identified for discontinuousstate con-straints(speedlimits)anddynamicequations(gradients)with bigjumps.

5. Knowncomputationalissuesofthepseudospectralmethod re-gardingsingularsolutionsaresolved byahybridapproach us-ingPontryagin’sMaximumPrinciple.

Thepaperisstructuredasfollows.Section 2definestheMTTC andEETC problemsandderivestheir necessary optimality condi-tions byapplicationof thePMP. Thisreflects thetraditional indi-rect optimaltrain control solutionapproach. Then,the numerical pseudospectral methodis introducedin Section 3.Section 4 pro-posesthecomputationalevaluationframeworkthatcombinesPMP and the pseudospectral method to verify and validate the solu-tions obtained from the pseudospectral methods. Section 5 ap-plies the computational framework to the MTTC andEETC prob-lems for various structured scenarios and identifies the incon-stistenciesofthe approximatedpseudospectralsolutions withthe PMP.Section6thendiscussesthenumericalchallengesofthe dis-cretized pseudospectral solutions and proposes a postprocessing step toobtain feasible continuoussolutions by exploiting knowl-edgefromthePMP.Section7endsthepaperwiththeconclusions.

2. Theoptimaltraincontrolproblemsandnecessaryoptimality conditions

Thissectiongivestheproblemformulationsoftheoptimal con-trolproblemforboth the EETC andMTTC, andderivesnecessary conditionsforoptimalitybyapplicationofthePMP(Lewis,Vrabie, &Syrmos,2012;Pontryaginetal.,1962;Ross,2015).First,theEETC problemisdiscussedandafterwardstheMTTCproblem.The prob-lemformulationsconsiderdistanceasindependentvariablerather thantime,becausediscontinuitypointsassociatedwithchangesin gradientsandspeedlimit arenaturallygivenintermsofdistance (Scheepmaker&Goverde,2015;Scheepmakeretal.,2017;Wang& Goverde, 2016a). To focus on the essence of the optimal control problems, we do not consider the effect of regenerative braking andwemodelthetrainwithoutlossofgeneralityasapointmass (Brünger&Dahlhaus,2014;Howlett&Pudney,1995).The applica-tionofthePMPtooptimaltraincontrolproblemsisnotnew,butis includedheretomakethepaperself-contained.Thederived equa-tionsfromthenecessaryoptimalityconditionswillbeusedinthe restofthepaper.Otherpapersthatalsogive in-depthderivations oftheoptimalityconditionsbyapplyingPMPconsiderslightly dif-ferent problem formulations, such as energy and time as state variables (Khmelnitsky, 2000), normalized control variables (Liu & Golovitcher, 2003),and adifferent objectivefunction including regenerativebraking (Albrecht et al., 2016a; 2016b). Howlett and Pudney(1995)considersvariousproblemformulationswith exten-sivederivations. Obviouslytheoptimalcontrolstructureisalways thesamebutthealgebraicexpressionsoftheoptimalityconditions depend on the specific problem formulation,and these are used laterinthepaperinconnectiontothepseudospectralsolutions.

Inourproblemformulationwemodelatypicalmaximum trac-tionforceasfunctionofspeedusingapurecontrolconstraintand amixedstate-controlconstraint.Inaddition,staticspeedlimitsare

Fig. 1. Typical traction force-speed diagram with a constant and hyperbolic part.

formulatedaspurestateconstraints.Theproblemformulationthus includesall typesofalgebraicpathconstraints(state,control,and mixed constraints),which will be usedlater to analyse their im-pactonthepseudospectralsolutionmethod.

2.1. Energy-efficienttraincontrol

We consider the problemof finding theoptimal control fora trainrunbetweentwostopsinagivenscheduledtimeT suchthat thetotaltractionenergy(J [m2_/s2_{])isminimized.Thiscanbe}

for-mulatedasthefollowingconstrainedoptimalcontrolproblem:

MinimizeJ= sf

s0

u+

(

s

)

ds, (1)

subjecttotheconstraints

˙ t

(

s

)

=

_v

₍

1 s

)

(2) ˙

v

(

s

)

=u

(

s

)

− r

(

v

)

− g

(

s

)

v

(

s

)

(3) u+

(

s

)

v

(

s

)

≤ pmax (4) 0≤

v

(

s

)

≤

v

max

(

s

)

(5) − umin≤ u

(

s

)

≤ umax (6) t

(

s0

)

=0,t

(

sf

)

=T,

v

(

s0

)

=0,

v

(

sf

)

=0,s0=0,sf=S, (7)

where the independent variableis distance s[m], thestate vari-ables are time t [s] and speed

v

[m/s], t˙=dt/ds and

v

˙=d

v

/ds

denote the derivatives of the state variables with respect to the independent variable s, and the control variable u [m/s2_{] is the}

mass-specificappliedforceu

(

s

)

₌F

(

s

)

_/

(

ρ

m

)

_,i.e.,theappliedforce divided by totalmassincludinga rotating-massfactor

ρ

[-], con-sisting of a traction part u+

(

s

)

and a braking part u−

(

s

)

. The control is bounded between a maximum specific braking rate −umin=−Fmin/

(

ρ

m

)

andamaximumspecifictractionforceumax=

Fmax/

(

ρ

m

)

. Moreover, the mass-specific traction power p=u+

v

[m2_/s3_{]islimitedbyamaximummass-specificpowerp}

max.Hence,

themaximumtractionforce isafunctionofspeedconsistingofa constant and hyperbolic partwhich is illustrated inthe traction-speed diagram of Fig. 1. It follows that the control variable is boundedby u∈U

(

v

)

=[−umin,min

(

umax,pmax/

v

)

].Weassumeda

constantmaximumbrakingforceasthiswastheonlybrakingdata

available. Note that based on the definition of the control vari-able, traction and braking control cannot be used at the same time. We use the notationu+

(

s

)

=max

(

u

(

s

)

,0

)

≥ 0 and u−

(

s

)

= min

(

u

(

s

)

_,0

)

_{≤ 0}sothatu

(

s

)

₌u+

(

s

)

₊u−

(

s

)

.Theresistanceforces consistofamass-specifictrainresistancer

(

v

)

=R

(

v

)

/

(

ρ

m

)

[m/s2_]

and a mass-specific line resistance g

(

s

)

=G

(

s

)

/

(

ρ

m

)

[m/s2_{]. The}

trainresistance isdefinedby theDavisequation r

(

v

)

=r0+r1

v

+

r2

v

2, withnon-negative coefficients r0,r1≥ 0 andr2>0 (Davis,

1926). The line resistance g

(

s

)

is defined as the specific gravity forceduetotrackgradients.Itisassumedthattrackshave piece-wiseconstantgradients.Notethatonuphillslopesg

(

s

)

>0andon downhillslopes g

(

s

)

<0.Inaddition,lineresistance mayhavean additionalnonnegativetermduetocurvaturewhichalsodepends onthelocationofthecurves.Finally,thespeedisboundedabove byaspeedlimit

v

max

(

s

)

,whichisassumedpiecewiseconstant.

Next,we derive thenecessary optimality conditionsusing the PMP.Forthis,definetheHamiltonianH [m/s2_]as

H

(

t,

v

,

λ

1,

λ

2,u,s

)

=−u++

λ

1

v

+

λ

2

(

u− r

(

v

)

− g

(

s

))

v

, (8)

where

λ

1[m2/s3]and

λ

2[m/s]arethecostatevariables,whichare

alsofunctionsoftheindependentvariables.Notethatthe Hamil-tonian depends on distance s via the line resistance g

(

s

)

, which is apiecewise constant function of distance.Forsections of con-stantgradienttheHamiltonianisindependentonsandthus con-stantover distance,i.e.

∂

H/

∂

s=0,although theHamiltonianmay have jumps at the points where the gradient changes and thus become piecewise constant. Also jumps in the costate may lead to jumps inthe Hamiltonian.Totake intoaccount the additional (mixed) state andcontrol path constraints (4)–(6), we define the augmentedHamiltonianH¯ [m/s2_{]astheLagrangianofthe}

Hamil-tonianasfollows: ¯

H

(

t,

v

,

λ

1,

λ

2,

μ

,u,s

)

=H+

μ

1

(

umax− u

)

+

μ

2

(

u+umin

)

+

μ

3

(

pmax− u+

v

)

+

μ

4

(

v

max−

v

)

, (9)

where

μ

1[-],

μ

2[-],

μ

3[s/m]and

μ

4[1/s]aretheLagrange

multi-pliers,with

μ

_{i≥ 0}(i_{= 1,. . .,}4).Thecostates

λ

1and

λ

2satisfythe

differentialequations

λ

˙1

(

s

)

=−

∂

H¯/

∂

t and

λ

˙2

(

s

)

=−

∂

H¯/

∂v

,which

gives ˙

λ

1

(

s

)

=0 (10) ˙

λ

2

(

s

)

=

λ

1+

v

λ

2 r

(

v

)

+

λ

2

(

u− r

(

v

)

− g

(

s

))

v

2 +

μ

3u++

μ

4. (11)

Fromthefirstequationfollowsimmediatelythat

λ

1

(

s

)

≡

λ

1is

con-stant. Note that these dynamicequations (10),(11) together with thestateequations(2),(3)defineaBVPwithfourdifferential equa-tions in the states and costates and four fixed (begin and final) endpointsforthe statesgivenin(7). Therefore,theendpointsfor thecostatesarefree.

AccordingtothePMPtheoptimalcontrolmaximizesthe Hamil-tonian (Pontryaginet al., 1962). Therefore, theoptimal control is definedby ˆ u

(

s

)

=argmax u∈U H

(

ˆ t

(

s

)

,

v

ˆ

(

s

)

,ˆ

_λ

₁

₍

_s

₎

_,

_λ

ˆ₂

₍

_s

₎

_,u,s

₎

_{, (12)} where

(

tˆ_,

v

ˆ

)

and

(

λ

ˆ1,

λ

ˆ2

)

arethestateandcostatetrajectories

cor-respondingtotheoptimalcontroltrajectoryuˆ.Moreover,the max-imizedHamiltonianisindependentonuandthusconstantwhenit isalsoindependentons.InourcasetheHamiltoniandependson distancesthrough thepiecewise constantgradient g

(

s

)

,andalso thecontrol,stateandcostatetrajectoriesmaydependonsthrough thestate constraint(5)ifthespeed limit

v

max

(

s

)

isnot constant.

Hence, on intervalswith constant g

(

s

)

≡ gr and

v

max

(

s

)

≡

v

max,r,

thereexistsaconstant

ϕ

suchthat

NotethattheHamiltonian(8)ispiecewiselinearinthecontrolu,

andcanbesplitintotractionandbrakingpartsas

H

(

t,

v

,

λ

1,

λ

2,u,s

)

=



₍

_λ2 v − 1

)

u+λ1−λ2 r(v)−λ2g(s) v ifu≥ 0 λ2 vu+λ1−λ2 r(v)−λ2g(s) v ifu<0. (14) The optimal control u that maximizes the Hamiltonian thus de-pends onthe signofuandtherelative valueofspeed

v

andthe costate

λ

2,whichgivesfivecases:

1. If

λ

2>

v

thenumustbemaximal.

2. If

λ

2=

v

thenu∈[0,umax] isundetermined(1stsingular

solu-tion).

3. If0<

λ

2<

v

thenu=0.

4. If

λ

2=0 then [−umin,0] is undetermined(2nd singular

solu-tion).

5. If

λ

2<0thenumustbeminimal.

Tofindafull characterizationoftheoptimalcontroland asso-ciated state andcostate trajectories,we apply the Karush–Kuhn– Tucker(KKT)conditionsontheaugmentedHamiltonian(9).These conditionsconsistofthestationaryandcomplementaryconditions associated to H¯ (Bertsekas, 1999). The complementary slackness conditionsonthepathconstraintsare

μ

_{i≥ 0,}i₌1_,...,4_,and

μ

1

(

umax− u

)

=0,

μ

2

(

u+umin

)

=0,

μ

3

(

pmax− u+

v

)

=0,

μ

4

(

v

max−

v

)

=0. (15)

Note that

μ

2=0 if u≥ 0, and likewise

μ

1=

μ

3=0 if u<0.

The stationaryconditionstatesthat

∂

H¯/

∂

u=0,wherethepartial derivative ofthe augmented Hamiltoniancan be split in two ac-cording to the signof the control variable anddoes not exist at

u=0duetothediscontinuityofH atu=0.Thestationary condi-tioncanthenbederivedas



−1+λ2 v −

μ

1−

μ

3

v

=0 ifu>0 λ2 v +

μ

2=0 ifu<0, (16) anditisundefinedforu=0.

BasedonthePMPandKKTconditionswecannowcharacterize each ofthe five driving regimes dependingon the costate

λ

2 as

follows:

1. Maximum acceleration (MA). The first case is

λ

2>

v

with u

maximal.Inthiscase,the stationarycondition(16)gives

μ

1+

μ

3

v

=−1+

λ

2/

v

>0andthus

μ

1>0or

μ

3>0.Therefore,

ei-theru=umaxoru=pmax/

v

.Notethatu=umaxandu=pmax/

v

cannot holdboth at the same time except at one point

v

1=

pmax/umax since the hyperbolic function is decreasing in

v

.

Hence, u=umax

(

v

)

=min

(

umax,pmax/

v

)

, whichspecifies

max-imumtraction,and

μ

1 and

μ

3 canbe expressedintermsof

v

and

λ

2 using thestationarycondition(16)witheither

μ

3=0

or

μ

1=0.Inaddition,duringthismaximalaccelerationregime

the speed bound

v

=

v

max cannot be maintained except at a

single point, so

μ

4=0. This completely defines the costate

equation(11).

2. Cruisingbypartialtraction(CR1).Thesecondcaseisthe singu-larsolutionwith

λ

2=

v

andu∈[0,umax].Inthiscase,the

sta-tionary condition(16)reads

μ

1+

μ

3

v

=0and therefore

μ

1=

μ

3=0. Moreover, since

λ

2

(

s

)

=

v

(

s

)

on a nontrivial interval,

we must have

λ

˙2

(

s

)

=

v

˙

(

s

)

. From (3) and (11), we can then

derive

v

2

₍

μ

4+r

(

v

))

+

λ

1=0. Thisequation hasa unique

so-lution

v

c which can be proved as follows. If we view the

left-hand side as a function h

(

v

)

, and recalling that r

(

v

)

is a quadratic function of speed with in particular r

(

v

)

=

r1+2r2

v

≥ 0,r

(

v

)

=2r2>0, andthirdderivative r(3)

(

v

)

=0,

then h

(

v

)

=2

v

μ

4+2

v

r

(

v

)

+

v

2r

(

v

)

≥ 0 and h

(

v

)

=2

μ

4+

2r

(

v

)

+4

v

r

(

v

)

>0.Hence,h

(

v

)

isanincreasing convex func-tion,whichhasauniquerooth

(

v

c

)

=0for

v

c>0,ifany.Now

assumethat

v

c<

v

max andso

μ

4=0.Thenthecruising speed

mustsatisfy

v

2

cr

(

v

c

)

+

λ

1=0, (17)

whichhasaunique solution

v

c if

λ

1<0.Moreover,the

maxi-mizedHamiltonianvaluecondition(13)appliedtothiscasefor anintervalwithconstantg

(

s

)

≡ g thengives

v

cr

(

v

c

)

+r

(

v

c

)

+g+

ϕ

=0 (18)

with

ϕ

piecewiseconstantnegativeincongruencetog

(

s

)

,and thus the Hamiltonian is negative(H=

ϕ

<0). Ifon the other hand

v

(

s

)

=

v

maxonanontrivialintervalthen

μ

4isdetermined

as

μ

4=−

λ

1/

v

2max− r

(

v

max

)

whichgivescruisingatthe

maxi-malspeed

v

max asthe unique solution.Note that also inthis

case

λ

1<0musthold.Inconclusion, thiscasecorresponds to

aconstant cruisingspeed

v

(

s

)

=min

(

v

c,

v

max

)

with

v

c the

so-lutionto(17).Thisimpliesthat

v

˙

(

s

)

=0andthusu

(

s

)

=r

(

v

)

+

g

(

s

)

_,i.e., thetractionforce equals thetotal resistanceforce to maintainthe cruisingspeed. Thiscasethusreducesto finding asolution

v

cand

λ

1 to(17).

3. Coasting (CO). In the third case

λ

2∈

(

0,

v

)

and u=0. Since

u

(

s

)

=0 the complementary slackness conditions give

μ

1=

μ

2=

μ

3=0andalso

μ

4=0sincespeedchangeswithout

trac-tion effortdue to the resistance, and, therefore,will not stay atthespeedlimitfora nontrivialinterval.Thiscompletely de-finesthecostateequation(11).Notethattherarecaseofr

(

v

)

= −g

(

s

)

forsomeintervalreducestocase2or4withcruisingat zerotractionu

(

s

)

₌0.

4. Cruising by partial braking (CR2). The fourth case is the sec-ondsingularsolution with

λ

2=0 andu∈[−umin,0].The

sta-tionarycondition(16)nowgives

μ

2=0,andalso

μ

1=

μ

3=0

sinceu<0.Forthenontrivalcasethat

λ

2=0onsomeinterval,

we must have

λ

˙2=0. Thus from (11) we get

λ

1/

v

2+

μ

4=0.

If

μ

4=0 then also

λ

1=0 and by (8) also

ϕ

=0. However,

both

λ

1 and

ϕ

should be negative otherwise other regimes

would be prohibited, so

μ

4>0. Therefore,

v

(

s

)

=

v

max with

μ

4=−

λ

1/

v

max2 . Hence, partial braking is optimal only when

cruising at the speed limit with u

(

s

)

=r

(

v

max

)

+g

(

s

)

, which

can occur only on downward slopes with −umin− r

(

v

max

)

<

g

(

s

)

<−r

(

v

max

)

.

5. Maximumbraking(MB)Inthefinalcase

λ

2<0withuminimal.

The stationary condition(16) now gives

μ

2=−

λ

2/

v

>0 and,

therefore,u=−umin.Thespecialcasethat

v

=

v

maxovera

non-trivialinterval,i.e., g

(

s

)

_{=−umin− r}

(

v

max

)

ispartofthefourth

case(with full braking). Hence, the costate equation (11) ap-plieswith

μ

4=0,andalso

μ

3=0sinceu<0.

Theabove analysisthusleadstothefollowing optimalcontrol structure: ˆ u

(

s

)

=

⎧

⎪

⎪

⎪

⎪

⎨

⎪

⎪

⎪

⎪

⎩

umax

(

v

(

s

))

if

λ

2

(

s

)

>

v

(

s

)

(MA) r

(

v

(

s

))

₊g

(

s

)

_∈[0_,umax] if

λ

2

(

s

)

=

v

(

s

)

(CR1) 0 0<

λ

2

(

s

)

<

v

(

s

)

(CO)

r

(

v

max

(

s

))

+g

(

s

)

∈[−umin,0] if

λ

2

(

s

)

=0 (CR2)

−umin if

λ

2

(

s

)

<0 (MB).

(19) The optimal energy-efficient train control thus consists of maxi-mumtraction forceumax

(

v

)

formaximumacceleration(MA),

par-tialtraction to counterbalancetheresistance forcesand cruiseat anoptimalcruisingspeedmin

(

v

c,

v

max

)

(CR1),zerotractionu=0

forcoasting(CO) (i.e., rollingwithout usingthe engine orbrakes ofthe train), partialbraking tocruise atthe maximumspeed on a sufficiently downward slope (CR2), and full braking force for maximumbraking(MB). Thekey challenge isto findtheoptimal

switchingpointsbetweenthedifferentdrivingregimesandthe or-deroftheregimestodeterminetheoptimaldrivingstrategy.Fora simpleflattrack(zerogradient),nospeedlimit,andsufficient run-ningtimesupplement,theoptimaldrivingstrategyisgivenbythe sequence MA-CR1-CO-MBorMA-CO-MBifthespeedlimit cannot be reachedwithinthegiventime.Notethatrunningtime supple-ments are the extrarunning timesabove thetechnicalminimum running time included in the timetable to deal with variations in the running time orto recover fromsmall delays.Ifthe train runs punctual,they canbe used forenergy-efficienttrain driving (Scheepmaker&Goverde,2015).

Note that the dynamic equations (2) and (3) are not defined at zero speed. Inpractice, we use slightlypositive speeds atthe endpointstoavoidthedivideby zeroissueandcorrectthetotal distance and time accordingly. From Pontryagins Maximum Prin-ciple it is known that a train starts with maximum acceleration and ends with maximum braking, so this control can be simply extrapolated from/to zero. For our case studies the error is less than 25 cm or 1s atthe start and endpoint. In this paper,we ignored thiserror,whichiswellwithin thepracticalstopping ac-curacy. Likewise, for theminimum-time train control ofthe next subsection.

2.2. Minimum-timetraincontrol

TheaimoftheMTTCproblemistominimizethetotalrunning time ofa train, sothe train should arrive asearly aspossible at the next station atthe minimizedtime t

(

sf

)

. Thisobjective J [s]

canbeformulatedas

MinimizeJ=t

(

sf

)

, (20)

subjecttotheconstraints(2)–(6)andtheendpointconditions

t

(

s0

)

=0,

v

(

s0

)

=0,

v

(

sf

)

=0,s0=0,sf=S, (21)

whilethefinaltimet

(

sf

)

isfree.

The PMP analysis is similar to Section 2.1 and, therefore, we focusonthemainpointshereonly. TheHamiltonianH [s/m]and augmentedHamiltonianH¯ [s/m]aredefinedas

H

(

t,

v

,

λ

1,

λ

2,u,s

)

=

λ

1

v

+

λ

2

(

u− r

(

v

)

− g

(

s

))

v

, (22) ¯ H

(

t,

v

,

λ

1,

λ

2,

μ

,u,s

)

=H+

μ

1

(

umax− u

)

+

μ

2

(

u+umin

)

+

μ

3

(

pmax− u+

v

)

+

μ

4

(

v

max−

v

)

, (23)

where

λ

1 [-]and

λ

2 [s2/m]arethecostate functionsofthe

inde-pendent variable s, and

μ

1 [s3/m2],

μ

2 [s3/m2],

μ

3 [s4/m3] and

μ

4 [s2/m2],

μ

i≥ 0

(

i=1...,4

)

, are the Lagrange multipliers

as-sociated to the pathconstraints. Again the Hamiltonian is piece-wise constant withpossible jumps at the distances s where the gradient g

(

s

)

or speed limit

v

max

(

s

)

change, and(13) also holds

for thisHamiltonian.Thecostates

λ

1 and

λ

2 satisfy the

differen-tial equations

λ

˙1

(

s

)

=−

∂

H¯/

∂

t and

λ

˙2

(

s

)

=−

∂

H¯/

∂v

, which gives

˙

λ

1=0andthus

λ

1

(

s

)

≡

λ

1 asintheEETCcase,and

˙

λ

2

(

s

)

=

λ

1+

v

λ

2

r

(

v

)

+

λ

2

(

u− r

(

v

)

− g

(

s

))

v

2 +

μ

3u++

μ

4. (24)

However,differentfromtheEETCcase,thefinalconditionsforthe statesare notall fixednow,since thefinaltime t

(

sf

)

isfree.The

transversality conditions now specify a fixed value for the final timeoftheassociatedcostate

λ

1.Thiscanbefoundbyconsidering

theendpointLagrangianE¯[s]definedas ¯

E

(

t

(

sf

)

,

v

(

sf

)

,sf

)

=−t

(

sf

)

+

γ

1

v

(

sf

)

+

γ

2

(

sf− S

)

, (25)

where

γ

1and

γ

2areLagrangemultipliers,withthecomplementary

slacknessconditions

γ

1

v

(

sf

)

=0and

γ

2

(

sf− S

)

=0.The

transver-salityconditionsstatethat

λ

1

(

sf

)

=

∂

E¯/

∂

t

(

sf

)

=−1andtherefore

wenowobtain

λ

1≡ −1.Likewise,

λ

2

(

sf

)

=

∂

E¯/

∂v

(

sf

)

=

γ

1andfor

thevalueofthemaximizedHamiltonianattheendpoints_f,weget

H[s_f]=−

∂

E¯/

∂

s_f=−

γ

2,whichhoweverdoesnotgiveanynew

in-formation.

The Hamiltonian is again linear in the control u with coeffi-cient

λ

2/

v

. The optimal control u that maximizes this

Hamilto-nian therefore depends on the sign ofthe costate

λ

2, and must

be maximal for

λ

2>0,isundetermined for

λ

2=0,andis

mini-malfor

λ

2<0.TheKKTconditionsfortheaugmentedHamiltonian

(23)givethesamecomplementaryslacknessconditions(15)onthe pathconstraintswith

μ

i≥ 0,i=1,...,4,astheEETC case,while

thestationaryconditionswith

∂

H¯_/

∂

u₌0nowbecome



_λ2

v −

μ

1−

μ

3

v

=0 ifu>0

λ2

v +

μ

2=0 ifu<0,

(26) and is undefined for u₌0. Like the EETC case, we have

μ

2=0

foru≥ 0,and

μ

1=

μ

3=0foru<0.Wecannowcharacterizethe

threedrivingregimesdependingonthesignofthecostate

λ

2

us-ingtheresultsofthePMPandKKTconditions.

1. Maximum acceleration (MA). In the first case,

λ

2>0 and u>

0. Then (26) gives

μ

1+

μ

3

v

=λ_v2 >0 and thus

μ

1>0 or

μ

3>0. So either u=umax or u=pmax/

v

, and therefore u=

umax

(

v

)

=min

(

umax,pmax/

v

)

, which specifies that maximum

tractionneedstobeapplied.Inaddition,duringamaximal ac-celerationregime the speed bound

v

=

v

max cannot be

main-tainedexceptatasinglepoint,so

μ

4=0.Themultiplier

μ

3 is

eitherzeroorequalto

μ

3=

λ

2/

v

2using(26)with

μ

1=0.This

completelydefinesthecostateequation(24).

2. Cruising(CR).Thesecondcaseisthesingularsolutionwith

λ

2=

0andu_∈[_−umin,umax].With

λ

2=0,(26)specifies

μ

1=

μ

2=

μ

3=0forbothu>0andu<0.TheHamiltonian(22)now

be-comesH=

λ

1/

v

(

s

)

=−1/

v

(

s

)

.Moreover,

λ

˙2=0onanontrivial

interval andthen

μ

4=1/

v

2 using (24). Inparticular, this

im-pliesthat

μ

4>0andtherefore,

v

=

v

max.Hence,thiscase

cor-respondstocruisingatthemaximumspeed

v

max withcontrol

u

(

s

)

=r

(

v

max

)

+g

(

s

)

∈[−umin,umax]tomaintainthemaximum

speed. The control can be anything from full to partial trac-tionor braking dependingon the value of theresistance and gradient at

v

max. Moreover, the Hamiltonian value is fixed at

H_=−1/

v

max<0.

3. Maximumbraking(MB)Inthefinalcase,

λ

2<0andu<0.Now

(26)gives

μ

2=−

λ

2/

v

>0andthereforeu=−umin.Thespecial

casethat

v

=

v

max over anontrivialinterval ispartofthe2nd

case(with full braking).Hence,

μ

4=0, andalso

μ

3=0since

u<0, which completelydefinesthe costateequation (24). Fi-nally,from(22)followsthat H=

ϕ

<0sincebothcostates are smallerthanzero.

Theabove analysisthus leadstothe optimalcontrol structure

ˆ u

(

s

)

=

⎧

⎨

⎩

umax

(

v

(

s

))

if

λ

2

(

s

)

>0 (MA) r

(

v

max

(

s

))

+g

(

s

)

if

λ

2

(

s

)

=0 (CR) −umin if

λ

2

(

s

)

<0 (MB). (27)

So,thetractionandbrakingeffortsareashighaspossibletokeep thetrainatitsmaximumspeedsaslongaspossible,whichindeed generatestheminimalrunningtime.Notimeiswastedoncoasting inthiscase.

3. Pseudospectralmethod

In thissection, the train trajectory optimization isformulated asa multiple-phase optimalcontrol problem(Betts, 2010; Darby et al., 2011a) andthen discretized according to the Radau pseu-dospectralmethod(Garg,2011;Gargetal.,2009;Raoetal.,2010).

Fig. 2. Illustration of phases in the multiple-phase optimal control problem.

The Radaupseudospecralmethod isimplementedinthe MATLAB toolbox GPOPS, which we use as solver for the MTTC and EETC problems.

Amultiple-phaseoptimalcontrolproblemisonewherethe tra-jectoryconsistsofacollectionofphases.Simplyspeaking,aphase isanysegmentofthecompletetrajectory.Ingeneral,anyparticular phaseofanoptimalcontrolproblemhasacostfunctional,dynamic constraints, path constraints, and boundary conditions. However, the models used toquantitatively describe the trajectory maybe different in different phases of the trajectory. The complete tra-jectory is then obtained by properly linking adjacent phases via linkageconditions.Similarly,thetotalcostfunctionalisthesumof thecostfunctionalsovereachphase.Theoptimaltrajectoryisthen found by minimizingthe totalcost functionalsubjecttothe con-straints withineach phaseandthelinkageconstraints connecting adjacentphases.

The multiple-phase optimal control problemdescribed in the previous section is now transcribed to a discrete NLPvia an ex-tension of the single-phase Radau pseudospectral method. Pseu-dospectral methods are orthogonal collocation methods using globalpoynomials. Acollocation methodcan beused forthe nu-merical solution of ordinary differential equations and integral equations.Theideaistochooseafinite-dimensionalspaceof can-didate solutions (here polynomials) and a number of points in the domain (called collocation points), and to select that solu-tion which satisfies the givenequation at the collocation points. TheRadaupseudospectralmethodusestheLegendre–Gauss–Radau (LGR)pointsplusadditionallythefinalpoint.TheresultingNLPcan thenbesolvedbyoneofthemanywelldevelopednonlinear opti-mizationalgorithms.

3.1. Multiple-phaseoptimalcontrolproblemformulation

The optimal train control problems can be formulated as finite-horizon multiple-phase optimal control problems (Wang & Goverde,2016a).Therunningsection fromthe departurepointto thearrival pointisthendividedintoa finitenumberofsegments by the critical points of speed limits, gradients and curves, and eachofthesesegmentsiscalleda“phase”.Fig.2givesan illustra-tion,wheretherunningsectionisdividedintosixphases.Within each phase, the gradient and speed limit are constant, but their valuesmaybedifferentoverthevariousphases.

Assume that the running section is divided into phases [s(₀r),s(_fr)], with s(₀r) and s(_fr), s(₀r)<s(_fr), denoting the initial and

terminal location of phase r_∈

{

1_,...,R

}

. The train departs from theinitial points(₀1) ofphase1,andarrives attheterminal point

s(_fR) ofphase R.Within phaser,the state vector consistsof time andspeedx(r)

₍

_s

₎

_=[t(r)

₍

_s

₎

_,

_v

(r)

₍

_s

₎

_]_{,andthe controlvector isthe}

mass-specificappliedforceu(r)

₍

_s

₎

_{.Moreover,denotetheboundary}

pointsofaphaser asx(₀r)₌x(r)

₍

_s(r) 0

)

andx (r) f =x( r)

₍

_s(r) f

)

.

Themultiple-phaseoptimalcontrolproblemcannowbe formu-latedasminimizingthecostfunctional

MinimizeJ= R  r=1 E(r)

_x(r) 0 ,s( r) 0 ,x( r) f ,s( r) f

+  s(r) f s(r) 0 F(r)

x(r)

(

s

)

,u(r)

(

s

)

,s

ds, (28) subjecttothedynamicconstraints

˙

x(r)

₍

_s

₎

_=f(r)

_x(r)

₍

_s

₎

_,u(r)

₍

_s

₎

_,s

_{, r=1,...,R, (29)}

thepathconstraints

h(r)

_x(r)

₍

_s

₎

_,u(r)

₍

_s

₎

_,s

_{≥ 0, r=1,...,R, (30)}

thelinkageconditions

l(r)

(

x(_fr),s(_fr),x(₀r+1),s(₀r+1)

)

=0, r=1,...,R− 1, (31) andtheboundary(orendpoint)conditions

e

x₀(1),s(₀1),x(_fR),s(_fR)

=0. (32) Note that the path constraints are inequalities over the phases andthe linkage conditions are equality constraints on the phase boundaries.

FortheEETCproblemthecostfunctionalisgivenbyE=0and

F(r)

₍

_u(r)

₎

_=u+(r)_withu+(r)_=max

₍

_u(r)_,0

₎

_{,andfortheMTTC}

prob-lemtheseareE(r)

₍

_x(r)

0 ,x (r) f

)

=t (r) f − t (r)

0 andF(r)=0.Thedynamic

andpathconstraintsaregivenforbothproblemsby(2)–(6),i.e.,

f(r)

_x(r)

₍

_s

₎

_,u(r)

₍

_s

₎

_,s

₌



1

v

(r)

₍

_s

₎

, u(r)

₍

_s

₎

_{− r}(r)

₍

_v

(r)

₍

_s

₎₎

_{− g}(r)

v

(r)

₍

_s

₎



 and h(r)

_x(r)

₍

_s

₎

_,u(r)

₍

_s

₎

_,s

₌

⎡

⎢

⎢

⎢

⎢

⎣

pmax− u+(r)

(

s

)

v

(r)

(

s

)

v

(r)

₍

_s

₎

v

(r) max−

v

(r)

(

s

)

u(r)

₍

_s

₎

_+u min umax− u(r)

(

s

)

⎤

⎥

⎥

⎥

⎥

⎦

,