Mission Planning of GEO Active Debris Removal Based on Revolver Mode

(emission planning of active debris removal (ADR) of revolver mode on geosynchronous orbit (GEO) is studied in this paper. It is assumed that there are one service satellite, one space depot, and some pieces of space debris in the ADR mission. (e service satellite firstly rendezvouses with the debris and then releases the thruster deorbit kits (TDKs), which are carried with the satellite, to push the debris to the graveyard orbit. Space depot will provide replenishment for the service satellite. (e purpose of this mission planning is to optimize the ADR sequence of the service satellite, which represents the chronological order, in which the service satellite approaches different debris. In this paper, the mission cost will be stated firstly, and then a mathematical optimization model is proposed. ADR sequence and orbital transfer time are used as designed variables, whereas the fuel consumption in the whole mission is regarded as objective for optimizing, and a specific number of TDKs is also a new constraint. (en, two-level optimization is proposed to solve the mission planning problem, which is low-level for finding optimal transfer orbit using accelerated particle swarm optimization (APSO) algorithm and up-level for finding best mission sequence using immune genetic (IGA) algorithm. Numerical simulations are carried out to demonstrate the effectiveness of the model and the optimization method. Results show that TDK number influences the fuel consumption through impacting the replenishing frequency and TDK redundancy. To reduce fuel consumption, the TDK number should be optimized and designed with suitable replenishing frequency and minimum TDK redundancy.


Introduction
Due to the human space activities and active spacecraft increasing, the quantities of failed/decayed satellites and space debris grow at the same time and become the major source of the space pollutions [1]. According to the reports from Space Surveillance Net (SSN), more than 18,000 onorbit space targets have been detected and catalogued, and 90% of them are space debris [2]. Moreover, the geosynchronous orbit (GEO) is an important orbital resource but crowded by debris with different shape and size, which will threaten the safety of on-orbit active spacecraft. By 2016, there are 1,484 objects on GEO, and only 31.7% of them are active satellites. According to European Space Agency (ESA), the collision probability between the debris and satellites will reach 3.7% till 2030 if the effective debris removal measure is not taken [3]. In recent decades, the active debris removal (ADR) has been come out and regarded as an effective measure for this terrible situation. ere are some methods to realize ADR mission: space robot arm catching [4], space rope net catching [5], increasing the resistance [6], and the method of laser cleaning [7]. e ADR requires that one service satellite rendezvouses with different on-orbit debris and remove them [8,9]. Obviously, according to the background of ADR mission, it is a multispacecraft rendezvous problem scientifically. e mission planning method can be raised as an important research topic for this scientific problem, which could reduce the fuel consumption and also find the rational sequence for multispacecraft rendezvous in the ADR mission.
Research on mission planning in ADR task has been investigated for several years. Braun [10] analyzed different types of propulsive systems carried by service satellite and obtained the optimal mission sequence by exhaustive method. Considering the fuel and time constraints, Madakat [11] optimized the mission planning for debris removal on low earth orbit (LEO). Olympio [12] studied the optimal sequence for removing debris on sun-synchronous orbit for low-thrust orbit maneuver, and the results showed that the proposed method can reduce fuel consumption effectively under the influence of J2 perturbation. Yu [13,14] studied the debris removal sequence on GEO and LEO with multiple service satellites. e previous literatures studied the mission planning in different debris removal scenarios with the orbit at different heights. Mostly, they focused on the phoenix mode [15] for removing debris; that is, the service satellite dragged debris to graveyard orbit firstly and then maneuvered to next debris and continued this action sequence. In this paper, the revolver mode [15] is utilized to realize ADR, the mission scenario of which can be shown in Figure 1. e scenario includes one service satellite, one space depot, and some space debris on GEO. Different from phoenix mode, the revolver mode avoids the service satellite transfer between graveyard orbit and GEO-debris' orbit. e service satellite carries several thruster deorbiting kits (TDKs) [16] and attaches them on the surface of debris. en, the TDK will open the propeller and push the debris to the graveyard orbit. In this way, the service satellite only needs to maneuver among the target debris, and the fuel consumption can be effectively lessened.
Mission planning for revolver mode brings new problem: the TDK is a new constraint. If fuel or TDK is used, the ADR mission will fail. is paper quotes the concept of space depot in ADR mission to provide replenishment for service satellite. e orbital depot, which is used to store TDK and fuel, is defined as a system that enables refueling of spacecraft elements in space [17]. When depleted of the TDKs, service satellite will return to the space depot for replenishment. One of the most prominent counterarguments to the propellant depot is centered on the thermal management of cryogenic fuels, and numerous studies and experiments have been done to solve it [18]. Recently, the propellant depot in space has proved practicable. Since the launch in December 2018, NASA's Robotic Refueling Mission 3 (RRM3) demonstrated the first ever long-term storage of cryogenic fluid with zero boil off, having successfully stored cryogenic fluid for four months [19]. Under the circumstances of revolver mode, the number of TDKs and the sequence of the service's maneuver should be optimized obviously. e purpose of this paper is to find the optimal sequence in the revolver mode for ADR mission and analyze the influence of the number of TDKs on the mission planning.
So far, several reported works have focused on the optimization of mission planning. Some approaches are the nested loops, including the outer-loop to optimize the mission sequence and inner loop to find the optimal transfer orbit. Wall and Conway [20] used genetic algorithm (GA) to solve the double outer-loop and inner-loop problem for ADR mission, and it was proved to be efficient. Ross and D'Souza [21] introduced an HOC method to solve the space mission planning problem and proposed a formalism that can free mission planners to focus on high-level decision making by automating and optimizing the details of the inner loops. Du [22] used the MIGA-SQP algorithm to tackle the nested loops in optimal scheduling of multispacecraft refueling. Motivated by Du's work of nested loops optimization, this paper integrates accelerated particle swarm optimization (APSO) algorithm with immune genetic algorithm (IGA) for mission planning problem. e IGA algorithm and APSO algorithm are used to solve outer loop and inner loop, respectively. In the previous study, PSO was a stochastic population-based optimization method proposed by Kennedy and Eberhart [23]. PSO is a promising evolutionary algorithm with a strong global search capability and simple effective operators [24]. e particle position represents the value of designed variables and particle velocity means the changing rate of variables during optimization process. e optimal solution is found through the moving of particle. Moreover, APSO algorithm aims to fully utilize the global exploration of PSO and the rapid local exploitation of the gradient-based method to speed up the global search process [25]. e IGA is the genetic algorithm inspired by the biological immune system. It can effectively deal with mission planning problem [26]. Compared with genetic algorithms, IGA has good global searching capability through inhibition of repetitive antibodies. Meanwhile, good antibodies are preserved by memory cells to obtain good convergence [27,28].
In this paper, we aim to investigate the mission planning problem considering the constrain of fuel consumption and number of TDKs. Compared to existing results, the major contributions and differences of this study are summarized as follows: (1) the mission planning for ADR mission in revolver mode is studied, and the mathematical model considering the number constraint of TDKs is modelled; (2) two-level optimization is organized to solve mission planning problem, low-level for finding optimal transfer orbit using APSO and up-level for finding best removing sequence using IGA. is paper is organized as follows. e mission scenario is described firstly, and mission costs are stated in Section 2.
en, ADR optimization method is organized to seek the optimal solution in Section 3. Two-level model came out for optimization, including APSO algorithm for low-level optimization and IGA for the up-level. Following that, numerical simulations are illustrated in Section 4 to demonstrate the efficiency of the proposed solution method. Influence on fuel consumption of different quantities of TDKs carried by service satellite is discussed in this part as well. Finally, some appropriate conclusions are presented in Section 5.

Mission Scenario.
ere are three different types of spacecraft in revolver mode, including service satellite, space depot, and space debris. As shown in Figure 1, ignoring the process of launching, the whole process for ADR mission can be divided into four steps. Firstly, the service satellite maneuvers to the position of the first debris. en, the satellite extends the robotic arms to grab the debris and installs TDK on the surface of debris. Next, two robotic arms are retracted, one for grabbing space debris and the other for attaching TDK on debris. Finally, TDK pushes the debris to graveyard orbit and service satellite maneuvers to the next debris until all pieces of debris are removed. When the fuels of the service or the number of TDKs are insufficient, service satellite will move to space depot for replenishment. In the mission, space depot is located on fixed location in orbit and does not maneuver due to its huge mass and limited fuel.
Summing up, the background of revolver mode studied in this paper is that substantial pieces of debris existed on GEO and one service satellite running among them. e service satellite needs to rendezvous with each piece of debris by sequences and attach TDKs to remove them. e removing sequence and transfer orbit need to be planned and optimized to reduce the fuel consumption.

Two-Impulse Rendezvous Model.
e two-body model and Earth J2000 coordinate are utilized to depict the relative and absolute motion of service satellite around the Earth. It should be noted that the third-body perturbation and solar pressure have significant influence on the inclination and eccentricity of GEO target, respectively, the fact of which has been shown in relevant references [29,30]. In the analysis of two references [31,32], scholars pointed out that, in a short period (within one year), the impact of perturbation on GEO target is very small and can be ignored, and the optimal rendezvous order does not change. erefore, the impact of atmospheric drag can be neglected, while since the ADR task is set within 5 days, it is assumed that the influence of other perturbations is ignored, including third-body perturbation and solar radiation pressure. Hence, it may be reasonable to neglect perturbations for simplifying the orbital model. Furthermore, the complexity of the problem will increase if the uncertainty of debris is taken into consideration. erefore, the perturbations and position error of debris are ignored in this paper.
Two-impulse rendezvous maneuvering model is adopted. Given the orbital elements of the two targetsE 1 , E 2 and the orbital transfer time t, the transfer orbit between the two targets can be obtained by solving the loopy Lambert problem [33]. en, the impulse required for the orbital transfer is

Mission Costs.
In the revolver mode, fuel, time, and number of TDKs are the main costs. Suppose that the quantity of debris to be cleaned is n and the TDKs are labelled as numbers 1, 2, 3, . . . , n, respectively. DefineΔT i , ΔV i , andΔM i as the transfer time, impulse, and fuel consumption referring to i-th debris, respectively. During the process of ADR mission, two cases will occur: e first case is that fuels and quantities of TDKs are sufficient, and the service maneuvers to next debris directly. Fuel consumption for removing i-th debris can be calculated by the Tsiolkovsky formula: where M i is the initial mass of the service before removing i-th debris. I sp denotes specific impulse, and g is the gravitational constant. e second case is that fuel or TDKs are insufficient. Service satellite should be transferred to space depot for replenishment and then maneuver to the next debris. Fuel consumption is obtained by where ΔM s i−1 and ΔM i s represent the fuel consumption of transferring to depot and to i-th debris, respectively. e reference [16] gives full analysis of how much fuel is required to place an object with the TDK into the graveyard orbit and weight; therefore, the fuel consumption constraint of TDK meets the result of the argument in this paper.
Denote the time for debris removing and replenishing as T g and T f , respectively. e time cost can be given as where q ∈ 0, 1 { } is the binary decision parameter, while q � 1 represents that service satellite transfers to space depot for replenishment, and q � 0 represents satellite transfers to the next debris directly. Define k as the quantity of TDKs the service satellite carries. Suppose that each debris should be removed by using one TDK. erefore, the TDK is reduced by n after removing all the debris.
In summary, the fuel, time, and TDK consumption can be calculated as follows.

ADR Optimization Method
On the basis of the maneuvering strategy developed in Section 2, the mission planning problem for ADR can be treated as an optimization problem for the purpose of minimizing the total fuel consumption.

Design Variables.
For mission planning problem, it is necessary to determine the ADR sequence [34]. e design variables consist of removal order X, binary decision variable S, and time distribution T. e X and S are integers, and T is a continuous number, respectively. e variables are of the form For example, the removing order X � [2, 5, 4, 3, 1] denotes that 5 pieces of debris need to be cleaned, and service satellite removes the debris as the order 2 ⟶ 5 ⟶ 4 ⟶ 3 ⟶ 1. S � [0, 0, 1, 0, 1]means that service satellite transfers to space depot for replenishment before removing the 3-th and 5-th debris. In T, [Δt 1 , Δt 2 , . . . Δt m ]indicate the time period for each maneuver, and [tw 1 , tw 2 , . . . , tw m ]are the waiting time before orbit transfer.
In the ADR mission, the objective is to minimize the total fuel consumption, MinΔM.

Constraints
(1) Fuel Constraint. e service satellite cannot maneuver if the remaining fuel is insufficient, and the removing task will be interrupted. e constraint for fuel can be stated as where M i f is remaining fuel before removing the i-th debris.
(2) TDK Number Constraint. TDK is being consumed during the mission. erefore, the number of TDKs is required to be sufficient and can be stated as where k left,i is the number of left TDKs the service satellite carries before removing the i-th debris.
(3) Time Distribution Constraint. e ADR mission usually requires timeliness. e constraint is e integrated penalty function is designed to deal with the constraints of fuel and number of TDKs, which is constructed as where the penalty function φ(·) is defined as e penalty function is adopted as objective for optimization to solve the problem of satellite layout optimization designed by Chen [25]. a, b, c are the given constants. In this paper, a, b, c are 1, 2, 0 respectively. M i f − ΔM i+1 i < 0 or k left,i < 0 means that the constraint is not satisfied, and the penalty function G(X, S) will add to the objective function.

Two-Level Optimization Model.
Method for the ADR problem can be divided to the one-level and two-level optimization. e one-level optimization means that all designed variables are optimized directly. e two-level optimization includes the up-level to optimize the mission sequence and the low-level to find the optimal transfer orbit. As shown in Figure 2, low-level optimization is nested in the up-level one; in other words, the up-level objective function is calculated by optimizing the low-level one. e up-level planning optimizes the sequence, while the low-level planning optimizes the object-to-object transfer impulse. Zhou [35] has demonstrated that the results of two-level optimization are better than the ones of one-level optimization. e reason is that the solution space of the up-level optimization is much smaller than that of the one-level optimization. In low-level optimization, APSO is used, which integrates a local acceleration mechanism based on gradient search in iterations, to solve optimal transfer impulse, and it can get faster calculation and maintain global search capability. In up-level optimization, IGA is an improved algorithm that inherits the genetic algorithm and integrates immune algorithm into it. IGA has the ability to generate diverse antibodies, self-regulating mechanisms, and immune memory. ese mechanisms can get global search capability without falling into the local optimum. rough antigen memory, the search speed can be accelerated, and search capability of the genetic algorithm can be improved.
In this way, IGA is more suitable for sequence optimization search problem of ADR. e two-level optimization model can be given as

Low-Level Optimization.
e low-level optimization is to determine the optimal transfer time for service satellite, the variable of which is T � [Δt 1 , Δt 2 , . . . Δt m , tw 1 , tw 2 , . . . , tw m ]. e objective for optimization is to minimize fuel consumption, shown in equation (5).
In this article, the novel accelerated particle swarm optimization (APSO) algorithm is adopted in the low-level optimization [25]. It aims to fully utilize the global exploration of PSO and the rapid local exploitation of the gradient-based method to speed up the global search process. e traditional PSO algorithm has problems such as slow convergence speed, low computational efficiency, and inability to maintain the global search ability in high-dimensional problems.
is paper uses PSO as basic framework and integrates a gradient-based search in its iterations, and in this way, APSO algorithm is proposed. When the number of iterations reaches a certain number T acc , the current position of each particle is used as an initial point for the sequential quadratic programming (SQP) algorithm to exploit its local best and update its position until all the particles in the swarm have been accelerated. e process of APSO is shown in the right part of Figure 3. In this way, APSO can get faster calculation speed when optimizing object-to-object impulse.

Up-Level Optimization.
e immune genetic algorithm (IGA) is adopted for up-level optimization. Compared with traditional GA, IGA is an improved genetic algorithm based on biological immune mechanism. e objective function of the problem corresponds to the antigen of the invading organism, and the solution corresponds to antibody produced by the immune system. e algorithm has the following features: (1) diverse antibodies: in IGA algorithm, the ability of the immune system to generate a large number of antibodies to resist various antigens is simulated.
is mechanism is used to improve the global search ability of the genetic algorithm rather than falling into the local optimum. (2) Self-regulation mechanism: the immune system has a mechanism to maintain immune balance. By inhibiting and promoting antibodies, it can self-regulate to produce an appropriate amount of antibodies. Using this function, it can improve local search ability of genetic algorithm. (3) Immune memory cell function: some of the cells that produce affinity will be preserved as memory cells. For similar antigens that will invade in the future, the corresponding memory cells will be quickly stimulated to produce a large amount of antibodies. Using this memory function, it can speed up search speed and improve the overall search ability of genetic algorithm.
In IGA, antigen stands for objective function, and antibody represents designed variables. e affinity between the antigen and the antibody indicates the matching degree between objective function and optimal variables. Each generation selects excellent antibodies by their affinity as memory cells. e affinity can be calculated as follows: where fit and con represent the objective function and antibody concentration, respectively, and p is the evaluation parameter. From equation (13), the calculation for affinity is not only based on the objective function, but also based on the antibody concentration. e high concentration will cause a decrease in affinity. is operation will keep individual varieties and avoid premature convergence.
Code, select, crossover, and mutation are used in IGA. Each antibody is coded asC � [X, S]. e roulette wheel method is used in selection [36]. e one-point method is utilized in crossover. Select one position k in antibody and two different antibodies randomly. e mutation means, select two genes in variable X and exchange them, and then select one gene in variable S and change it to another value. e process of IGA is shown in the left part of Figure 3. In summary, the procedure of two-level optimization is summarized as follows: Step 1: set parameters (shown in Tables 1 and 2) to initialize APSO and IGA

Mathematical Problems in Engineering
Step 2: calculate the fitness function using APSO in lowlevel optimization Step 3: calculate the affinity based on the fitness value and antibody concentration Step 4: select the memory cells and operate the process of selection, crossover, and mutation in IGA of up-level optimization Step 5: repeat the Steps 2-5 until reaching maximum iterations e whole process of two-level optimization is shown in Figure 3.

Simulation and Analysis
Select 9 pieces of space debris on near GEO orbit. e orbital elements of debris are shown in Table 3. To avoid the singularity in classical orbital elements, the eccentricity is initialized as 0.001 as near-circle orbit. e semimajor axis and argument of perigee are initialized as 42166. 3    depot and service satellite are initially located in the near GEO orbit, where the semimajor axis, eccentricity, inclination, right ascension of ascending node (RAAN), argument of latitude, and true anomaly are 42166.3 km, 0.001, 0.001°, 0°, 0°, and 0°, respectively. e constants used in this simulation are I sp � 300 s and g � 9.8. Time period for ADR mission is t j � 1000 s. e refueling time is t f � 8000 s. e maximum time t max for each maneuver is 5 orbital periods. e parameters of APSO and IGA are provided in Tables 1 and 2.

Case 1.
In this case, the efficiency of the two-level optimization is tested by comparing with exhaustive search (ES) method, that is, the enumeration of possible solutions and the choice of the best one. e consumption time for ES grows exponentially with the instance size, so the serial numbers 1-6 of objects in Table 2 are selected to be removed. In that way, there are 2 6 · 6! � 46080 different sequences in total, which can be the search space of simulation. e mass of propellant and TDK is 400 kg and 50 kg, respectively, the consumption of which meets the result of the argument in reference [16]. Set the structural mass of service satellite as 400 kg. Each satellite carries 4 TDKs. erefore, the total mass of service satellite is 1000 kg. e fuel consumption for all sequence is shown in Figure 4. e extreme data with fuel consumption, which is greater than 2000 kg, is removed from figure. Each point in the figure stands for one possible sequence. According to the frequency of replenishment, the mission sequence can be divided into 7 categories, which, respectively, means that the service satellite is replenished 0-6 times. e minimum fuel cost is 443.22 kg, as the optimal mission sequence isX � [1, 3, 6, 5, 4, 2], S � [0, 0, 0, 0, 1, 0]. It takes about 49827.5 s for a 3.6 GHZ Intel Core personal computer with a memory of 16 GB to complete this calculation. Using the proposed two-level optimization method, the convergence curve is shown in Figure 5.
e algorithm converges after about 25 generations. e optimal mission sequence is X � [1,3,6,5,4,2], S � [0, 0, 0, 0, 1, 0]and the minimum fuel cost is 443.22 kg, which is consistent with result of ES method. However, in two-level optimization, the calculation time in the same computer is 4326.9 s. Clearly, these two methods (ES and two-level optimization method) can both find the optimal solution, but the time consumed by two-level optimization is reduced by 91.3%.
From Figure 4, the fuel cost of each category is concentrated in a certain range. For example, 84.65 percent (6616 of 7815 in total) of the data with twice replenishment distribute in the range of fuel consumption 600-1200 kg. Data analysis of the fuel consumption with different sequence category is shown in Figure 6. e data of zero replenishing frequency is removed, because the TDK constraint is not satisfied. In general, when the initial mass is the same, higher replenishing frequency means larger fuel consumption.

Case 2.
In the revolver mode for debris removing, the quantity of TDKs is an important part for the whole mission. e different TDK number constraints in mission planning are discussed in case 2. Keep structural mass and TDK mass the same as those in case 1. All nine objects in Table 3 are selected to be removed. e number of TDKs carried by service satellite is 2, 3, 4, 5, 6, 7, 8, and 9, respectively, and the total mass changes to 900 kg, 950 kg, 1000 kg, 1050 kg, 1100 kg, 1150 kg, 1200 kg, and 1250 kg with the change of the TDKs number. Table 4 addresses the optimized mission sequence with different number of TDKs. It is supposed that each debris should be removed by one attached TDK. e minimum fuel consumption for different TDK numbers is 1014.1 kg, 854.4 kg, 904.9 kg, 726.2 kg, 934.1 kg, 981.7 kg, 986.8 kg, and 1066.4 kg, respectively. As the quantity of total pieces of space debris is 9, for the TDK number with 2, 3, and 4, service satellite maneuvers to space depot more than once for TDK replenishment, because the TDKs the service satellite carries are not enough for removing all pieces of debris.
is situation increases the fuel consumption as shown in case 1, in which higher frequency of replenishment means larger fuel consumption. For the number of 5, service satellite needs to maneuver to depot once, and the fuel cost is small. When the TDK number increases to 6, 7, 8, and 9, the service satellite maneuvers to space depot twice due to the lack of propellant, and the fuel consumption is large. It can be analyzed that the TDK number influences the replenishing frequency through the constraints of TDKs number and fuel. Fewer TDKs the service satellite carries will cause high replenishing frequency due to constraint of TDK, and larger quantity will have high replenishing frequency as well, due to the constraint of fuel.
In this case, a new concept is introduced: TDK redundancy. e TDK redundancy indicates the total remaining quantity of TDK before replenishment and after the whole mission.
is index reflects the usage rate of TDK.  Obviously, the initial TDK number and replenishing frequency determine the TDK redundancy. As shown in Table 4, for the TDK numbers of 7, 8, and 9, when the replenishing frequency is the same, higher TDK number causes larger TDK redundancy and further leads to larger fuel cost. e reason is that, during the orbital transfer, service satellite will have higher mass (redundant TDKs).
is indicates that the lower replenishing frequency does not always mean fewer fuel cost, as the mass may be large in orbit transfer. For the TDK numbers as 2 and 3, the replenishing frequency is high, and the fuel consumption is large. For the TDK numbers as 4, 5, 7, 8, 9, 10, and 11, the replenishing frequency is the same, and the larger TDK redundancy leads to more fuel consumption.
To sum up, TDK is the main constraint of mission. e TDK number influences replenishing frequency and TDK redundancy and then impacts the fuel consumption. In general, high replenishing frequency means large fuel consumption. But in special conditions, the lowest replenishing frequency may not mean minimum fuel cost, because more TDKs carried by service satellite will lower the frequency of replenishment, while, in turn, they may cause higher mass in orbit transfer. When the replenishing frequency is the same, high TDK redundancy causes more fuel cost. e two indexes (replenishing frequency and TDK redundancy) almost determine the fuel cost in case the target debris is specific. erefore, for the revolver mode, the TDK number the service satellite carries should be designed   beforehand with suitable replenishing frequency and minimum TDK redundancy. is is, in fact, take full advantage of TDK carried by service satellite and reduce orbital maneuver with small satellite mass.

Conclusions
e ADR mission is an effective solution for GEO debris removal. In this paper, the mission planning for ADR mission of revolver mode is studied. A mathematical optimization model for revolver mode is proposed, considering TDK number, mission sequence, and orbital transfer time as constraints of mission. Two-level optimization is organized to solve this problem, including low-level for finding optimal transfer orbit using APSO and up-level for finding best removing sequence using IGA. Compared with ES method, the simulations demonstrate the efficiency of the two-level optimization.
Experimental results indicate that the proposed method could solve the mission planning problem in ADR successfully. TDK is the main constraint of mission. e number of TDKs carried by the service satellite influences the fuel consumption by impacting the replenishing frequency and TDK redundancy. More TDKs carried in service satellite will lower the frequency of replenishment while in turn causing large TDK redundancy and higher mass in orbit transfer. Fewer TDKs carried by service satellite will cause high replenishing frequency. When the initial mass is the same, higher frequency of replenishment means larger fuel consumption, and when the replenishing frequency is the same, higher TDK redundancy causes more fuel cost. To reduce fuel consumption, the quantities of TDKs carried by service satellite should be optimized and designed with suitable replenishing frequency and minimum TDK redundancy.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.