Neighboring Optimal Guidance and Constrained Attitude Control for Accurate Orbit Injection

Accurate orbit injection represents a crucial issue in several mission scenarios, e.g., for spacecraft orbiting the Earth or for payload release from the upper stage of an ascent vehicle. This work considers a new guidance and control architecture based on the combined use of (i) the variable-time-domain neighboring optimal guidance technique (VTD-NOG), and (ii) the constrained proportional-derivative (CPD) algorithm for attitude control. More specifically, VTD-NOG & CPD is applied to two distinct injection maneuvers: (a) Hohmann-like finite-thrust transfer from a low Earth orbit to a geostationary orbit, and (b) orbit injection of the upper stage of a launch vehicle. Nonnominal flight conditions are modeled by assuming errors on the initial position, velocity, attitude, and attitude rate, as well as actuation deviations. Extensive Monte Carlo campaigns prove effectiveness and accuracy of the guidance and control methodology at hand, in the presence of realistic deviations from nominal flight conditions.


Introduction
Precise orbit insertion maneuvers represent a crucial issue in several mission scenarios, e.g., for spacecraft orbiting the Earth or for payload release from the upper stage of an ascent vehicle. In fact, the dynamical conditions at injection affect the subsequent phases of spaceflight, because corrective maneuvers may be needed if orbit insertion occurs with unsatisfactory accuracy.
In the scientific literature, only a limited number of works dealt with the joint application of guidance and control (G&C) algorithms to aerospace vehicles. In Ref. [1] proportional-derivative (PD) control is employed for both guidance and control algorithms. Guidance and control based on Nonlinear Dynamic Inversion is studied by Marcos et al. [2], and a comparison between Dynamic Inversion and State Dependent Riccati Equation approaches is presented in Ref. [3]. An integrated G&C method is proposed by Tian et al. [4], while the use of G&C based on sliding-mode is investigated in Ref. [5].
This work considers a recently-introduced guidance and control architecture based on the combined use of (i) the variable-time-domain neighboring optimal guidance technique (VTD-NOG), and (ii) the constrained proportionalderivative (CPD) algorithm for attitude control, with the final aim of achieving precise orbit insertion, in the presence of nonnominal flight conditions.
The Variable-Time-Domain Neighboring Optimal Guidance (VTD-NOG) [6,7] belongs to the class of implicit guidance approaches [8], relies on the second-order sufficient conditions for optimality, and aims at finding the corrective control actions in the neighborhood of the reference trajectory. This must be an optimal path that satisfies the second-order sufficient conditions for optimality. Adoption of a normalized time scale allows overcoming the main difficulties encountered in former formulations of neighboring optimal guidance [8][9][10], which are the occurrence of singularities for the gain matrices and the challenging implementation of the updating law for the time-to-go.
VTD-NOG identifies the trajectory corrections by assuming a thrust direction always aligned with the longitudinal axis of the space vehicle. As a result, VTD-NOG iteratively generates the desired attitude, which can be eventually discontinuous across subsequent guidance intervals. This circumstance implies that the actual 1 3 orientation, which is subject to the space vehicle attitude dynamics, does not coincide with the desired orientation. Hence, the attitude control system must be capable of maintaining the actual orientation sufficiently close to the desired one. This work considers a control algorithm inspired by a proportional-derivative law (PD) [11]. However, a proper saturation action might be necessary to avoid excessive rates for the thrust deflection and the roll torque. The Constrained Proportional Derivative (CPD) algorithm introduces such an appropriate saturation, with the final aim of maintaining these rates within acceptable limits. In past researches VTD-NOG & CPD was successfully applied to lunar ascent [12]. Moreover, closely related G&C algorithms were adopted for low-thrust orbit transfers [13,14].
In this research, VTD-NOG & CPD is applied to two distinct orbit injection maneuvers: (a) Hohmann-like finitethrust transfer from a low Earth orbit to a geostationary orbit, and (b) orbit injection of the upper stage of a launch vehicle. Nonnominal flight conditions are modeled by assuming errors on the initial position, velocity, attitude, and attitude rates, as well as actuation deviations. All of these perturbations are simulated stochastically, in the context of extensive Monte Carlo (MC) campaigns. The ultimate purpose is in demonstrating that VTD-NOG & CPD indeed represents an effective methodology for spacecraft guidance and control, capable of guaranteeing accurate orbit injection, in the presence of deviations from nominal flight conditions.

Spacecraft Dynamics
This paper addresses the problem of accurate orbit injection in two distinct mission scenarios. However, in both cases a maneuvering spacecraft subject to the gravitational attraction of a single body is considered, and the related equations of motion coincide as a result. This section describes the equations that govern both the trajectory and the attitude of the spacecraft of interest.

Trajectory Equations
The spacecraft motion takes place about the Earth, and the dynamics of its mass center is investigated under the following assumptions: (a) the initial and final trajectories are coplanar, (b) the Earth has spherical mass distribution, and (c) the propulsive thrust is continuous.
Assumption (b) implies that the gravitational attraction is directed radially. In light of assumption (a), an optimal trajectory coplanar with the two terminal paths can be reasonably conjectured to outperform any hypothetical alternative three-dimensional trajectory. In fact, due to symmetry of the gravitational field, out-of-plane thrusting has the only effect of rotating the instantaneous velocity, and would imply a useless waste of propellant. This means that the optimal thrust direction is assumed coplanar with the two terminal paths for the entire time of flight, and the optimal powered arc preceding orbit injection is sought in the same plane as a result. However, the trajectory equations (as well as the necessary conditions for optimality) are three-dimensional, because perturbed paths are no longer planar.
The spacecraft motion can be described in a convenient inertial reference frame, associated with the right-hand sequence of unit vectors ĉ 1 ,ĉ 2 ,ĉ 3 , and with origin located at the center of the attracting body. The two terminal trajectories lie on the ĉ 1 ,ĉ 2 -plane (cf. Fig. 1a). The timevarying position can be identified by the following three variables: radius r, right ascension , and declination , portrayed in Fig. 1a. The spacecraft velocity can be projected into the rotating frame r,t,n , where r is aligned with the position vector r and t is parallel to the ĉ 1 ,ĉ 2 -plane (and in the direction of the spacecraft motion, cf. Fig. 1a). The related components are denoted with v r , v t , v n and termed, respectively, radial, transverse, and normal velocity component. The state vector x (with components denoted with x k (k = 1, … , 6) ) of the spacecraft includes the variables associated with the position and velocity vectors and is given by The spacecraft is controlled through the thrust direction, defined by the in-plane angle and the out-of-plane angle , both illustrated in Fig. 1b (in which T is aligned with the thrust direction). Thus, the control vector u is The equations of motion, also termed state equations hence forward, govern the spacecraft trajectory, and involve the state vector x and the control vector u: where (T∕m) is specified and = 398, 600.4 km 3 s 2 is the Earth gravitational parameter. Equations (3) through (6) can be written in the general compact form: In general, the boundary conditions are problem-dependent, and are written in compact form as where subscripts i and f are associated with the initial and final value of the state, t f is the final time, and ̃ is the parameter vector, which contains all the time-independent unknown quantities.
the two inertial frames, i.e., â 1â2â3 The spacecraft instantaneous orientation is associated with the body frame x b ,ŷ b ,ẑ b , whose origin is in the instantaneous center of mass of the vehicle, its axes coincide with the principal axes of inertia, and x b is aligned with the longitudinal axis. The body frame is obtained from the inertial frame â 1 ,â 2 ,â 3 through three elementary counterclockwise rotations, about axes 3, 2, and 1: The angles Ψ,Θ, and Φ are the yaw, pitch, and roll angle, respectively. Thrust vector control (TVC) is used to control pitch and yaw motion of the spacecraft, and side jet system (SJS) is used to control roll motion. The attitude kinematics equations are given by [15] where P, Q, R are the body coordinates of the spacecraft angular velocity. The attitude dynamics equations are given by [15] where I x , I y , and I z are the principal moments of inertia, M x is the torque generated by the side jet system, T is the thrust magnitude, and l is the distance between the center of mass and the swivel point of the TVC, which lies on the longitudinal axis. Variables Δ y and Δ z denote the thrust deflection angles portrayed in Fig. 1c.
The electro-hydraulic servoactuator that controls the engine deflection angles is modeled by the following two first-order systems: [11] In Eq. (16) Δ yc is the commanded Δ y and represents one of the three control inputs for the attitude control system, while the actual angle Δ y (appearing in Eq. (14)) is obtained by saturating Δ y to its maximum value Δ y : Similar considerations apply to Δ zc , Δ z , and Δ z . The actuator of the side jet system is model by the following first order system: for which similar considerations apply. In particular, M xc is the commanded M x , and represents the third control input. Moreover, the saturation value for M x is denoted with M x .

Nominal Trajectory
This study addresses real-time closed-loop guidance and attitude control of a space vehicle during the last powered arc preceding orbit injection. The state components related to the motion of the mass center are collected in x (cf. Eq. (1)), whereas the control vector u includes the thrust direction angles (cf. Eq. (2)). The formal and numerical methodology for the derivation of the nominal trajectory is identical for the two orbit injection maneuvers that are being addressed, and is described in the following.

Formulation of the problem
With regard to the applications of interest, the nominal initial conditions on some state variables are specified = 0, = 0, and v n = 0 , whereas the initial values of the remaining components of x are either all specified (A) or expressed in terms of true anomaly at ignition, f 0 (B): where subscript 0 denotes a specified value, whereas p and e are, respectively, the semilatus rectum and eccentricity of the Keplerian path prior to ignition. The desired final conditions are associated with orbit injection into a circular orbit of specified radius R f and coplanar with the initial trajectory, thus As continuous thrust is employed in the powered arc, propellant consumption is minimized through minimization of the time of flight. Therefore, if the ignition time t 0 is set to 0, the objective function is J = t f .
The problem at hand can be reformulated using the dimensionless normalized time defined as Let the dot denote the derivative with respect to hence forward. Equation (7) is rewritten as In the last relation, the p-dimensional vector a contains the unknown time-independent parameters of the problem. In summary, the optimal control problem consists of identifying a feasible solution that minimizes the objective functional J, through selection of the optimal control law * ( ) and the optimal parameter vector * .

Necessary conditions for optimality
To state the necessary conditions for optimality, which are extensively employed in this research, a Hamiltonian H and a function of the boundary conditions Φ are introduced as [16] where the time-varying, n-dimensional costate vector ( ) and the time-independent, q-dimensional vector are the adjoint variables conjugate to the state Eqs. (22) and to the conditions (8), respectively.
In the presence of an optimal (locally minimizing) solution, the following conditions hold: For the Hamiltonian of Eq. (23) the Pontryagin minimum principle (Eq. (24)) yields the control variables as functions of the adjoint variables and the state variables: where the superscript * denotes the optimal value of the respective variable. Equations (25) are the adjoint (or costate) equations, together with the related boundary conditions; Eq. (26) is equivalent to p algebraic scalar equations. If the control u is unconstrained, then Eq. (24) implies that H is stationary with respect to u along the optimal path, i.e., H * = 0 T . Equations (24) through (26) are well established in optimal control theory (and are proven, for instance, in Ref. [16]), and allow translating the optimal control problem into a two-point boundary-value problem. It is straightforward to demonstrate that the condition (26) is equivalent to As previously remarked, planar trajectories leading to injection can be reasonably assumed to outperform any alternative three-dimensional path. As a result, the problem of determining the minimum-fuel path can be simplified by assuming that at any time the out-of-plane variables equal 0, i.e., = 0,v n = 0 , * 3 = 0, and * 6 = 0 , leading to * = 0.

Indirect heuristic method
Unlike alternative heuristic approaches, the indirect heuristic method (IHM) does not assume any functional structure for the control variables. Instead, the analytical necessary conditions for optimality are derived and employed, for the primary purpose of expressing the control variables as functions of the adjoint variables. These are in turn subject to the costate equations (and the related boundary conditions), which are to be integrated together with the dynamics equations. As a result, a reduced parameter set (mainly composed of the unknown initial values of the adjoint variables) suffices to transcribe the optimal control problem into a parameter optimization problem. This is solved using a heuristic technique (PSO [17] in this work), which employs a population of individuals that evolve toward the optimal solution. Each individual is associated with a particular choice of the parameter set. The optimal control variables are determined without any restriction, because no particular representation is assumed. Furthermore, satisfaction of all the necessary conditions guarantees (at least the local) optimality of the solution. The methodology at hand is thus capable of circumventing some disadvantages of using heuristic techniques, while retaining the main advantage, which is the absence of any starting guess. On the other hand, considerable analytical developments are needed for the application of this technique, and the size of the problem doubles as the adjoint equations must also be integrated for each individual. IHM has already been successfully applied to several space trajectory optimization problems. [18,19]. For minimum-time orbit injection maneuvers the parameter set, associated with each individual, includes the initial values of the adjoint variables, the final time, and the remaining time-independent parameters collected in a. The state Eqs. (22) are integrated together with the adjoint Eqs. (25), while the Pontryagin minimum principle is employed with the intent of supplying the optimal control as a function of the adjoint variables. For each individual, the final condition violations are evaluated. Further details on the specific implementation of IHM are provided in Refs. [18] and [19].

Variable-Time-Domain Neighboring Optimal Guidance
The Variable-Time-Domain Neighboring optimal guidance (VTD-NOG) uses the optimal trajectory as the reference path, with the final intent of determining the control correction at each sampling time t k k=0,…,n S , with t 0 = 0 . These are the times at which the displacement between the actual trajectory, associated with , and the nominal trajectory, corresponding to * , is evaluated, to yield The total number of sampling times, n S , is unspecified, whereas the actual time interval between two successive sampling times is given and denoted with Δt S ,

Time-to-Go Updating Law and Termination Criterion
The fundamental principle that underlies the VTD-NOG scheme consists in finding the control correction ( ) in the generic interval k , k+1 such that the second differential of J is minimized [6,7], while holding the first-order expansions of the state equations, the related final conditions, and the parameter condition. Minimizing the second differential of J is equivalent to solving the accessory optimization problem, defined in the interval k , 1 . The solution of the same problem in the overall interval [0, 1] leads to deriving all the relations reported in Ref. [16]. This means that the latter relations need to be extended to the generic interval k , 1 .
For the sake of brevity, only a few fundamental equations that form VTD-NOG are reported. Full details are contained in Refs. [6] and [7]. The following equation yields the control correction in k , k+1 : where the vectors reported as subscripts denote partial derivatives. Moreover, the iterative correction of the flight time, dt (k) f , derives again from generalization of the accessory optimization problem, defined and solved in each interval Because the actual sampling interval Δt S is specified, the general formula for k is The overall number of intervals n S is found at the first occurrence of the following condition: The introduction of the normalized time is extremely useful, because all the gain matrices are defined in the normalized interval [0,1] and cannot become singular. Moreover, the limiting values k are dynamically calculated at each sampling time using Eq. (31). Also the termination criterion has a consistent definition, and corresponds to the upper bound of the interval [0,1], to which is constrained.

Gain Matrices
The definition of a neighboring optimal path requires the numerical backward integration of the sweep Eqs. [6,7,16], for the purpose of obtaining several time-varying gain matrices, which are essential for real-time neighboring optimal guidance. These matrices are evaluated along the optimal trajectory, which represents the reference solution. For the sake of brevity, only an outline of the steps needed to get the gain matrices is presented in the following. Full details can be found again in Refs. [6] and [7].
As a first step, matrices A, B, C, D, E, and F are derived analytically [7], and evaluated along the optimal path. Then, the following sweep equations for Ŝ , R, Q, n, and are integrated backward in time: where represents a constant, auxiliary matrix [6,7], and The displacement can be written in terms of the preceding matrices, as

Preliminary Offline Computations and VTD-NOG & CPD Structure
The implementation of VTD-NOG requires several preliminary computations that can be completed offline, and yield several nominal quantities, stored in the onboard computer. First of all, the optimal trajectory is to be determined, together with the related state, costate, and control variables, which are assumed as the nominal ones. These are obtained in the time domain and are represented as sequences of equally-spaced values, e.g., * i = * i . The successive step is the analytical der ivation of the matr ices , , , H , H , H , H , H , H , H , H , H , H , H , H , which are evaluated along the optimal path. Then, also the matrices A, B, C, D, E, and F are introduced and evaluated. All of these quantities are interpolated, to be available at arbitrary values of in the interval [0, 1] . Subsequently, the two-step backward integration of the sweep equations is performed, and yields the gain matrices Ŝ , R, m, Q, n, and , using also the analytic expressions of W, U, and V (written in terms of R, m, Q, n, and ). Finally, all the remaining nominal quantities are interpolated, to be available again at arbitrary values of .
Using the nominal quantities computed offline, at each time k the VTD-NOG algorithm determines the time of flight t (k) f , the value k+1 , and the control correction ( ) . VTD-NOG ends when Eq. (31) is satisfied. Figure 2 portrays a block diagram that illustrates the sampled-data feedback structure of the VTD-NOG algorithm, in which the control and flight time corrections definitely depend on the state displacement (evaluated at specified discrete times) through the time-varying gain matrices, which are computed offline and stored onboard. The attitude control loop (encircled by the dotted line) is being described in detail in the following.

Constrained Proportional-Derivative Attitude Control
The attitude control system is designed for the purpose of guaranteeing the correct spacecraft orientation, on the basis of the corrected control u yielded iteratively by VTD-NOG. The output of the attitude control system is represented by the actual control a (cf. Fig. 2).

Commanded Attitude
With reference to Fig. 2, VTD-NOG yields the corrected control u, i.e., the thrust direction identified by the angles and , under the approximating assumption that this direction is always aligned with the longitudinal axis of the spacecraft. This means that these two corrected angles represent the commanded values, denoted with c and c , that the attitude control system must pursue through the use of thrust vectoring.
The commanded yaw, pitch, and roll angles are denoted with Ψ c ,Θ c , and Φ c , and correspond to a thrust direction T c that points toward x b . This alignment condition leads to On the other hand, the commanded thrust direction is identified through the c and c : Identification of Eqs. (38) and (39) leads to obtaining Ψ c and Θ c as functions of c , c , , and . The commanded roll angle Φ c is independent of T c , and can be selected on the basis of practical requirements.

Control Algorithm
A baseline attitude control action for such spacecraft is given by the following PD control: (38) T c = cos c sin c cos c cos c sin c rtn T = cos c sin c cos c cos c sin c ⋅ To analyze convergence properties achieved by the considered controller, first it is worth noticing that most of times the commanded angles Φ c , Θ c , and Ψ c can be modeled as constant, since the guidance command usually changes slowly compared to attitude maneuvers. Next, it will be shown that if Φ c , Θ c , and Ψ c are constant, then the proposed PD control guarantees local convergence to the desired attitude. In fact, linearizing Eqs. (10) through (15) about Φ = Θ = Ψ = 0 , P = Q = R = 0 , M x = 0, Δ y = Δ z = 0 , the following relations are obtained: Linearization of the actuator's Eqs. (16) and (18) about Note that equations relative to each single axis are decoupled from the others. Then, it is easy to obtain that the linearized closed-loop system given by Eqs. (40)-(44) possesses the following property. If Φ c , Θ c , and Ψ c are constant and k dx > k px x > 0 , k dy > k py y > 0 , k dz > k pz z > 0 , then The considered PD control can lead to excessive rates for the thrust deflection angles Δ y and Δ z , and roll control torque M x . In fact, high values for the gains k px , k dx , k py , k dy , k pz , and k dz might be required to obtain a fast response of the attitude control system in comparison with the guidance command. Then, high gains can in turn lead to high amplitudes for the rates of M x , Δ y , and Δ z . If the rates are too high, then clearly they become physically infeasible. The latter issue is here tackled using Constrained Proportional and Derivative (CPD) control, which is described by the following equations: which replace Eqs. (40) and (41). In Eqs. (45)-(46) M xc > 0 , Δ yc > 0 , and Δ zc > 0 are additional design parameters. Even though the formal proof is omitted for the sake of conciseness [12], the use of CPD control in Eqs. (45)-(46) ensures that (44) Thus, under those conditions also CPD control achieves local convergence to the desired attitude.

Determination of Control Gains
The goal of the current subsection is presenting a method for determining at least first guess values for the gains k px , k dx , k py , k dy , k pz , and k dz . The method is here illustrated only for gains k py and k dy , because it can be easily extended to the other gains. Neglect dynamics of the actuator in the first relation of Eq. (44) and consider the linearized closed-loop system in Eqs. (41) and (42) ceeding by trial-and-error. Because G y k py = 2 ny and G y k dy = 2 y ny , it is easy to verify that the specifications y ≥ y and ny ≥ ny are fulfilled for all G y ≤ G y ≤ G y by setting k py = ny 2 G y and k dy = 2 y ny G y .

Actual Attitude and Thrust Angles
Using the relation â 1â2â3 the actual spacecraft orientation is identified through the following equation: Moreover, the actual thrust direction T a depends on the two angles Δ y and Δ z : T a = cos Δ z cos Δ y sin Δ z cos Δ z sin Δ y x bŷbẑb T .

Leo-to-Geo Orbit Transfer
As a first application, this research considers the injection maneuver that concludes a high-thrust orbit transfer, from a low Earth orbit (LEO) to a geostationary orbit (GEO), in the presence of nonnominal flight conditions. As continuous, constant thrust is employed, the nominal thrust acceleration is T∕m = n 0 c c − n 0 t , where c is the (constant) effective exhaust velocity of the propulsive system, n 0 is the initial thrust acceleration (at t 0 = 0 ), and t is the actual time ( c = 3 km∕sec and n 0 = 1.112 g). The initial LEO is assumed to be circular, at altitude of 400 km, and coplanar with the final geostationary orbit. The optimal transfer is composed of three arcs: (a) a short-duration powered phase, (b) a coast (ballistic) arc, and (c) a final thrust arc, which concludes at orbit injection.
This work considers the problem of minimizing the time of flight of the last powered phase, through selection of the ignition time (along the coast arc) and determination of the optimal thrust direction time history. The Keplerian ellipse preceding the last ignition is the one that follows the first velocity change of the LEO-to-GEO Hohmann transfer. This means that p = 24471 km and e = 0.723 . The initial conditions at ignition are written in terms of f 0 (cf. Equation (19) case (B)), whose optimal value, found through IHM, is 179.9 deg. The final conditions are given by Eqs. (20), where R f = 42164 km . The minimum time needed for orbit injection is 105.9 s, and is obtained together with the optimal control time history, through the use of IHM. (51) (52) T a = cos a sin a cos a cos a sin a rtn T .
The initial spacecraft mass m 0 equals 2400 kg. Further characteristics are the maximal deflection angles Δ y and Δ z (both set to 3 deg), the maximal torque generated by the side jet system M x (set to 60 Nm), the time-varying distance l , and inertia moments I x , I y , and I z : where l 0 = 0.25 m , ̇l = 2.4 ⋅ 10 −3 m∕sec , I x0 = 1200 kg m 2 , ̇I x = −4.36 kg m 2 sec , I y0 = I z0 = 800 kg m 2 , a n d ̇I y =̇I z = −2.91 kg m 2 sec .
T h e v a l u e s This alignment condition was assumed for the derivation of the optimal ascent path. However, the actual spacecraft dynamics is driven by a thrust direction not exactly aligned with the longitudinal axis, due to the use of thrust vectoring for attitude control. This circumstance is apparent also by inspection of Fig. 2, which illustrates clearly that the corrected control u does not coincide with the actual control a , which affects the real dynamics of the center of mass. As a first step, VTD-NOG & CPD has been tested to evaluate these deviations, exclusively related to the alignment assumption. Table 6 reports the related results (obtained in a single simulation), i.e., the final displacements from the nominal final altitude, declination, and velocity components, and testifies to the excellent accuracy of VTD-NOG & CPD in this context. However, perturbations can exist that affect the overall spacecraft dynamics. Monte Carlo (MC) campaigns are usually run, with the intent of obtaining some useful statistical information on the accuracy of the guidance and control algorithm of interest, in the presence of the existing perturbations, which are simulated stochastically.
More specifically, for the initial conditions an error on the initial radius and on the initial declination are assumed, with Gaussian distribution, zero mean value and standard deviation r ( ) 0 = 5 km (for r 0 ) and ( ) 0 = 6.8 ⋅ 10 −3 deg (for 0 ) corresponding to 5 km in direction normal to the plane of the nominal trajectory. Errors on the initial velocity with Gaussian distribution, zero mean and standard deviation of 30 m/s are introduced as well. These errors are consistent with some values found in the scientific literature [13,14,20]. Moreover, errors on the initial attitude angles and rates are introduced. All these displacements have Gaussian distribution and zero mean. Their standard deviation equals 10 deg for the initial attitude angles and 10 deg/sec for the initial attitude rates.
Moreover, propulsive perturbations are included in simulations. In fact, usually the thrust magnitude (and the related acceleration, as a result) exhibits small fluctuations. This timevarying behavior is modeled through a trigonometric series with random coefficients: where n (p) 0 denotes the perturbed value of n 0 , whereas the coefficients ã k k=1,…,10 have a random Gaussian distribution centered around the zero and a standard deviation equal to 0.02.
Furthermore, actuation errors may affect the spacecraft performance. They are modeled through Gaussian noise terms, with zero mean. A term with standard deviation 0.06 Nm is added to the roll control torque M x . Terms with standard deviation of 0.05 deg are added to the thrust angles Δ y and Δ z in Eq. (16). All the variables subject to uncertainties are reported in Tables 3, 4, 5, 6.
A single MC campaign is performed, including 100 numerical simulations, which assume all the previously described perturbations. Figures 3, 4 through 5 portray the time histories of the relevant trajectory variables, M x , the engine deflection angles, as well as the commanded and actual attitude angles. Moreover, two statistical quantities are evaluated, i.e., the mean value and the standard deviation for all of the outputs of interest. The symbols ____ Δ and ( ) will denote the mean error (with respect to the nominal value) and standard deviation of henceforth. Table 7 reports the statistics on the errors at injection and the time of flight. Inspection of this table reveals that VTD-NOG & CPD guarantees orbit injection with excellent accuracy. Furthermore, the average time of flight is very close to the nominal value, and the corresponding standard deviation is modest. Table 3 LEO-to-GEO transfer: standard deviations on the initial position and velocity variables  Table 4 LEO-to-GEO transfer: standard deviations on the initial attitude and angular rates

Upper Stage Orbit Injection
As a second application, the terminal ascent path of the upper stage of the Scout launch vehicle is considered. The actual time histories of the rocket thrust magnitude and mass are employed [21], and this implies that the thrust acceleration time history is specified as well.
The overall ascent trajectory that ends with circular orbit injection at 800 km of altitude is addressed in Ref. [22]. In this, research, only the minimum-time path of the upper stage is needed as the reference trajectory, and is found by where R E is the Earth radius. The final conditions at orbit injection are given by Eqs. (20), where R f = R E + 800 km . The minimum time needed for orbit injection is 25.8 s.
Further characteristics of the upper stage are the maximal deflection angles Δ y and Δ z (both set to 3 deg), the maximal torque generated by the side jet system M x (set to 200 Nm), the time-varying distance l , and inertia moments I x , I y , and I z , whose time evolution is given by Eq.  Tables 8 and 9. Moreover, the following values are selected for VTD-NOG & CPD. The sampling interval Δt S is set to 0.2 s, and the CPD gains are determined using the same steps adopted in the preceding section. Because T min = 23.7 N is the minimum thrust [20], the lower bounds are r 0 ) and ( ) 0 = 0.018 deg (for 0 ) corresponding to 2 km in direction normal to the plane of the nominal trajectory. Errors on the initial velocity with Gaussian distribution, zero mean and standard deviation of 10 m/sec are introduced as well. Moreover, errors on the initial attitude angles and rates are introduced, as well as actuation errors. These are all modeled as in the previous section, except for the error on M x , whose standard deviation is set to 0.2 Nm. All the variables subject to uncertainties are reported in Tables 10,11,12. A single MC campaign is performed, including 100 numerical simulations, which assume all the previously described perturbations. Figure 6 portrays the relevant state variables, i.e., altitude, declination, and the three velocity components, together with the roll torque M x . It is remarkable that both the radial and normal velocity components exhibit remarkable deviations from the respective nominal values. Nevertheless, the final convergence to the desired conditions is achieved to a great accuracy. Torque M x reaches the maximal value only in the first part of the injection arc (Table 13). Figure 7 depicts the nozzle deflection angles, which assume modest values and never saturate, unlike what occurs for the preceding case, focused on orbit injection at GEO. Finally, Fig. 8 illustrates the time histories of the commanded and actual attitude angles, showing the correct tracking of the former by the latter ones. At the end of the MC campaign, two statistical quantities are evaluated, i.e., the mean value and the standard deviation for all of the outputs of interest. Table 14 reports the statistics on the errors at injection and the time of flight. Inspection of this table reveals that VTD-NOG & CPD guarantees orbit injection with excellent accuracy. Furthermore, the average time of flight is very close to the nominal value, and the corresponding standard deviation is modest.

Conclusion
This work considers a new guidance and control architecture based on the combined use of (i) the variable-time-domain neighboring optimal guidance technique (VTD-NOG), and (ii) the constrained proportional-derivative (CPD) algorithm for attitude control, with the final aim of achieving precise orbit insertion, in the presence of nonnominal flight conditions. VTD-NOG relies on the second-order sufficient conditions for optimality, and finds the corrective control actions in the neighborhood of the reference trajectory. This 3.5e −6 0.95 0.60 0.83 0.32 must be an optimal path that satisfies the second-order sufficient conditions for optimality. A fundamental original feature of VTD-NOG is the use of a normalized time scale as the domain in which the nominal trajectory and the related vectors and matrices are defined. The attitude control system employs thrust vectoring and a side jet system. The CPD attitude control law is proportional-derivative-like, and guarantees local convergence of the actual attitude toward the desired attitude prescribed by VTD-NOG. Moreover, a proper saturation action is introduced, to avoid excessive rates for the thrust deflection and the roll torque. VTD-NOG & CPD is applied to two distinct orbit injection maneuvers: (a) Hohmann-like finite-thrust transfer from a low Earth orbit to a geostationary orbit and (b) orbit injection of the upper stage of a launch vehicle. Nonnominal flight conditions are modeled by assuming errors on the initial position, velocity, and attitude, as well as actuation deviations. Extensive Monte Carlo campaigns prove effectiveness and accuracy of the guidance and control methodology at hand, in the presence of realistic deviations from nominal flight conditions.
Funding Open access funding provided by Università degli Studi di Roma La Sapienza within the CRUI-CARE Agreement.