A constellation optimization design method for providing navigation services for earth-moon space

Abstract

The application discloses a constellation optimization design method for providing navigation services for the earth-moon space, and belongs to the field of space constellation optimization. The application is realized by the following method: a nominal orbit of each orbit type is designed under the CRTBP model to generate a reference working orbit library; a plurality of performance index models are established, and then comprehensive scoring is carried out based on the weighted average method; the earth-moon navigation constellation scheme optimization algorithm takes the number of satellites in the constellation as the only input, and the maximum number of satellite deployment under each orbit type is constrained; the simulation start and end time is given, and the weight distribution of each index is carried out; then, the dynamic programming method is used to generate all constellation combination schemes; the index calculation and comprehensive evaluation sorting are carried out for each scheme; and finally, the optimal scheme is obtained. The orbit number combination under the optimal scheme is taken as the optimization variable, and the continuous amplitude optimization problem is simplified into a discrete number integer programming problem. After the orbit amplitude combination is determined, the orbit initial phase is further optimized to obtain the optimal constellation configuration.

Description

Technical Field

[0001] This invention relates to a constellation optimization design method for providing navigation services in the Earth-Moon space, applicable to universal constellation optimization design for providing navigation services to users on the ground, lunar surface, and Earth-Moon space, and belongs to the field of aerospace constellation optimization technology. Background Technology

[0002] With the increasing number of deep space exploration missions, especially the growing demand for lunar exploration and development, reliable and real-time communication and navigation services are crucial for the precise execution of space missions, whether manned lunar landings, unmanned lunar surface missions, or deep space exploration missions. Existing GNSS navigation satellite systems are insufficient to provide comprehensive navigation and communication services for lunisolar space missions; therefore, it is necessary to establish a lunisolar navigation constellation with full lunisolar coverage. Constellation configuration optimization design is key to lunisolar navigation constellation design. Traditional methods are often applied to lunisolar constellation design problems with limited mission scenarios and specific orbit types, often lacking universal applicability and difficult to widely apply to lunisolar constellation design problems with multiple service targets and multiple orbit types. Therefore, it is necessary to propose a more universal lunisolar constellation design method that is applicable to multiple service targets and consists of multiple orbits.

[0003] Currently, among the design methods for Earth-Moon constellations, prior art 1, based on numerical extension, has designed constellations that take the periodic orbit of the Earth-Moon system's L1 point as the working orbit. However, this method is not applicable to Earth-Moon constellation designs that include other orbit types.

[0004] Prior technology 2 uses a hybrid integer genetic algorithm based on orbital stability and synchronization constraints to solve the lunar communication constellation design problem for LPOs, but it cannot be applied to constellation design problems for other service objects, and the constellation evaluation index considered is relatively simple.

[0005] Prior technology 3 uses the hypervolume ant colony optimization algorithm to solve the multi-objective optimization problem to solve the design problem of lunar communication and navigation constellations in frozen lunar orbits. However, this method is not applicable to more general constellation design problems and it is difficult to guarantee the convergence of other types of constellation design problems.

[0006] In summary, existing technologies use relatively few orbit types when designing constellations, do not comprehensively consider multiple orbit types and their combinations to analyze the navigation performance of target areas, and often only use coverage as an evaluation indicator, resulting in a relatively simplistic evaluation system. Furthermore, the aforementioned constellation design methods are highly targeted and only effective for constellation design problems in specific scenarios. The feasibility of extending these methods to constellation design in different scenarios is unknown, and they lack universality. Summary of the Invention

[0007] The purpose of this invention is to address the problems of limited orbit types, singular performance evaluation indicators, and poor universality of constellation design methods in existing technologies. It provides a constellation optimization design method for providing navigation services in the Earth-Moon space. This method first generates 17 major three-body orbit families in the Earth-Moon space and reference working orbits for each family based on three-body orbit design theory, thus creating a corresponding constellation working orbit library and reference working orbit library. The constellation service target area is defined and meshed using cubic or spherical space methods. Using Earth-Moon shadows, coverage, and navigation accuracy factors as constellation design constraints, calculation models for each constraint indicator are constructed, and statistical processing is performed on each indicator under simulation scenarios. A comprehensive constellation evaluation model is established based on weighted averages. The constellation scheme is optimized using the reference working orbit library combined with dynamic programming and ergonomic methods according to the comprehensive evaluation model. After the constellation scheme is determined, a two-layer amplitude-phase optimization strategy based on the working orbit library is further adopted, and heuristic optimization methods such as genetic algorithms are used to solve for the constellation with the best performance in the current scenario. This method provides an effective tool for the construction of Earth-Moon space navigation systems, and the optimization results provide a reliable reference for the deployment of Earth-Moon space constellations.

[0008] The objective of this invention is achieved through the following technical solution.

[0009] This invention discloses a constellation optimization design method for providing navigation services in the Earth-Moon space. It is a method for designing Earth-Moon constellations based on the generation of a working orbit library and multi-index modeling. The method involves two processes: constellation scheme determination and constellation optimization, resulting in the optimal constellation for a specific objective. First, a working orbit library is established, containing 17 types of orbits, including DROs and periodic orbits near five translation points. Nominal orbits for each orbit type are designed under the CRTBP model, generating a reference working orbit library. Multiple performance index models are established, including Earth-Moon shadow, coverage, and DOP value performance indicators, and a comprehensive score is performed based on a weighted average method. The Earth-Moon navigation constellation scheme optimization algorithm uses the number of satellites within the constellation as the sole input, constrains the maximum number of satellites deployed under each orbit type, provides a simulation start and end time, and assigns weights to each index. Then, dynamic programming is used to generate all constellation combination schemes, and index calculations and comprehensive evaluations are performed for each scheme to obtain the optimal scheme. Based on the optimal scheme, the orbital number combinations under the selected orbital family of the optimal scheme are used as optimization variables, simplifying the continuous amplitude optimization problem into a discrete integer numbering programming problem. Here, the orbit number is equivalent to the amplitude. After determining the combination of orbital amplitudes, the initial phase of the orbit is further optimized. This invention realizes the design of a Moon-Earth constellation for a customized target object, comprehensively considering multiple constellation performance indicators, maximizing constellation performance, and providing the optimal constellation configuration.

[0010] This invention discloses a constellation optimization design method for providing navigation services in the Earth-Moon space, comprising the following steps:

[0011] Step 1: Generate Constellation Working Orbit Library B and Reference Working Orbit Library B r .

[0012] The dynamic model of the spacecraft under CRTBP is as follows:

[0013]

[0014] in

[0015] Based on the circularly restricted three-body theory, N is generated through three-body orbit design and numerical extension methods. t Orbital families of a type. The set of all orbits under each orbital type is the orbital family of that type, and the set of orbital families formed by all orbital types is the working orbital library B of the constellation.

[0016] The orbital families under each type contain C1, C2, ..., C F (F = 1, 2, ..., N) t There are ) orbits, where F represents the index of the orbit type. The orbit families under each orbit type are sorted and numbered according to the magnitude of the amplitude along the X-axis in the Earth-Moon rotation system. From the orbit family of the Fth orbit type, take the number C. F The track with a length of 2 is called the reference working track for the current track type. The collection of reference working tracks for all track types is called the reference working track library B. r .

[0017] Step 2: Determine the target area for navigation services provided by the Earth-Moon constellation and perform gridding to obtain the set of spatial points N. tar .

[0018] When the target region is a spatial target region, it is defined using a cubic space approach, with the cubic space center coordinates P given. c0 =(x c0 ,y c0 ,z c0 The set of coordinates P of the target point in cubic space:

[0019] P c ={(x c ,y c ,z c )|x c =X0+n x ·dX+x c0 ,y c =Y0+n y ·dY+yc0 ,z c =Z0+n z ·dZ+z c0} (2)

[0020] in

[0021] This indicates rounding down. Cubic space is defined by giving coordinate ranges parallel to the three axes in a lunar rotation system. The X-axis defines the starting point X0 and the ending point X... f The Y direction is defined as starting point Y0 and ending point Y0. f The Z direction is defined as the starting point Z0 and the ending point Z0. f The spacing in the X direction is dX, the spacing in the Y direction is dY, and the spacing in the Z direction is dZ.

[0022] Based on equation (2), equidistant division is performed to achieve cubic spatial gridding, then we have

[0023] N tar =P c (3)

[0024] When the target area is a ground target area or a lunar target area, a spherical space method is used to define the ground target area and the lunar target area, generating a set of coordinates P of the target points in the spherical space. s

[0025] P s ={(α,λ,R)α=α0+n α ·dα,λ=λ0+n λ ·dλ} (4)

[0026] in

[0027] Spherical space is defined within the Earth-Moon rotation system and described by longitude, latitude, and spherical radius. Longitude is defined with the X-axis of the Earth-Moon rotation system as the starting point, with counterclockwise rotation being positive, and the longitude range α∈[0,360]°. Latitude is defined with the lunar orbit plane as the starting point, with the North Pole as positive, and the latitude range λ∈[-90,90]°. When defining the target region of spherical space, the spherical radius R is given, the starting point of longitude is α0, and the ending point of longitude is α... f The latitude starting point is λ0, and the latitude ending point is λ. f The longitude interval is dα, and the latitude interval is dλ.

[0028] Based on equation (4), the spherical space is divided at equal intervals to achieve spherical spatial gridding.

[0029] To simplify calculations, the spherical coordinates are transformed to a rectangular coordinate system. Meanwhile, the coordinates of the sphere's center are given as P. s0 =(x s0 ,ys0 ,z s0 The set of coordinates of the target points in any spherical space under a Cartesian coordinate system is updated as follows:

[0030] P s ′={(x s ,y s ,z s )|x s =Rcosλcosα+x s0 ,y s =Rcosλsinα+y s0 ,z s =Rsinλ+z s0} (5)

[0031] Then there is

[0032] N tar =P s ′ (6)

[0033] Step 3: Establish the criteria for determining the Earth-Moon shadow and coverage, establish a calculation model for the navigation accuracy factor, and obtain the Earth-Moon shadow, coverage, and navigation accuracy factor based on the calculation model for the navigation accuracy factor.

[0034] Step 3.1: Establish criteria for determining the Earth-Moon shadow.

[0035] Let r b r s and r are the position vectors of the obscuring celestial body, the Sun, and the satellite, respectively. s =r s -r, s = r b -r represents the position vector of the sun and the satellite relative to the celestial body being obscured. R s R and R are the radii of the Sun and the occluding object, respectively. The shadow axis is the straight line connecting the Sun's center of mass and the occluding object's center of mass. The shadow determination angle is the angle between the line connecting the shadow cone vertex V1 and the satellite's position and the shadow axis. The opening direction of the angle is consistent with the penumbra angle f1, denoted as θ. s1 is the geocentric vector pointing from the Earth's center to the shadow cone vertex V1. The above angles and vectors are determined by the following formula:

[0036]

[0037] The penumbra cone angle f1 refers to the semi-cone angle of the penumbra cone with vertex V1. f1 is expressed as follows:

[0038]

[0039] The shadow determination criteria are as follows:

[0040]

[0041] When equation (9) is satisfied, the shaded state e = 1; otherwise, e = 0.

[0042] Step 3.2: Establish coverage determination criteria.

[0043] Geometric coverage ensures that the line connecting the satellite and the target is not obstructed by celestial bodies. Based on the positions of the celestial bodies and the target relative to the satellite, the visible area is divided into Zone I and Zone II. The angle δ between the celestial body, satellite, and target is calculated.

[0044]

[0045] Where, r sat r represents the position coordinates of the constellation satellites. tar Let r be the position coordinates of the target object. c The position vector of the celestial body being obscured.

[0046] For geometrically visible objects r that satisfy the condition of equation (10) tar Further judgment is made. Based on the cone-shaped region formed by the tangency of the satellite's line of sight with the celestial body, region II is divided into three sub-regions: II-1, II-2, and II-3. The value of r is determined according to equation (11). tar Location:

[0047]

[0048] Where d is the distance from the celestial body to the line connecting the observation position and the target position. The radius of the celestial body is R. c .

[0049] When r tar When located in zones II-2 and II-3, further determine r tar Visibility:

[0050]

[0051] Among them, the satellite-target-celestial body angle δ′. Based on the above judgment criteria, when the satellite and the target are visible, the coverage state c = 1; otherwise, c = 0.

[0052] Step 3.3: Establish a navigation accuracy factor calculation model.

[0053] When calculating the navigation accuracy factor (DOP), the number of satellites participating in the positioning solution at a certain moment in the constellation is n. t The position vector r of the target object d =(x d ,y d ,z d The position vector of the i-th satellite is... Calculate the pseudorange measurement residual matrix A:

[0054]

[0055] in i = 1, 2, ..., n t The covariance matrix Q is calculated using the least squares method:

[0056]

[0057] σ represents covariance.

[0058] The calculation model for the navigation accuracy factor is as follows:

[0059]

[0060] The Earth-Moon shadow, coverage, and navigation accuracy factor are obtained based on the calculation model of the navigation accuracy factor.

[0061] Step 4: Given a time series T, perform statistical processing on the results of Earth-Moon shadow, coverage, and navigation accuracy factors obtained in Step 3, and establish a calculation model for the comprehensive constellation evaluation index using the weighted average method.

[0062] Step 4.1, the time series T is represented as:

[0063] T = [t1 t2 ... t] M (16)

[0064] M represents the number of moments in the time series.

[0065] Based on the shadow situation determined in step 3.1, the n×M dimensional shadow matrix E is represented as follows:

[0066]

[0067] Where e ij Let represent the shadow state of the i-th satellite at time j. Let n represent the number of satellites in the constellation. The shadow vector E of the constellation is obtained by summing the rows of the shadow matrix. c :

[0068]

[0069] Then with E c The maximum value in the shaded index S e ,Right now

[0070] S e =max[E c (19)

[0071] Step 4.2: Based on the coverage determined in Step 3.2, the coverage matrix C is represented as:

[0072]

[0073] N is the set from step two. tar The number of elements.

[0074] N c To satisfy the requirement of having enough spatial points for coverage multiplicity, where Cover represents the coverage multiplicity threshold, we have:

[0075]

[0076] Therefore, the coverage ratio vector of the constellation is obtained as follows:

[0077]

[0078] With C r The average value is the coverage index S c ,Right now

[0079]

[0080] Step 4.3: Based on the DOP value results calculated in Step 3.3, the DOP value matrix D is represented as follows:

[0081]

[0082] N d To satisfy the requirement of having enough spatial points to cover the multiplicity of DOP, where DOPT represents the DOP threshold, we have:

[0083]

[0084] The average of the spatial proportions of DOP values ​​within the target area that satisfy DOPT over time is the DOP value index S. d

[0085]

[0086] Step 4.4, based on S in step 4.1 e S in step 4.2 c And S in step 4.3 d A calculation model for the comprehensive evaluation index of constellations is established using the weighted average method.

[0087] Given a weight vector W = [w1 w2 w3], where the internal elements represent shadow weight, overlay weight, and DOP weight respectively. Define a performance vector S = [S e S c S d The comprehensive index of the constellation is

[0088]

[0089] Step 5: Use dynamic programming algorithm to generate all constellation combinations I that satisfy the basic constraints. c Based on the reference work track library B in step one r The set of spatial points N in step two tar Then, using the comprehensive index from step four, all constellation combinations are evaluated, and the constellation scheme with the highest evaluation score is selected as the optimal constellation scheme Cons0.

[0090] Generate all constellation combinations I that satisfy the basic constraints using a dynamic programming algorithm. c , represented as

[0091]

[0092] a lF This indicates the number of satellites deployed for orbit type F in the l-th constellation combination.

[0093] to I c Evaluate by row-by-row, and denote the l-th row as A. l =[a l1 a l2 ...a lNt ] represents the l-th constellation scheme, obtaining A l The non-zero elements in the vector form a vector in accordance with From B r Obtain the orbit of the l-th constellation and perform orbital integration over time series T. Let N... tar To achieve the goal, the l-th constellation is evaluated using equation (27) in step four to obtain S. l .

[0094] Complete I c After traversing row by row, we obtain the pair N. tar Cons0 is the optimal constellation scheme for providing navigation services.

[0095]

[0096] b represents I c The row number where the highest score in the constellation assessment is located.

[0097] Step Six: Construct amplitude optimization and phase optimization problems based on Cons0 obtained in Step Five. Based on the constellation orbit library B generated in Step One, and combined with the spatial point set N given in Step Two... tar Using the comprehensive index from step four as the optimization objective, the genetic algorithm is used to solve the amplitude optimization problem and the phase optimization problem in sequence, thereby obtaining the optimal amplitude constellation Cons1 and the optimal phase constellation Cons2, thus realizing the optimization of the navigation constellation for the custom target region in Earth-Moon space.

[0098] Construct and solve the amplitude optimization problem PRO1:

[0099] Obtain the index vector of non-zero elements in Cons0 M0 represents the number of non-zero elements in Cons0, v o This represents the index position of the o-th non-zero element in Cons0. The v-th orbital type is retrieved from the constellation orbital library B, and the v-th type orbital library contains C orbitals. v The orbital numbers are 1 to C. v Let the index be f. v f v ∈[1,C v And f v ∈N + .

[0100] Therefore, the optimization variable is

[0101]

[0102] Upper and lower bound constraints are

[0103]

[0104] The optimization goal is

[0105]

[0106] Construct and solve the phase optimization problem PRO2:

[0107] The phase optimization process uses the initial phase angle combination of each orbit within the constellation as the optimization variable, therefore we have

[0108]

[0109] Upper and lower bound constraints are

[0110]

[0111] The optimization objective is the same as equation (32).

[0112] The optimal constellation Cons1 for amplitude is obtained by solving PRO1 using a genetic algorithm. Based on Cons1, the optimal constellation Cons2 for phase is obtained by solving PRO2 using a genetic algorithm. This completes the optimization design of the navigation constellation for the custom target region in Earth-Moon space.

[0113] It also includes step seven, which is the optimal navigation constellation Cons2 for the target area in the Earth-Moon space under the constellation size constraint obtained in step six. This constellation can provide high-precision and high-efficiency positioning and navigation services to spacecraft users in this area, thereby improving the mission execution efficiency of spacecraft users in the Earth-Moon space.

[0114] Beneficial effects:

[0115] 1. The present invention discloses a constellation optimization design method for providing navigation services in the Earth-Moon space. By establishing a target area definition model, the service area of ​​the Earth-Moon navigation constellation is customized, and the navigation constellation is optimized based on this model, thereby enabling the method to have a wide range of applications.

[0116] 2. The present invention discloses a constellation optimization design method for providing navigation services in the Earth-Moon space. Based on the weighted average method, a comprehensive evaluation index calculation model for the constellation is established considering three performance indicators: Earth-Moon shadow, coverage, and accuracy factor. Using this model to evaluate satellite resource allocation can improve the reliability of the constellation system.

[0117] 3. The present invention discloses a constellation optimization design method for providing navigation services in the Earth-Moon space. It constructs a constellation optimal scheme generation algorithm based on dynamic programming. The algorithm takes constellation size as the main constraint and the maximum number of satellites deployed under each orbit type as the secondary constraint. It can handle not only single-orbit single-satellite deployment schemes, but also single-orbit multi-satellite deployment schemes, and has broad application potential.

[0118] 4. The constellation optimization design method for providing navigation services in the Earth-Moon space disclosed in this invention combines a comprehensive evaluation model with a global optimization algorithm. It uses a working orbit library to optimize orbit amplitude through integer programming, simplifying the continuous amplitude optimization problem into a discrete numbered integer programming problem. This avoids orbit design and correction during the optimization process, significantly improves optimization efficiency, and realizes the optimal design of navigation constellations for custom target areas in the Earth-Moon space. Attached Figure Description

[0119] Figure 1 This is an overall flowchart of an implementation case of a constellation optimization design method for providing navigation services in the Earth-Moon space disclosed in this invention.

[0120] Figure 2 These are schematic diagrams of 17 reference working tracks for embodiments of the present invention.

[0121] Figure 3 It is the cubic space target point in the embodiment of this invention.

[0122] Figure 4 This is a schematic diagram of the orbit of the optimal constellation scheme in an embodiment of the present invention.

[0123] Figure 5 This is the curve showing the spatial proportion of the optimal constellation whose DOP value meets the threshold in the embodiment of this invention.

[0124] Figure 6This is the spatial proportion change curve of the coverage multiple of the optimal constellation in the embodiment of the present invention, which meets the threshold. Detailed Implementation

[0125] To better illustrate the purpose and advantages of the present invention, the invention will be further described below in conjunction with the accompanying drawings and examples.

[0126] Example 1:

[0127] Centered on the Earth's center, a target cubic space of 500,000 kilometers is defined as the coverage area within a rotating frame. The Earth-Moon constellation is selected to consist of 4 satellites, with a maximum deployment of 1 satellite for each orbital type. The simulation period is from January 1, 2030, 00:00:00 to February 1, 2030, 00:00:00, a total duration of 31 days. This example uses the method proposed in this invention to present the constellation scheme with the optimal overall performance under the above scenario.

[0128] like Figure 1 As shown in the figure, the constellation optimization design method for providing navigation services in the Earth-Moon space disclosed in this embodiment has the following specific implementation steps:

[0129] Step 1: Generate Constellation Working Orbit Library B and Reference Working Orbit Library B r .

[0130] The dynamic model of the spacecraft under CRTBP is as follows:

[0131]

[0132] in

[0133] Based on the circularly restricted three-body theory, N is generated through three-body orbit design and numerical extension methods. t Orbital families of a type. The set of all orbits under each orbital type is the orbital family of that type. The set of orbital families formed by all orbital types in the table below is the working orbital library B of the constellation.

[0134] The orbital families under each type contain C1, C2, ..., C F (F = 1, 2, ..., N) t There are ) orbits, where F represents the index of the orbit type. The orbit families under each orbit type are sorted and numbered according to the magnitude of the amplitude along the X-axis in the Earth-Moon rotation system. From the orbit family of the Fth orbit type, take the number C. F The track with a length of 2 is called the reference working track for the current track type. The collection of reference working tracks for all track types is called the reference working track library B. r .

[0135] In this embodiment, 17 orbit types are considered as candidate orbits for constellation operation, namely N t =17. The number of tracks included in each type of track family in the working track library B is shown in Table 1.

[0136] Table 1. Number of tracks included in each track family of the working track library

[0137]

[0138]

[0139] Take the number C from each orbital family. F / 2 track composition reference working track library B r As shown in Table 2. See the schematic diagram of the reference working track. Figure 2 .

[0140] Table 2 Reference Work Track Library

[0141]

[0142] Step 2: Determine the target area for navigation services provided by the Earth-Moon constellation and perform gridding to obtain the set of spatial points N. tar .

[0143] When the target region is a spatial target region, it is defined using a cubic space approach, with the cubic space center coordinates P given. c0 =(x c0 ,y c0 ,z c0 The set of coordinates P of the target point in cubic space:

[0144] P c ={(x c ,y c ,z c )x c =X0+n x ·dX+x c0 ,y c =Y0+n y ·dY+y c0 ,z c =Z0+n z ·dZ+z c0} (36)

[0145] in

[0146] This indicates rounding down. Cubic space is defined by giving coordinate ranges parallel to the three axes in a lunar rotation system. The X-axis defines the starting point X0 and the ending point X... fThe Y direction is defined as starting point Y0 and ending point Y0. f The Z direction is defined as the starting point Z0 and the ending point Z0. f The spacing in the X direction is dX, the spacing in the Y direction is dY, and the spacing in the Z direction is dZ.

[0147] Based on equation (2), equidistant division is performed to achieve cubic spatial gridding, then we have

[0148] N tar =P c (37)

[0149] This embodiment uses a spatial region as the coverage target, which is defined as P centered on the Earth's center. c0 = (0,0,0), a cubic space region of 100×100×100 million kilometers in the Earth-Moon rotating system, i.e., X0 = -50, X f =50, Y0=-50, Y f =50, Z0=-50, Z f =50. Given a three-axis interval of 50,000 kilometers, i.e. dX=5, dY=5, dZ=5. The spatial point set generated using equation (36) has a total of 9261 elements. The distribution of the cubic spatial target points is shown in [reference needed]. Figure 3 .

[0150] Step 3: Establish Earth-Moon shadow and coverage determination criteria, and establish a navigation accuracy factor calculation model.

[0151] Step 3.1: Establish criteria for determining the Earth-Moon shadow;

[0152] Let r b r s and r are the position vectors of the obscuring celestial body, the Sun, and the satellite, respectively. s =r s -r, s = r b -r represents the position vector of the sun and the satellite relative to the celestial body being obscured. R s R and R are the radii of the Sun and the occluding object, respectively. The shadow axis is the straight line connecting the Sun's center of mass and the occluding object's center of mass. The shadow determination angle is the angle between the line connecting the shadow cone vertex V1 and the satellite's position and the shadow axis. The opening direction of the angle is consistent with the penumbra angle f1, denoted as θ. s1 is the geocentric vector pointing from the Earth's center to the shadow cone vertex V1. The above angles and vectors are determined by the following formula:

[0153]

[0154] The penumbra cone angle f1 refers to the semi-cone angle of the penumbra cone with vertex V1. f1 is expressed as follows:

[0155]

[0156] The shadow determination criteria are as follows:

[0157]

[0158] When equation (9) is satisfied, the shaded state e = 1; otherwise, e = 0.

[0159] Step 3.2: Establish coverage determination criteria;

[0160] Geometric coverage ensures that the line connecting the satellite and the target is not obstructed by celestial bodies. Based on the positions of the celestial bodies and the target relative to the satellite, the visible area is divided into Zone I and Zone II. The angle δ between the celestial body, satellite, and target is calculated.

[0161]

[0162] Where, r sat r represents the position coordinates of the constellation satellites. tar Let r be the position coordinates of the target object. c The position vector of the celestial body being obscured.

[0163] For geometrically visible objects r that satisfy the condition of equation (10) tar Further judgment is made. Based on the cone-shaped region formed by the tangency of the satellite's line of sight with the celestial body, region II is divided into three sub-regions: II-1, II-2, and II-3. The value of r is determined according to equation (11). tar Location:

[0164]

[0165] Where d is the distance from the celestial body to the line connecting the observation position and the target position. The radius of the celestial body is R. c .

[0166] When r tar When located in zones II-2 and II-3, further determine r tar Visibility:

[0167]

[0168] Among them, the satellite-target-celestial body angle δ′. Based on the above judgment criteria, when the satellite and the target are visible, the coverage state c = 1; otherwise, c = 0.

[0169] Step 3.3: Establish a navigation accuracy factor calculation model;

[0170] When calculating the navigation accuracy factor (DOP), the number of satellites participating in the positioning solution at a certain moment in the constellation is n. t The position vector r of the target object d =(x d ,y d ,zd The position vector of the i-th satellite is... Calculate the pseudorange measurement residual matrix A:

[0171]

[0172] in This embodiment is designed for a four-star constellation, therefore n t =4. The covariance matrix Q is calculated using the least squares method:

[0173]

[0174] σ represents covariance.

[0175] The calculation model for the navigation accuracy factor is as follows:

[0176]

[0177] Step 4: Given a time series T, perform statistical processing on the results of Earth-Moon shadow, coverage, and navigation accuracy factors obtained in Step 3, and establish a calculation model for the comprehensive constellation evaluation index using the weighted average method.

[0178] The time series T is represented as:

[0179] T = [t1 t2 ... t] M (47)

[0180] M represents the number of moments in the time series.

[0181] In this embodiment, the simulation start and end times are first converted into normalized time, and a time series T is generated with a normalized time step of 0.03. T has a total of 238 time points, i.e., M = 238.

[0182] Step 4.1: Based on the shadow situation determined in Step 3.1, the n×M dimensional shadow matrix E is obtained as follows:

[0183]

[0184] Where e ij This represents the shadow state of the i-th satellite at time j. n represents the number of satellites in the constellation. In this example, it is a four-star constellation, so n = 4. The shadow vector E of the constellation is obtained by summing the rows of the shadow matrix. c :

[0185]

[0186] Then with E c The maximum value in the shaded index S e ,Right now

[0187] Se =max[E c (50)

[0188] Step 4.2: Based on the coverage determined in Step 3.2, the coverage matrix C is represented as:

[0189]

[0190] N is the set from step two. tar The number of elements. In this embodiment, according to step two, N = 9261.

[0191] N c To satisfy the requirement of having enough spatial points for coverage multiplicity, where Cover represents the coverage multiplicity threshold, we have:

[0192]

[0193] In this embodiment, Cover = 4.

[0194] Therefore, the coverage ratio vector of the constellation is obtained as follows:

[0195]

[0196] With C r The average value is the coverage index S c ,Right now

[0197]

[0198] Step 4.3: Based on the DOP value results calculated in Step 3.3, the DOP value matrix D is represented as follows:

[0199]

[0200] N d To satisfy the requirement of having enough spatial points to cover the multiplicity of DOP, where DOPT represents the DOP threshold, we have:

[0201]

[0202] In this embodiment, DOPT = 20.

[0203] The average of the spatial proportions of DOP values ​​within the target area that satisfy DOPT over time is the DOP value index S. d

[0204]

[0205] Step 4.4, based on S in step 4.1 e S in step 4.2 c And S in step 4.3d A calculation model for the comprehensive evaluation index of constellations is established using the weighted average method.

[0206] Given a weight vector W = [w1 w2 w3], where the internal elements represent shadow weight, overlay weight, and DOP weight respectively. Define a performance vector S = [S e S c S d The comprehensive index of the constellation is

[0207]

[0208] In this embodiment, there is no index preference, so let W = [1 1 1].

[0209] Step 5: Use dynamic programming algorithm to generate all constellation combinations I that satisfy the basic constraints. c Based on the reference work track library B in step one r The set of spatial points N in step two tar Using equation (27) from step four, all constellation combinations are evaluated, and the constellation scheme with the highest evaluation score is taken as the optimal constellation scheme Cons0.

[0210] Generate all constellation combinations I that satisfy the basic constraints using a dynamic programming algorithm. c , represented as

[0211]

[0212] a lF This indicates the number of satellites deployed for orbit type F in the l-th constellation combination.

[0213] to I c Evaluate by row-by-row, and take the l-th row as... Let A represent the l-th constellation scheme. l The non-zero elements in the vector form a vector in accordance with From B r Obtain the orbit of the l-th constellation and perform orbital integration over time series T. Let N... tar To achieve the goal, the l-th constellation is evaluated using equation (27) in step four to obtain S. l .

[0214] Complete I c After traversing row by row, we obtain the pair N. tar Cons0 is the optimal constellation scheme for providing navigation services.

[0215]

[0216] b represents I cThe row number where the highest score in the constellation assessment is located.

[0217] In this embodiment, the Earth-Moon constellation is selected to have 4 satellites, with a maximum deployment of 1 satellite per orbital type. A dynamic programming algorithm is used to generate 2380 constellation combinations. After evaluating the constellation's performance and traversing all possible configurations, the optimal constellation schemes for providing navigation services to a specified space target area are determined to be the L3 South Halo, L3 Vertical, L4 Vertical, and L2 Vertical configurations. Their orbital parameters are shown in Table 3, and the orbital diagrams are as follows. Figure 4 As shown.

[0218] Table 3. Orbital parameters of the optimal constellation scheme

[0219]

[0220] The combined Earth-Moon shadow score is 100, the coverage score is 99.99, and the DOP value score is 68.1497. The weight coefficients for all three indicators are 1, resulting in a comprehensive score of 89.3815. The evaluation scores for the optimal constellation scheme are shown in Table 4.

[0221] Table 4 Evaluation Scores of the Optimal Constellation Scheme

[0222]

[0223] Step Six: Construct amplitude optimization and phase optimization problems based on Cons0 obtained in Step Five. Based on the constellation orbit library B generated in Step One, and combined with the spatial point set N given in Step Two... tar Using equation (27) from step four as the optimization objective, the genetic algorithm is used to solve the amplitude optimization problem and the phase optimization problem in turn.

[0224] Construct and solve the amplitude optimization problem PRO1:

[0225] Obtain the index vector of non-zero elements in Cons0 M0 represents the number of non-zero elements in Cons0, v o This represents the index position of the o-th non-zero element in Cons0. The v-th orbital type is retrieved from the constellation orbital library B, and the v-th type orbital library contains C orbitals. v The orbital numbers are 1 to C. v Let the index be f. v f v ∈[1,C v And f v ∈N + .

[0226] Therefore, the optimization variable is

[0227]

[0228] Upper and lower bound constraints are

[0229]

[0230] The optimization goal is

[0231]

[0232] The optimal amplitude constellation Cons1 is obtained by solving PRO1 using a genetic algorithm.

[0233] In this embodiment, when solving PRO1, the genetic algorithm is set with a population size of 5, an iteration count of 50, a crossover rate of 0.9, a mutation rate of 50%, a selection pressure of 0.9, and a mutation probability of 50%. The optimal variables obtained are [1017 8717 41045227], and the optimal performance is 93.3178. The orbital parameters corresponding to Cons1 are shown in Table 5.

[0234] Table 5. Orbital parameters of Cons1

[0235]

[0236] Construct and solve the phase optimization problem PRO2:

[0237] The phase optimization process uses the initial phase angle combination of each orbit within the constellation as the optimization variable, therefore we have

[0238]

[0239] Upper and lower bound constraints are

[0240]

[0241] The optimization objective is the same as equation (32).

[0242] Based on Cons1, a genetic algorithm is used to solve PRO2 to obtain the optimal constellation Cons2.

[0243] In this embodiment, when solving PRO2, the genetic algorithm is set with a population size of 5, an iteration count of 50, a crossover rate of 0.9, a mutation rate of 50%, a selection pressure of 0.9, and a mutation probability of 50%. The optimal variables obtained are [1017 8717 41045227], and the optimal performance is 93.3178. The orbital parameters corresponding to Cons1 are shown in Table 6. The spatial proportion of the optimal constellation whose DOP value meets the threshold is shown in the curve. Figure 5 The curve showing the change in the spatial proportion that meets the coverage threshold is shown below. Figure 6 .

[0244] Table 6. Orbital parameters of Cons1

[0245]

[0246] This completes the optimization design of the navigation constellation for the custom target area in the Earth-Moon space.

[0247] Step 7: Based on the constellation size constraint obtained in Step 6, the optimal navigation constellation Cons2 for the target region in the Earth-Moon space can provide high-precision and high-efficiency positioning and navigation services to spacecraft users in this region, thereby improving the mission execution efficiency of spacecraft users in the Earth-Moon space.

[0248] The above detailed description further illustrates the purpose, technical solution, and beneficial effects of the invention. It should be understood that the above description is only a specific embodiment of the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A constellation optimization design method for providing navigation services in the Earth-Moon space, characterized in that, Includes the following steps: Step 1: Generate Constellation Working Orbit Library B and Reference Working Orbit Library B r ; Step 2: Determine the target area for navigation services provided by the Earth-Moon constellation and perform gridding to obtain the set of spatial points N. tar ; Step 3: Establish the criteria for determining the Earth-Moon shadow and coverage, establish a calculation model for the navigation accuracy factor, and obtain the Earth-Moon shadow, coverage, and navigation accuracy factor based on the calculation model for the navigation accuracy factor; Step 4: Given a time series T, perform statistical processing on the results of Earth-Moon shadow, coverage, and navigation accuracy factors obtained in Step 3, and establish a calculation model for the constellation comprehensive evaluation index using the weighted average method. Step 5: Use dynamic programming algorithm to generate all constellation combinations I that satisfy the basic constraints. c Based on the reference work track library B in step one r The set of spatial points N in step two tar Using the comprehensive index from step four, all constellation combinations are evaluated, and the constellation scheme with the highest evaluation score is taken as the optimal constellation scheme Cons0. Step Six: Construct amplitude optimization and phase optimization problems based on Cons0 obtained in Step Five; based on the constellation working orbit library B generated in Step One, combined with the spatial point set N given in Step Two. tar Using the comprehensive index from step four as the optimization objective, the genetic algorithm is used to solve the amplitude optimization problem and the phase optimization problem in sequence, thereby obtaining the optimal amplitude constellation Cons1 and the optimal phase constellation Cons2, thus realizing the optimization of the navigation constellation for the custom target region in Earth-Moon space.

2. The constellation optimization design method for providing navigation services in the Earth-Moon space as described in claim 1, characterized in that, The implementation method for step one is as follows: Based on the circularly restricted three-body theory, N is generated through three-body orbit design and numerical extension methods. t Orbital family of a type; the set of all orbits under each orbital type is the orbital family of that type, and the set of orbital families formed by all orbital types is the working orbital library B of the constellation; The orbital families under each type contain C1, C2, ..., C F (F = 1, 2, ..., N) t There are ) orbits, where F represents the index of the orbit type; sort and number the orbit families under each orbit type according to the magnitude of the amplitude in the X-axis direction under the Earth-Moon rotation system; take the number C from the orbit family of the Fth orbit type. F The track with a value of / 2 is called the reference working track for the current track type; the collection of reference working tracks for all track types is called the reference working track library B. r .

3. The constellation optimization design method for providing navigation services in the Earth-Moon space as described in claim 2, characterized in that, The second step is implemented as follows: When the target region is a spatial target region, it is defined using a cubic space approach, with the cubic space center coordinates P given. c0 =(x c0 ,y c0 ,z c0 The set of coordinates P of the target point in cubic space: P c ={(x c ,y c ,z c )x c =X0+n x ·dX+x c0 ,y c =Y0+n y ·dY+y c0 ,z c =Z0+n z ·dZ+z c0 } (1) in This indicates rounding down; cubic space is defined by giving coordinate ranges parallel to the three axes in the Earth-Moon rotation system, with the X-direction defining the starting point X0 and the ending point X... f The Y direction is defined as starting point Y0 and ending point Y0. f The Z direction is defined as the starting point Z0 and the ending point Z0. f The interval in the X direction is dX, the interval in the Y direction is dY, and the interval in the Z direction is dZ. Based on equation (1), equidistant division is performed to achieve cubic spatial gridding, then we have N tar =P c (2) When the target area is a ground target area or a lunar target area, a spherical space method is used to define the ground target area and the lunar target area, generating a set of coordinates P of the target points in the spherical space. s P s ={(α,λ,R)α=α0+n α ·dα,λ=λ0+n λ ·dλ} (3) in Spherical space is defined within the Earth-Moon rotation system and described by longitude, latitude, and spherical radius. Longitude is defined with the X-axis of the Earth-Moon rotation system as the starting point, with counterclockwise as positive, and the longitude range α∈[0,360]°. Latitude is defined with the ecliptic plane as the starting point, with the North Pole as positive, and the latitude range λ∈[-90,90]°. When defining the target region of spherical space, the spherical radius R is given, the starting point of longitude is α0, and the ending point of longitude is α. f The latitude starting point is λ0, and the latitude ending point is λ. f The longitude interval is dα, and the latitude interval is dλ. Based on the equidistant division according to equation (3), a spherical spatial grid is achieved; Transform the spherical coordinates to a rectangular coordinate system; simultaneously, give the coordinates of the sphere's center as P. s0 =(x s0 ,y s0 ,z s0 The set of coordinates of the target points in any spherical space under a Cartesian coordinate system is updated as follows: P s ′={(x s ,y s ,z s )|x s =Rcosλcosα+x s0 ,y s =Rcosλsinα+y s0 ,z s =Rsinλ+z s0 } (4) Then there is N tar =P s ′ (5)。 4. The constellation optimization design method for providing navigation services in the Earth-Moon space as described in claim 3, characterized in that, The method for implementing step three is as follows: Step 3.1: Establish criteria for determining the Earth-Moon shadow; Let r b r s r and r are the position vectors of the obscuring celestial body, the sun, and the satellite, respectively; s s =r s -r, s = r b -r are the position vectors of the Sun and satellite relative to the occluded celestial bodies; R s R and R are the radii of the Sun and the occluding object, respectively; the shadow axis is the straight line connecting the Sun's center of mass and the occluding object's center of mass; the shadow determination angle is the angle between the line connecting the shadow cone vertex V1 and the satellite's position and the shadow axis, with the opening direction of the angle being consistent with the penumbra angle f1, denoted as θ; s1 is the geocentric vector pointing to the shadow cone vertex V1; the above angles and vectors are determined by the following formula: The penumbra cone angle f1 refers to the semi-cone angle of the penumbra cone with V1 as its vertex; f1 is expressed as follows: The shadow determination criteria are as follows: When equation (8) is satisfied, the shaded state e = 1; otherwise, e = 0. Step 3.2: Establish coverage determination criteria; Geometric coverage ensures that the line connecting the satellite and the target is not obstructed by celestial bodies. Based on the positions of the celestial bodies and the target relative to the satellite, the visible area is divided into Zone I and Zone II. The angle δ between the celestial body, satellite, and target is calculated. Where, r sat r represents the position coordinates of the constellation satellites. tar Let r be the position coordinates of the target object. c The position vector of the celestial body being obscured; For geometrically visible objects r that satisfy the condition of equation (9) tar Further judgment is made; based on the cone-shaped region formed by the tangency of the satellite line of sight with the celestial body, region II is divided into three sub-regions: II-1, II-2, and II-3; r is judged according to equation (10). tar Location: Where d is the distance from the celestial body to the line connecting the observation position and the target position; the radius of the celestial body is R. c ; When r tar When located in zones II-2 and II-3, further determine r tar Visibility: Among them, the satellite-target-celestial body angle δ′; based on the above judgment criteria, when the satellite and the target are visible, the coverage state c = 1, otherwise c = 0; Step 3.3: Establish a navigation accuracy factor calculation model; When calculating the navigation accuracy factor (DOP), the number of satellites participating in the positioning solution at a certain moment in the constellation is n. t The position vector r of the target object d =(x d ,y d ,z d The position vector of the i-th satellite is... i = 1, 2, ..., nt; Calculate the pseudorange measurement residual matrix A: in i = 1, 2, ..., n t The covariance matrix Q is calculated using the least squares method. σ represents covariance; The calculation model for the navigation accuracy factor is as follows: .

5. A constellation optimization design method for providing navigation services in the Earth-Moon space as described in claim 4, characterized in that, Step four is implemented as follows: Step 4.1, the time series T is represented as: T=[t1 t2 ... t M ] (15) M represents the number of moments in the time series; Based on the shadow situation determined in step 3.1, the n×M dimensional shadow matrix E is represented as follows: Where e ij The shadow vector E represents the shadow state of the i-th satellite at time j; n represents the number of satellites in the constellation; the shadow matrix is ​​summed row by row to obtain the shadow vector E of the constellation. c : Then with E c The maximum value in the shaded index S e ,Right now S e =max[E c ] (18) Step 4.2: Based on the coverage determined in Step 3.2, the coverage matrix C is represented as: N is the set from step two. tar The number of elements; N c To satisfy the requirement of having enough spatial points for coverage multiplicity, where Cover represents the coverage multiplicity threshold, we have: Therefore, the coverage ratio vector of the constellation is obtained as follows: With C r The average value is the coverage index S c ,Right now Step 4.3: Based on the DOP value results calculated in Step 3.3, the DOP value matrix D is represented as follows: N d To satisfy the requirement of having enough spatial points to cover the multiplicity of DOP, where DOPT represents the DOP threshold, we have: The average of the spatial proportions of DOP values ​​within the target area that satisfy DOPT over time is the DOP value index S. d Step 4.4, based on S in step 4.1 e S in step 4.2 c And S in step 4.3 d A calculation model for the comprehensive evaluation index of constellations is established using the weighted average method; Given a weight vector W = [w1 w2 w3], where the internal elements represent shadow weight, overlay weight, and DOP weight respectively; define a performance vector S = [S e S c S d The comprehensive indicators of a zodiac sign are: 。 6. A constellation optimization design method for providing navigation services in Earth-Moon space as described in claim 5, characterized in that, Step five is implemented as follows: Generate all constellation combinations I that satisfy the basic constraints using a dynamic programming algorithm. c , represented as a lF This indicates the number of satellites deployed for orbital type F in the l-th constellation combination; to I c Evaluate by row-by-row, and take the l-th row as... Let A represent the l-th constellation scheme. l The non-zero elements in the vector form a vector in accordance with From B r Obtain the orbit of the l-th constellation and perform orbital integration over time series T; with N tar To achieve the goal, the l-th constellation is evaluated using equation (26) in step four to obtain S. l ; Complete I c After traversing row by row, we obtain the pair N. tar Cons0 is the optimal constellation scheme for providing navigation services. b represents I c The row number where the highest score in the constellation assessment is located.

7. A constellation optimization design method for providing navigation services in Earth-Moon space as described in claim 6, characterized in that, Step six is ​​implemented as follows: Construct and solve the amplitude optimization problem PRO1: Obtain the index vector of non-zero elements in Cons0 M0 represents the number of non-zero elements in Cons0, v o This indicates the index position of the o-th non-zero element in Cons0; retrieve the orbital library of the v-th orbital type from the constellation orbital library B, where the v-th type orbital library contains C orbitals. v The orbital numbers are 1 to C. v Let the index be f. v f v ∈[1,C v And f v ∈N + ; Therefore, the optimization variable is Upper and lower bound constraints are The optimization goal is Construct and solve the phase optimization problem PRO2: The phase optimization process uses the initial phase angle combination of each orbit within the constellation as the optimization variable, therefore we have Upper and lower bound constraints are The optimization objective is the same as in equation (31); The optimal constellation Cons1 with amplitude is obtained by solving PRO1 using a genetic algorithm; the optimal constellation Cons2 with phase is obtained by solving PRO2 using a genetic algorithm based on Cons1; thus, the navigation constellation optimization design for the custom target region in Earth-Moon space is completed.

8. A constellation optimization design method for providing navigation services in Earth-Moon space as described in claims 1, 2, 3, 4, 5, 6, or 7, characterized in that, It also includes step seven, which is the optimal navigation constellation Cons2 for the target area in the Earth-Moon space under the constellation size constraint obtained in step six. This constellation can provide high-precision and high-efficiency positioning and navigation services to spacecraft users in this area, thereby improving the mission execution efficiency of spacecraft users in the Earth-Moon space.

Images