Low-Thrust Orbital Transfers in the Two-Body Problem

Abstract

Low-thrust transfers between given orbits within the two-body problem are considered; the thrust is assumed power limited. A simple method for obtaining the transfer trajectories based on the linearization of the motion near reference orbits is suggested. Required calculation accuracy can be reached by means of use of a proper number of the reference orbits. The method may be used in the case of a large number of the orbits around the attracting center; no averaging is necessary in this case. The suggested method also is applicable to the cases of partly given final orbit and if there are constraints on the thrust direction. The method gives an optimal solution to the linearized problem which is not optimal for the original nonlinear problem; the difference between the optimal solutions to the original and linearized problems is estimated using a numerical example. Also examples illustrating the method capacities are given.

1. Introduction

In this paper spiral transfers with power-limited low thrust between two given orbits are considered. The transfer trajectory may include a large number of orbits around the attracting center if the transfer is performed in a strong gravity field; this situation takes place, for example, near Earth. This makes optimization of the transfer more difficult.

There are various mathematical methods for optimization of multirevolution orbital transfers [1–7], most of them are based on averaging of motion. However, the methods have some defects: they are complicated [3, 5] or have a limited application, that is, are applicable only to the circular or neighboring orbits [1, 6, 7] or only to the planar or coaxial orbits [2, 4, 6, 7].

This paper suggests a simple mathematical method for calculation of the orbital transfers within the two-body problem. (A paper describing a generalization of the method to any force field will be published in the Cosmic Research journal.) The method is based on a linearization of motion near some reference orbits; in this respect the method is similar to the Method of Transporting Trajectory (MTT) for solving the two-point boundary value problem (TPBVP) suggested in [8] and developed in [9–12]. The method suggested here makes it possible obtaining transfer trajectories of the following types:(i)transfer between a given position in the initial orbit and an optimally determined position in the final orbit (this case takes place, for example, if the spacecraft launches from a given orbit at a given time and phasing of the final orbit is not necessary);(ii)transfer between the initial orbit and a given position in the final orbit with obtaining an optimal position of the launch (this situation takes place, for instance, in the transfer to a given point of geostationary orbit if the time of launch from the initial orbit may be selected arbitrarily);(iii)transfer between two given orbits with obtaining optimal launch and arrival positions (this is a classical case of the orbital transfer).

The method suggested in this paper does not need any averaging of motion. Required calculation accuracy is reached by means of use of a proper number of the reference orbits: the bigger is the number, the shorter are intervals of linearization and the higher is accuracy. The method does not put any limits on the number of orbits around the attracting center and on the number of the reference orbits. The suggested method is applicable also to the case of a partly given final orbit (e.g., only semimajor axis and eccentricity or only orbital energy may be given) and to the case of constraints imposed on the thrust direction.

However, the method gives an optimal solution to the linearized problem and this solution is not optimal for the initial nonlinear problem. Note that this nonoptimality is not caused by the linearization errors, but is a principal sequence of the substitution of the original problem for the linearized one. A comparison of the optimal solutions to the initial and linearized problems made for a numerical example is given in the paper. Other numerical examples illustrate the method capacities.

2. Statement of the Problem

The spacecraft motion is described by the following:𝐱̇=𝐟(𝐱,𝑡)+𝐠, (2.1) where 𝐱={𝐫,𝐯} is the spacecraft state vector, 𝐟=𝐟(𝐱,𝑡)={𝐯,𝐟_𝑣}, 𝐟_𝑣=𝐟_𝑣(𝐱,𝑡) is the acceleration caused by the external forces,𝐠={𝟎,𝛼}. (2.2) In the considered case of the two-body problem,𝐟_𝑣𝜇=−𝑟³𝐫. (2.3) Let us assume for simplicity that the electric power of the propulsion system is constant. This case takes place, for example, if a nuclear power source is used or if the transfer is performed near a planet with nearly circular orbit. Then the performance index for the power-limited thrust is1𝐽=2^𝑇₀𝛼²𝑑𝑡, (2.4) where 𝛼=|𝛼|. Minimal value of the functional (2.4) gives minimal propellant consumption [13].

Let 𝐪_𝑖,𝐪_𝑓 be 5-dimensional vectors of orbital elements defining the initial and final orbits. The problem is to find transfer trajectory between the initial and final orbits in given time 𝑇 in which minimal value of the performance index (2.4) is reached.

3. Transfer between Given State and Given Orbit

The boundary values for the problem formulated in Section 2 may be written as𝐱₀=𝐲_𝑖0,𝐱_𝑇=𝐲_𝑓𝑇, (3.1) where 𝐲_𝑖,𝐲_𝑓 are state vectors of the initial and final orbits. The case will be considered in this section when vector 𝐲_𝑖0 is given and vector 𝐲_𝑓𝑇 is not given.

3.1. Transfer between Neighboring Orbits

First let us assume that the initial and final orbits are close to each other. Then the equation of motion (2.1) can be linearized near the initial orbit as follows (Figure 1):𝝃̇=𝐅𝝃+𝐠, (3.2) where𝜉=𝜉(𝑡)=𝐱(𝑡)−𝐲_𝑖(𝑡), (3.3) is state vector of the linearized motion and𝐅=𝜕𝐟𝜕𝐲_𝑖. (3.4) Vector f and matrix F in (3.2) and (3.4) are calculated in the initial orbit, this is why matrix F is a function of time even if vector f does not depend on time explicitly. The Hamiltonian for (3.2) and performance index (2.4) is [14]𝛼𝐻=−²2+𝐩^𝑡𝐅𝝃+𝐏^𝑡_𝑣𝜶+𝑝_𝑡, (3.5) where 𝐩={𝐩_𝑟,𝐩_𝑣} is a vector of the costate variables corresponding to (3.2), 𝐩_𝑣 is Lawden’s primer vector, 𝐩_𝑡 is costate variable corresponding to the additional equation 𝑡̇=1 making the system autonomous. Vector 𝐩 satisfies the following costate variational equation:𝐩̇^𝑡=−𝜕𝐻𝜕𝜉^𝑡=−𝐩^𝑡𝐅, (3.6) variable 𝑝_𝑡 satisfies the following equation:𝑝_𝑡̇=−𝜕𝐻𝜕𝑡=−𝐩^̀‌o𝐅̇𝝃. (3.7) Let us consider sixth-order matrixΨ=Ψ(𝑡) which is a general solution to (3.6) with initial value Ψ(0)=𝐈. Matrix Ψ is costate transition matrix given by𝚿=𝜕𝐲_𝑖0/𝜕𝐲_𝑖. (3.8) Dividing matrix Ψ into two 6×3-dimensional submatrices as follows:𝚿=𝚿_𝑟𝚿_𝑣, (3.9) the costate variables may be represented as𝐩=𝚿^𝑡𝛽,𝐩_𝑣=𝚿^𝑡_𝑣𝛽, (3.10) where 𝛽 is a constant 6-dimensional vector. The Hamiltonian (3.5) reaches its maximum if𝛼=𝐩_𝑣=𝚿^𝑡_𝑣𝛽. (3.11) Solution to (3.2) is given by the Cauchy formula𝜉=𝝃(𝑡)=𝝃₀+^𝑡₀Φ(𝑡,𝜏)𝐠𝑑𝜏, (3.12) where𝚽𝑡₁,𝑡₂=𝜕𝐲_𝑖𝑡₁𝜕𝐲_𝑖𝑡₂ (3.13) is the state transition matrix. Due to (3.1) and (3.3) 𝜉₀=0 in (3.12). Using (2.2), (3.9), and (3.11) and relations𝚽(𝑡,𝜏)=𝚽(𝑡,0)Φ(0,𝜏)=𝚽(𝑡,0)𝚽⁻¹(𝜏,0),𝚽=𝚽(𝑡,0)=𝚿⁻¹, (3.14) the solution (3.12) can be represented as𝜉=𝚽𝐒𝜷, (3.15) where𝐒=𝐒(𝑡)=^𝑡₀𝚿_𝑣𝚿^𝑡_𝑣𝑑𝑡 (3.16) is matrix of sixth order (also see [9]). As follows from (3.3), (3.11), and (3.15), in order to find optimal thrust vector and the state vector of the transfer trajectory, it is sufficient to obtain vector 𝜷. Nonspecified state vector 𝐲_𝑓𝑇 may be found using transversality condition which in the considered case is𝐩_𝑇=𝐩(𝑇)=𝜕𝐪_𝑓𝜕𝐲_𝑓𝑇^𝑡𝝈, (3.17) where 𝝈 is an arbitrary constant 5-dimensional vector. Due to the closeness of the initial and final orbits, the approximate equality 𝜕𝐪_𝑓/𝜕𝐲_𝑓≈𝐔 is fulfilled, where𝐔=𝐔(𝑡)=𝜕𝐪_𝑖𝜕𝐲_𝑖(𝑡), (3.18) is a 5×6-dimensional matrix. Then the condition (3.17) takes the form𝐩_𝑇=𝐔^𝑡_𝑇𝝈. (3.19) Equations (3.10) and (3.19) and the last relation of (3.14) give𝜷=𝐔_𝑇𝚽_𝑇^𝑡𝝈. (3.20) The closeness of the initial and final orbits makes it possible to use the following approximate linear relation:Δ𝐪=𝐪_𝑓−𝐪_𝑖=𝐔_𝑇𝝃_𝑇. (3.21) Equations (3.15), (3.20), and (3.21) giveΔ𝐪=𝐖𝝈, (3.22) where𝐖=𝐔_𝑇𝚽_𝑇𝐒_𝑇𝐔_𝑇𝚽_𝑇^𝑡 (3.23) is 5-dimensional matrix. Using (3.8), (3.9), (3.13), (3.14), (3.16), and (3.18), we obtain𝐔_𝑇𝚽_𝑇=𝜕𝐪_𝑖𝜕𝐲_𝑖𝑇𝜕𝐲_𝑖𝑇𝜕𝐲_𝑖0=𝜕𝐪_𝑖𝜕𝐲_𝑖0=𝐔₀, (3.24)𝐔_𝑇𝚽_𝑇𝐒_𝑇𝐔_𝑇𝚽_𝑇^𝑡=𝜕𝐪_𝑖𝜕𝐲_𝑖0^𝑇₀𝜕𝐲_𝑖0𝜕𝐯_𝑖𝜕𝐲_𝑖0𝜕𝐯_𝑖^𝑡𝑑𝑡𝜕𝐪_𝑖𝜕𝐲_𝑖0^𝑡, (3.25) where 𝐯_𝑖=𝐯_𝑖(𝑡) is velocity vector in the initial orbit. (3.23), (3.25) give𝐖=𝐔₀𝐒_𝑇𝐔^𝑡₀=^𝑇₀𝐐𝐐^𝑡𝑑𝑡, (3.26) where𝐐=𝐐(𝑡)=−𝜕𝐪_𝑖𝜕𝐯_𝑖 (3.27) is a 5×3-dimensional matrix; as is seen from (3.18) and (3.27), matrix 𝐐 is a part of matrix 𝐔. Using (3.20), (3.22), and (3.24), vector 𝜷 can be found as follows:𝜷=𝐔^𝑡₀𝐖⁻¹Δ𝐪. (3.28) Now due to (3.3), (3.11), (3.15), and (3.28), the optimal thrust vector and the spacecraft state vector are𝛼=𝚿^𝑡_𝑣𝐔^𝑡₀𝐖⁻¹Δ𝐪, (3.29)𝐱(𝑡)=𝐲_𝑖(𝑡)+𝚽𝐒𝐔^𝑡₀𝐖⁻¹Δ𝐪. (3.30) Equations (3.1) and (3.30) give the state vector of entry into the final orbit𝐲_𝑓𝑇=𝐲_𝑖𝑇+𝚽_𝑇𝐒_𝑇𝐔^𝑡₀Δ𝐪. (3.31) Substituting (3.29) into (2.4) and taking into account definitions (3.16), (3.26), the minimal value of the performance index can be found as follows:1𝐽=2Δ𝐪^t𝐖⁻¹Δ𝐪. (3.32) Equations (3.29), (3.30), and (3.32) give the solution to the problem.

Figure 1 

Transfer between given orbits.

Note that matrices 𝐔, 𝐖, Φ, Ψ are calculated in [9, 15, 16] in an explicit form which makes the suggested method analytical.

3.2. Transfer between Arbitrary Orbits

Let us consider transfer between arbitrary orbits and divide the time interval 𝐓 into 𝐧 subintervals defined by instants 𝐭_𝟎=𝐭_𝟏,…,𝐭_𝐧−𝟏,𝐭_𝐧=𝐓. Also we assume that 𝐧−𝟏 intermediate reference orbits between the initial and final orbits are specified somehow and 𝐪_𝟏,…,𝐪_𝐧−𝟏 are 5-dimentional vectors of elements of the reference orbits (see Figure 2). These vectors may be given, for instance, in the following way:𝐪_𝑗=𝐪_𝑗+𝑗𝑛𝐪_𝑓−𝐪_𝑖,𝑗=1,…,𝑛−1. (3.33) Let us define vectorsΔ𝐪_𝑗=𝐪_𝑗−𝐪_𝑗−1,𝑗=1,…,𝑛, (3.34) where 𝐪₀=𝐪_𝑖,𝐪_𝑛=𝐪_𝑓 (here subscript “0” is the number of the reference orbit). Equations (3.21) and (3.34) give the following equality:^𝑛_𝑗=1Δ𝐪_𝑗=Δ𝐪. (3.35) Let us divide the transfer trajectory into 𝑛 arcs corresponding to the time subintervals and assume that the jth arc begins in the j–1st reference orbit and ends in the jth one. Assuming number 𝑛 big enough to make the j–1st and jth reference orbits close to each other for all 𝑗=1,…,𝑛, the results of the Section 3.1 may be applied to each pair of the reference orbits. Now the problem is to find the reference orbits providing optimality of the transfer. Due to (3.32) the performance index for the jth transfer arc is𝐽_𝑗=12Δ𝐪^𝑡_𝑗𝐖_𝑗⁻¹Δ𝐪_𝑗,𝑗=1,…,𝑛, (3.36) where similarly to (3.26), (3.27) the following definitions are used:𝐖_𝑗=^𝑡𝑗_𝑡𝑗−1𝐐_𝑗𝐐^𝑡_𝑗𝑑𝑡,𝐐_𝑗=𝜕𝐪_𝑗−1𝜕𝐯_𝑗−1𝑗=1,…,𝑛, (3.37) where (41) 𝐯_𝑗 is velocity vector of the jth reference orbit. Performance index of the whole problem is𝐽=^𝑛_𝑗=1𝐽_𝑗. (3.38) In order to find the transfer trajectory that gives minimum value for 𝐽, it is sufficient to find intermediate reference orbits providing a minimum for (3.38). Thus, function (3.38) should be minimized with respect to vectors Δ𝐪_𝑗,𝑗=1,…,𝑛 taking into account (3.35). Let us introduce the helping function𝐿=𝐽−𝜆^𝑡𝑛_𝑗=1Δ𝐪_𝑗−Δ𝐪, (3.39)

Figure 2 

Transfer between arbitrary orbits.

where 𝜆 is a Lagrange multiplier. Necessary conditions of a minimum for the functional (3.38) are𝜕𝐿𝜕Δ𝐪_𝑗^𝑡=𝐖_𝑗⁻¹Δ𝐪_𝑗−𝜆=𝟎𝑗=1,…,𝑛. (3.40) Thus,Δ𝐪_𝑗=𝐖_𝑗𝝀𝑗=1,…,𝑛, (3.41) and (3.35) and (3.41) give𝝀=𝐖⁻¹Δ𝐪, (3.42) where𝐖=^𝑛_𝑗=1𝐖_𝑗. (3.43) Now new values of vectors Δ𝐪_𝑗 can be found from (3.41) and (3.42) as follows:Δ𝐪_𝑗=𝐖_𝑗𝐖⁻¹Δ𝐪𝑗=1,…,𝑛. (3.44) New reference orbits are defined by the new values of the orbital elements given by𝐪_𝑗+1=𝐪_𝑗+Δ𝐪_𝑗𝑗=0,…,𝑛−1. (3.45) Assuming optimal locations of the reference orbits found and applying (3.29) and (3.30) to the jth arc of the transfer trajectory, we obtain optimal thrust vector and the trajectory state vector in the time interval 𝑡_𝑗−1≤𝑡≤𝑡_𝑗(𝑗=1,…,𝑛) as follows:𝜶(𝑡)=𝚿^𝑡_𝑗𝑣𝐔^𝑡_𝑗𝐖⁻¹Δ𝐪, (3.46)𝐱(𝑡)=𝐲_𝑗−1(𝑡)+𝚽_𝑗𝐒_𝑗𝐔^𝑡_𝑗0𝐖⁻¹Δ𝐪, (3.47) where𝐒_𝑗=𝐒_𝑗(𝑡)=^𝑡_𝑡𝑗−1𝚿_𝑗𝑣𝚿^𝑡_𝑗𝑣𝑑𝑡,𝐔_𝑗=𝜕𝐪_𝑗−1𝜕𝐲_𝑗−1, (3.48) matrix 𝐖 is given by (3.43), matrices Φ=Φ_𝑗(𝑡,0),Ψ_𝑗𝑣=Ψ_𝑗𝑣(𝑡) are calculated in the 𝑗−1st reference orbit.

Multiplying (3.40) by Δ𝐪_𝑗 and taking into account (3.35), (3.36), and (3.38), we obtain1𝐽=2𝜆^𝑡Δ𝐪. (3.49) Due to (3.42) and (3.49) minimal value of the performance index in the considered case is given by (3.32) with matrix W defined by (3.43).

3.3. Calculation Procedure

Let 𝐲⁰_𝑗,𝐲¹_𝑗 be the state vectors of the jth reference orbit at the beginning and at the end of the jth time subinterval (i.e., at times 𝑡_𝑗−1, 𝑡_𝑗 resp.). Then, the solution to the problem considered here may be obtained by means of the following iterative calculation procedure.(1)𝑛−1 intermediate reference orbits are specified somehow, for example, using (3.33). The launch position in the initial orbit is specified in the case considered here (i.e., state 𝐲⁰₀=𝐲_𝑖(0) is given) and the respective initial state vector of the transfer trajectory is 𝐱₀=𝐲⁰₀.(2)Vector 𝐲⁰_𝑗 is calculated for 𝑗=1 using the following equation, that is, similar to (3.30):𝐲⁰_𝑗=𝐲¹_𝑗−1+𝚽¹_𝑗𝐒¹_𝑗𝐔^𝑡_𝑗0𝐖_𝑗⁻¹Δ𝐪_𝑗, (3.50) where 𝐒¹_𝑗=𝐒_𝑗(𝑡_𝑗),𝐔_𝑗0=𝐔_𝑗(0) matrix Ψ_𝑗𝑣=Ψ_𝑣(𝑡_𝑗) is calculated in the 𝑗−1st reference orbit.(3)Step  2 is repeated for 𝑗=2,…,𝑛−1 and for 𝑗=𝑛 (3.50) gives state vector 𝐲_𝑓𝑇 in the final orbit.(4)New vectors Δ𝐪_𝑗 are calculated using (3.44) and new reference orbits with elements given by (3.45) are found.(5)Performance index is calculated using (3.32) and steps  2–4 are repeated until decrement Δ𝐽 of the performance index gets smaller than a given parameter 𝜀>0. As soon as |Δ𝐽|<𝜀, the thrust vector 𝜶 and the state vector 𝐱 of the transfer trajectory may be calculated at each time subinterval using (3.46) and (3.47).

The suggested method is approximate, although any desired accuracy can be reached by means of selecting an appropriate amount 𝑛 of subintervals. It can be shown that if 𝑛→∞ then the solution converges to a limit which is an accurate solution to the linearized problem.

4. Other Types of the Transfer

4.1. Transfer between Given Orbit and Given State Vector

Now let us consider the case when the launch position can be selected in the initial orbit arbitrarily and the position of the entry into the final orbit is given (i.e., vector y_𝑖0 is not given and vector 𝐲_𝑓𝑇 is given in (3.1)). In this case the method described in Section 3 should be applied in the backward direction with decreasing time, that is, vector 𝐲⁰_𝑛−1 of the start from 𝑛−1st reference orbit can be found for a given state vector 𝐲_𝑓𝑇 of the arrival to the final orbit, and so forth, until vector 𝐲_𝑖0 is found. The equations necessary to solve the problem considered here can be easily derived from the equations given in Section 3.

4.2. Transfer between Two Given Orbits

Let us assume that neither launch position in the initial orbit nor arrival position in the final orbit are given and should be determined in an optimal way during solving the transfer problem. This case may be solved using the suggested method in the following way.

A first guess for the launch position should be given somehow. This position defines the vector 𝐲_𝑖0 and in the first iteration of the calculation procedure described in Section 3.3 the final state vector 𝐲_𝑓𝑇 may be found. In the second iteration of the calculation procedure for this state vector a new value of the vector 𝐲_𝑖0 may be found as described in Section 4.1, and so forth, that is, odd iterations of the calculation procedure use the case described in Section 3 and even iterations use the case described in Section 4.1.

5. Partly Given Final Orbit

The suggested method can be also used in the case of partly given elements of the final orbit, that is, if vector 𝐪_𝑓 has dimension 𝑚<5. For instance, only energy of the final orbit (𝑚=1) or semimajor axis and eccentricity (𝑚=2) may be given. In this case the method described above is applied with 𝑚-dimensional vectors 𝐪_𝑖,𝐪_𝑓,𝐪_𝑗(𝑗=1,…,𝑛−1), 𝑚×6-dimensional matrix 𝐔_𝑗, and 𝑚-order matrices 𝐖_𝑗,𝐒_𝑗. Nongiven orbital elements are determined using transversality condition (3.17) and the respective conditions for vectors 𝐪_𝑗.

6. Constraint Imposed on the Thrust Direction

Here the case when a constraint is imposed on the thrust direction is considered. This constraint may be caused by specific features of the spacecraft design or of its attitude control system. Let us assume the constraint given by𝐁𝜶=𝟎, (6.1) where 𝐁=𝐁(𝐱,𝑡) is a matrix of dimension 1×3 (i.e., B is a row) or 2×3 (in this case the thrust direction is given and the problem is to find optimal thrust value). In this case, as is shown in [17, 18], the suggested optimization method is also applicable with matrices 𝐖_𝑗,𝐒_𝑗 and vector 𝜶_𝑗 from (41), (3.48), and (3.46) replaced by𝐖_𝑗=^𝑡𝑗_𝑡𝑗−1𝐐_𝑗𝐏𝐐^𝑡_𝑗𝑑𝑡,𝐒_𝑗=𝐒_𝑗(𝑡)=^𝑡_𝑡𝑗−1𝚿_𝑗𝑣𝐏𝚿^𝑡_𝑗𝑣𝑑𝑡,𝜶_𝑗=𝐏𝚿^𝑡_𝑣𝐔^𝑡_𝑗𝐖_𝑗⁻¹Δ𝐪_𝑗, (6.2) where third-order matrix𝐏=𝐈−𝐁^𝑡(𝐁𝐁^𝑡)⁻¹𝐁, (6.3) is a projector onto the constraining set given by (6.1). Note that rank of matrix P is less than 3. This is why matrices 𝐖_𝑗 given by (6.2) may be singular, while in the case of no constraint on the thrust direction matrices 𝐖_𝑗 given by (41) are nonsingular in any interval of integration [9]. As is shown in [18] nonsingularity of all matrices 𝐖_𝑗 given by (6.2) is a sufficient condition of feasibility of the transfer with constraint (6.1).

7. On Optimality of the Method

One of the necessary conditions of optimality of the thrust vector is constancy of the Hamiltonian in the whole time interval of the transfer [14]. It can be shown that the Hamiltonian of the linearized problem given by (3.5) is constant in the solution given by the suggested method, although the Hamiltonian of the original nonlinear problem calculated in the solution of the linear problem is not constant. Thus, the method described in the paper gives the solution which is not optimal in the original problem.

A comparison of the solutions to the problem of the orbital transfer in the original formulation and linearized one was performed. Transfer between coplanar circular orbits of radii 𝑎_𝑖 = 1 and 𝑎_𝑓 = 4 with gravitational parameter of the primary body 𝜇=1 was considered. Solution to the original problem was provided by Petukhov. He used the mathematical developed by his method for solving two-point boundary value problem [19]; optimal solution to the orbital transfer problem was obtained by Petukhov by means of variation of the arrival point in the final orbit. Results of the comparison are given in Table 1.

Table 1 

Comparison of the solutions to the original and linearized problems.

As is seen in Table 1, the difference between the values of the performance index for the optimal solutions to the original and linearized problems is small and practically does not depend on the transfer duration. Although in more complicated transfers the difference may be bigger.

8. Illustrative Examples

This section demonstrates potentialities of the suggested method by means of examples of transfers in the Earth’s sphere of influence. Orbits are given by the orbital elements𝐪=𝑟_𝜋,𝑟_𝛼,𝑖,Ω,𝜔, (8.1)

where 𝑟_𝜋,𝑟_𝛼 are radii of perigee and apogee in thousands of kilometers (Mm), i is the orbital inclination, 𝛀 is the longitude of the ascending node, 𝜔is the argument of the perigee. Angular elements are given in degrees. Direction of the thrust vector is given in the examples by angles between the projection of the thrust vector onto the instantaneous orbital plane (angle 𝜑) and between the thrust vector and orbital plane (angle Ψ). In all examples given below the optimal positions of departure from the initial orbit and arrival to the final orbit were obtained (i.e., the classical case of orbital transfer described in Section 4.2 was considered). The changes ΔV of the velocity by means of the thrust also are given in the examples. These changes were calculated using numerical integration of the magnitude of the thrust vector given by (3.46) or (6.2). The number of the time subintervals was taken 𝑛=5000 in all examples.

8.1. Transfer between Elliptical Orbits with High Mutual Inclination

Transfer between two orbits given by 𝐪_𝑖={7,30,50,80,−60},  𝐪_𝑓={40,80,80,−80,70} with 𝑇=400 hour is considered here. The transfer trajectory is shown in Figure 3 in projections onto the equator plane xy and the polar plane 𝑥𝑧. Dashed lines show node lines of the initial and final orbits. The jet acceleration value divided by 𝑔=9.8066 m/s² is shown in Figure 4; Figure 5 shows angles 𝜑, Ψ. Performance index and total ΔV for the transfer are 𝐽=44.42m²s⁻³, ΔV = 10.05 km/s.

Figure 3 

Transfer between two elliptic orbits with high mutual inclination.

Figure 4 

Jet acceleration for transfer between two elliptic orbits with a high mutual inclination.

Figure 5 

Angles of the thrust direction transfer between two elliptic orbits with a high mutual inclination.

8.2. Transfer to an Orbit with Given Perigee and Apogee Radii

Transfer in time 𝑇=400hrs to a partly given orbit, namely, to an orbit with given only perigee and apogee radii, is considered here. Optimal transfer is planar in this case. Only perigee and apogee radii of the initial orbit are specified, because the other initial orbital elements are not important and may be taken equal to zero. Elements of the initial and final orbits are taken as follows: 𝐪_𝑖={7,20,0,0,0},𝐪_𝑓={40,80,0,0,0}. Transfer orbit is shown in Figure 6; as is seen in this figure optimal is coaxial collocation of the final orbit with respect to the initial orbit. Figure 7 shows respective propulsion acceleration value divided by g, and angle φ is shown in Figure 8 (Ψ=0 in the planar transfer). Performance index and total ΔV are J = 3.01 m²s⁻³, ΔV = 2.82 km/s.

Figure 6 

Transfer trajectory to an orbit with given perigee and apogee radii.

Figure 7 

Jet acceleration for transfer to orbit with given perigee and apogee distances.

Figure 8 

Angle 𝜑 for transfer to orbit with given perigee and apogee distances.

8.3. Transfer to an Orbit with Given Energy

A transfer to an orbit with a given energy, namely, to a hyperbolic orbit given only by 𝐶₃ = 1 km²/s², is considered here. Initial orbit with elements 𝐪_𝑖={7,20,0,0,0} is taken, the transfer duration is 𝑇=1000hrs. The transfer trajectory (which is obviously planar in this case) and the jet acceleration are shown in Figures 7 and 8. The dot at the end of the spiral trajectory in Figure 7 marks the entry into the hyperbolic orbit. Corresponding performance index and total Δ𝑉 are 𝐽=4.85 m² s⁻³, Δ𝑉=4.93 km/s. The thrust direction is tangential in this case; this confirms nonoptimality of the solution to the linearized problem of the orbital transfer, because in the original nonlinear problem the optimal thrust is not tangential (see, e.g., [20]).

8.4. Constrained Thrust Direction

A transfer between the orbits given by the elements 𝐪_𝑖={7,20,0,0,0},𝐪_𝑓={50,100,60,60,60}, in a time 𝑇=400 hours, is considered here. The transfer trajectory is shown in Figure 9. The performance index and total ΔV are 𝐽=10.83 m²s⁻³, ΔV=5.39 km/s.

Figure 9 

Transfer trajectory to the hyperbolic orbit with given energy.

Now let us consider the following constraint on the thrust direction: the thrust is always orthogonal to the spacecraft position vector, that is, 𝐁=𝐫^𝑡 in (6.1). The projector (6.3) in this case is 𝐏=𝐈−𝐫𝐫^𝑡/𝑟². The transfer trajectory for the constrained thrust direction visually does not differ from the one for the unconstrained direction shown in Figure 9. Performance index and total ΔV in the case of the constrained thrust direction are 𝐽=11.44 m²s ⁻³, ΔV = 5.52 km/s.

Acceleration value versus time for the unconstrained and the constrained thrust direction is shown in Figure 10, the transfer trajectories are shown in Figure 11, the Jet acceleration for the transfers are shown in Figure 12 and the angles of the thrust direction for the transfers are shown in Figure 13.

Figure 10 

Jet acceleration for transfer to hyperbolic orbit with given energy.

(a)

(b)

Figure 11 

Transfer trajectories without (a) and with (b) constraint on the thrust direction.

(a)

(b)

Figure 12 

Jet acceleration for transfer without (a) and with (b) constraint on the thrust direction.

(a)

(b)

Figure 13 

Angles of the thrust direction for transfer without (a) and with (b) constraint on the thrust direction.

9. Conclusion

The mathematical method for calculation of low-thrust orbital transfers presented in this paper has two essential disadvantages, as follows.(1)The method is applicable to the power-limited thrust, while the existing thrusters are close to ones with constant or given exhaust velocity.(2)The original nonlinear problem is replaced in the method by a linearized problem solution to which it is not optimal for the original problem.

These disadvantages are compensated by simplicity of the method and its analytical form at each iteration, and also by the wide applicability of the method: it works well in the case of a big difference between the initial and final orbits, for a very high number of orbits around the attracting center; also it is applicable for different transfer types (such as point-to-orbit, orbit-to-point, and orbit-to-orbit transfers), in the cases of partly given final orbit and of a constraint imposed on the thrust direction. Despite the fact that the method is based on linearization of motion, any necessary accuracy of calculations may be reached by means of augmentation of the number 𝑛 of the reference orbits.

The suggested method may be used at early phases of the mission design when a high optimization accuracy is not needed and at the same time massive calculations are necessary for selection of a best mission scheme.

Nomenclature

r, v:	Position and velocity vectors
𝑡:	Current time
𝑡=0, 𝑡=𝑇:	Initial and final instants of the transfer
I:	Unit matrix
α:	Jet acceleration vector (thrust vector)
α:	|α|
𝜇:	Gravitational parameter of the attracting center
𝜑:	Angle between the projection of the thrust vector onto the orbital plane and the velocity vector
ψ:	Angle between the thrust vector and the orbital plane
Subscripts “0” and “𝑇”:	Values of the parameters at the time instants 0 and T (if another meaning of the “0” subscript is not stipulated), superscript “𝑡” denotes transposition.

Acknowledgment

The authors are grateful to the Brazilian São Paulo Research Foundation (FAPESP) for financial support of this study.