Orbit-injection Strategy to a Multistage Launch Vehicle Based on Six-Degree of Freedom and Coupled Axes Attitude Control System

ABSTRACT

Understanding of various aerodynamic factors involved in flight trajectories is fundamental to design launch vehicles. First and foremost, computer simulation is an efficient way of predicting its behavior in the movement across the atmosphere. Considering that the available Brazilian version of Analysis, Simulation and Trajectory Optimization Software for Space Applications (Astos) does not simulate a controlled vehicle in six degrees of freedom (DoF), the aim of this article is to complement the Astos outcomes, particularly evaluating the trajectory of a controlled launch vehicle from liftoff to orbit injection, considering the model of rigid body dynamics with a six DoF. This approach carried out with an in-house developed simulator called Scott that simulated a multistage launcher with three flight configurations. In the Scott computer program, a launcher was modeled with differential equations in six DoF, coupled axes attitude control system, and aerodynamic coefficients that changed as a function of Mach number. These features improved the results generated by Astos software for the same configurations and the same initial conditions. Additionally, the results provided by Scott were close to actual vehicle in terms of attitude change and Mach number reached.

Keywords
Launch vehicle; Attitude control; Vehicle trajectory

INTRODUCTION

Over the last few years, there has been a great evolution in space transport systems such as composite structural materials, more efficient engines, more accurate control systems, as well as reusable vehicles. Furthermore, a growing number of private companies, namely SpaceX, Boeing, and Lockheed Martin, have developed their own launchers and vehicles in order to control this market. Consequently, costs have decreased and vehicles have become more reliable too. In this scenery, not only have embedded systems evolved, but also computer simulation has become more efficient in previewing flight performance and also preventing catastrophic events.

The classical space transportation system simulation commences with mission definition that includes vehicle configuration, payload, launch site, and orbit definition. In general, the first simulation takes place in a three-degrees of freedom (DoF) computer program that needs propulsion features, aerodynamic coefficients, and mass properties. As the simulation result, a feasible trajectory and a vehicle reference attitude are achieved, being enough to preliminary trajectory analysis. On the other hand, a three-DoF model does not represent an actual vehicle motion; for this reason, a six-DoF simulation could provide more reliable outcomes.

In order to get the aerodynamic coefficients, the Missile Datcom computer program could be used for determining the aerodynamic characteristics of missile-types bodies of revolution as launch vehicles. Additionally, there is a great deal of computer codes for space vehicle design and analysis, including Satellite Tool Kit (STK) developed by Analytical Graphics Inc. (AGI) and Analysis, Simulation and Trajectory Optimization Software for Space Applications (Astos) made by the German company Astos Solutions GmbH. Furthermore, Astos computer program provides this reference based on optimized vehicle and flight path.

Taking into account that the available Brazilian version of Astos does not simulate a six-DoF vehicle with attitude control, the aim of this article is to complement the Astos outcomes, particularly evaluating the trajectory of a controlled launch vehicle from liftoff to orbit injection, considering the model of rigid body dynamics with a six-DoF. This approach carried out with an in-house developed computer program called Scott that simulated a multistage launcher with three flight configurations. In the Scott computer program a launcher was modeled with differential equations in six DoF, coupled axes attitude control system, and aerodynamic coefficients that changed as a function of Mach number. These features improved the results generated by Astos software for the same configurations and the same initial conditions. Additionally, the results provided by Scott were close to actual vehicle in terms of attitude change and Mach number reached. To solve this gap, a sophisticated approach is necessary, including a six-DoF vehicle dynamic and coupled axes control system. In other words, translational and rotational dynamics, with control in terms of pitch and yaw attitudes, and an interaction among rotation axes. Undoubtedly, these features could improve the reference outcomes.

For this reason, a simulation environment was created with Earth gravity and atmosphere models to get gravity and aerodynamic forces as well as atmosphere density and temperature at each point of the flight trajectory. Furthermore, the simulator was assembled in a six-DoF with orientation described by quaternions, a control system of coupling axes, time-variable mass properties, as well as variable aerodynamic properties in function of the Mach number.

Then a controlled three-stage launch vehicle simulation was carried out in Scott and the results compared with Astos. Particularly, the control system followed the attitude and trajectory references with significative attitude differences, which are more coherent with an actual flight. In spite of the qualities presented in the Astos computer program, it is crystal clear that rotational dynamics and others features present in the Scott results in a more realistic simulation, and its outcomes could be useful in a refined launch vehicle trajectory reference.

METHODOLOGY

The shape and gravitation of the Earth

The launch vehicle simulation, from surface to orbit injection, requires a high precision shape of the Earth, gravity, and navigation models. Taking into account that Earth is not a perfect sphere, but a slightly oblate spheroid with an irregular surface, as a result an accurate shape modelling is necessary. For that, the World Geodetic System 1984 (WGS-84) model was adopted, as shown by Stevens and Lewis (2003). Firstly, the authors define the geoid as an equipotential surface of the Earth’s gravity field that coincides with the undisturbed mean sea level. Due to the uneven mass geoid distribution, a waving surface is generated. In addition to that, the equatorial radius of the Earth is roughly 21 km wider than the polar radius, which fits an ellipse. Therefore, the rotation about minor axis produces an Earth’s oblate spheroidal model. Both geoid and oblate spheroidal models, as well as geocentric attitude (vertical local) and geodetic altitude, are illustrated in Fig. 1.

Figure 1
The geoid, spheroid, and height definitions.

Secondly, the angle λ_gd between the spheroid normal and the equator plane is called the geodetic latitude, and h_gd represents the geodetic altitude, as presented in Fig. 2.

Figure 2
WGS-84 oblate spheroidal model of the Earth.

Besides, the angle λ_gc means the geocentric latitude, r_gc is the distance from Earth center to P, and the longitude value is the same to geocentric (lon_gc) or geodetic (lon_gd) representation.

Thirdly, the gravity model is presented with respect to gradients of the potential function (V), which currently uses Earth Gravity Model 86 (EGM86) coefficients. The potential function is shown in Eq. 1:

$$V(r_{gc},{\lambda }_{gc})=\frac{GM}{r_{gc}}\left[1-0.5J_2{\left(\frac{a}{r_{gc}}\right)}^2(3 {\sin }^{2 }{\lambda }_{gc}-1)\right]$$ (1)

where GM is gravitational constant given by GM = 3986004.418 x 10⁸ m³/s², r is the geocentric radius in meters, J₂ is the second Earth spherical harmonic, given by J₂ = 1.0826267 x 10^-3, a is semimajor axes of earth ellipse given by a = 6378137.0 m, and λ_gc that is geocentric latitude.

The Eq. 2 presents the gradients of the potential function calculated in geocentric coordinates and transformed into Earth Centered Earth Fixed reference system (ECEF):

$$G=\frac{-GM}{r_{gc}^2}\left[\begin{array}{l}\left[1+1.5J_2{\left(\frac{a}{r_{gc}}\right)}^2(1-5{\sin }^2 {\lambda }_{gc})\right]\frac{r_x}{r_{gc}}\\ \left[1+1.5J_2{\left(\frac{a}{r_{gc}}\right)}^2(1-5{\sin }^2 {\lambda }_{gc})\right]\frac{r_y}{r_{gc}}\\ \left[1+1.5{/}_2{\left(\frac{a}{r_{gc}}\right)}^2(3-5{\sin }^2 {\lambda }_{gc})\right]\frac{r_z}{r_{gc}}\end{array}\right]$$ (2)

where r_x, r_y and r_z are the r_gc components into ECEF and the r_gc module. Additionally, the geocentric latitude is evaluated as λ_gc = asin r_z/r_gc.

The gravitational acceleration is determined by gravitational attraction minus centripetal acceleration, which in terms of vector equation is given by Eq. 3:

$$g=G-{\omega }_{e/i}x\left({\omega }_{e/i}x{r}_{gc}\right)$$ (3)

where ω_e/i is the fixed Earth rotation rate vector.

Finally, the navigation model presented by Stevens and Lewis (2003) allows calculating geodetic longitude as Eq. 4,

$${\mathrm{lon}}_{gd}=\mathrm{atan} \left(\frac{{\mathrm{r}}_{\mathrm{y}}}{{\mathrm{r}}_{\mathrm{x}}}\right)$$ (4)

while geodetic latitude (λ_gd) and height (h_gd) are evaluated following the iterative algorithm, from Eqs. 5-8:

$${\lambda }_{gd}=\mathrm{atan} \frac{r_z}{\sqrt{(r_x^2+r_y^2)}\left(1-\frac{N{e}^2}{N+h_{gd}}\right)}$$ (5)

$$N=\frac{a}{{(1-e^2{\sin }^2 {\lambda }_{gd})}^{0.5}}$$ (6)

$$N+h_{gd}=\frac{\sqrt{(r_x^2+r_y^2)}}{\cos {\lambda }_{gd}}$$ (7)

$$h_{gd}=(N+h_{gd})-N$$ (8)

where e is the Earth eccentricity given by e = 0.08181919.

To sum up, the gravity forces are evaluated in geocentric coordinates (λ_gc and lon_gc) and consider the distance from Earth center (r_gc) and, for navigation purposes, the geodetic coordinates (λ_gd and lon_gd) and heigh (h_gd) are suitable.

Equations of motion

The dynamic equations describe the six-DoF vehicle movement through the space. These equations are based on initial and time update states, as well as mass center inertial position, quaternion orientation, mass center inertial velocity, angular velocity body, and total mass body.

In order to completely describe the six-DoF vehicle movement, a set of 13 differential equations are necessary, namely: four rotational kinematic equations, three translational kinematic equations, three rotational dynamic equations, and three translational dynamic equations. The set of equations governing the motion of the vehicle is (Eqs. 9-12)

$$\dot{q}=\frac{1}{2}{q}_{b/i}{\omega }_{b/i}^b$$ (9)

$${\dot{p}}_{CM/O}=L_{BI}{V}_{CM/e}+{\omega }_{e/i}{p}_{CM/O}$$ (10)

$${\dot{\omega }}_{b/t}={\left(J^b\right)}^{-1}(M_{A,T}-{\mathrm{\Omega }}_{b/l}{J}^b{\omega }_{b/i})$$ (11)

$${\dot{v}}_{CM/e}=\frac{1}{m}{F}_{A,T}-({\omega }_{b/i}+{\omega }_{e/l})v_{CM/e}+L_{BI}g$$ (12)

The previous matrix movement equations and the variable mass represent the state vector X (Eq. 13):

$$X^T=\left[{\left(p_{CM/o}\right)}^T{\left(q_{b/i}\right)}^T{\left(v_{CM/e}\right)}^T{\left({\omega }_{b/i}\right)}^{T}m\right]$$ (13)

Aerodynamic models

According to Anderson (2010) and Sirious (2004), the aerodynamic forces over the “axisymmetric vehicle” are defined as (Eqs.14-16):

$$F_{AX}=-C_A{p}_{\text{dyn}}{S}_{\text{ref}}$$ (14)

$$F_{AY}=-C_{Na}\beta p_{dyn}{S}_{\text{ref}}$$ (15)

$$F_{AZ}=-C_{N\alpha }\alpha p_{dyn}{S}_{ref}$$ (16)

where F_AX is axial force, F_AY is lateral force, and F_AZ is normal force. In terms of aerodynamic derivatives, C_A is the axial force coefficient and C_Nα is the normal force coefficient in pitch, P_dyn is dynamic pressure and S_ref is the reference area. Finally, α is the angle of attack and β is the sideslip angle.

Additionally, the aerodynamic roll, pitch, and yaw moments are given as follows (Eqs. 17-19):

$$M_{AX}=\left(C_{i\delta }\delta -C_{lp}\frac{p{L}_{ref}}{2V}\right)pdin{S}_{ref}{L}_{ref}-y_{cm}{F}_{AZ}+z_{cm}{F}_{AY}$$ (17)

$$M_{AY}=\left[-C_{N\alpha }\left(\frac{x_{cp}-x_{cm}}{L_{ref}}\right)\alpha -C_{mq}q\frac{L_{ref}}{2V}\right]pdin{S}_{ref}{L}_{ref}$$ (18)

$$M_{Az}=\left[-C_{N\alpha }\left(\frac{x_{cp}-x_{cm}}{L_{ref}}\right)\beta -C_{mq}r\frac{L_{ref}}{2V}\right]pdin{S}_{ref}{L}_{ref}$$ (19)

where C_lδ is the roll-driving derivative, C_lp is the roll-damping derivative, C_mq is the pitching moment coefficient, L_ref is the reference length, x_cm, y_cm, and z_cm are center of gravity position, x_cp is the center of pressure position, and δ is the fin deflection.

Aerodynamics effects over the vehicle structure

There are two main aerodynamic effects over the vehicle structure: vibration level peak and aerodynamic loads. The former happens in the transonic flight region, and it is defined by Anderson (2010) as a range of relative speed that produces subsonic and supersonic flow over the body. In general, this region begins at Mach 0.8, with the first shock waves, and ends at Mach 1.2, when the flow is fully supersonic. The flight into this region causes severe unsteadiness and vibration over the vehicle as result of asymmetry and instability in the flow around the body.

The latter aspect is dynamic pressure, which is explained by Anderson (2010) and Gunstom (2004) as a pressure at body stagnation point due to relative movement inside the fluid. In particular, the aerodynamic load value is proportional to dynamic pressure. In other words, in the point of maximum dynamic pressure, the vehicle is undergoing the maximum mechanical stress and, for this reason, it becomes a significant factor in the vehicle design. According to Anderson (2010), the equation of dynamic pressure, for incompressible flow, is given by Eq. 20,

$$p_{dyn}=\frac{1}{2}\rho V^2$$ (20)

where P_dyn is the dynamic pressure, ρ is air density, and V is the relative speed between fluid and the body.

Considering the mechanical stress over the vehicle, there are three points of interest: dynamic pressure maximum value (MaxQ), exactly time which MaxQ happens (MaxQ time), and total time above 90% of MaxQ (HighQ region).

Control system model

It is suggested by Kadan (2019) a way to control a vehicle in terms of pitch and yaw, regardless of the roll attitude, by resolving attitude errors in order to generate signals for control system. The concept is comparing a general-desired-stored trajectory, in terms of θ_d, ψ_d, and ϕ_d to θ_a, ψ_a, and ϕ_a that are provided by inertial sensors. This difference results in θ_ε, ψ_ε, and ϕ_ε, which could be applied to thrust vector actuators (TVA) or reaction control system (RCS) control equations. Considering the elementary rotation matrices (Eqs. 21-23):

$$L_x\left(\varphi \right)=\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \left(\varphi \right)& \sin \left(\varphi \right)\\ 0& -\sin \left(\varphi \right)& \cos \left(\varphi \right)\end{array}\right]$$ (21)

$$L_y\left(\theta \right)=\left[\begin{array}{ccc}\cos \left(\theta \right)& 0& -\sin \left(\theta \right)\\ 0& 1& 0\\ \sin \left(\theta \right)& 0& \cos \left(\theta \right)\end{array}\right]$$ (22)

$$L_z\left(\psi \right)=\left[\begin{array}{ccc}\cos \left(\psi \right)& \sin \left(\psi \right)& 0\\ -\sin \left(\psi \right)& \cos \left(\psi \right)& 0\\ 0& 0& 1\end{array}\right]$$ (23)

where θ, ψ, and ϕ could be written in terms of actual (a), desirable (d), or error (e) attitudes.

A desirable body orientation with respect to actual body orientation, in 231-rotation sequence, is given by Eq. 24:

$$L_{bd}=L_x\left({\varphi }_d\right)L_z\left({\psi }_d\right)L_y\left({\theta }_d\right)L_y{\left({\theta }_a\right)}^T{L}_z{\left({\psi }_a\right)}^T{L}_x{\left({\varphi }_a\right)}^T{L}_{ba}$$ (24)

Additionally, for small angles, the angular rotation about the actual body axes required to align them with

$$L_{bd}=\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& {\varphi }_e\\ 0& -{\varphi }_e& 1\end{array}\right]\left[\begin{array}{ccc}1& {\psi }_e& 0\\ -{\psi }_e& 1& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{ccc}1& 0& -{\theta }_e\\ 0& 1& 0\\ {\theta }_e& 0& 1\end{array}\right]L_{ba}$$ (25)

From Eqs. 24 and 25, as well as considering small angles, hence Eqs. 26-28:

$${\theta }_e=({\theta }_d-{\theta }_a)\cos \left({\psi }_d\right)\cos \left({\varphi }_a\right)+({\psi }_d-{\psi }_a)\sin \left({\varphi }_a\right)$$ (26)

$${\psi }_e=-({\theta }_d-{\theta }_a)\cos \left({\psi }_d\right)\sin \left({\varphi }_a\right)+({\psi }_d-{\psi }_a)\cos \left({\varphi }_a\right)$$ (27)

$${\varphi }_e=({\varphi }_d-{\varphi }_a)+({\theta }_d-{\theta }_a)\sin \left({\psi }_d\right)$$ (28)

Taking into account a six-DoF and coupled axes dynamics, the Eqs. 26-28 provide the error data in the body reference that are useful in the proportional-integral-derivative (PID) control system. Considering the thrust vector control (TVC) in terms of pitch and yaw direction, the displacement angle nozzle is modeled by the Eqs. 29 and 30:

$${\beta }_{\text{pitch}}=k_p^{pitch}{\theta }_e-k_d^{pitch}q+k_i^{pitch}\int {\theta }_{e}dt$$ (29)

$${\beta }_{\text{yaw}}=k_p^{\text{yaw}}{\psi }_e-k_d^{\text{yaw}}r+k_i^{\text{yaw}}\int {\psi }_{e}dt$$ (30)

Especially in ballistic flight phase, the same approach can be used to point the vehicle into correct attitude. Eqs. 31 and 32 provide a final effect over the vehicle likewise the RCS.

$$M_y=k_p^{pitch}{\theta }_e-k_d^{pitch}q+k_i^{pitch}\int {\theta }_{e}dt$$ (31)

$$M_z=k_p^{\text{yaw}}{\psi }_e-k_d^{\text{yaw}}r+k_i^{\text{yaw}}\int {\psi }_{e}dt$$ (32)

In the previous equations, k_p, k_d, and k_i are the proportional, derivative, and integral gains, respectively.

Gain schedule calculation

According to Campos (2004) arguments, an analytical method (based on control theory provided by Nise (2012) and Ogata (1997) could be used to calculate a set of gains that are able to follow the reference, as well as stabilize the system during the flight. By hypothesis, the author considered that the roll angle and velocity are small, which implies that the pitch and yaw axes are uncoupled. This assumption is a particular case of control system showed early. Then, the Eqs. 33-35 evaluated the pitch-plane gains:

$$k_p=\frac{{\omega }^2+2\xi \omega p_0+M_{\alpha }}{-M_{\beta z}}$$ (33)

$$k_d=\frac{2\xi \omega +p_0}{-M_{\beta z}}$$ (34)

$$k_i=\frac{{\omega }^2{p}_0}{-M_{\beta z}}$$ (35)

where ω = 3.1464 rad/s, ξ = 0.8645, and P₀ = 0.4519 rad/s.

Furthermore, the parameters M_α and M_βz are given by Eqs. 36 and 37:

$$M_{\alpha }=\frac{C_{N\alpha }{p}_{dyn}{S}_{ref}{I}_a}{I_{yy}}$$ (36)

$$M_{\beta z}=\frac{T{\mathrm{I}}_{\mathrm{c}}}{I_{yy}}$$ (37)

where C_Nα is the normal force coefficient in pitch, P_dyn is dynamic pressure, S_ref is the reference area, I_α is the difference between center of gravity and center of pressure positions, T is the thrust module, I_c is the distance from center of gravity and the thrust force application point, and finally I_yy is the moment of inertia about lateral body axis.

Aerodynamic coefficients calculation

All the aerodynamic coefficients were calculated with Missile Datcom computer program (Blake 1998). Some of first flight configuration coefficients – C_A, C_Nα, C_mq and X_cp – are presented below. As shown in the Fig. 3d, by hypothesis, the center of gravity position is constant at 10 m from rocket tip. As a result, the vehicle is aerodynamic-stable until Mach 3.5 and is unstable above this value. Regardless of this aerodynamic instability, the control system must be able to keep it on the reference trajectory and attitude.

Figure 3
Aerodynamic coefficients for the first flight configuration: (a) axial force coefficient, (b) normal force coefficient in pitch, (c) pitching moment coefficient, and (d) center of pressure and center of gravity position.

Scott development

As stated in Aguiar (2017), Scott computer program was developed to simulate a multistage launch vehicle, considering translational and rotational dynamics, with coupling-axes-attitude control system. Especially for this work, Scott was configurated to simulate a three-stage launch vehicle, with a fairing and a payload. Undoubtedly, the ability to follow an attitude and trajectory reference is essential for this simulation. In addition to that, the use of a coupling-axes-attitude control improves the reliability of the system and makes a control-roll bay unnecessary under certain rotation values that does not affect the inertial unit accuracy.

Case study

Astos simulation

According to Astos Solutions website, the newer version of computer programmer is the 9.23. Furthermore, the current version and the version available for this work (e.g., version 7.1.0) are encompassed reduced and full three DoF aerodynamics models, as well as full six DoF aerodynamic model. Although the computer program supports six DoF model, just reduced and full three DoF are controlled. As described in Astos 9 Model Reference (Astos Solutions GmbH 2019) and Astos 7 Model Library (Astos Solutions GmbH 2013), full three DoF is related to full three DoF Euler angles that is included yaw, pitch, and roll angles in a vehicle without rotational symmetry. On the other hand, a reduced three DoF could be applied in a rotational symmetry vehicle, which is controlled in terms of yaw and pitch attitudes. In this case, an artificial roll angle is computed in order to guide the air-path axis system. For these reasons, the reduced three DoF aerodynamic and Euler angles will be used in this work.

Vehicle description

The vehicle showed in Fig. 4 is composed by three stages, and a fairing and a payload. Each stage has a solid propellant engine, equipped with a TVC and a set of TVA that moves each nozzle in pitch and yaw directions. Additionally, on the top of the rocket, a fairing gives an aerodynamic shape to the body and protects a payload during the atmospheric flight. The vehicle´s longitudinal dimension is 20 m and its diameter is about 1.5 m.

Figure 4
Complete vehicle: three stages, payload, and fairing.

During the powered flight, the TVC deflects each nozzle and provides an attitude control in order to follow a preload reference. Furthermore, a complex coupled axes control system was employed in the launcher. Regardless of the roll position, this control system holds the rocket in the desired attitude and trajectory. On the other hand, during ballistic flight phase, a RCS maneuvers the vehicle into the correct attitude.

Following a traditional approach, the vehicle was designed in three-stages and a payload tandem arrangement, which generated three flight configurations:

• Configuration 1 has three stages, besides a fairing and a payload. The first stage is active for 80 seconds and a PID control system keeps the vehicle in the reference trajectory during the propelled flight. After the burnout time, the first stage is jettisoned from the whole vehicle.

• Configuration 2 is composed of second and third stages, as well as a fairing and a payload. As in the first configuration, a PID control system holds the rocket in the reference trajectory for 73 seconds until the separation of second stage and jettisoning of fairing.

• Configuration 3 has the third stage and a payload. The aim of this configuration is to put the payload in the correct altitude and attitude to inject it in the target orbit. Thus, it is necessary a ballistic flight to reach the desired altitude and a pointing maneuver to grasp the correct attitude before starting the third stage engine.

Table 1 summarizes the flight configurations and the payload. Particularly, the table provides data about structure and propellant mass, burn time, length of each configuration, as well as center of gravity and pressure with respect to the tip of the vehicle.

Table 1
Flight configurations and payload specifications.

vehicle properties	Configuration 1	Configuration 2	Configuration 3	Payload
Burn time (sec)	80	73	100	-
Total mass (kg)	30,600	14,100	3,100	450
Active stage propellant mass (kg)	15,000	10,000	2,000	-
Active stage structural mass (kg)	1,500	1,000	650 (structural mass + fairing)	450
Total length (m)	20	10	5	1.12
Center of pressure position (m)	-13 (min) -6 (max)	-3.04 (min) -0.91 (max)	-0.69 (min) -1.93 (max)	-
Center of gravity position (m)	- 10	- 5	- 4	- 0.10

• Source: Elaborated by the authors.

In addition to vehicle description, there are eight flight phases that must be described. Firstly, by hypothesis, the vehicle body axes and the flight plane are aligned up before the launch time. Therefore, an initial roll maneuver is not necessary. Secondly, after the engine ignition, the vehicle performs a vertical rise until clear the launch pad. Thirdly, to steer the rocket’s longitudinal axis toward downrange direction, the pitch-over one and two maneuvers were performed, followed by first and second stages burn, under gravity turn law control and constant yaw. Consequently, the vehicle attitude tended to local horizontal. Finally, tipper actuators guarantee the third stage precision start and the payload injection into the planned orbit. Table 2 summarizes the trajectory sequence of events.

Table 2
Trajectory flight phases.

Trajectory phase	Initial time (sec)	Final time (sec)	Control law
Vertical climb	0.0	4.3	Constant pitch Constant yaw
Pitch-over 1	4.3	7.7	Linear pitch Constant yaw
Pitch-over 2	7.7	11.2	Constant pitch Constant yaw
First stage burn	11.2	80.0	Gravity turn Constant yaw
Second stage burn	80.0	153	Gravity turn Linear pitch Constant yaw
Third stage pointing maneuver	153	230	Linear pitch Linear yaw
Third stage injection payload	230	330	Constant pitch Constant yaw
Third stage deorbit	330	350	None

• Source: Elaborated by the authors.

RESULTS AND DISCUSSION

Taking into account the proposed method, this section illustrates the results achieved for the three flight configurations. These outcomes were focused on control system performance, structural loads analyses, altitude reached, as well as jettisoned stages impact point prediction. In some of them, there is a comparison between Scott results and Astos reference; in others, Scott simulator provided values regardless of reference, such as dynamic pressure and ground impact point.

First configuration results and discussions

The ability to follow a reference trajectory is a key factor to inject the satellite into a precise orbit. In order to track that reference, a PID control system and a set of gains are crucial for a real launch operation. Especially for the first flight configuration, a dedicated PID control system and a proper set of gains were added to hold the vehicle close to desired pitch and yaw attitudes desired. Therefore, a couple of pitch attitude as function of flight time were generated and are shown in Figs. 5 and 6.

Figure 5
Pitch attitude reached by PID control in comparison to pitch reference provided by Astos during first configuration flight.

Figure 6
Pitch attitude reached by PID control in comparison to pitch reference provided by Astos from liftoff to 20 seconds.

Firstly, from lift-off to 7.4 seconds, the vehicle carried out a vertical climb in order to avoid launch tower strike, as well as prevent launch pad damage due to high temperature gases. Secondly, the vehicle performed a gradually pitch-over maneuver that has been ranging its attitude to -11 degrees. Finally, a gravity turn maneuver shifted the pitch attitude to -52.8 degrees. Taking into account that Astos reference did not consider the rotational dynamics, it is expected a time delay during pitch-over maneuver, as shown below. Particularly, Scott simulator provided a smoothing change attitude when compared with Astos.

As seen in Fig. 7, the controller tracked properly the yaw reference despite the fact that there is a peak of lateral deviation in the order of tenth of degrees. This peak took place at the vertical climb phase (between lift-off and 7.4 seconds) and suggests an interaction between pitching and yawing attitude change. Considering the tenth of degrees in lateral deviation, as well as the low vehicle velocity, it is reasonable to conclude that there was little mechanical stress over the structure in this part of the flight.

Figure 7
Yaw attitude reached by PID control in comparison to yaw reference provided by Astos.

Since the control system was working properly and the vehicle overpassed the launch tower, the following aspect to be considered is the stress over the structure caused by, mainly, vibration and aerodynamic loads. According to Vos and Farokhi (2015), the vibration is presented all over the flight and its peak value happens during the transonic phase. As seen in Fig. 8, the transition from Mach 0.8 to Mach 1.2 (transonic phase) happened between 27.3 and 38.5 seconds, still in a dense atmosphere. Hence, the maximum value of vibration took place in this time range, as result the vehicle structure and control system must be designed to endure this effect. Particularly, this event could be considered as first design point.

Figure 8
Mach number as function of flight time. The graph presents transonic region, since Mach 0.8 until Mach 1.2, which happens between 27.3 and 38.5 seconds. At 80 seconds, the first stage separation in Mach 4.0 occurs.

In addition to transonic phase, the dynamic pressure could be considered the second design point on account of mechanical stress over the structure. Especially, above 90% of the simulated maximum dynamic pressure, the vehicle comes into the highest mechanical stress region (HighQ region). For this vehicle, the HighQ region began at 8.75 km and finished at 15.45 km. As seen in Fig. 9, the maximum dynamic pressure was 60.1κN/m², at altitude of 12.5 km, while it was travelling at Mach 2.1. Despite of shape and performance differences, the maximum dynamic pressure takes place near to 12 km for the most of launcher vehicles as Space Shuttle, that was close to 11 km, and SpaceX Falcon, that happened at 12.5 km.

Figure 9
Dynamic pressure as a function of altitude. Between 8.75 and 15.45 km, the vehicle achieves more than 90% of MaxQ. In terms of range, HighQ extends from 54.1 until 60.1 κN/m². This is the part of flight that had major mechanical stress over the structure.

Considering the vibration peak caused by transonic flight and the maximum level of mechanical stress due to dynamic pressure, a time separation between these events is mandatory in order to avoid the overlap effects. As shown in Fig. 10, the transonic region happened from 27.5 to 38 seconds, and the HighQ region took place between 43.8 and 59.2 seconds after the liftoff. In fact, this time separation supports a future structural project, and optimization of payload mass or altitude orbit injection.

Figure 10
Dynamic pressure and Mach number as a function of time. It is possible to see that transonic region had occurred between 27.5 and 38 seconds, while HighQ region happened from 43.8 to 59.2 seconds.

Another point to be considered in the first configuration is the vehicle’s altitude as function of flight time. As shown in Fig. 11, from the lift off until burn out time, the vehicle ascended to 29 km above the surface. At 80 seconds, three events happened, namely: engine burnt out, first stage separation, and second stage ignition. This sequence is called “cold separation”.

Figure 11
Altitude in function of flight time. At 80 seconds, the separation of the first stage happened near to 29 km above the surface.

The simulation of the first stage had been continuing during a non-controlled ballistic flight up to impact with the ground. As seen in Fig. 12, the jettisoned stage reached the apogee point close to 56 km above the surface and 93 km from the launch pad. In addition to that, a first impact point approximation took place at 28° azimuth and 185 km from launch point.

Concluding this section, it may be argued that the vehicle’s control system tracked, satisfactorily, the attitude reference. Additionally, the six DoF attitude control implemented in Scott simulator results in a feasible attitude response when compared to Astos. Another point explored was the qualitative stress over the structure in terms of vibration and dynamic pressure. The results showed that the peak of vibration was associated with transonic flight, and the maximum aerodynamic loads were related to HighQ region. Furthermore, the graphic of Mach number and aerodynamic loads as function of flight time showed that vibration peak and MaxQ happened in different flight times and altitude. Finally, the flight trajectory was propagated and generated information about apogee altitude, as well as the first stage ground impact estimation.

Figure 12
Altitude as function of horizontal range. The first stage separation happened at 24.7 km from the launch pad, it reached 56 km above the surface (apogee point) and the impact point was 185 km from the launch pad.

Second configuration results and discussions

As presented before, the second flight configuration is composed of second and third stages, as well as a fairing and a payload. This configuration performed a powered and controlled flight during 73 seconds until second stage burnout, which happened at 153 seconds of flight time. In addition to that, a PID control system managed the vehicle attitude in order to track a preload reference, as done in configuration 1.

Considering that there were not abrupt changes in pitch and yaw, the controller kept vehicle’s attitude very close to reference, as expected. The Fig. 13 shows that the pitch attitude coincided with reference, during the completely second configuration flight. The same is true for yaw attitude reference that changed just only a tenth of degree, hence the control system held the lateral reference without significant deviation, as seen in Fig. 14. These outcomes imply that the PID control system and its gains set were suitable to conduct the second configuration in a correct trajectory.

Figure 13
Comparison between reference pitch attitude provided by Astos and reached in configuration 2.

Figure 14
Comparison between reference yaw attitude provided by Astos and the one obtained in configuration 2.

Another point to be considered is the Mach number as function of flight time. As shown in Fig. 15, from lift off to second stage separation in 137 seconds, the Mach number reached is close to that predicted by Astos. However, between 146 and 153 seconds we have a drop of two Mach number units. This drop happened at the same time the velocity increased, as shown in Fig. 16.

Figure 15
Mach number as function of flight time.

Taking into account that the Mach number was evaluated as a ratio between vehicle velocity and local sound speed, the first question is whether the velocity provided by Astos and reached by Scott are close. Both velocities are illustrated in Fig. 16 and their outcomes are very similar. Hence, there is not enough difference to produce that Mach number deviation.

Since the Astos and Scott velocities had similar behavior, the second question is about how the Mach number evaluation is biased by local velocity and local sound speed deviation. The former is dependent for the thrust and the latter is straightly proportional to square root of the local temperature. As seen in Fig. 17, the temperature at 153 seconds (281 km) is almost 1.67 times at 146 seconds (100 km). On the other hand, Fig. 16 shows that the velocity at 153 seconds is 1.14 times more than the velocity at 146 seconds. Hence, the Mach number at 153 seconds should be almost 80% of Mach number at 146 seconds.

Figure 16
Comparison between velocity provided by Astos and the velocity achieved in the simulation.

Figure 17
Atmosphere temperature as a function of flight time.

Finally, as presented in Fig. 18, the second stage separation happened at 112 km height and the simulation of the second stage had been continuing during a non-controlled ballistic flight up to impact on the ground. Due to kinetic and potential energy in the second stage separation point, the jettisoned stage reached 377 km high, at 1,309 km from the launch pad, splashed down at 28° azimuth and 2,644 km from the launch point, shown in Fig. 19.

Figure 18
Altitude as a function of flight time. In 80 seconds, separation of the first stage and ignition of second stage happened. In addition to that, the vehicle had risen about 83 km before the burn out and second stage separation.

Figure 19
Altitude as a function of horizontal range. After second stage separation, which occurred at 112 km high, it reached apogee point 377 km above of surface and the impact point was 2,644 km.

In conclusion, for the second configuration flight, the control system commanded properly pitch and yaw vehicle attitudes, apart from holding it very close to trajectory reference. Additionally, the Mach number outcomes presented a value reduced between 146 and 153 seconds of flight. In spite of increase value of velocity, the Mach number rate reduces due to atmospheric temperature changes. There is a strong relationship between sound velocity and atmospheric temperature. Similarly, for the first configuration, the flight trajectory was propagated and generated data about apogee altitude and second stage ground impact.

Third configuration results and discussions

The third configuration flight was divided in ballistic and propelled phases. The former took place from 153 to 230 seconds and the latter burnt for 100 seconds, summed 177 seconds of flight time. In the ballistic flight phase, the vehicle performed a pointing maneuver in order to get a local horizontal and ignited the third stage engine. This powered phase was controlled by a PID control system, in pitch and yaw axes, that tracked a reference trajectory as the early configurations.

As shown in Fig. 20, the pointing maneuver began at 153 seconds, pitched the attitude from -68 to -113 degrees, and there was a time delay between Astos reference and Scott outcomes. As earlier explained, Astos did not consider rotational dynamics and this fact results in time delay whenever an abrupt attitude change happens. Considering that the pitch angle was measured from launch system reference, its large distance from launch point, and the Earth curvature, the local horizontal was measured as -113 degrees rather than -90 degrees.

Figure 20
Third flight configuration pitch attitude as a function of flight time.

Not only pitch attitude changes showed a time delay, but also yaw attitude presented the same behavior, as can be seen in Fig. 21. The graphic shows a time delay caused by the rotational dynamics present in the Scott modeling. After pointing maneuver, both curves had some close values until third stage separation.

Figure 21
Third flight configuration yaw attitude as a function of flight time.

Taking into account the proposed sequence of events, after the second stage separation, the vehicle performed a ballistic flight. In this phase, a pointing maneuver happened as well as a velocity reduction. At 230 seconds, the third stage engine was ignited, increasing the vehicle velocity until 7.24 km/s, and injecting the payload in a circular orbit. At this time, a payload separation and the third stage deorbit engine ignition took place. These outcomes, summarized in Fig. 22, supported the trajectory of a launch vehicle evaluation proposed in this paper.

Figure 22
Third flight configuration velocity as a function of flight time.

Another point to be explored is the orbit altitude reached. As shown in Fig. 23, after second stage separation, the vehicle rose until 233 km before the third stage ignition. In the 100 seconds of burning time, it got almost 70 km and the payload was injected in a 300 km high orbit.

Figure 23
Third flight configuration altitude as a function of flight time.

As both, payload and third stage, were injected into orbit, a deorbit maneuver was necessary in order to bring the stage back to the Earth. As presented in Fig. 23, a 2.44 km/s delta V impulse was applied, and the stage reentered the Earth’s atmosphere. In the Fig. 24, a flight path propagation shows the ground impact at 2,670 km from the launch point, with a vertical descent due to atmospheric drag force. In addition to that, this is just 25 km from second stage impact point, but not at the same time. Despite the time difference between second and third stages, the delta V impulse turned both descent trajectories very close, as predicted by the model.

Figure 24
Third jettisoned stage altitude as a function of horizontal range.

Concluding this section, it may be argued that the payload was injected at a near circular orbit at 300 km above de Earth surface. In addition to that, an efficient deorbit maneuver brought the third stage back in a position near to second stage impact point. In spite of the qualities presented in the Astos computer program, it is crystal clear that rotational dynamic results in a more realistic simulation, and its outcomes could be useful in a refined launch vehicle simulation.

CONCLUSIONS

This work introduces the Scott trajectory simulation tool and present their complementary outcomes in an Astos simulation in terms of attitude changes control. To demonstrate these improved results a three-stages launch vehicle flight was simulated in Scott. This approach included a six-DoF control system, as well as advanced Earth and aerodynamic models. Despite of easier software implementation, a decoupled axes strategy needs a roll control bay in order to keep a low roll rate. Consequently, this system has some disadvantages, such as addition weight and possible mechanical failures.

To alleviate this problem, a coupled axes control strategy was investigated and implemented. This control-system operated in order to keep pitch and yaw attitudes close to reference, regardless of the roll position. As discussed, the couple between axes were modeled inside the error equations, with respect to the body axes frame. This model allows to manage the body attitude by using a reaction control system and/or thrust vector actuator.

In addition to that, a reference attitude and trajectory was created with use of Astos computer program. Furthermore, this reference was optimized in terms of payload mass and orbit altitude. During the Scott simulation, this reference was followed by PID control system. In fact, the Scott software simulated an optimized vehicle following a reference trajectory. However, the outcomes achieved suggested that Scott could provide similar results with a non-optimized path.

The simulations demonstrated that the control system followed satisfactorily the reference trajectory proposed and injected a payload into pre-specified orbit. During the flight, the vehicle was exposed to dynamic pressure and crossed the sound barrier. As result, the time and altitude of transonic flight, maximum dynamic pressure, and their quantitative values were gathered. The outcomes showed the transonic flight and the peak of dynamic pressure happened in separated time. As result of six Dof simulation, Scott computer program provided smoothing and feasible changes attitude, which complemented the Astos results. Finally, the findings that were presented may support a preliminary design of new launch vehicles apart from help in flight safety analyzes.

Future work involves a comparison between Scott’s results and a real rocket launch in order to validate the outcomes reached in this paper.

ACKNOWLEDGMENTS

Not applicable.

• Peer Review History: Single Blind Peer Review.

• FUNDING

  Not applicable.

DATA AVAILABILITY STATEMENT

The data will be available upon request.

REFERENCES

• Astos Solutions GmbH] (2013) Astos Model Library. Unterkirnach: Astos Solutions GmbH.
• Astos Solutions GmbH] (2019) Astos Model Reference. Unterkirnach: Astos Solutions GmbH.
• Aguiar DF (2017) Mission space simulation: launch, orbit injection and re-entry (master’s thesis), São José dos Campos: Instituto Nacional de Pesquisas Espaciais. In Portuguese.
• Anderson J (2010) Fundamentals of aerodynamics. 5th ed. New York: McGraw-Hill.
• Blake WB (1998) Missile DATCOM: User’s Manual – 1997 FORTRAN 90 Version, Air Force Research Laboratories Document AFRL-VA-WP-TR-1998-3009. Columbus: Air Force Research Lab.
• Campos DC (2004) Methodology for evaluating a gain schedule control attitude for a launcher vehicle (master’s thesis). São José dos Campos: National Institute for Space Research. In Portuguese.
• Gunstom B (2004) The Cambridge Aerospace Dictionary. New York: Cambridge University.
• Kadan NV (2019) Practical design of flight control systems for launch vehicles and missiles. New Delhi: Allied Publishers.
• Nise S (2012) Engenharia de sistemas de controle. 6th ed. Rio de Janeiro: LTC.
• Ogata K (1997) Modern control engineering. 3rd ed. New Jersey: Prentice Hall.
• Sirious GM (2004) Missile guidance and control systems. New York: Springer.
• Stevens BJ, Lewis FL (2003) Aircraft control and simulation. 2nd ed. New Jersey: John Wiley & Sons.
• Vos R, Farokhi S (2015) Introduction to transonic aerodynamics. New York: Springer.

Edited by

Section editor: Ana Morais https://orcid.org/0000-0003-0810-179X

Publication Dates

• Publication in this collection
  18 Mar 2024
• Date of issue
  2024

History

• Received
  17 Aug 2023
• Accepted
  07 Dec 2023