21st IFAC Symposium on Automatic Control in Aerospace 21st IFAC Symposium on Automatic Control in Aerospace August 27-30, 2019. Cranfield, UK Control 21st Symposium on 21st IFAC IFAC Symposium on Automatic Automatic Control in in Aerospace Aerospace August 27-30, 2019. Cranfield, UK Available at www.sciencedirect.com August 27-30, 2019. Cranfield, Cranfield, UK Control in online 21st IFAC Symposium on Automatic Aerospace August 27-30, 2019. UK August 27-30, 2019. Cranfield, UK
ScienceDirect
IFAC PapersOnLine 52-12 (2019) 274–279
Ascent Ascent Trajectory Trajectory Optimization Optimization of of Launch Launch Ascent Trajectory Optimization of Launch Vehicles with Air-Breathing Propulsion Ascent Trajectory Optimization of Launch Vehicles with Propulsion Vehicles with Air-Breathing Air-Breathing Propulsion Vehicles with Air-Breathing Propulsion∗∗∗ ∗ ∗∗ Vijith Mukundan ∗ Arnab Maity ∗∗ Shashi Ranjan Kumar ∗∗∗
Vijith Mukundan ∗∗ Arnab Maity Shashi Ranjan Kumar ∗∗ ∗∗∗ ∗∗ ∗∗∗∗ Vijith Maity Shashi U Rajeev ∗∗∗∗ Vijith Mukundan Mukundan ∗ Arnab Arnab Shashi Ranjan Ranjan Kumar Kumar ∗∗∗ UP PMaity Rajeev ∗∗ ∗∗∗∗ Vijith Mukundan Arnab Shashi Ranjan Kumar ∗∗∗ U Rajeev UP PMaity Rajeev ∗∗∗∗ ∗∗∗∗ ∗ U P Rajeev ∗ Vikram Sarabhai Space Centre, Indian Space Research Organization, ∗ Vikram Sarabhai Space Centre, Indian Space Research Organization, ∗ VikramTrivandrum Sarabhai Space Space Centre, Indian Space Research Organization, Organization, 695022, India (e-mail:
[email protected]) Vikram Sarabhai Centre, Indian Space Research 695022, India (e-mail:
[email protected]) ∗ ∗∗ Trivandrum Vikram Sarabhai of Space Centre, Indian Space Research Organization, Trivandrum 695022, India (e-mail:
[email protected]) Aerospace Engineering at Indian Institute of ∗∗ Department Trivandrum 695022, India (e-mail:
[email protected]) Department of Aerospace Engineering at Indian Institute of ∗∗ ∗∗ Department Trivandrum 695022, India (e-mail:
[email protected]) of Aerospace Engineering at Indian Institute of Technology, Mumbai 400076, India (e-mail:
[email protected]) Department of Aerospace Engineering at Indian Institute of Technology, Mumbai 400076, India (e-mail:
[email protected]) ∗∗ ∗∗∗ Department of Aerospace Engineering at Indian Institute of Technology, Mumbai 400076, India (e-mail:
[email protected]) Department of Aerospace Engineering at Indian Institute of ∗∗∗ Technology, Mumbai 400076, India (e-mail:
[email protected]) Department of Aerospace Engineering at Indian Institute ∗∗∗ ∗∗∗ Technology, Mumbai 400076, India (e-mail:
[email protected]) Department of Aerospace Aerospace Engineering at Indian Institute Institute of of Technology, Mumbai 400076, India (e-mail:
[email protected]) Department of Engineering at Indian of Technology, Mumbai 400076, India (e-mail: ∗∗∗∗∗∗∗ Department of Aerospace Engineering
[email protected]) Indian Institute of Technology, Mumbai 400076, India (e-mail:
[email protected]) Vikram Sarabhai Space Centre, Indian Space Research ∗∗∗∗ Technology, Mumbai 400076, India (e-mail:
[email protected]) Vikram Sarabhai Space Centre, Indian Space Research ∗∗∗∗ ∗∗∗∗ Vikram Technology, Mumbai 400076, India (e-mail:
[email protected]) Sarabhai Space Centre, Indian Space Research Organization, Trivandrum 695022, India Vikram Sarabhai Space Centre, 695022, Indian Space Organization, Trivandrum India Research ∗∗∗∗ Vikram Sarabhai Space Centre, 695022, Indian Space Organization, Trivandrum 695022, India Research Organization, Trivandrum India Organization, Trivandrum 695022, India Abstract: This This paper paper aims to to investigate optimization optimization of of a highly constrained constrained three-dimensional three-dimensional Abstract: Abstract: This paper aims aims to investigate investigate optimization of aaa highly highly constrained three-dimensional trajectory of a launch vehicle stage powered by air-breathing engine, with oblate gravitational Abstract: This paper aims to investigate optimization of highly constrained three-dimensional trajectory of a launch vehicle stage powered by air-breathing engine, with oblate gravitational Abstract: This paper aims to investigate optimization of aoptimize, highly constrained three-dimensional trajectory of a launch vehicle stage powered by air-breathing engine, with oblate gravitational field. The trajectories of such vehicles are difficult to owing to strong trajectory of a launch vehicle stage powered by air-breathing engine, with oblate gravitational field. The trajectories of such vehicles are difficult to optimize, owing to strong nonlinear nonlinear trajectory of a launch vehicle stage powered by air-breathing engine, with oblate gravitational field. The trajectories of vehicles are difficult to owing to strong nonlinear coupling between the thrust, thrust, aerodynamics and the vehicle’s states, and stringent constraints on field. Thebetween trajectories of such such vehicles areand difficult to optimize, optimize, owing to strong nonlinear coupling the aerodynamics the vehicle’s states, and stringent constraints on field. The trajectories of such vehicles are difficult to optimize, owing to strong nonlinear coupling between the thrust, aerodynamics and the vehicle’s states, and stringent constraints on structural load and aerodynamic heating. The objective of trajectory optimization is to impart coupling theaerodynamic thrust, aerodynamics andobjective the vehicle’s states, andoptimization stringent constraints on structuralbetween load and heating. The of trajectory is to impart coupling the thrust, aerodynamics and the vehicle’s states, andoptimization stringent constraints on structural load and aerodynamic heating.and The objective of trajectory trajectory is to to impart impart maximumbetween velocity at desired altitude direction, with minimum fuel using structural load and aerodynamic heating. The objective of is maximum velocity at aa desired altitude and direction, with minimumoptimization fuel consumption, consumption, using structural load and aerodynamic heating. The objective of trajectory optimization is to impart maximum velocity at a desired altitude and direction, with minimum fuel consumption, using angle of of attack attack andatengine engine throttle as control control inputs. with Path minimum constraints onconsumption, dynamic pressure, pressure, maximum velocity a desired altitude and direction, fuel using angle and throttle as inputs. Path constraints on dynamic maximum velocity atengine aheat desired altitude and direction, with using angle of attack and throttle as control inputs. Path constraints on dynamic pressure, bending moment, and flux as well well as control constraints need to to be befuel satisfied. Trajectory is angle of moment, attack and engine throttle as as control inputs. Path minimum constraints onconsumption, dynamic pressure, bending and heat flux as control constraints need satisfied. Trajectory is angle of moment, attack and engine throttle as as control inputs. Path constraints on dynamic pressure, bending moment, and heat flux flux as well as control constraints need to to be satisfied. satisfied. Trajectory is optimized using the approach of discretization of control variables, which converts this optimal bending and heat as well control constraints need be Trajectory is optimized using the approach of discretization of control variables, which converts this optimal bending moment, and fluxof well as control constraints needThis to beNLP satisfied. Trajectory is optimized using approach discretization of control variables, which converts this optimal control problem problem to aa heat nonlinear programming problem. problem solved optimized using the the ofas discretization control variables, which converts thisis control to approach nonlinear programming of(NLP) (NLP) problem. This NLP problem isoptimal solved optimized using the approach of discretization of control variables, which converts this optimal control problem to a nonlinear programming (NLP) problem. This NLP problem is solved using a metaheuristic method (harmony search optimization), and the sensitivity of the solution control problem to amethod nonlinear programming (NLP) problem. NLP problem solved using a metaheuristic (harmony search optimization), andThis the sensitivity of theissolution control problem to amethod nonlinear programming (NLP) problem. This NLPplanning problem issolution solved using a (harmony search optimization), and the of with respect respect to discretization discretization is also also analyzed. The versatility of trajectory trajectory method is using a metaheuristic metaheuristic method (harmony search The optimization), and the sensitivity sensitivity of the the solution with to is analyzed. versatility of planning method is using a metaheuristic method (harmony search optimization), and theobjectives sensitivity of the solution with respect tobydiscretization discretization is also also analyzed. analyzed. The versatility of trajectory trajectory planning method is demonstrated solving the problem with different initial conditions, and constraints. with respect to is The versatility of planning method is demonstrated by solving the problem with different initial conditions, objectives and constraints. with respect toby is also analyzed. Theinitial versatility of trajectory planning method is demonstrated solving with conditions, objectives and demonstrated bydiscretization solving the the problem problem with different different initial conditions, objectives and constraints. constraints. Copyright © 2019.byThe Authors. Published by with Elsevier Ltd. Allinitial rightsconditions, reserved. objectives and constraints. demonstrated solving the problem different Keywords: Trajectory, Trajectory, launch launch vehicle, vehicle, optimization, optimization, constraints, constraints, harmony harmony search, search, nonlinear nonlinear Keywords: Keywords: Trajectory, Trajectory, launch vehicle, vehicle, optimization, optimization, constraints, constraints, harmony harmony search, search, nonlinear nonlinear programming problem Keywords: launch programming problem Keywords: Trajectory, programming problem launch vehicle, optimization, constraints, harmony search, nonlinear programming problem programming problem 1. INTRODUCTION INTRODUCTION different approach approach of of trajectory trajectory design design than than that that using using 1. different 1. different approach of trajectory design than that using chemical rockets. Launch vehicles using chemical rockets 1. INTRODUCTION INTRODUCTION different approach of trajectory design than that using chemical rockets. Launch vehicles using chemical rockets 1. launch INTRODUCTION different approach of trajectory design thantime thatrockets using rockets. Launch vehicles using chemical Traditionally, space space vehicles have have been been propelled propelled chemical have trajectory trajectory designed to spend minimal within chemical rockets.designed Launch vehicles using chemical Traditionally, launch vehicles have to spend minimal time rockets within Traditionally, space launch vehicles have been propelled chemical rockets. Launch vehicles using chemical rockets have trajectory designed to spend minimal time within by rocket propulsion systems which use fuel and oxidizer the sensible atmosphere to reduce the aerodynamic loads Traditionally, space launch havefuel been havesensible trajectory designedtotoreduce spendthe minimal time within by rocket propulsion systemsvehicles which use andpropelled oxidizer the atmosphere aerodynamic loads Traditionally, space vehicle. launch vehicles havefuel been propelled by rocket propulsion systems which use and oxidizer have trajectory designed to spend minimal time within the sensible atmosphere to reduce the aerodynamic loads stored in the launch The necessity of carrying both experienced by the vehicle. On the other hand, launch by rocket propulsion systems which use fuel and oxidizer the sensible atmosphere to reduce theother aerodynamic loads stored in the launch vehicle. The necessity of carrying both experienced by the vehicle. On the hand, launch by propulsion systems use reason fuel and oxidizer stored in the launch vehicle. necessity of both the sensible to propulsion reduce theother aerodynamic loads by the On fuelrocket and oxidizer onboard, is The thewhich primary behind the experienced vehicle usingatmosphere air-breathing have tohand, dwelllaunch in the stored inoxidizer the launch vehicle. The necessity of carrying carrying both experienced by the vehicle. vehicle. On the the other hand, launch fuel and onboard, is the primary reason behind the vehicle using air-breathing propulsion have to dwell in the stored inoxidizer the launch vehicle. The necessity of leading carrying both fuel and onboard, is the primary reason behind the experienced by the vehicle. Onlonger the other hand, launch vehicle using air-breathing propulsion have to dwell in the huge mass and size of space launch vehicles, to low atmosphere for a considerably period in order to fuel and oxidizer onboard, is the primary reason behind the vehicle using air-breathing propulsion have to dwell in the huge mass and size of space launch vehicles, leading to low atmosphere for a considerably longer period in order to fuel and oxidizer onboard, islaunch the primary reason behind the atmosphere huge mass and space launch vehicles, leading to vehicle air-breathing propulsion to dwell in the for aa considerably longer period in order to payload fraction ofof space systems. Air-breathing generateusing thrust and optimally accelerate to high high velocities. huge mass and size size ofspace space launch vehicles, leading to low low atmosphere for considerably longer have period invelocities. order to payload fraction of launch systems. Air-breathing generate thrust and optimally accelerate to huge mass and sizeof space launchusing vehicles, leading tooxylow atmosphere payload fraction space launch systems. Air-breathing for and a considerably longer period invelocities. order to thrust accelerate to propulsion systems like scramjet atmospheric Trajectory design ofoptimally air-breathing launch vehicles pose payload fraction ofoflike space launch systems. Air-breathing generate thrust and optimally accelerate to high high velocities. propulsion systems scramjet using atmospheric oxy- generate Trajectory design of air-breathing launch vehicles pose payload fraction of space launch systems. Air-breathing propulsion systems like scramjet using atmospheric oxygenerate thrust and optimally accelerate to high velocities. Trajectory design of air-breathing launch vehicles pose gen as oxidizer, removes the need of carrying oxidizer the challenge of accelerating to high velocity within the propulsion systems like scramjet using atmospheric oxyTrajectory design of air-breathing launch vehicles pose gen as oxidizer, removes the need of carrying oxidizer the challenge of accelerating to high velocity within the propulsion systems like reduction scramjet using oxygen as removes the of carrying Trajectory design of air-breathing launch challenge of to velocity within the onboard. The resulting in gross gross lift-off oxidizer mass of the atmosphere, and simultaneously keep the vehicles thermal and gen as oxidizer, oxidizer, removes the need need of atmospheric carrying oxidizer the challenge of accelerating accelerating to high high velocity withinpose the onboard. The resulting reduction in lift-off mass of atmosphere, and simultaneously keep the thermal and gen as oxidizer, removes the need of carrying oxidizer onboard. The resulting reduction in gross lift-off mass of the challenge of accelerating to high velocity within the atmosphere, and simultaneously keep the thermal and the vehicle increases the payload capability of the launch aerodynamic loads on the vehicle under check, giving rise onboard. resulting in gross lift-off of aerodynamic atmosphere, and simultaneously keep the thermal and the vehicleThe increases thereduction payload capability of themass launch loads on the vehicle under check, giving rise onboard. The resulting reduction in gross lift-off mass the vehicle increases the payload capability of the launch atmosphere, and simultaneously keep the thermal and aerodynamic loads on the vehicle under check, giving rise vehicle. This can potentially lower the operational cost of to a highly constrained trajectory optimization problem. the vehicle increases the payload capability of the launch aerodynamic loads on the vehicle under check, giving rise vehicle. This can potentially lower the operational cost of to a highly constrained trajectory optimization problem. the vehicle the of thelaunches launch vehicle. This can lower the cost of aerodynamic loads on the under check, giving rise aa highly constrained trajectory optimization problem. space transportation and payload making frequent space Trajectory optimization ofvehicle air-breathing launch vehicle vehicle. Thisincreases can potentially potentially lowercapability the operational operational cost of to to highly optimization constrained trajectory optimization problem. space transportation and making frequent space launches Trajectory of air-breathing launch vehicle vehicle. This can potentially lower the operational cost of space transportation and making frequent space launches to a highly constrained trajectory optimization problem. Trajectory optimization of air-breathing launch vehicle economically more viable. Increase in payload fraction by also encounters encounters strong of nonlinear coupling between the space transportation and making space launches Trajectory optimization air-breathing launch vehicle economically more viable. Increasefrequent in payload fraction by also aa strong nonlinear coupling between the space and frequent space launches economically more viable. Increase in fraction by Trajectory optimization of air-breathing launch vehicle encounters aa by strong nonlinear coupling between the using transportation scramjet engine hasmaking been demonstrated demonstrated by Smart thrust generated engine, the aerodynamic aerodynamic forces on economically more viable. Increase in payload payload fraction by also also encounters strong nonlinear coupling between the using scramjet engine has been by Smart thrust generated by engine, the forces on economically more viable. Increase in payload fraction bya thrust using scramjet engine has been demonstrated by Smart also encounters astates strong nonlinear coupling between the generated by engine, the aerodynamic forces on and Tetlow (2009) in their trajectory calculations for vehicle, and the and control of vehicle. Such a couusing scramjet engine has been demonstrated by Smart thrust generated by engine, the aerodynamic forces on and Tetlow (2009) in their trajectory calculations for a vehicle, and the states and control of vehicle. Such a couusing scramjet engine has been demonstrated by Smart and Tetlow (2009) in trajectory calculations for thrust by engine, the aerodynamic forces on and states and control of vehicle. aa coumultistage rocket-scramjet-rocket system. The system system de-aa vehicle, pling is is generated absent launch vehicle using propulsion, and Tetlow (2009) in their their trajectory calculations fordevehicle, and the thein control of chemical vehicle. Such Such coumultistage rocket-scramjet-rocket system. The pling absent instates launchand vehicle using chemical propulsion, and Tetlow (2009) in their trajectory calculations for a multistage rocket-scramjet-rocket system. The system devehicle, and the states and control of vehicle. Such a coupling is absent in launch vehicle using chemical propulsion, signed to deliver a payload of 100 kg to 200 km equatorial whereisthe the thrust produced by engine engine is fairly independent independent multistage rocket-scramjet-rocket system. de- where pling absent in launch vehicle using is chemical propulsion, signed to deliver a payload of 100 kg to 200The kmsystem equatorial thrust produced by fairly multistage rocket-scramjet-rocket system. de- where signed to aa payload of kg to km equatorial pling absent launch vehicle using is chemical propulsion, thrust produced by fairly orbit showed showed payload fraction of 1.5%, which is better better of the theisthe states ofinvehicle. vehicle. signed to deliver deliver payload of 100 100of kg to 200 200The kmsystem equatorial where the thrust produced by engine engine is fairly independent independent orbit aa payload fraction 1.5%, which is of states of signed to deliver a payload of 100of kg1.5%, to 200 km equatorial orbit showed aarocket payload fraction which is better where the thrust produced by engine is fairly independent of the states of vehicle. than existing based systems with payload fraction orbit showed payload fraction of 1.5%, which is better of the states of vehicle. than existing rocket based systems with payload fraction The problem of trajectory optimization of air-breathing orbit showed arocket payload fraction of 1.5%, which better The than existing based systems with payload of theproblem states ofofvehicle. trajectory optimization of air-breathing of 0.9% 0.9% for similar similar payloads. Flaherty et al. al. (2010)isfraction include than existing rocket based systems with payload fraction of for payloads. Flaherty et (2010) include The problem trajectory optimization of launch vehiclesof has been addressed addressed by Corban Corban et al. al. (1991); (1991); The problem ofhas trajectory optimization of air-breathing air-breathing than existing rocket based systems with payload fraction of 0.9% for similar payloads. Flaherty et al. (2010) include launch vehicles been by et flexibility in mission and launch window as benefits of air of 0.9% forinsimilar payloads. Flaherty et al. flexibility mission and launch window as (2010) benefitsinclude of air launch The problem of trajectory optimization of air-breathing vehicles has been addressed by Corban et al. (1991); Lu (1993); Brinda et al. (2006); Murillo and Lu (2010); launch vehicles has been addressed by Corban et al. (1991); of 0.9% forin similar payloads. Flaherty et al. (2010) flexibility mission and launch window as benefits of (1993); Brinda et al. (2006); Murillo and Lu (2010); breathing propulsion based space launch systems. flexibility in mission and launch window as benefitsinclude of air air Lu breathing propulsion based space launch systems. launch vehicles has been addressed by Corban et al. (1991); Lu (1993); Brinda et al. (2006); Murillo and Lu (2010); Forbes-Spyratos et al. (2018). Corban et al. (1991) deLu (1993); Brinda et al. (2006); Murillo and Lu (2010); flexibility in mission and launch window as benefits of air breathing propulsion based space systems. et al. (2018). Corban et al. (1991) debreathing propulsion based space launch launch systems. system Forbes-Spyratos Lu (1993); Brinda et al.(2018). (2006); Murilloet and Lu (2010); et al. Corban al. (1991) deThe dependency dependency of the the air-breathing propulsion scribed onboard trajectory optimization of air-breathing Forbes-Spyratos et al. (2018). Corban et al. (1991) debreathing propulsion based space launch systems. system Forbes-Spyratos The of air-breathing propulsion scribed onboard trajectory optimization of air-breathing The dependency air-breathing system Forbes-Spyratos al. (2018). Corban etof (1991) detrajectory optimization air-breathing on the the atmosphereof tothe generate thrust, propulsion results in in aa system rather scribed The dependency ofto the air-breathing propulsion scribed onboard onboard et trajectory optimization ofal. air-breathing on atmosphere generate thrust, results rather The dependency of the air-breathing propulsion system on on the the atmosphere atmosphere to to generate generate thrust, thrust, results results in in aa rather rather scribed onboard trajectory optimization of air-breathing on the atmosphere generate thrust, results a rather Copyright © 2019 IFAC 2405-8963 Copyright © to 2019. The Authors. Published by in Elsevier Ltd. 274 All rights reserved. Copyright © 2019 IFAC 274 Copyright 2019 IFAC 274 Peer review© of International Federation of Automatic Copyright ©under 2019 responsibility IFAC 274 Control. 10.1016/j.ifacol.2019.11.255 Copyright © 2019 IFAC 274
2019 IFAC ACA August 27-30, 2019. Cranfield, UK
Vijith Mukundan et al. / IFAC PapersOnLine 52-12 (2019) 274–279
single stage to orbit vehicle using energy state approximation. Lu (1993) optimized the two dimensional trajectory of an air-breathing aerospace plane. An inverse dynamics approach was used wherein the optimal output was obtained and the corresponding control was derived by inverting the system dynamics. Trajectory optimization of air-breathing single stage to orbit launch vehicle was also reported by Brinda et al. (2006). The problem was converted to NLP, which was solved using sequential quadratic programming. Murillo and Lu (2010) generated fuel optimal ascent trajectories for air-breathing launch vehicle. The problem was posed as a two-point boundary value problem and was solved using finite difference method for finding the optimum angle of attack control to achieve the required target states with minimum fuel consumption. The optimum angle of attack profile was obtained without considering the dynamic pressure constraint. The engine throttle command was then separately used to satisfy only one path constraint, i.e. dynamic pressure. Forbes-Spyratos et al. (2018) described trajectory optimization of a three-stage launch vehicle, which has scramjet propulsion in second stage. Here too, only dynamic pressure constraint was considered without addressing the other path constraints. The study presented in this paper aims to explore optimization of a highly constrained three dimensional trajectory of an air-breathing launch vehicle. The objective is to deliver the upper stage with maximum velocity at a desired direction and altitude. The trajectory has to be fuel optimal and satisfy the path constraints of dynamic pressure, bending load and heat flux, as well as the constraints on the control variables, the angle of attack and the engine throttle. Trajectory optimization is attempted using a direct method. The control variables are discretized with respect to time. Values of control variables are initialized at the time nodes, selecting them from their respective allowable ranges. The control variables in between the nodes are interpolated using piecewise cubic Hermite interpolation. Using the control profile generated, the trajectory is propagated in time and the optimality of the trajectory is checked based on evaluation of an objective function. The proposed method converts the optimal control problem to NLP, in which the control variables at the nodes act as design variables deciding the optimality of the trajectory. The above NLP is solved using a gradient free, metaheuristic optimization algorithm, called harmony search optimization (HSO). This algorithm proposed by Geem et al. (2001) is inspired by the way the musicians achieve the best harmony. Similar to genetic algorithm, HSO also converges to the optimal solution by evolution of the candidates from an initial population. This paper is organized as follows. In Section II, the trajectory optimization problem is formulated. The trajectory optimization procedure is presented in Section III. In Section IV, the results of trajectory optimization are discussed. Conclusion and scope of future work are given in Section V. 2. PROBLEM FORMULATION In this section, the problem formulation will be presented. A three dimensional point mass dynamics of launch vehicle 275
275
is considered over rotating Earth with assumption of oblate gravitational field, as given by Tewari (2007). The states of the vehicle are given by radial distance (r), latitude (φ), longitude (λ), velocity (V ), flight path angle (γ) and heading angle (ψ). The dynamic model for such vehicle can be given by r˙ = V sin γ (1) V cos γ cos ψ φ˙ = (2) r V cos γ sin ψ (3) λ˙ = r cos φ T cos α − D V˙ = − gr sin γ + gφ cos γ cos ψ m 2 (4) + ωE r cos φ(sin γ cos φ − sin φ cos γ cos ψ) T sin α + L gr cos γ V cos µ − γ˙ = cos γ + r mV V gφ sin γ cos ψ + 2ωE sin ψ cos φ − V ω 2 r cos φ(cos ψ sin γ sin φ + cos γ cos φ) + E (5) V T sin α + L gφ sin ψ V sin µ − ψ˙ = tan φ sin ψ cos γ + r mV cos γ V cos γ 2 rωE sin ψ cos φ sin φ + 2ωE (sin φ − cos ψ cos φ tan γ) + V cos γ (6) where, m, T , L, and D denote mass, thrust, lift, and drag, respectively. µ is bank angle and α is angle of attack. Thrust is computed using air-breathing engine propulsion model, by Chuang and Morimoto (1997). This model gives thrust and specific impulse of the air-breathing engine as a function of Mach number (M ), angle of attack (α) and altitude (h). Thrust is throttleable. The coefficient of maximum thrust (CT max ) is given by the following relation 0.4736M 1.5 + 1.6947M −2 (M < 4) CT max (M, α) = 15(α + 5)0.25 M e C (M ≥ 4) M 1.15 (7) where, 2 35 M 0.08 MC = − (α + 5) − 0.6 (8) 200 M CT max is dimensionless and α is in degrees in the above relation. The thrust (T , in Newtons) and specific impulse (Isp , in seconds) of air-breathing engine are given by T = sQCT max Ae f (0 ≤ s ≤ 1) (9) 4500 − 10(h − 20) (M < 4) Isp (h, M ) = −245M + 5480 − 10(h − 20) (M ≥ 4) (10) where s, Ae , and Q are throttle (dimensionless, varies from 0 to 1), inlet area of air-breathing engine (in m2 ) and dynamic pressure (in Pascals), respectively. Inlet area is considered as 0.61575 m2 . f is a dimensionless factor for scaling the thrust given by the model to a suitable value for the vehicle considered in the current study and is taken as 40. h is altitude in kilometers. The mass flow rate (m, ˙ in kg/s) is given by T (11) m ˙ = Isp g0
2019 IFAC ACA 276 August 27-30, 2019. Cranfield, UK
Vijith Mukundan et al. / IFAC PapersOnLine 52-12 (2019) 274–279
Lift and drag are obtained from the aerodynamic coefficients, stored in table as a function of Mach number and angle of attack, for a generic launch vehicle. The Earth’s rotation speed is assume to be ωE = 7.29211 × 10−5 rad/s. The terms gr and gφ are components of oblate Earth gravity and are given by 2 µE 3 RE 2 (3sin φ − 1) (12) g r = 2 1 − J2 r 2 r 2 3µE RE sin φ cos φ (13) g φ = − 2 J2 r r
where, µE = 3.98600448 × 1014 m3 /s2 is Earth’s gravitation constant, J2 = 1.0826368×10−3 m3 /s2 is second zonal harmonic and RE = 6378137 m is the equatorial radius of Earth. The state dynamics can be represented as x˙ = f (x, α, s) (14) where the state vector, x, is given by x = [r φ λ V γ ψ]T (15) and the angle of attack (α) and the throttle (s) are used as the control variables.
The objective of trajectory optimization is to reach maximum velocity at the desired altitude, flight path angle, and heading angle while minimizing fuel consumption, without violating the path constraints. Fuel consumption minimization can be achieved by maximizing the final mass. The desired terminal conditions for the above optimization problem are h(tf ) = h∗ = 30 km (16) γ(tf ) = γ ∗ = 25◦ (17) (18) ψ(tf ) = ψ ∗ = 0◦ where the superscript ‘∗’ denotes the desired value of the variables. These constraints act as terminal constraints for the current problem. The path constraints can be written as αmin ≤ α ≤ αmax (19) (20) smin ≤ s ≤ smax Q ≤ Qmax (21) (22) |Qα| ≤ Qαmax q˙ ≤ q˙max (23) ◦ ◦ where αmin and αmax are -5 and 5 , and smin and smax are 0.1 and 0.8, respectively. The dynamic pressure, Q, with the atmospheric density denoted by ρ, is given by 1 (24) Q = ρV 2 2 The maximum allowed value of Q is 75 kPa. The term Qα, the product of dynamic pressure and angle of attack, is an indicator of the bending moment experienced by the vehicle, and it is limited to 7000 Pa-rad. The stagnation point heat flux (q, ˙ in W/cm2 ) due to aerodynamic heating is given by Scott et al. (1985) as, 3.05 ρ V q˙ = 18300 (25) rn 10000 where, rn denotes the nose radius and is taken as 0.5 m. The limit on q˙ is 30 W/cm2 . 276
Path constraints of dynamic pressure, bending moment and heat flux constraints are included as the part of objective using exterior penalty functions. To satisfy the terminal constraint on altitude, the state dynamics are propagated till altitude h∗ is reached. Terminal constraints on flight path angle (γ) and heading angle (ψ) are considered in the objective function. The objective function to be minimized is given as J = − wV V (tf ) − wm m(tf ) + wγ (γ(tf ) − γ ∗ )2 + wψ (ψ(tf ) − ψ ∗ )2 + w1 G1 + w2 G2 + w3 G3 (26) where G1 = max(0, g1 ), G2 = max(0, g2 ), G3 = max(0, g3 ) and the terms g1 , g2 , g3 are defined as below (27) g1 (x(t), α(t), s(t)) = Q − Qmax ≤ 0 (28) g2 (x(t), α(t), s(t)) = |Qα| − Qαmax ≤ 0 g3 (x(t), α(t), s(t)) = q˙ − q˙max ≤ 0 (29) In equation (26), velocity is in m/s, mass is in kg, flight path angle and heading angle are in degrees, dynamic pressure is in Pascals, bending moment indicator is in Pa-rad and heat flux is in W/cm2 . Note that the terms wV , wm , wγ and wψ , are the weights associated with the objective and the desired terminal constraints, and is given by (30) wV = 0.1, wm = 0.01, wγ = 100, wψ = 10 and terms w1 , w2 and w3 are weights associated with path constraints, and is given by (31) wk = (0.5i)2 (i = iteration number) for k = 1, 2 and 3 . 3. TRAJECTORY OPTIMIZATION PROCEDURE The trajectory optimization aims to obtain the control profiles, that is, time variation of angle of attack and throttle, such that the launch vehicle reaches the desired terminal conditions, with maximum velocity and minimum fuel consumption, satisfying the path constraints. This is an optimal control problem, which is converted to a nonlinear programming problem by discretizing the control variables. 3.1 Discretization of control variable Control variables, angle of attack (α) and throttle (s), are discretized at a sequence of the time nodes. (32) tα = {0, tα2 , tα3 , . . . , tαi , . . . , tαN } (33) ts = {0, ts2 , ts3 , . . . , tsi , . . . , tsN } The selection of time nodes is an engineering judgement. Increasing the number of nodes improves the accuracy of the solution but also results in high computational load. Thus there is a trade-off between the two. To begin with, the time nodes can be considered at a larger interval and then finer time intervals can be taken till the solutions converge. Values of control variables, angle of attack and throttle, at the time nodes are selected between their allowable limits, such that the constraints on the control variables, given by equations (19) and (20), are also satisfied. The angle of attack and throttle at the nodes are as given below. (34) αc = {α1 , α2 , α3 , . . . , αN } sc = {s1 , s2 , s3 , . . . , sN } (35)
2019 IFAC ACA August 27-30, 2019. Cranfield, UK
Vijith Mukundan et al. / IFAC PapersOnLine 52-12 (2019) 274–279
Piecewise cubic Hermite interpolation is used to compute the values of the control variables between two nodes. This generates the control variable profiles, as shown in Figure 1. Using these profiles, the trajectory of launch vehicle is
277
Step 2. Improvise a new harmony: A new solution candidate [x1 , x2 , x3 , ...xN ] from the harmony memory is generated, as follows: (1) Define Harmony Memory Consideration Rate, which is the probability of selecting a component (xj ) of the new solution candidate from the current harmony memory. If the component is not selected from the current memory, then it is generated randomly. (2) If xj is selected from the harmony memory, it is obtained from the j th dimension of a randomly selected solution candidate from the harmony memory. (3) Define Pitch Adjustment Rate (PAR) and Bandwidth (BW). PAR is probability that xj , obtained from harmony memory, will be mutated. BW defines the extent of mutation. Step 3. Update the harmony memory: If the new solution candidate obtained in Step 2 gives a better objective function than that of the worst solution candidate of harmony memory, then it is used to replace that worst solution candidate in the harmony memory. Otherwise, the new candidate is discarded.
Fig. 1. Angle of Attack and throttle Profile propagated in time. The nature of the trajectory depends on the control variable profiles, decided by the selected values of control variables at the nodes. Hence, the angle of attack and throttle at the nodes are the design variables which decide the optimality of the trajectory, evaluated by the objective function given by equation (26). Obtaining the design variables to minimize the objective function (26) is a nonlinear programming problem, which is solved using a optimization technique, harmony search optimization, described below. 3.2 Harmony Search Optimization Harmony Search Optimization is a metaheuristic optimization algorithm, proposed by Geem et al. (2001). This algorithm imitates the process of musician improving to achieve the best harmony aesthetically. The aesthetics of harmony is determined by the correct sound of each instrument, the same way as in optimization a objective function is determined by the set of design variables. HSO algorithm consists of four steps which are given as follows. Step 1. Initialize harmony memory (HM) array: The harmony memory is initialized by generating a certain number of randomly generated solution candidates of the optimization problem. For a N-dimensional optimization problem, harmony memory array is given as x11 x12 x13 · · · x1N x21 x22 x23 · · · x1N HM = . (36) . .. .. .. .. .. . . . S S S S xHM xHM xHM · · · xHM 1 2 3 N
where [xi1 , xi2 , xi3 , ...xiN ] is a solution candidate and HMS denotes harmony memory size.
277
Step 4. Check the termination criterion: If a predefined termination criterion (like maximum number of iterations or minimum improvement in objective function) is satisfied, the process is stopped, or else Step 2 and Step 3 are repeated. HSO algorithm is implemented to obtain the angle of attack and throttle at the nodes, so as to minimize the objective function. The versatility of the trajectory optimization procedure is demonstrated by solving the problem with different objectives and constraints. 4. RESULTS AND DISCUSSION This section describes the results of trajectory optimization and is followed by the pertinent discussion. As discussed in previous section, the objective of the airbreathing stage is to impart maximum velocity to the upper stage in the desired direction and altitude with minimum fuel consumption, simultaneously satisfying the path constraints. The initial and final conditions considered for the problem are as given in Table 1. HSO algorithm described in Section 3 is implemented with following parameters of algorithm, HMS = 50, HMCR = 0.75, PAR = 0.5. These were selected based on a parametric study to ensure that the design space is searched extensively and premature convergence to a local optimum is prevented. In order to assess the sensitivity and convergence of solution, with respect to discretization, the above described methodology is implemented by discretizing control variables at different uniform time intervals (∆t) of 20 s, 15 s, 10 s and 5 s. Figure 2 shows the variation of control variables, angle of attack and throttle, for all the four cases. It is seen that the angle of attack and throttle variations have converged for the cases ∆t = 10 s and ∆t = 5 s. Same can be stated about the variation of altitude and velocity, shown in Figure 3. Table 2 shows the maximum final velocity and mass for the four cases. The improvement in final velocity and mass achieved by
2019 IFAC ACA 278 August 27-30, 2019. Cranfield, UK
Vijith Mukundan et al. / IFAC PapersOnLine 52-12 (2019) 274–279
Table 1. Initial and Final States
35
Initial condition
Final condition
Parameter
at t = 0
at t = tf
r
RE + 11 km
RE + 30 km
φ
0◦
free
λ
0◦
free
V
614.3 m/s
To Maximize
γ
35◦
25◦
ψ
0◦
0◦
γ
35◦
25◦
m
33 T
To Maximize
Altitude(km)
State/
∆t = 20 s ∆t = 15 s ∆t = 10 s ∆t = 5 s
30 25 20 15 10 0
5
10
15
20 time(s)
trajectory optimization decreases as we go for discretization at finer time intervals. This means that the solution becomes less sensitive to discretization time interval below a certain limit. Hence, for a given problem making finer discretizations after attaining convergence is not required and thus heavy computational loads can be prevented. In the problem under study, the convergence has been attained for ∆t = 10 s and the same has been chosen for further investigation of the problem.
Velocity(km/s)
3
25
30
35
40
35
40
∆t = 20 s ∆t = 15 s ∆t = 10 s ∆t = 5 s
2
1
0 0
5
10
15
20 time(s)
25
30
Fig. 3. Altitude and velocity variation
Angle of attack (deg)
10 ∆t = 20 s ∆t = 15 s ∆t = 10 s ∆t = 5s
8 6 4
Fig. 4. Path constraints: dynamic pressure, bending moment indicator and heat flux
2 0 0
5
10
15
20 time(s)
25
30
35
40 30 Altitude (km)
1
Throttle
0.8 0.6
∆t = 20 s ∆t = 15 s ∆t = 10 s ∆t = 5 s
0.4 0.2 0 0
20
heat flux boundary 10 0 0
5
10
15
20 time(s)
25
30
35
m(tf ) (in T)
20
1919.05
33.80
15
2045.17
38.19
10
2048.54
38.53
5
2049.23
38.54
500
region violating constraints
1000 1500 Velocity (m/s)
2000
Fig. 5. Altitude- velocity trajectory plot and maximum dynamic pressure boundary
Table 2. Final velocity and mass for cases of ∆t V (tf ) (in m/s)
dynamic pressure boundary
40
Fig. 2. Angle of attack and throttle variation
Case (∆t in s)
region within constraints
Fig. 4 shows the variation of path constraints, dynamic pressure, bending moment indicator and heat flux. It can be seen that all the path constraints are satisfied throughout the trajectory. The optimal trajectory is flown close to the maximum dynamic pressure limit. This is 278
because the thrust of air-breathing engine increases with dynamic pressure, as seen in Eq. (9). In altitude-velocity trajectory plot shown in Fig. 5, it can be seen that the trajectories of all the four cases are seen to be grazing the maximum dynamic pressure boundary. The vehicle is able to maintain the trajectory around the boundary without crossing the boundary and violating the constraints, by using the throttle control shown in Fig. 2. Engine is throttled down when the vehicle has to be flown around the dynamic pressure boundary. When the heat flux constraint is about to get activated, the engine is throttled up to increase the thrust. To study the effect of initial altitude on the trajectory, the trajectory optimization was solved with different initial altitudes. The objective was modified to meet a required
2019 IFAC ACA August 27-30, 2019. Cranfield, UK
Vijith Mukundan et al. / IFAC PapersOnLine 52-12 (2019) 274–279
variables, angle of attack and engine throttle, with respect to time. The optimal control problem was converted to a nonlinear programming problem, which is solved by a metaheuristic optimization algorithm, Harmony Search. The sensitivity of solution with respect to discretization was studied and it was found that, the solution converges as the discretization becomes finer. The solution becomes less sensitive to the discretization time interval below a certain limit. This is an important result based on which excessively fine discretization and increasing computational load can be avoided. The trajectory was found to be grazing the dynamic pressure and heat flux constraints, in order to meet the objectives of maximizing velocity and minimizing fuel consumption. On studying the effect of initial conditions on the trajectory, it was observed that for the air-breathing launch vehicle, meeting the required velocity becomes infeasible if the trajectory is initiated at high altitude. Future work will explore development of schemes suitable for real-time trajectory optimization and onboard implementation. These will cater the requirement of closed loop guidance to meet the desired mission in the presence of vehicle parameter and environmental uncertainties.
35 Altitude(km)
30 25 h0 = 11 km h0 = 12 km h0 = 13 km h0 = 14 km
20 15 10 0
5
10
15
20 time(s)
25
30
35
40
Fig. 6. Altitude variation for different initial altitudes
Angle of attack (deg)
6 4 h0 = 11 km h0 = 12 km h0 = 13 km h0 = 14 km
2 0 -2 0
5
10
15
20 time(s)
25
30
35
40
REFERENCES
1
Throttle
0.8 h0 = 11 km h0 = 12 km h0 = 13 km h0 = 14 km
0.6 0.4 0.2 0 0
5
10
15
20 time(s)
279
25
30
35
40
Fig. 7. Angle of attack and throttle variation for different initial altitudes final velocity, instead of maximizing velocity. In these cases, the objective function is expressed as below J = − wV (V (tf ) − V ∗ )2 − wm m(tf ) + wγ (γ(tf ) − γ ∗ )2 + wψ (ψ(tf ) − ψ ∗ )2 + w1 G1 + w2 G2 + w3 G3 (37) ∗ where, V is assumed to be equal to 2048.54 m/s. Fig. 6 shows the trajectories reaching the desired final altitude. Fig. 7 shows the angle of attack and throttle profile, for these trajectories. It is observed that as initial altitude increases, the angle of attack initially becomes more negative and vehicle flying a shallow trajectory to increase the thrust generated so that the required final velocity is met. Fig. 7 also shows that throttle increases for higher initial altitudes. The throttle reached its upper limit, for initial altitude of 14 km, beyond which it becomes infeasible for the vehicle to meet the final velocity requirement at the desired altitude and flight path angle. This defines the upper bound of altitude at which the air-breathing stage can start and meet the required terminal conditions. 5. CONCLUSIONS AND FUTURE WORK This paper presented the ascent trajectory optimization of air-breathing launch vehicle by discretizing the control 279
Brinda, V., Dasgupta, S., and Lal, M. (2006). Trajectory optimization and guidance of an air breathing hypersonic vehicle. In 14th AIAA/AHI Space Planes and Hypersonic Systems and Technologies Conference, 7997. Chuang, C.H. and Morimoto, H. (1997). Periodic optimal cruise for a hypersonic vehicle with constraints. Journal of Spacecraft and Rockets, 34(2), 165–171. Corban, J., Calise, A., and Flandro, G. (1991). Rapid near-optimal aerospace plane trajectory generation and guidance. Journal of guidance, control, and dynamics, 14(6), 1181–1190. Flaherty, K.W., Andrews, K.M., and Liston, G.W. (2010). Operability benefits of airbreathing hypersonic propulsion for flexible access to space. Journal of Spacecraft and Rockets, 47(2), 280–287. Forbes-Spyratos, S.O., Kearney, M.P., Smart, M.K., and Jahn, I.H. (2018). Trajectory design of a rocket– scramjet–rocket multistage launch system. Journal of Spacecraft and Rockets, 56(1), 53–67. Geem, Z.W., Kim, J.H., and Loganathan, G.V. (2001). A new heuristic optimization algorithm: harmony search. simulation, 76(2), 60–68. Lu, P. (1993). Inverse dynamics approach to trajectory optimization for an aerospace plane. Journal of Guidance, Control, and Dynamics, 16(4), 726–732. Murillo, O. and Lu, P. (2010). Fast ascent trajectory optimization for hypersonic air-breathing vehicles. In AIAA Guidance, Navigation, and Control Conference, 8173. Scott, C.D., Ried, R.C., Maraia, R., Li, C., and Derry, S. (1985). An aotv aeroheating and thermal protection study. Thermal design of aeroassisted orbital transfer vehicles, 96, 198–229. Smart, M.K. and Tetlow, M.R. (2009). Orbital delivery of small payloads using hypersonic airbreathing propulsion. Journal of Spacecraft and Rockets, 46(1), 117–125. Tewari, A. (2007). Atmospheric and space flight dynamics. Springer.