## 1 Introduction  
With the deployment of huge satellite groups, constellation networking has ushered in new challenges. In particular, orbit-control technology and orbit-maintenance technology have become the constraints of large-scale constellation networking control. Determining how to efficiently, economically, and stably control large-scale constellation initialization and constellation configuration maintenance in low orbit has become a challenging problem that needs to be focused on and solved in the field of large-scale constellation orbit control. Constellation orbit control involves aspects such as control technology, control methods, and optimal control. Among them, the selection of the orbit operation control strategy will directly affect the constellation networking efficiency and the constellation operation cost, and further affect the geometric integrity, structural stability, and service availability of the constellation system. Large-scale constellation networking in low orbit is mainly based on constellation service characteristics and satellite orbit characteristics, which are crucial for completing the establishment and maintenance of the overall constellation geometry. The main tasks of a constellation network are the deployment and maintenance of the constellation orbital plane and the satellite phase in the orbital plane. In terms of drag, [PERSON] et al. proposed a phase space method for accomplishing differential drag maneuvers, based on the relative phase drift rate and collision avoidance of the satellite constellation. In this method, orbit modification was the only means of collision avoidance of the satellite constellation [1]. [PERSON] et al. proposed a differential drag control method for a large fleet of propulsionless satellites deployed in the same orbit and used the simulated annealing algorithm to minimize the required phasing time given the available control authority [2]. [PERSON] et al. proposed a station-keeping strategy to compensate for perturbations due to the atmospheric drag, in which the errors in the orbital measurements and control as well as in the estimation of the semimajor axis decay rate were taken into account [3]. [PERSON] proposed a new linear quadratic controller for satellite constellations, which helped to increase the formation-keeping accuracy for the low-altitude satellite constellations affected by atmospheric drag [4]. In terms of low thrust, [PERSON] et al. applied Hamiltonian theory to transform the low-thrust station-keeping optimization problem into a two-point boundary value problem, which was used for the station-keeping strategy of satellites with low thrust [5]. [PERSON] et al. established a long-term perturbation model for relative motion based on a Hamiltonian model. By analyzing the long-term evolution process of the BDS medium Earth orbit (MEO) constellation over 10 years, the validity of the proposed analytical perturbation compensation calculation method was verified [6]. In terms of perturbations, [PERSON] et al. employed perturbation analysis to formulate the relative perturbation motion equations for constellation orbits, elucidating the orbital perturbation principles and the overall drift characteristics of relative perturbation motion. In addition, they established the governing equations for both absolute and relative motion in the context of constellation orbital perturbation compensation control [7; 8]. [PERSON] et al. conducted an in-depth analysis of the orbital evolution process influenced by Earth's oblateness perturbation and satellite orbit deviations. Their findings indicated that Earth's oblateness perturbation resulted in the general drift of the constellation structure, while satellite orbit deviations played a critical role in determining the stability of the constellation geometry [9; 10]. [PERSON] et al. accounted for three-body and \(J_{2}\) perturbations based on the integral over time of the undesired perturbation forces.